domain: phys sdv

text:
Thomas Dupic

Meriem Bensouda Koraichi

Anastasia Minervina

Mikhail Pogorelyy

Thierry Mora

Aleksandra M Walczak

Immune Fingerprinting through Repertoire Similarity

Immune repertoires provide a unique fingerprint reflecting the immune history of individuals, with potential applications in precision medicine. However, the question of how personal that information is and how it can be used to identify individuals has not been explored. Here, we show that individuals can be uniquely identified from repertoires of just a few thousands lymphocytes. We present "Immprint," a classifier using an information-theoretic measure of repertoire similarity to distinguish pairs of repertoire samples coming from the same versus different individuals. Using published T-cell receptor repertoires and statistical modeling, we tested its ability to identify individuals with great accuracy, including identical twins, by computing false positive and false negative rates < 10 -6 from samples composed of 10,000 T-cells. We verified through longitudinal datasets and simulations that the method is robust to acute infections and the passage of time. These results emphasize the private and personal nature of repertoire data.

INTRODUCTION

Personalized medicine is a frequent promise of nextgeneration sequencing. These high-throughput and lowcost sequencing technologies hold the potential of tailored treatment for each individual. However, progress comes with privacy concerns. Genome sequences cannot be anonymized: a genetic fingerprint is in itself enough to fully identify an individual, with the rare exception of monozygotic twins. The privacy risks brought by these pseudonymized genomes have been highlighted by multiple studies [START_REF] Homer | Resolving Individuals Contributing Trace Amounts of DNA to Highly Complex Mixtures Using High-Density SNP Genotyping Microarrays[END_REF][START_REF] Naveed | Privacy in the Genomic Era[END_REF][START_REF] Sweeney | Identifying Participants in the Personal Genome Project by Name[END_REF], and the approach is now routinely used by law enforcement. Sequencing experiments that focus on a limited number of expressed genes should be less prone to these concerns. However, as we will show, Band T-cell receptor (BCR and TCR) genes are an exception to this rule. BCR and TCR are randomly generated through somatic recombination [START_REF] Hozumi | Evidence for somatic rearrangement of immunoglobulin genes coding for variable and constant regions[END_REF], and the fate of each B-or T-cell clone depends on the environment and immune history. The immune T-cell repertoire, defined as the set of TCR expressed in an individual, has been hailed a faithful, personalized medical record, and repertoire sequencing (RepSeq) as a potential tool of choice in personalized medicine [START_REF] Robins | Immunosequencing: applications of immune repertoire deep sequencing[END_REF][START_REF] Attaf | αβ T cell receptors as predictors of health and disease[END_REF][START_REF] Woodsworth | Sequence analysis of T-cell repertoires in health and disease[END_REF][START_REF] Bradley | Using T Cell Receptor Repertoires to Understand the Principles of Adap-tive Immune Recognition[END_REF][START_REF] Davis | Recent progress in the analysis of αβ T cell and B cell receptor repertoires[END_REF]. In this report, we describe how, from small quantities of blood (blood spot or heel prick), one can extract enough information to uniquely identify an individual, providing an immune fingerprint. The "Immprint" classifier analyzes this immune fingerprint to decide whether two samples were sampled from the same individual. * Corresponding authors. These authors contributed equally.

RESULTS

Given two samples of peripheral blood containing respectively M 1 and M 2 T-cells, we want to distinguish between two hypothetical scenarios: either the two samples come from the same individual ("autologous" scenario), or they were obtained from two different individuals ("heterologous" scenario), see Fig. 1a.

TCR are formed by two protein chains α and β. They each present a region of high somatic variability, labeled CDR3α and CDR3β, randomly generated during the recombination process. These regions are coded by short sequences (around 50 nucleotides), which are captured by RepSeq experiments. The two chains are usually not sequenced together so that the pairing information between α and β is lost. Most experiments focus on the β chain, and when not otherwise specified, the term "receptor sequence" in this paper will refer only to the nucleotide sequence of the TRB gene coding for this β chain (which include CDR3β). Similarly, as most cells expressing the same beta chain are clonally related, we will be using the terms "clone" and "clonotype" to refer to set of cells with the same nucleotide TRB sequence, even if they were produced in separate generation events and are not real biological clones (since we have no means of distinguishing the two cases). CDR3β sequences are very diverse, with more than 10 40 possible sequences [START_REF] Mora | Quantifying lymphocyte receptor diversity[END_REF]. For comparison, the TCRβ repertoire of a given individual is composed of 10 8 to 10 10 unique clonotypes [START_REF] Qi | Diversity and clonal selection in the human T-cell repertoire[END_REF][START_REF] Lythe | How many TCR clonotypes does a body maintain?[END_REF]. As a result, most of the sequences found in a repertoire are "private". cleotide receptor sequences shared between the two samples, which we call S. Samples coming from the same individual should have more receptors in common because T-cells are organized in clones of cells carrying the same TCR. By contrast, S should be low in pairs of samples from different individuals, in which sharing is due to rare convergent recombinations. Appropriately setting a threshold to jointly minimize the rates of false positives and false negatives (Fig. 1b), we can use S as a classifier to distinguish autologous from heterologous samples. S is not normalized for sequencing depth and values of S should not be compared between samples of different size.

