So begins another season

Microbe and virii have an important advantage over homo sapiens: quick generation times.

The adaptive immune system

Fortunately, the jawed vertebrates also have an adaptive immune system unmoored from the germline.

Somatic recombination and diversity

• In T-cells, TCR is generated by splicing together gene segments and non-templated nucleotides.
  Purely combinatorial calculations show \(>10^{18}\) possible \(\alpha\)-\(\beta\) pairs. However, gene usage and cell maturation biases mean the effective sample size much reduced.

• B cells also use combinatorial splicing, as well as somatic mutations.

Simultaneous repertoire and functional expression

• Long history of methods for repertoire sequencing
• Lu, McDavid… (2018) used a well-based method to recover B-cell immunoglobulin region and gene expression from single cells.
• 10X Genomics Chromium has “beads” with forest of bead- and molecule-specific oligos.
• Cell lysis and cDNA template happens in aqueous droplets

Single cell expression

Vignettes

In situ synovial B cells from Rheumatoid Arthritis

• 4 patients with RA and large numbers of B cells infiltrating the synovium (knees)
• 1 patient with matched blood

Autoreactive B Cell Responses in RA

• 6 RA patients, B cells selected by tetramer for response to rheumatoid factor (Rf) or citrullinated protein (Cit) or background (Tet-)

Pre- and postnatal T cell landscape in human infants

• T cells from cord blood, 8 infants, from very premature to full-term.
• Sorted into CD31+/- subsets. CD31+: recent thymic emigrants, more quiescent, IL8 secreting, rhinovirus protection.

[TB]CR clustering

• 15 thousand cells \(\rightarrow\) over 100 million pairwise distances. So don’t calculate a distance matrix directly.
• Use greedy clustering (CD-HIT; Li and Godzik [2006]) to define equivalence classes, then cluster within each set to yield “isomorphs.”
• It’s tempting to call them “clones,” but that’s generally too stringent of a term.

Statistical models for single cell [TB]CR data

• Interested in some statistic \(Y_i\), design \(x_i\), and a replicate index \(R(i)\). For instance, \(Y_i\) is the CDR3 length, \(x_i\) is a tetramer specificity and \(R(i)\) indexes subjects.
• Given some parametric assumptions on \(Y_i\), we might estimate and test parameters \(\beta\) with (generalized-)linear mixed models \[ E(Y_i | Z_{R(i)}) = {\mathbf{x}}_i \beta + Z_{R(i)}, \] \(Z \sim N(0, \tau^2)\).

CDR3 length

Here \(Y_i\) is Gaussian.

Class switching rates

Here \(Y_i\) is Bernoulli.

Expanded isomorphs

Tests for isomorph expansion

Omnibus models via permutation

• \(Y_i\) is discrete, and \(R_i\) is nested within \(x_i\), so we can write \(X_r\). Let \(\mathcal{I}_r\) be the cell indices in replicate \(R\). We are interested in some function \(T_r = T( \{ Y_i: i \in \mathcal{I}_r \} )\), a function that depends only on the CDR3 properties of a replicate. For example, \(T_r\) counts the average multiplicity of isomorphs within a replicate, and we might test that CD31+ replicates have higher average multiplicity.

• Consider a test of independence between a statistic : \[ P(T_r | X_r) = P(T_r) \]

• If we permute cells across replicates, we break link between \(Y_i\) and \(R_i\), then we will also have broken link between \(T_r\) and \(X_r\)

Example

\(T_r\) is the number of isomorphs present at \(>1\) copies in each replicate. Test if \(T(\text{CD31}^-) - T(\text{CD31}^+)>0\):

chain	observed	expected(se)	p.value
TRA	197	121 (0.6)	<0.01
TRB	230	65 (0.5)	<0.01

Bulk [TB]CR

Pros:

• Cheaper
• Much greater depth

Cons:

• Lose cellular-level resolution, instead associate a sample-frequency \(F_{is}\) with each sequence \(s\) and replicate \(i\). Can still estimate \(E(F_{is}) = x_{is} \beta\).
  
• Best case scenario is that \(F_{is} = \frac{ \text{# in } i \text{ with } s}{\text{# in } i}\). Implicitly must condition on the (effective) number of cells sequenced.
• No obvious way to derive permutation tests
• No pairing (without excessive cleverness)

Combining single immune cell repertoire sequencing and functional expression

1. Conditional tests – differential expression, between isomorphs
2. Quality control – how to tell cells, from doublets, from droplets containing only ambient RNA
3. Multiview clustering

Conditional tests

• What is an isomorph? Use 5’ expression to find latent subpopulations (cell types), then test for isomorph enrichment.
  
• What does an isomorph do? Test for differential expression between isomorphs (maybe within subtype).

Ex: synovial B cells

• If you aren’t careful, 5’ expression confounded with isomorphs.
• Might be a feature, rather than a bug: BCR inference possible from expression alone? (With some annoying pipeline to write…)

Synovial B cells

Synovial B cell isomorphs

Joint QC

• Repertoire pool provides additional information. Rather than hard thresholds, well-calibrated probabilities would be best.
• DNA-barcoded antibodies (Citeseq/Reapseq/cell hashing) probably offer a better “multiview” way forward
• But:

Multiplets and non-canonical [TB]CR receptors

Conclusions

• The adaptive immune system plays a large, pathogenic and heterogenous role in autoimmune disease
• The adaptive immune system (eg CAR T cells) is poised to play a large, beneficial role in cancer treatment
• Multidomain data complexity feels like \(\mathcal{O}( NP^D)\), but modest attempts at integration will bare fruit.
• Immune cell receptor data is fun

Acknowledgments

Daniel Lu
Bill Robinson

Kristin Scheible
Jennifer Anolik
Nathan Laniewski
UR Genomics Core

Software availability

https://github.com/amcdavid/CellaRepertorium

Appendix

10X protocol

1. First strand synthesis with GEM Bead which has cell barcode + UMI
2. PCR + target enrichment (using TCR/BCR constant regions)
3. 5" expression from chopping off 3’ end and performing end repair

From 10X “Assay Scheme and Configuration” technical note

10X bioinformatic pipeline

• Read trimming: trim adapters and primers from ends
• Read filtering: keep reads with at least some alignment to VDJ segments
• Assembly: de novo from each cell into contigs
• Contig productivity: Full length = both V and J. Productive: recognizable start codon, CDR3 is in-frame, no non-sense mutations.
• Contig filtering. When there are >2 contigs with distinct CDR3, they must be supported by >2 UMI and 20% of the top contig.