Citation

If you use immLynx in your research, please cite:

• Belean, A. et al. Single-cell sequencing of non-infectious and viral uveitis defines distinct mechanisms in ocular inflammation. Frontiers in Immunology 17, 1916535 (2026). doi:10.3389/fimmu.2026.1916535

immLynx is an interface to existing Python tools, so please also cite the original publications for each tool you use. The full immLynx citation, including a BibTeX entry, is available with citation("immLynx").

Python Dependencies

immLynx uses basilisk to manage Python dependencies automatically. The first time you call a function that invokes Python, basilisk will create an isolated conda environment containing all required packages (tcrdist3, olga, soNNia, clusTCR, metaclonotypist, pyrepseq, and PyTorch with the ESM-2 model). This initial setup may take several minutes but only needs to occur once per installation. Subsequent calls reuse the existing environment with minimal overhead.

A few tools pin dependency versions that cannot coexist with that main environment, so they get their own. exportToScanpy() uses a scanpy/scirpy environment, and runDeepTCR() uses a TensorFlow 2.12 environment. Each is built on first use of the function that needs it, so you only pay for the ones you actually call.

Two functions use their own environments rather than that one. exportToScanpy() needs the scanpy/scirpy stack, and runScXpand() needs Python 3.11 with a current PyTorch, which conflicts with the pinned versions above. Each builds on first use. The scXpand environment is the largest, at several gigabytes.

Loading the Package and Example Data

We begin by loading immLynx along with scran and scater, which provide utilities for single-cell analysis and visualization. The immLynx_example dataset is a SingleCellExperiment object containing scRNA-seq data from multiple patients, with TCR repertoire information added via scRepertoire.

library(immLynx)
library(scran)
library(scater)

data("immLynx_example")

immLynx_example
#> class: SingleCellExperiment 
#> dim: 2000 500 
#> metadata(0):
#> assays(2): counts logcounts
#> rownames(2000): TRBV7-2 HSPA1B ... AC119396.1 BHLHE40
#> rowData names(0):
#> colnames(500): P17B_AAACCTGCAACACGCC-1 P17B_AAACGGGTCTCCAGGG-1 ...
#>   P20L_TCAACGATCCCTCAGT-1 P20L_TCTGGAAGTGGAAAGA-1
#> colData names(16): orig.ident nCount_RNA ... ident clusters
#> reducedDimNames(1): UMAP
#> mainExpName: RNA
#> altExpNames(0):

Extracting TCR Data

Before running downstream analyses, it is often useful to extract TCR data into a data.frame format for inspection and quality control. The extractTCRdata() function wraps immApex::getIR() and supports both long format (one row per chain) and wide format (one row per cell with columns for each chain). The validateTCRdata() function checks for common issues such as missing CDR3 sequences, non-standard amino acids, and inconsistent chain annotations.

tcr_data <- extractTCRdata(immLynx_example,
                           chains = "TRB")
head(tcr_data)
#>             cdr3_aa        v    d       j     c                 barcode chain
#> 1     CAISEQGKGELFF TRBV10-3 <NA> TRBJ2-2 TRBC2 P17B_AAACCTGCAACACGCC-1   TRB
#> 2 CASSVRRERANTGELFF    TRBV9 <NA> TRBJ2-2 TRBC2 P17B_AAACGGGTCTCCAGGG-1   TRB
#> 4 CASSVRRERANTGELFF    TRBV9 <NA> TRBJ2-2 TRBC2 P17B_AACTTTCTCAGCGACC-1   TRB
#> 6   CASSSRQDSTDTQYF   TRBV27 <NA> TRBJ2-3 TRBC2 P17B_AATCGGTCAAGCCCAC-1   TRB
#> 8 CASSVRRERANTGELFF    TRBV9 <NA> TRBJ2-2 TRBC2 P17B_ACGAGCCTCATTTGGG-1   TRB
#> 9  CASSPRAGTPNTEAFF  TRBV4-1 <NA> TRBJ1-1 TRBC1 P17B_ACGCCAGGTTCATGGT-1   TRB

tcr_wide <- extractTCRdata(immLynx_example,
                           chains = "both",
                           format = "wide")