The S score can be improved upon by exploiting the fact that some receptors are much more likely than others to be generated during V(D)J-recombination, with variations in generation probability (P gen , [START_REF] Murugan | Statistical inference of the generation probability of Tcell receptors from sequence repertoires[END_REF][START_REF] Marcou | Highthroughput immune repertoire analysis with IGoR[END_REF][START_REF] Sethna | Population variability in the generation and thymic selection of T-cell repertoires[END_REF]) spanning 15 orders of magnitude. Public sequences (with high P gen ) are likely to be found in multiple individuals [START_REF] Venturi | The molecular basis for public T-cell responses?[END_REF], while rare sequences (low P gen ) are unlikely to be shared by different individuals, and thus provide strong evidence for the autologous scenario when found in both samples. To account for this information, we define the score:

I = shared s [ln (1/P gen (s)) -γ] , (1)
which accounts for Shannon's "surprise" ln(1/P gen )-a measure of unexpectedness-associated with each shared sequence s, so that rare shared sequences count more than public ones. The constant γ depends on the repertoire's clonal structure and is set to 12 in the following (see Methods for an information-theoretic derivation). P gen is computed using models previously trained on data from multiple individuals [START_REF] Marcou | Highthroughput immune repertoire analysis with IGoR[END_REF]. Small differences reported between the P gen of distinct individuals justify the use of a universal model [START_REF] Sethna | Population variability in the generation and thymic selection of T-cell repertoires[END_REF].

Measuring error rates

We tested the classifiers based on the S and I scores on TCRβ RepSeq datasets from 656 individuals [START_REF] Emerson | Immunosequencing identifies signatures of cytomegalovirus exposure history and HLAmediated effects on the T cell repertoire[END_REF] labeled according to their cytomegalovirus (CMV) serological status. Sequences were downsampled to mimic experiments where M 1 = M 2 = M cells were sampled and their receptors sequenced. The frequency of a particular clonotype in the sample (the proportion of cells expressing a particular beta chain) was estimated using the read counts of unique TCRβ sequences, and the mean values of S and I computed with a procedure designed to correct for the limited diversity of the sampled repertoire relative to the full repertoire, see Methods . Similar results may be obtained when M 1 and M 2 are different (see Methods). The clones most often shared between two autologous samples are also the most clonally expandedand hence are probably antigen experienced. We verified that the sequence statistics of those expanded clonotypes did not differ from generic ones (Fig. S1).

In Fig. 1c, we plot the mean value of S (over many draws of M = 5000 cells) for each individual (autologous scenario, in blue) and between pairs of different individuals (heterologous scenario, in yellow). The two scenarios are clearly discernable under both scores. This result holds for pairs of monozygotic twins obtained from a distinct dataset [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] (green dots), consistent with previous reports that twins differ almost as much in their repertoires as unrelated individuals [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF][START_REF] Zvyagin | Distinctive properties of identical twins' TCR repertoires revealed by high-throughput sequencing[END_REF][START_REF] Tanno | Determinants governing T cell receptor α/β-chain pairing in repertoire formation of identical twins[END_REF]. Heterologous scores (yellow dots) vary little, and may be bounded from above by a theoretical prediction (dashed line) based on a model of recombination [START_REF] Elhanati | Predicting the spectrum of TCR repertoire sharing with a data-driven model of recombination[END_REF] (see Methods). On the other hand, autologous scores (blue dots) show several orders of magnitude of variability across individuals. These variations stem from the clonal structure of the repertoire, and correlates with measures of diversity (Fig. S2), which is known to vary a lot between individuals and correlates with age [START_REF] Britanova | Dynamics of Individual T Cell Repertoires: From Cord Blood to Centenarians[END_REF], serological status, and infectious disease history [START_REF] Sylwester | Broadly targeted human cytomegalovirus-specific CD4+ and CD8+ T cells dominate the memory compartments of exposed subjects[END_REF][START_REF] Khan | Cytomegalovirus Seropositivity Drives the CD8 T Cell Repertoire Toward Greater Clonality in Healthy Elderly Individuals[END_REF]. To explore the worst case scenario of discriminability, hereafter we will focus on the individual with the lowest autologous S found in the dataset.

The sampling process introduces an additional source of variability within each individual. Two samples of blood from the same individual do not contain the exact same receptors, and the values of S and I is expected to vary between replicates. Examples of these variations are shown in Fig. 1d. The blue (respectively yellow) curves correspond to the sample distributions in the autologous (heterologous) scenario for different individuals (pairs of individuals). The distribution of S is well-approximated by a Poisson distribution, while I follows approximately a compound distribution of a normal and Poisson distributions (see Methods for details). Armed with these statistical models of variations, we can predict upper bounds for the false negative and false positive rates. As seen from the detection error trade-off (DET) graph Fig. 1e, the Immprint classifier performs very well for a few thousand receptors with an advantage for I.

