Python bio data visualization

Data Visualization for Bioinformatics

When to Use

• Plotting differential expression results as a volcano plot or MA plot.
• Building an expression heatmap/clustermap from a genes x samples matrix.
• Laying out multi-panel figures (QC dashboards, timecourse + enrichment + MA plot).
• Comparing expression across groups with box/violin/strip/KDE plots.
• Debugging a blank saved figure, a broken log-axis, or a clipped legend.

Version Compatibility

matplotlib >= 3.8, seaborn >= 0.13, pandas >= 2.0, scipy >= 1.11, Python >= 3.10. seaborn 0.13 renamed the old palette without hue pattern (now warns/deprecated) — always pass hue alongside palette.

Prerequisites

• pip install matplotlib seaborn pandas numpy scipy
• Familiarity with pandas DataFrames (see python-bio-pandas) and NumPy arrays (see python-bio-numpy).
• A tidy (long-form) DataFrame or a genes x samples numeric matrix to plot.

Critical Gotchas

• savefig before plt.show(): show() clears the figure; if you call it first, savefig writes a blank file.
• Log axes and zeros: plt.yscale('log') fails on/distorts zeros. Log-transform the data before plotting, e.g. np.log2(counts + 1), instead of a log-scale axis on raw counts.
• Seaborn expects long-form data: if data is wide (one column per sample), melt first: pd.melt(df, id_vars=['gene'], value_vars=samples).
• OO interface for multi-panel: use fig, ax = plt.subplots() + ax.plot(). plt.plot() acts on the current active axes and silently breaks inside loops over subplots.
• Colormaps: use RdBu_r/coolwarm for diverging fold-change data, viridis/cividis for sequential expression levels, and seaborn's colorblind palette for categorical groups. Never jet (not perceptually uniform, misleads readers).

Multi-Panel Figures with GridSpec

Goal: lay out panels of unequal size (e.g. one wide timecourse plot on top, two square plots below) in a single figure. Approach: use fig.add_gridspec and slice it with gs[row, col]; a full-row slice gs[0, :] spans all columns.

import matplotlib.pyplot as plt
import numpy as np

fig = plt.figure(figsize=(12, 7))
gs = fig.add_gridspec(2, 2, height_ratios=[1, 1], hspace=0.35, wspace=0.3)

ax_top = fig.add_subplot(gs[0, :])   # spans both columns
ax_bl = fig.add_subplot(gs[1, 0])
ax_br = fig.add_subplot(gs[1, 1])

# Top: expression timecourse for a few genes after DNA damage
t = np.linspace(0, 48, 100)
for name, phase in [('p53', 0), ('MDM2', 4), ('CDKN1A', 8)]:
    signal = 2 * np.sin(2 * np.pi * t / 24 - phase) * np.exp(-t / 80) + 5
    ax_top.plot(t, signal, linewidth=2, label=name)
ax_top.set_xlabel('Time (hours)')
ax_top.set_ylabel('Expression (log2)')
ax_top.legend()

# Bottom-left: pathway enrichment
pathways = ['Apoptosis', 'Cell Cycle', 'DNA Repair', 'p53 Signaling', 'Autophagy']
neg_log_p = [8.5, 7.2, 6.1, 10.3, 3.2]
ax_bl.barh(pathways, neg_log_p, color=plt.cm.viridis(np.linspace(0.2, 0.8, 5)))
ax_bl.axvline(-np.log10(0.05), color='red', linestyle='--', label='p=0.05')
ax_bl.set_xlabel('-log$_{10}$(p-value)')
ax_bl.legend(fontsize=8)

plt.show()

Volcano and MA Plots

Goal: visualize DE results — significance vs effect size (volcano), or effect size vs expression level (MA). Approach: compute a boolean significance mask from fold-change and adjusted p-value thresholds, plot non-significant points first (gray, low alpha) then significant points on top (colored), and add threshold reference lines.

import numpy as np
import matplotlib.pyplot as plt


def plot_volcano(log2fc, neg_log10_padj, ax=None, fc_thresh=1.0, padj_thresh=0.05):
    """Volcano plot: log2 fold-change vs -log10(adjusted p-value).

    log2fc, neg_log10_padj: 1D array-likes of equal length, one row per gene.
    fc_thresh: |log2FC| cutoff for calling a gene significant.
    padj_thresh: adjusted p-value cutoff (converted internally to -log10).
    """
    log2fc = np.asarray(log2fc)
    neg_log10_padj = np.asarray(neg_log10_padj)
    if ax is None:
        _, ax = plt.subplots(figsize=(7, 5))
    sig = (np.abs(log2fc) > fc_thresh) & (neg_log10_padj > -np.log10(padj_thresh))
    ax.scatter(log2fc[~sig], neg_log10_padj[~sig], s=5, alpha=0.4, c='gray')
    ax.scatter(log2fc[sig], neg_log10_padj[sig], s=8, alpha=0.7, c='red')
    ax.axvline(-fc_thresh, color='gray', linestyle='--', linewidth=0.8)
    ax.axvline(fc_thresh, color='gray', linestyle='--', linewidth=0.8)
    ax.axhline(-np.log10(padj_thresh), color='gray', linestyle='--', linewidth=0.8)
    ax.set_xlabel('log$_2$ Fold Change')
    ax.set_ylabel('-log$_{10}$(padj)')
    return ax