validation <- validateTCRdata(tcr_data)
print(validation)
#> $valid
#> [1] TRUE
#> 
#> $errors
#> character(0)
#> 
#> $warnings
#> character(0)
#> 
#> $summary
#> $summary$n_rows
#> [1] 306
#> 
#> $summary$n_unique_barcodes
#> [1] 306
#> 
#> $summary$n_unique_cdr3
#> [1] 275
#> 
#> $summary$n_na_cdr3
#> [1] 0
#> 
#> $summary$chains
#> [1] "TRB"

Summarizing TCR Repertoire

The summarizeTCRrepertoire() function computes standard repertoire statistics including clonotype frequencies, diversity indices (Shannon entropy, Simpson index, clonality), CDR3 length distributions, and V/J gene usage. These metrics provide a high-level overview of repertoire structure before detailed analysis.

summary <- summarizeTCRrepertoire(immLynx_example,
                                  chains = "TRB")
print(summary)
#> === TCR Repertoire Summary ===
#> Chain(s): TRB 
#> 
#> --- Basic Statistics ---
#> Total cells with TCR: 306 
#> Unique clonotypes: 275 
#> Clonotype ratio: 0.8987 
#> 
#> --- Diversity Metrics ---
#> Shannon entropy: 5.511 
#> Clonality: 0.0188 
#> Simpson index: 0.0057 
#> Inverse Simpson: 174.69 
#> 
#> --- Top Clonotypes ---
#>            cdr3_aa count  proportion
#>  CASSVRRERANTGELFF    14 0.045751634
#>      CASSPTVAGEQFF     3 0.009803922
#>   CASSQAPFSTSGELFF     3 0.009803922
#>   CASSQDRTGLDYEQYF     3 0.009803922
#>   CASSRLRTGALYEQYF     3 0.009803922
#> 
#> --- CDR3 Length Distribution ---
#> Mean: 14.8 
#> Median: 15 
#> Range: 10 - 21

TCR Clustering with clusTCR

clusTCR groups TCR sequences by CDR3 similarity using the Markov Cluster (MCL) algorithm. The inflation parameter controls cluster granularity: higher values produce more, smaller clusters, while lower values yield fewer, larger clusters. Cluster assignments are stored directly in the colData of the SingleCellExperiment object, making them immediately available for downstream visualization and statistical testing.

sce <- runClustTCR(
  immLynx_example,
  chains = "TRB",
  method = "mcl",
  inflation = 2.0
)

table(sce$clustcr_TRB)

Near-Neighbor Search with symdel

Symmetric deletion lookup finds every pair of CDR3 sequences within a given edit distance without ever building a full pairwise distance matrix. It indexes each sequence by its set of single-character deletions, so two sequences within one edit share at least one deletion variant and collide in the index. Cost scales with the number of sequences rather than with their square, which makes it practical on repertoires where a full distance matrix would not fit in memory.

runSymdelNeighbors() searches over unique sequences and maps results back onto cells. By default it adds a symdel_degree column to colData giving the number of neighbors each cell’s CDR3 has. High degree marks sequences sitting in dense similarity neighborhoods, the signature of convergent recombination and antigen-driven expansion.

sce <- runSymdelNeighbors(
  immLynx_example,
  chains = "TRB",
  max_edits = 1
)

table(sce$symdel_degree, useNA = "ifany")

Setting return_object = FALSE returns the neighbor pairs themselves, which is the more useful form when building a similarity graph.

edges <- runSymdelNeighbors(
  immLynx_example,
  chains = "TRB",
  max_edits = 2,
  return_object = FALSE
)

head(edges)

Note that this is the same algorithm benchmarked in the XT-neighbor work. immLynx uses the CPU implementation from pyrepseq, so no GPU is required.

TCR Distance Calculations with tcrdist3

tcrdist3 computes pairwise distances between TCR sequences using a substitution-matrix-based metric that accounts for CDR3 amino acid sequence, V gene, and J gene identity. The resulting distance matrix can be used for hierarchical clustering, nearest-neighbor analysis, or as input to dimensionality reduction methods. Note that tcrdist3 deduplicates sequences internally, so the distance matrix dimensions correspond to the number of unique TCR sequences rather than the total number of cells.