With 10, 000 cells, corresponding to ∼ 10 µL of blood, Immprint may simultaneously achieve a false positive rate of < 10 -16 and false negative rate of < 10 -6 , allowing for the near-certain identification of an individual based on the I score in pairwise comparisons against the world population ∼ 10 10 . When a large reference repertoire has been collected (M 1 = 1, 000, 000, corresponding to ∼ 1mL of blood), an individual can be identified with just 100 cells (Fig. S3).

The AUROC estimator (Area Under the Curve of the Receiver Operating Characteristic), a typical measure of a binary classifier performance, can be used to score the quality of the classifier with a number between 0.5 (chance) and 1 (perfect classification). The I score outperforms the S score (Fig. 1f), particularly above moderate sample sizes (M ≈ 5000). Both scores can be readily generalized to the case of paired receptors TCRαβ, when the pairing of the two chains is available (through singlecell sequencing [START_REF] Dash | Paired analysis of TCRα and TCRβ chains at the single-cell level in mice[END_REF][START_REF] Redmond | Single-cell TCRseq: Paired recovery of entire T-cell alpha and beta chain transcripts in T-cell receptors from single-cell RNAseq[END_REF][START_REF] Grigaityte | Single-cell sequencing reveals αβ chain pairing shapes the T cell repertoire[END_REF] or computational pairing [START_REF] Howie | High-throughput pairing of T cell receptor a and b sequences[END_REF]), using P gen (α, β) = P gen (α) × P gen (β) [START_REF] Dupic | Genesis of the Aβ T-cell receptor[END_REF] for the gener-

FIG. 1. a)

The two samples A and B can either originate from the same individual (autologous) or two different individuals (heterologous). In both scenarios, sequences can be shared between the two samples, but their quantity and quality vary. b) Schematic representation of the distribution of the S or I scores for multiple pairs of samples extracted from the same individual (in the autologous scenario) or the same pair of individuals (heterologous). The dashed vertical line represents the threshold value. c) Expected value of S and I for different pairs of samples, sampled from the same individual (in blue) or different ones (yellow). Green dots represent samples extracted from pairs of identical twins. The dashed lines represents the theoretical upper bound for heterologous repertoires (see Methods) for both S and I (γ = 12). d) Sampling distributions of S for 6 different patients (autologous case, each blue curve is one patient) or 6 different pairs of patients (heterologous, each yellow or green curve is a pair of patients) for M = 5000. The y-axis scale on the left is adapted to the heterologous distributions while the scale on the right corresponds to the (much wider) autologous ones. The 3 sampling distributions in green correspond to a pair of samples extracted from identical twins. e) Detection Error Trade-off (DET) graph for both summary statistics and different sample sizes M . I (γ = 12) outperforms S in all scenarios. f) AUROC (Area Under Receiver Operating Characteristic), as a function of M . The AUROC is a traditional measure of the quality of a binary classifier (a score closer to one indicates a better classifier). The results are shown for S and I both in the default case (only the β chain considered) or for the full (α-β) receptor.

ation probability of the full TCR. Because coincidental sharing of both chains is substantially rarer than with the β chain alone, using the paired chain information greatly improves the classifier.

While this paper focuses on T-cells and TCR sequences, the structure of the B-cell receptors (BCR) repertoire is very similar to the TCR repertoire and we expect to find qualitatively similar results. As an example we use the dataset obtained in Ref. [START_REF] Briney | Commonality despite exceptional diversity in the baseline human antibody repertoire[END_REF] to measure S and I for IGH chains (forming half of the BCR receptor), Fig S4 . We see that 5000 IgG+ memory B-cells are enough to reliably identify the individuals in the study. However, B-cells are 5 times less common in peripheral blood than T-cells, and somatic hypermutation tends to distort the statistics of receptors, reducing the reliability of our classifier. Hence, for practical applications, T-cells are a better means of identification.

Evolution with time

The previous results used samples obtained at the same time. However, immune repertoires are not static: interaction with pathogens and natural aging modify their composition. The evolution of clonal frequencies will decrease Immprint's reliability with time, especially if the individual has experienced immune challenges in the meantime.

To study the effect of short-term infections, we analyzed an experiment where 6 individuals were vaccinated with the yellow fever vaccine, which is regarded as a good model of acute infection, and their immune system was monitored regularly through blood draws [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF]. We observe an only moderate drop in I caused by vaccination (Fig. 2 a). This is consistent with the fact that infections lead to the strong expansion of only a limited number of clones, while the rest of the immune system stays stable [START_REF] Dewitt | Dynamics of the Cytotoxic T Cell Response to a Model of Acute Viral Infection[END_REF][START_REF] Wolf | Identifying and Tracking Low-Frequency Virus-Specific TCR Clonotypes Using High-Throughput Sequencing[END_REF][START_REF] Qi | Diversification of the antigen-specific T cell receptor repertoire after varicella zoster vaccination[END_REF][START_REF] Sycheva | Quantitative profiling reveals minor changes of T cell receptor repertoire in response to subunit inactivated influenza vaccine[END_REF]. While other types of infections, auto-immune diseases, and cancers may affect Immprint in more substantial ways, our result suggests that it is relatively robust to changes in condition.

