# Python Bio Data Visualization

> Build volcano/MA plots, clustermap heatmaps, and multi-panel GridSpec figures with matplotlib/seaborn. Use when plotting DE results, expression data, or QC distributions, or fixing savefig, log-axis, colormap bugs.

- Skill: `pavel-kravchenko/python-bio-data-visualization` (Agent Skill)
- Install (CLI): `npx skillmds@latest add pavel-kravchenko/python-bio-data-visualization`
- Raw SKILL.md: https://api.skillmd.com/api/skills/pavel-kravchenko/python-bio-data-visualization/raw
- Safety review: pending
- Works with: Claude Code, Claude.ai, OpenAI Codex
- Category: Coding & Dev Tools
- Author: pavel-kravchenko (https://skillmd.com/u/pavel-kravchenko)
- Updated: 2026-09-17
- Page: https://skillmd.com/skills/pavel-kravchenko/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.

```python
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.

```python
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`:

```r
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.

```python
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.

```python
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

```python
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.