Bio Tcr Bcr Analysis Repertoire Visualization

Version Compatibility

Reference examples tested with: matplotlib 3.8+, seaborn 0.13+, pandas 2.2+, numpy 1.26+; R circlize 0.4+, iNEXT 3.0+; VDJtools 1.2.1+

Before using code patterns, verify installed versions match. If versions differ:

• Python: pip show <package> then help(module.function) to check signatures
• R: packageVersion('<pkg>') then ?function_name to verify parameters
• CLI: <tool> --version then <tool> --help to confirm flags

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Note: matplotlib 3.9 removed plt.cm.get_cmap(name, N); use plt.get_cmap(name) and sample it, or matplotlib.colormaps[name].resampled(N). seaborn 0.13 deprecated bare palette= without hue=; assign hue= and legend=False.

Repertoire Visualization

"Plot my TCR/BCR repertoire" -> Render V-J usage, CDR3-length spectratype, clonal-space, clonal tracking, diversity/rarefaction, overlap, and similarity-network figures, and read each one correctly.

• CLI: vdjtools PlotFancyVJUsage, vdjtools RarefactionPlot (delegates plotting to R)
• R: circlize::chordDiagram (V-J), iNEXT::iNEXT + ggiNEXT (rarefaction/extrapolation)
• Python: matplotlib/seaborn for bespoke figures; numpy multinomial resampling for rarefaction; networkx + rapidfuzz for similarity networks

The governing principle: every figure inherits two upstream choices

A repertoire figure is only as valid as the numbers behind it, and those numbers depend on two decisions made before any plot is drawn:

1. The clonotype definition. A clonotype = CDR3 (nucleotide OR amino acid) + optionally V and J. Amino-acid CDR3 collapses convergent recombination (many nt rearrangements -> one aa clonotype), inflating apparent sharing and deflating richness; nt CDR3 is more conservative (Venturi et al. 2006 PNAS 103:18691). Every count on every axis shifts with this choice, so it must be stated in the figure and held constant across all samples in a comparison.
2. The sequencing depth. Richness, Shannon, clonality, and set-overlap (Jaccard, shared-clonotype counts) are strictly increasing functions of reads sequenced - the rare-clone tail never saturates. Diversity, rarefaction, and overlap figures are comparable across samples ONLY after depth normalization: downsample all samples to a common depth, or read rarefaction curves at a shared x. Comparing raw values across libraries of unequal depth measures sequencing effort, not biology (Chao et al. 2014 Ecol Monogr 84:45; Greiff et al. 2015 Genome Med 7:49). This caveat governs the diversity, rarefaction, and overlap recipes below.

Choosing the figure and reading it

V-J usage: chord/circos and heatmap

The chord/circos shows which V and J segments pair, weighted by clonotype count, within one sample.

VDJtools route (delegates plotting to R; requires RInstall once):

vdjtools PlotFancyVJUsage -m metadata.txt output_dir/

R with circlize - a V-by-J count matrix becomes a chord diagram:

library(circlize)

plot_vj_chord <- function(clone_df) {
    vj_matrix <- table(clone_df$v_gene, clone_df$j_gene)
    chordDiagram(vj_matrix, transparency = 0.5, annotationTrack = c('grid', 'name'))
}

Python heatmap alternative (samples x V gene, or V x J for one sample) avoids a chord dependency and reads more quantitatively for many segments:

import seaborn as sns
import matplotlib.pyplot as plt

def plot_vj_heatmap(clone_df):
    vj = clone_df.pivot_table(index='v_gene', columns='j_gene', values='frequency', aggfunc='sum', fill_value=0)
    fig, ax = plt.subplots(figsize=(8, 6))
    sns.heatmap(vj, cmap='viridis', ax=ax)
    ax.set_title('V-J pairing frequency')
    return fig

CDR3 spectratype: length distribution reveals clonal structure

The spectratype is a histogram of CDR3 length. A roughly Gaussian, bell-shaped distribution indicates a diverse, polyclonal (naive-like) repertoire; skew or discrete spikes at particular lengths indicate clonal expansion(s). The classic immunoscope spectratype (and VDJtools CalcSpectratype) bins nucleotide length, where in-frame clones sit 3 nt apart and the periodicity is part of the readout; amino-acid length is a common simplification - state which is plotted. Weighting by read/UMI frequency shows expansions; weighting by unique clonotypes shows the underlying diversity - the two views can look opposite, so label which is plotted.

def plot_spectratype(clone_df, length_col='cdr3_length'):
    fig, ax = plt.subplots(figsize=(9, 5))
    bins = range(clone_df[length_col].min(), clone_df[length_col].max() + 2)
    ax.hist(clone_df[length_col], bins=bins, weights=clone_df['frequency'], color='steelblue')
    ax.set_xlabel('CDR3 length (aa)')
    ax.set_ylabel('Frequency')  # Gaussian = polyclonal; spikes = clonal expansion
    return fig

Clonal-space / proportion: how much repertoire the big clones hold

Bin clones into frequency strata (Rare / Small / Medium / Large / Hyperexpanded) and stack their summed frequency. Because strata are defined on frequency, this view is comparatively depth-robust and shows clonal-space homeostasis at a glance. A treemap of top clones is an alternative when individual dominant clones matter.

import numpy as np

def plot_clonal_space(clone_df, sample_col='sample'):
    # Frequency-based strata (immunarch homeo convention); frequency makes this depth-robust
    edges = [0, 1e-4, 1e-3, 1e-2, 1e-1, 1.0]
    labels = ['Rare', 'Small', 'Medium', 'Large', 'Hyperexpanded']
    clone_df = clone_df.copy()
    clone_df['stratum'] = pd.cut(clone_df['frequency'], bins=edges, labels=labels)
    space = clone_df.groupby([sample_col, 'stratum'], observed=True)['frequency'].sum().unstack(fill_value=0)
    fig, ax = plt.subplots(figsize=(8, 5))
    space[labels].plot(kind='bar', stacked=True, ax=ax, colormap='viridis')
    ax.set_ylabel('Fraction of repertoire')
    return fig

Clonal tracking across timepoints

Track individual clone frequencies over an ordered sample set (vaccination, infection, therapy). A line plot of the top clones, or an alluvial/stream for the same data, shows expansion, contraction, and persistence. Downsample timepoints to a common depth before declaring contraction - a clone scored absent at one timepoint is frequently below the sampling floor, not truly gone.

def plot_clone_tracking(clone_df, top_n=10, clone_col='cdr3_aa', time_col='timepoint'):
    top = clone_df.groupby(clone_col)['frequency'].sum().nlargest(top_n).index
    fig, ax = plt.subplots(figsize=(9, 5))
    for clone in top:
        d = clone_df[clone_df[clone_col] == clone].sort_values(time_col)
        ax.plot(d[time_col], d['frequency'], marker='o', label=clone[:12])
    ax.set_xlabel('Timepoint'); ax.set_ylabel('Clone frequency')
    ax.legend(bbox_to_anchor=(1.02, 1), loc='upper left', fontsize=7)
    return fig