We then asked how stable Immprint is over long times. Addressing this issue is hampered by the lack of longitudinal datasets over long periods, so we turn to mathematical models [START_REF] Lythe | How many TCR clonotypes does a body maintain?[END_REF][START_REF] Borghans | Quantification of T-cell dynamics: From telomeres to DNA labeling[END_REF][START_REF] Thomas-Vaslin | Comprehensive Assessment and Mathematical Modeling of T Cell Population Dynamics and Homeostasis[END_REF][START_REF] Desponds | Fluctuating fitness shapes the clone-size distribution of immune repertoires[END_REF][START_REF] De Greef | The naive t-cell receptor repertoire has an extremely broad distribution of clone sizes[END_REF] to describe the dynamics of the repertoire. Following the model of fluctuating growth rate described in Ref. [START_REF] Desponds | Fluctuating fitness shapes the clone-size distribution of immune repertoires[END_REF], we define two typical evolutionary timescales for the immune system: τ , the typical turnover rate of T-cell clones, and θ, which represents the typical time for a clonotype to grow or shrink by a factor two as its growth rate fluctuates. The model predicts a power-law distribution for the clone-size distribution, with exponent -1 -τ /2θ. This exponent has been experimentally measured to be ≈ -2 [START_REF] Mora | Quantifying lymphocyte receptor diversity[END_REF], which leaves us with a single parameter τ , and θ = τ /2. An example of simulated evolution of Immprint with time is shown in Fig. 2 b. The highlighted histogram represents a data point at two years obtained from [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF]. While a fit is possible for this specific individual (τ = 0.66 in Fig. 2 b), the τ parameter is not universal, and we expect it to vary between individuals, especially as a function of age. In Fig. 2c we explore the effect on the stability of Immprint for a range of reasonable values for the clone turn-over rate τ , from 6 months to 10 years, encompassing both previous estimates of the parameter [START_REF] Desponds | Fluctuating fitness shapes the clone-size distribution of immune repertoires[END_REF] and measured turnover rates for different types of T-cells [START_REF] Boer | Quantifying T Lymphocyte Turnover[END_REF]. While Fig. 2 focuses on I, the behaviour is similar for S (see Fig. S5). We observe that under this model, for most individuals and bar exceptional events, Immprint should conserve its accuracy for years or even decades.

DISCUSSION

In summary, the T-cells present in small blood samples provide a somatic and long-lived barcode of human individuality, which is robust to immune challenges and stable over time. While the uniqueness of the repertoire was a well known fact, we demonstrated that the most common T-cells clones are still diverse enough to uniquely define an individual and frequent enough to be reliably sampled multiple times. Unlike genome sequencing, repertoire sequencing can discriminate monozygotic twins with the same accuracy as unrelated individuals. However, a person's unique immune fingerprint can be completely wiped out by a hematopoietic stem cell transplant [START_REF] Buhler | Genetic T-cell receptor diversity at 1 year following allogeneic hematopoietic stem cell transplantation[END_REF]. The different datasets used cover a range of different sequencing methods , however different approaches may lead to slightly different threshold choices. In particular, in practical implementations, sequencing depth is an important concern. One needs enough coverage to sequence TCRβ genes from as many as possible of the T-cells present in the sample, in order to measure a more precise a immune fingerprint. In addition, the specific calculations presented here only apply to peripheral blood cells. Specific cell types or cells extracted from tissue samples may have different clonal distributions and potentially different receptor statistics. For example the value of S in the autologous case varies between CD4+ and CD8+ T-cells (Fig. S6), different individuals remain distinguishable using each subset.

Immprint is implemented in a python package and webapp (see Methods) allowing the user to determine the autologous or heterologous origin of a pair of repertoires. Beyond identifying individuals, the tool could be used to check for contamination or labelling errors between samples containing TCR information. The repertoire information used by Immprint can be garnered not only from RepSeq experiments, but also from RNA-Seq experiments, which contain thousands of immune receptor transcripts [START_REF] Li | Ultrasensitive detection of TCR hypervariable-region sequences in solid-tissue RNA-seq data[END_REF][START_REF] Bolotin | Antigen receptor repertoire profiling from RNA-seq data[END_REF]. Relatively small samples of immune repertoires are enough to uniquely identify an individual even among twins, with potential forensics applications. At the same time, unlike genetic data from genomic or mRNA sequencing, Immprint provides no information [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] and the blue curves show theoretical estimates, fitted to match (τ = 0.66). c) Evolution of the (normalized) mean of I (M=5000) as a function of time for different values of the turnover rate τ . The dashed line represents the threshold value divided by the smallest value of It=0 (M=5000) in the data. The data points were obtained from the datasets [START_REF] Chu | Longitudinal immunosequencing in healthy people reveals persistent T cell receptors rich in highly public receptors[END_REF] (yellow), [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] (green) and [START_REF] Britanova | Dynamics of Individual T Cell Repertoires: From Cord Blood to Centenarians[END_REF] (orange). Different markers indicate different individuals.

