9: Data Analysis and Visualization with Python Libraries

Chapter 9: Data Analysis and Visualization with Python Libraries

Modern genomics relies on the Python scientific computing stack: NumPy for multi-dimensional matrix operations, Pandas for high-throughput DataFrame analysis, SciPy for statistical modeling, and Matplotlib / Seaborn for publication-grade visualizations.

1. NumPy Arrays & Matrix Operations

NumPy vectors and matrices execute array operations up to 100x faster than standard Python loops via optimized C-level SIMD instructions.

import numpy as np

# Calculating Position Frequency Matrix (PFM) for transcription factor binding
# Columns: Positions 1 to 6 | Rows: A, C, G, T
pfm = np.array([
    [12,  3,  0, 24, 18,  2],  # A
    [ 2, 18,  0,  1,  3,  1],  # C
    [ 8,  2, 24,  0,  2,  0],  # G
    [ 3,  2,  1,  0,  2, 22]   # T
])

# Normalize columns to calculate Position Weight Matrix (PWM) probabilities
column_sums = pfm.sum(axis=0)
pwm = pfm / column_sums
print("Position Weight Matrix (Probabilities):\n", np.round(pwm, 3))

2. Pandas DataFrames & Data Manipulation

Pandas DataFrames provide tabular manipulation designed for RNA-seq counts, variant annotations, and metadata tables.

import pandas as pd

# Load RNA-Seq differential expression table
df = pd.read_csv("differential_expression.tsv", sep="\t", index_col="gene_id")

# Filter for significantly upregulated genes (log2FC >= 2.0 and FDR < 0.01)
upregulated = df[(df["log2FoldChange"] >= 2.0) & (df["padj"] < 0.01)]
print(f"Upregulated Genes Count: {len(upregulated)}")
print(upregulated[["gene_name", "log2FoldChange", "padj"]].head())

3. SciPy for Statistical & Scientific Computing

SciPy provides statistical testing (Mann-Whitney U, t-test) and hierarchical clustering for biological matrices.

4. Matplotlib Library for Charts & Graphs

Generating a standard Volcano Plot displaying statistical significance vs biological effect size:

import matplotlib.pyplot as plt

plt.figure(figsize=(8, 6))
# Non-significant genes
plt.scatter(df["log2FoldChange"], -np.log10(df["padj"]), color="gray", alpha=0.4, s=15, label="NS")

# Significant genes
sig = (df["padj"] < 0.05) & (df["log2FoldChange"].abs() >= 1.5)
plt.scatter(df.loc[sig, "log2FoldChange"], -np.log10(df.loc[sig, "padj"]), color="crimson", s=25, label="Significant")

plt.axvline(1.5, color="black", linestyle="--", alpha=0.6)
plt.axvline(-1.5, color="black", linestyle="--", alpha=0.6)
plt.axhline(-np.log10(0.05), color="blue", linestyle="--", alpha=0.6)

plt.xlabel("log2 Fold Change")
plt.ylabel("-log10 Adjusted P-value")
plt.title("RNA-Seq Volcano Plot: Treated vs Control")
plt.legend()
plt.tight_layout()
plt.savefig("volcano_plot.pdf", dpi=300)

5. Seaborn for Clustermaps & Heatmaps

Seaborn provides two-way hierarchical clustering of gene expression matrices across tumor and healthy tissue replicates.

import seaborn as sns

# Hierarchically clustered expression heatmap
sns.clustermap(df.iloc[:50, :8], cmap="vlag", z_score=0, figsize=(10, 12), metric="correlation")
plt.savefig("clustered_heatmap.png", dpi=300)

6. Creating Interactive Visualizations with Bokeh

Interactive web graphics allowing researchers to hover over points, zoom into genomic regions, and inspect gene annotations.

7. Publication-Ready Figure Design

Best practices for Nature/Science/Cell figures: vector output (PDF/SVG), high DPI, accessible font families, and color-blind friendly palettes (viridis, cividis).