def plot_ma(base_mean, log2fc, sig_mask, ax=None):
    """MA plot: log10(mean expression) vs log2 fold-change, significant genes highlighted."""
    base_mean, log2fc, sig_mask = np.asarray(base_mean), np.asarray(log2fc), np.asarray(sig_mask)
    if ax is None:
        _, ax = plt.subplots(figsize=(7, 5))
    ax.scatter(np.log10(base_mean[~sig_mask]), log2fc[~sig_mask], s=10, alpha=0.5, c='gray')
    ax.scatter(np.log10(base_mean[sig_mask]), log2fc[sig_mask], s=15, alpha=0.8, c='red', label='Significant')
    ax.axhline(0, color='black', linewidth=0.8)
    ax.set_xlabel('log$_{10}$(Base Mean)')
    ax.set_ylabel('log$_2$(Fold Change)')
    ax.legend(fontsize=8)
    return ax

R practitioners typically make the same volcano plot from a DESeq2/edgeR results() table with EnhancedVolcano or plain ggplot2:

library(ggplot2)

# res: data.frame with columns log2FoldChange, padj (e.g. from DESeq2::results())
res$sig <- with(res, !is.na(padj) & padj < 0.05 & abs(log2FoldChange) > 1)

ggplot(res, aes(x = log2FoldChange, y = -log10(padj), color = sig)) +
  geom_point(alpha = 0.6, size = 1) +
  scale_color_manual(values = c(`FALSE` = "grey60", `TRUE` = "red")) +
  geom_vline(xintercept = c(-1, 1), linetype = "dashed", color = "grey40") +
  geom_hline(yintercept = -log10(0.05), linetype = "dashed", color = "grey40") +
  labs(x = expression(log[2]~"fold change"), y = expression(-log[10]~"(padj)")) +
  theme_bw()

Heatmap / Clustermap Pattern

Goal: show a genes x samples expression matrix with hierarchical clustering on both axes. Approach: z-score each row (per-gene) so the color scale reflects relative change across samples, then let sns.clustermap cluster and reorder rows/columns.

import seaborn as sns
from scipy.stats import zscore

# mat: genes (rows) x samples (columns) numeric expression matrix
mat_z = zscore(mat, axis=1)  # per-gene z-score across samples

g = sns.clustermap(
    mat_z,
    cmap='RdBu_r', center=0, vmin=-3, vmax=3,
    row_cluster=True, col_cluster=True,
    figsize=(10, 8), yticklabels=False,
)
# clustermap returns a ClusterGrid, not a Figure/Axes:
g.fig.savefig('heatmap.pdf', dpi=300, bbox_inches='tight')

Comparing Groups: Box / Violin / Strip / KDE

Goal: compare a gene's expression across categorical groups (cell type, condition). Approach: keep data in long form (one row per observation) and let seaborn's hue handle the grouping; split=True on violins compares exactly two hue levels side by side.

import seaborn as sns
import matplotlib.pyplot as plt

fig, axes = plt.subplots(1, 3, figsize=(15, 4))

sns.boxplot(data=expr_df, x='Cell_Type', y='TNF_alpha', hue='Condition',
            ax=axes[0], palette='Set2')
sns.violinplot(data=expr_df, x='Cell_Type', y='CD19', hue='Condition',
               split=True, ax=axes[1], palette='muted', inner='quartile')
sns.stripplot(data=expr_df, x='Cell_Type', y='CD3', hue='Condition',
              dodge=True, ax=axes[2], palette='deep', alpha=0.6, size=4)
plt.tight_layout()

Saving Figures

fig.savefig('figure.pdf', dpi=300, bbox_inches='tight')  # vector for publication
fig.savefig('figure.png', dpi=150, bbox_inches='tight')  # raster for display/slides
plt.close(fig)  # free memory when saving many figures in a loop

Colormap Quick Reference

Data type	Recommended colormap
Expression fold-change (diverging)	RdBu_r, coolwarm
Expression level (sequential)	viridis, YlOrRd
Categorical (cell types)	sns.color_palette('colorblind')
P-value / significance	plasma or custom threshold-based

Pitfalls

• Legend outside plot area gets clipped on save: use bbox_inches='tight' in savefig.
• Tick labels overlap: rotate with ax.set_xticklabels(labels, rotation=45, ha='right').
• seaborn changes global rcParams on import: import order matters; set custom plt.rcParams after import seaborn.
• clustermap returns a ClusterGrid, not a Figure: access the figure via g.fig, the heatmap axes via g.ax_heatmap.
• hue without matching palette length: seaborn cycles or drops colors silently if the palette has fewer entries than hue levels — pass an explicit palette sized to the number of groups.

See Also

• bio-data-visualization-volcano-customization — advanced volcano plot styling and gene labeling.
• bio-data-visualization-heatmaps-clustering — clustering algorithms and dendrogram options in depth.
• bio-data-visualization-multipanel-figures — complex multi-panel layout patterns.
• python-bio-pandas — reshaping data (wide-to-long) before plotting.