about kin relationships, very much like classical fingerprints, and avoids privacy concerns about disclosing genetic information shared with non consenting relatives.

METHODS

Datasets & Pre-processing

We use five independant RepSeq datasets in this study: (i) genomic DNA from Peripheral blood mononuclear cells (PBMCs) from 656 healthy donors [START_REF] Emerson | Immunosequencing identifies signatures of cytomegalovirus exposure history and HLAmediated effects on the T cell repertoire[END_REF]; (ii) cDNA of PBMCs sampled from three pairs of twins, before and after a yellow-fever vaccination [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF]; (iii), (iv) two longitudinal studies of healthy adults [START_REF] Britanova | Dynamics of Individual T Cell Repertoires: From Cord Blood to Centenarians[END_REF][START_REF] Chu | Longitudinal immunosequencing in healthy people reveals persistent T cell receptors rich in highly public receptors[END_REF];(v) cDNA dataset of IGH genes (B-cells) from 9 individuals, with multiple replicates [START_REF] Briney | Commonality despite exceptional diversity in the baseline human antibody repertoire[END_REF].

CDR3 nucleotide sequences were extracted with MIGEC [START_REF] Shugay | Towards error-free profiling of immune repertoires[END_REF] (for the second dataset) coupled with MiXCR [START_REF] Bolotin | MiXCR: Software for comprehensive adaptive immunity profiling[END_REF]. We also extract the frequency of reads from the three datasets. The non-productive sequences were discarded (outof-frame, non-functional V gene, or presence of a stop codon). The generation probability (P gen ) was computed using OLGA [START_REF] Sethna | OLGA: Fast computation of generation probabilities of B-and T-cell receptor amino acid sequences and motifs[END_REF], with the default TCRβ model. The frequency of each clone was estimated by summing the frequencies of all reads that shared the same nucleotide CDR3 sequence and identical V,J genes.