dist_results <- runTCRdist(
  immLynx_example,
  chains = "beta",
  organism = "human"
)

dim(dist_results$distances$pw_beta)

Generation Probability with OLGA

OLGA (Optimized Likelihood estimate of immunoGlobulin Amino-acid sequences) computes the generation probability (Pgen) of a TCR sequence—the probability that V(D)J recombination produces that specific CDR3 amino acid sequence. Low Pgen values indicate sequences that are unlikely to arise by chance, which may suggest antigen-driven selection. The runOLGA() function adds log10-transformed Pgen values to the colData of the input object.

sce <- runOLGA(
  immLynx_example,
  chains = "TRB",
  model = "humanTRB"
)

hist(log10(sce$olga_pgen_TRB), breaks = 50)

The generateOLGA() function can also produce synthetic TCR sequences sampled from the recombination model, which is useful for constructing background distributions or benchmarking analyses.

random_seqs <- generateOLGA(n = 100,
                            model = "humanTRB")
head(random_seqs)

Protein Language Model Embeddings

Protein language models such as ESM-2 learn dense vector representations of amino acid sequences from large-scale protein databases. When applied to CDR3 sequences, these embeddings capture biochemical and structural properties that complement sequence-based distance metrics. The runEmbeddings() function generates per-sequence embeddings and stores them as a dimensionality reduction in the SingleCellExperiment object, enabling visualization with standard tools such as UMAP or PCA. Multiple pooling strategies (mean, cls, max) are available to aggregate per-residue representations into a fixed-length vector.

sce <- runEmbeddings(
  immLynx_example,
  chains = "TRB",
  model_name = "facebook/esm2_t12_35M_UR50D",
  pool = "mean"
)

sce <- scater::runUMAP(sce,
                       dimred = "tcr_esm")
scater::plotReducedDim(sce, dimred = "UMAP")

Deep Generative Featurization with DeepTCR

DeepTCR learns a latent representation of CDR3 sequences with a variational autoencoder. Where ESM-2 brings general protein knowledge learned from UniProt, the DeepTCR VAE is trained on the repertoire in front of it, so its features describe the structure of that specific dataset. The two are complementary rather than interchangeable.

runDeepTCR() trains on unique sequences and writes the learned features to a dimensional reduction. DeepTCR prunes uninformative latent dimensions during training, so the result may be narrower than latent_dim requests.

sce <- runDeepTCR(
  immLynx_example,
  chains = "TRB",
  latent_dim = 64,
  verbose = FALSE
)

dim(SingleCellExperiment::reducedDim(sce, "tcr_deeptcr"))

This chunk is not evaluated when the vignette is built. DeepTCR pins TensorFlow 2.12 and its own numpy, pandas, and scipy versions, so it runs in a separate deepTCREnv environment that is a large first-use download. Set eval = TRUE to run it locally.

Passing use_genes = TRUE includes V and J gene usage alongside the CDR3 sequence, which helps when gene usage carries specificity signal that the CDR3 alone does not.

Metaclone Discovery with Metaclonotypist

Metaclonotypist identifies metaclones—groups of TCR clonotypes that share enough sequence similarity to potentially recognize the same antigen. This approach extends traditional exact-match clonotyping by allowing controlled levels of CDR3 edit distance and tcrdist-based similarity. The max_edits and max_dist parameters control the stringency of metaclone grouping.

sce <- runMetaclonotypist(
  immLynx_example,
  chains = "beta",
  method = "tcrdist",
  max_edits = 2,
  max_dist = 20
)

table(sce$metaclone)

Selection Inference with soNNia

soNNia uses a neural network to model post-thymic selection by comparing observed TCR sequences against a background distribution generated by V(D)J recombination. The model estimates a selection factor for each sequence, quantifying how much more (or less) likely a sequence is to appear in the observed repertoire relative to the recombination baseline. This requires a background sequence file, which can be generated using generateOLGA().

background <- generateOLGA(n = 10000, model = "humanTRB")
write.csv(background, "background.csv", row.names = FALSE)

sce <- runSoNNia(
  immLynx_example,
  chains = "TRB",
  background_file = "background.csv"
)