The preprocessing code is distributed on the Git repository associated with the paper. We also developed a command-line tool (https://github.com/statbiophys/immprint) that discriminates between sample origins, and a companion webapp (https://immprint.herokuapp.com).

Discrimination scores

To discriminate between the autologous and heterologous scenarios, we introduce a log-likelihood ratio test between the two possibilities:

I = s ln P (y 1 (s), y 2 (s)|autologous) P (y 1 (s), y 2 (s)|heterologous) , (2)
where y 1 (s) = 1 if the sequence s is found in sample 1, and 0 otherwise; likewise y 2 (s) = 1 if s is in sample 2. The sum runs over all potential sequences s, including unseen ones. To be present in a sample, a sequence s first has to be present in the repertoire. This occurs with probability 1 -(1 -p(s)) Nc , where N c is the total number of clonotypes in the repertoire, and p(s) is the probability of occurence of sequence s (resulting from generation and selection, see below). Second, it must be picked in a sample of size M , with probability 1 -(1 -f ) M ≈ M f (assuming M f 1) depending on its frequency f , which is distributed according to the clone size distribution ρ(f ). We checked that f (s) and P gen (s) were not correlated (Fig. S1) and that the effects of a shared infection between different individuals were limited to a handful of clones (Fig. S7). Then one can write

P (y 1 (s) = 1, y 2 (s) = 1 | autologous) ≈ 1 -e -Ncp(s) M 1 M 2 df ρ(f ) f 2 , (3)
P (y 1 (s) = 1, y 2 (s) = 0 | autologous) ≈ 1 -e -Ncp(s) M 1 N c and 1 ↔ 2, (4)
P (y 1 (s) = 0, y 2 (s) = 0 | autologous) ≈ 1 -1 -e -Ncp(s) M 1 + M 2 N c , (5)
where we've used df ρ(f ) f = 1/N c . For the heterologous case the probability factorizes as:

P (y 1 (s), y 2 (s) | heterologous) = P 1 (y 1 (s))P 2 (y 2 (s)), (6)
with

P a (y a (s) = 1) ≈ 1 -e -Ncp(s) M a N c , a = 1, 2. ( 7
)
Since only the term y 1 (s) = y 2 (s) = 1 (shared sequences) is different between the autologous and heterologous cases, we obtain:

I = shared s ln(N 2 c f 2 ) -ln 1 -e -Ncp(s) . ( 8
)
Further assuming N c p(s) 1, and p(s) = P gen (s)q -1 (where q accounts for selection [START_REF] Elhanati | Predicting the spectrum of TCR repertoire sharing with a data-driven model of recombination[END_REF] and P gen (s) is the probability of sequence generation [START_REF] Marcou | Highthroughput immune repertoire analysis with IGoR[END_REF]), the score simplifies to Eq. 1, with γ = -ln(qN c f 2 ) = ln(q -1 f / f 2 ). The factor γ depends on unknown parameters of the model, but can be estimated assuming a power-law for the clone size distribution [START_REF] Touzel | Inferring the immune response from repertoire sequencing[END_REF], ρ(f ) ∝ f -2 extending from f = 10 -11 to f = 0.01, and q = 0.01 [START_REF] Elhanati | Predicting the spectrum of TCR repertoire sharing with a data-driven model of recombination[END_REF], yielding γ ≈ 12.24. Alternatively we optimized γ to minimize the AUROC, yielding γ ≈ 15 (SI Fig. S8). Since performance degrades quickly for larger values, we conservatively set γ = 12.

Estimating mean scores from RepSeq datasets

The sampling of M cells from blood is simulated using large repertoire datasets. In a bulk repertoire sequencing dataset, the absolute number of cells for each clonotype (cells with a specific receptor) is unknown, but the fraction of each clonotype can be estimated using the proportion of reads that are associated with this specific receptor. To estimate the autologous S and I of two samples of size M 1 and M 2 in the absence of true replicates, we computed their expected values from a single dataset containing N reads, from which two random subsamples of sizes M 1 and M 2 were taken. The mean value of S is equal to

S = s (1 -(1 -f (s)) M1 )(1 -(1 -f (s)) M2
), where f (s) is the true (and unknown) frequency of sequence s. A naive estimate of mS may be obtained by repeatly resampling subsets of sizes M 1 and M 2 from the observed repertoire, calculate S for each draw, and average. One get the same result by replacing f (s) by fs = n(s)/N in the previous formula, where n(s) is the number of s reads in the full dataset, and N = s n(s). However, this naive estimate leads to a systematic overestimate of the sharing (visible when compared with biological replicates, see Fig. S9), simply because this procedure overestimates the probability of resampling rare sequences, in particular singletons whose true frequency may be much lower that 1/N . A similar bias occurs when computing I. To correct for this bias, we look for a function h(n) that satisfies for all f :

h(n) ≡ n N n f n (1 -f ) N -n h(n) = (1 -(1 -f ) M1 ) (1 -(1 -f ) M2 ), (9)
so that S and I can be well approximated by:

S ≈ s h(n(s)), (10)
I ≈ - s h(n(s)) [ln(1/P gen (s)) -γ] . (11)
Expanding the right-hand side of Eq. 9 into 4 terms, we find that h

(n) = 1 -g M1 (n) -g M2 (n) + g M1+M2 (n) satisfies Eq. 9 provided that: n N n f n (1 -f ) N -n g M (n) = (1 -f ) M . (12)
Under the change of variable x = f /(1 -f ), the expression becomes:

n N n x n g M (n) = (1 + x) N -M = n N -M n x n . ( 13
)
Identifying the polynomial coefficients in x n on both sides yields:

g M (n) = N -M n N n . (14)
These corrected estimates agree with the direct estimates using biological replicates (Fig. S9).

Similarly, S and I in heterologous samples can be estimated using:

S ≈ s [1 -g M1 (n(s))][1 -g M2 (n (s))], (15)
I ≈ s [1 -g M1 (n(s))][1 -g M2 (n (s))] [ln(1/P gen (s)) -γ] . (16)
where n(s) and n (s) are the empirical counts of sequence s in the two samples.

Theoretical upper bound on heterologous scores

When the two samples were extracted from two different individuals (heterologous scenario), we can use the universality of the recombination process to give upper bounds on the values of S and I. These bounds are represented by the dashed lines in Fig1c). If two samples of respectively M 1 and M 2 cells, containing T 1 ≤ M 1 and T 2 ≤ M 2 unique sequences are extracted from two different individuals, the number of shared sequences between them is given by [START_REF] Elhanati | Predicting the spectrum of TCR repertoire sharing with a data-driven model of recombination[END_REF]:

S heterologous ≤ s 1 -(1 -p(s)) T1 1 -(1 -p(s)) T2 (17) T 1 T 2 s p(s) 2 = T 1 T 2 p(s) ≤ M 1 M 2 p(s) . (18)
p(s) is the probability of finding a sequence s in the blood. Following [START_REF] Elhanati | Predicting the spectrum of TCR repertoire sharing with a data-driven model of recombination[END_REF], we make the approximation p(s) = P gen (s)q -1 , where the q = 0.01 factor is the probability that a generated sequences passes selection. Then p(s) can be estimated from the mean over generated sequences. Similarly, we can estimate I as

I heterologous T 1 T 2 s p(s) 2 [ln (1/P gen (s)) -γ] (19)
= -T 1 T 2 p(s)[γ + ln(qp(s))] ≤ -M 1 M 2 p(s)[γ + ln(qp(s))] , (20)
which is also estimated from the mean over generated sequences.

Error rate estimates

To make the quantitative predictions shown in Fig. 1, we need to constrain the tail behavior of the distributions of S and I, for the two scenarios.

The S statistic can be rewritten as a sum of Bernouilli variables over all possible sequences, each with a parameter corresponding to its probability of being present in both samples, either in the autologous or the heterologous case. Therefore S is a Poisson binomial distribution, a sum of independent Bernouilli variables with potentially different parameters. The variance and tails of that distribution are bounded by those of the Poisson distribution with the same mean, denoted by m a for the autologous case, and m h for the heterologous case (Fig. S10).

Thanks to that inequality, the rates of false negatives and false positives for a given threshold r are bounded by:

P (S < r|autologous) ≤ Q(r + 1, m a ), P (S > r|heterologous) ≤ 1 -Q(r + 1, m h ), (21)
where Q is the regularized gamma function, which appears in the cumulative distribution function of the Poisson distribution. The mean autologous score m a is estimated from experimental data: we use the smallest value of S in the Emerson dataset and Eq. 10. To compute m h , we use the semi-theoretical prediction made using the universality of the recombination process, Eq. 17.

Similarly, I can be viewed as a sum of S independent random variables, all following the distribution of ln(1/P gen )γ. However, this distribution differs in the two scenarios. Sequences shared between more than one donor have an higher probability of being generated, their ln(P gen ) distribution has higher mean and smaller variance (Fig. S11).

The sum is composed of a relatively large number of variables in most realistic scenarios. Hence, we rely on the central limit theorem to approximate it by a normal distribution, of mean and variance proportional to S. Explicitly:

P (I < r|autologous) = 1 2 ∞ S=0 (m a ) S e -ma S! 1 + erf r -S ln(1/P gen ) -γ 2SVar[ln(1/P gen ) -γ] , (22)
P (I > r|heterologous) = 1 2 ∞ S=0 (m h ) S e -m h S! 1 -erf r -S ln(1/P gen ) -γ shared 2SVar[ln(1/P gen ) -γ] shared . ( 23
)
The AUROC are computed based on these estimates, by numerically integrating the true positive rate P (S, I < r|heterologous) with respect to the false negative rate P (S, I < r|autologous) as the threshold r is varied.

Modeling the evolution of autologous scores

We use the model of Ref. [START_REF] Desponds | Fluctuating fitness shapes the clone-size distribution of immune repertoires[END_REF] to describe the dynamics of individual T-cell clone frequencies f under a fluctuating growth rate reflecting the changing state of the environment and the random nature of immune stimuli:

df dt = - 1 τ + 1 2θ + 1 √ θ η(t) f (t), (24)
where η(t) is a Gaussian white noise with η(t) = 0 and η(t)η(t ) = δ(t -t ).

With the change of variable x = ln(f ), these dynamics simplify to a simple Brownian motion in log-frequency: ∂ t x = -τ -1 + θ -1/2 η(t). In that equation, τ appears as the decay rate of the frequency, while θ is the timescale of the noise, interpreted as the typical time it takes for the frequency to rise or fall by a logarithmic unit owing to fluctuations. Considering a large population of clone, each with their independent frequency evolving according to Eq. 24, and a source term at small f corresponding to thymic exports, one can show that the steady-state probability density function of f follows a power-law [START_REF] Desponds | Fluctuating fitness shapes the clone-size distribution of immune repertoires[END_REF], ρ(f ) ∝ f -α , with exponent α = 1 + 2θ/τ . α was empirically found to be ≈ 2 in a wide variety of immune repertoires [START_REF] Mora | Quantifying lymphocyte receptor diversity[END_REF][START_REF] Touzel | Inferring the immune response from repertoire sequencing[END_REF][START_REF] Weinstein | High-throughput sequencing of the zebrafish antibody repertoire[END_REF][START_REF] Oakes | Quantitative characterization of the T cell receptor repertoire of naïve and memory subsets using an integrated experimental and computational pipeline which is robust, economical, and versatile[END_REF], implying 2θ ≈ τ . The turn-over time τ is unknown, and was varied from 1/2 year to 10 years in the simulations.

We the evolution of human TRB repertoires by starting with the empirical values of the frequencies of each observed clones, f (s, 0) = f (s, 0) = n(s, 0)/N from the analysed datasets. A sample of size M was randomly selected with respect to these frequencies, and the frequencies of the clones captured in that sample were then evolved with a time-step of 2 days using Euler-Maruyama's method, which is exact in the case of Brownian motion. Clones with frequencies falling below 10 -11 (corresponding to a single cell in the organism) were removed. At each time t > 0, we measured the mean value of S with the formula s (1 -(1 -f (s, t)) M ) where the sum runs over the sequences captured in the initial sample. [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] and the blue curves show theoretical estimates, fitted to match (τ = 0.66). c) Evolution of the (normalized) mean of S (M = 5000) as a function of time for different values of the turnover rate τ . The dashed line represents the threshold value divided by the smallest value of St=0 (M = 5000) in the data. The data points were obtained from the datasets [START_REF] Chu | Longitudinal immunosequencing in healthy people reveals persistent T cell receptors rich in highly public receptors[END_REF] (yellow), [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] (green) and [START_REF] Britanova | Dynamics of Individual T Cell Repertoires: From Cord Blood to Centenarians[END_REF] [START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] (6 individuals). In the autologous case (blue), there is a significant difference between the two cell types, caused by the differences between their clonal distributions. The top twenty clones, when shared between donors, have a smaller probability of being generated than expected by chance. This difference is likely to be driven by convergent selection against common pathogens, since CMV positive donors show a more prononced effect than CMV negative ones.

FIG. 2

 2 FIG.2. a) Evolution of I (M=5000) during vaccination, between a sample taken at day 0 (vaccination date) and at a later timepoint. Each color represents a different individual. Each pair timepoint/individual has two biological replicates. The dashed line represents the threshold value. b) Evolution of I between a sample taken at year 0 and a later timepoint. The red histogram corresponds to one of the individuals sampled in[START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] and the blue curves show theoretical estimates, fitted to match (τ = 0.66). c) Evolution of the (normalized) mean of I (M=5000) as a function of time for different values of the turnover rate τ . The dashed line represents the threshold value divided by the smallest value of It=0 (M=5000) in the data. The data points were obtained from the datasets[START_REF] Chu | Longitudinal immunosequencing in healthy people reveals persistent T cell receptors rich in highly public receptors[END_REF] (yellow),[START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] (green) and[START_REF] Britanova | Dynamics of Individual T Cell Repertoires: From Cord Blood to Centenarians[END_REF] (orange). Different markers indicate different individuals.

τ

  FIG. S2.Comparison between the mean of S (autologous case), and three common diversity measures: the number of unique sequences found in the dataset (top left), the Shannon index, -fs ln fs (top right), the Simpson index (bottom left), and the total number of reads in each datasets (bottom right). All the diversity measures show a strong correlation with S, but the correlation with the sequencing depth is low.

  FIG. S5. a) Evolution of S (M = 5000) during vaccination, between a sample taken at day 0 (vaccination date) and at a later timepoint. Each color represents a different individual. Each pair timepoint/individual has two biological replicates. The dashed line represents the threshold value. b) Evolution of S between a sample taken at year 0 and a later timepoint. The red histogram corresponds to one of the individuals sampled in[START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] and the blue curves show theoretical estimates, fitted to match (τ = 0.66). c) Evolution of the (normalized) mean of S (M = 5000) as a function of time for different values of the turnover rate τ . The dashed line represents the threshold value divided by the smallest value of St=0 (M = 5000) in the data. The data points were obtained from the datasets[START_REF] Chu | Longitudinal immunosequencing in healthy people reveals persistent T cell receptors rich in highly public receptors[END_REF] (yellow),[START_REF] Pogorelyy | Precise tracking of vaccineresponding T cell clones reveals convergent and personalized response in identical twins[END_REF] (green) and[START_REF] Britanova | Dynamics of Individual T Cell Repertoires: From Cord Blood to Centenarians[END_REF] (orange). Different markers indicate different individuals.

  FIG. S7. Left: Mean value of Pgen as a function of the rank of the clonotype, for generic sequences (blue) and sequences shared between more than two donors (orange). The mean stays flat indicating that the probability of being generated does not generally depend on the clonotype size. There is an exception (black rectangle), shown as a close-up on the right panel. The top twenty clones, when shared between donors, have a smaller probability of being generated than expected by chance. This difference is likely to be driven by convergent selection against common pathogens, since CMV positive donors show a more prononced effect than CMV negative ones.

= 15 )

 15 FIG. S8. Left panel: AUROC (Area Under Receiver Operating Characteristic) of I, as a function of γ (M = M1 = M2 = 5000). We observe an optimum near γ = 15. Right panel: AUROC as a function of M , for S, I(γ = 0), and I(γ = 15).

  FIG. S10. Comparison between the distribution of S obtained by computationally and repeatedly downsampling a single repertoire from Ref.[START_REF] Emerson | Immunosequencing identifies signatures of cytomegalovirus exposure history and HLAmediated effects on the T cell repertoire[END_REF] with M = 5, 000 (histogram), and a Poisson distribution of the same mean (full line).

  FIG. S11. Distribution of ln(Pgen) for generic sequences, and for sequences shared between heterologous samples.

ACKNOWLEDGMENTS

The study was supported by the European Research Council COG 724208. A.A.M and M.V.P. are supported by RSF-20-15-00351. The authors are grateful for the help of Natanael Spizak with the analysis of BCR repertoire datasets.