TCR Repertoires of Thymic Conventional and Regulatory T Cells: Identification and Characterization of Both Unique and Shared TCR Sequences Free

The ability to generate a rapid and sustained T cell response to external pathogens and transformed malignant cells is essential for the protection of the host. At the same time, deletion of autoreactive T cells during development and repression of excessive or autoreactive responses in peripheral tissues is essential to proper T cell protective function (1). This duality in the regulation of T cell function is accomplished by two types of T cells: conventional T cells that provide “helper” (CD4⁺) and “killer” (CD8⁺) functions and regulatory T cells (Treg) that suppress conventional T cell–dependent responses. Treg have been assigned to two subsets based on the origin of their generation: thymic Treg (tTreg) (or natural Treg) that develop in the thymus (2–5) and peripheral Treg generated in the periphery from thymic conventional CD4⁺ T cells (Tconv) under specific conditions (6–8). The development of tTreg appears to require both TCR signals and other factors, such as costimulatory signaling and cytokines, but the precise mechanisms of tTreg generation have not been fully elucidated.

A key factor in tTreg generation is the specificity of the TCR whose interaction with self-antigen/MHC plays a critical role in tTreg differentiation. Several studies using transgenic (Tg) mouse models suggest that the signal strength of TCR recognition of self-antigen/MHC ligand differs in tTreg and Tconv, with tTreg differentiation involving a higher level of signal strength (9, 10). Indeed, the disruption of normal self-antigen/MHC ligand expression in the thymus because of Aire deficiency causes a change in the fate of self-antigen–specific T cells from tTreg to Tconv (11, 12). Subsequent studies suggest that Ag presentation by different APCs (classical and plasmacytoid dendritic cells, cortical and medullary thymic epithelial cells, and B cells) at different thymic locations (cortex and medulla) influences the deletion of autoreactive thymocytes and the differentiation of tTreg (3, 13). In addition, factors such as cytokines (including IL-2 and TGF-β) (14–16) and costimulatory receptors (CD28) (17) have also been implicated in the development of tTreg. Collectively, it is clear that no single factor alone determines differentiation to the tTreg fate, but precisely how these factors act in combination remains to be determined.

Initial analyses of TCR sequences in Treg and Tconv of TCRα or TCRβ Tg mice reported that CDR3α and CDR3β sequence repertoires of Treg and Tconv are different either by exclusive appearance in only one of these lineages or by their relative abundance in Treg or Tconv when sequences were found in both cell types (18, 19). Subsequent studies using high-throughput sequencing generated larger numbers of TCR sequences in Treg and Tconv. Studies using TCRβ Tg mice to compare TCRα sequences between tTreg and Tconv reported little overlap of TCRα sequence between tTreg and effector T cells (20) or between tTreg and Tconv that recognize the same foreign Ag (21). Study of a TCRα Tg mouse to compare TCRβ sequences between Treg and Tconv from spleen and peripheral lymph nodes found that 12% of TCRβ sequences are shared by peripheral Treg and Tconv and are thus presumed to be derived from common progenitors (22). However, there has been no reported deep sequencing analysis examining endogenous TCRα and TCRβ of tTreg and Tconv from the thymus of non-TCR Tg mice. It is therefore unclear what degree of TCR sequence uniqueness and similarity exists overall between tTreg and Tconv or, importantly, whether there are general sequence features that distinguish αβ TCR of tTreg from those of Tconv.

In this study, we addressed the role of TCR sequence in determining whether T cells develop into tTreg or Tconv lineages. We report a comprehensive comparison of TCRα and TCRβ sequences of tTreg and Tconv using a Unique Molecular Identifier (UMI) methodology incorporating a 5′ single universal primer for PCR amplification of all V genes, significantly reducing PCR bias of amplification and sequencing errors affecting the quantitation of TCR frequency (23, 24). Comparison of TCRα and TCRβ sequences between tTreg and Tconv from two normal mouse strains revealed that, although many sequences were unique to either Treg or Tconv, a substantial proportion of TCRα (21–30%) and TCRβ (5–20%) sequences from tTreg were also found in Tconv. Analysis of a TCRβ Tg mouse line revealed an even higher proportion (71%) of TCRα sequences found in tTreg that were also found in Tconv. Interestingly, these shared TCRα clonotypes that were common to tTreg and Tconv were significantly more abundant than nonshared TCRα sequences of tTreg and Tconv. Finally, we used machine learning (ML) to develop an algorithm that was capable of distinguishing nonshared TCRα and TCRβ sequences expressed by tTreg from those of Tconv and, in addition, found that specific amino acid trimers were differentially expressed in either tTreg or Tconv. When we applied the same ML algorithm to an analysis of those TCR sequences that were shared by tTreg and Tconv, the vast majority of these sequences were classified as characteristic of Tconv and not tTreg. Taken together, our findings identify TCR sequence characteristics that bias to tTreg or Tconv fate, in addition to the presence of factors that can drive cells with an identical TCR sequence into either Tconv or tTreg lineages.

tTreg and Tconv were isolated from 4- to 8-wk-old mice of three strains, all on a C57BL/6 background: 1) Rag-GFP-Foxp3-RFP (25, 26): tTreg and Tconv were isolated from three individual mice based on GFP and Foxp3-RFP expression (Supplemental Fig. 1) and were sequenced independently; 2) Tcrd^CreERZsGreen-Foxp3-RFP: Tcrd^CreER knock-in with tamoxifen-induced ZsGreen reporter for TCRδ expression (three doses of 1 mg tamoxifen i.p. every other day, cells isolated 2 wk after last injection) (27): tTreg and Tconv cells were isolated from nine mice based on expressions of ZsGreen, Foxp3-RFP, and CD25 (Supplemental Fig. 1), and sequencing was performed on three pools (three mice pooled in one sample); and 3) Foxp3-GFP TCRα^+/− TCRβ-Tg mice carrying the DO11.10 TCRβ (28, 29): tTreg and Tconv samples from three individual mice were sorted based on Foxp3 reporter + and −, respectively, and individual samples were sequenced independently. Foxp3-GFP mice were provided by Vijay Kuchroo (30). All mice were maintained under specific pathogen free conditions at the animal facility of National Cancer Institute and Duke University. Animal procedures were reviewed and approved by National Institutes of Health or Duke Institutional Animal Care and Use Committee.

cDNA library construction was previously described (24). Briefly, total RNA was isolated from tTreg and Tconv using a QIAGEN RNeasy Micro Kit. Isolated total RNA (50–500 ng) was used for cDNA synthesis using TCRα and TCRβ C region–specific primers mTRAC1 (5′-GGCGTTGGTCTCTTTGAAG-3′) and mTRBC1 (5′-CACTTGTCCTCCTCTGAAAG-3′) (all oligos made by Eurofins USA), SMARTScribe Reverse Transcriptase (Takara Bio), and SmartN oligos (5′-AAGCAGUGGTAUCAACGCAGAGUNNNNUNNNNUNNNNUCTTrGrGrGrGp-3′) for template switching at the 5′ end to incorporate a UMI and M1SS sequence for PCR. The cDNA products were treated with uracil-DNA glycosylase (New England BioLabs) at 37°C for 30 min to remove SmartN oligos.

Three rounds of PCR using Super Fidelity Platinum Taq DNA Polymerase (Thermo Fisher Scientific) were performed to prepare libraries for sequencing. The first PCR (18–24 cycles) was used to enrich TCRs using M1SS (5′-AAGCAGTGGTATCAACGCA-3′, part of SmartN used as a 5′ PCR primer) and TCR C region primers (mTRAC2: 5′-CGGCACATTGATTTGGGAG-3′ and mTRBC2: 5′-TGTGGACCTCCTTGCCATTC-3′), and primers were removed by the QIAquick PCR Purification Kit (QIAGEN). The second PCR (20–32 cycles) was used to add an 8-bp sample barcode to each sample at 5′ end (P7M1S-n: 5′-CGTGTGCTCTTCCGATC(N)_1–2-NNNNNNNN(8 bp barcode)-CAGTGGTATCAACGCAGAG-3′) and internal C region primers (mTRAC3: 5′-AGGTTCTGGGTTCTGGATG-3′ and mTRBC3: 5′-GGTGGAGTCACATTTCTCAG-3′). PCR products were separated by 2% agarose (UltraPure; Thermo Fisher Scientific) gel electrophoresis, and DNA fragments (400–800 bp) were further purified by a QIAquick Gel Extraction Kit (QIAGEN). Purified DNAs of each sample were quantitated by a BioAnalyzer (Agilent Technologies) and combined for the third round of PCR (10 cycles), which incorporates the Illumina adaptor (P7: 5′-CAAGCAGAAGACGGCATACGAGATGTGACTGGAGTTCAGACGTGTGCTCTTCCGATC-3′, mP5TA: 5′-AATGATACGGCGACCACCGATCGTCGAGGTTCTGGGTTCTGGATG-3′ and mP5TB: 5′-AATGATACGGCGACCACCGATCGTCGGGTGGAGTCACATTTCTC-3′). Amplified DNA was separated by 2% agarose gel electrophoresis and further purified by QIAquick Gel Extraction and followed by the PCR purification kits. The amount of purified DNA was measured using a Qubit Fluorometer (Thermo Fisher Scientific). Fifty picomoles of DNA were used for sequencing on an Illumina HiSeq 2500 system. A modified paired end sequencing protocol was used: TCR-specific sequencing primers TRA (5′-TCGTCGAGGTTCTGGGTTCTGGATG-3′) and TRB (5′-TCGTCGGGTGGAGTCACATTTCTCAG-3′) were used for first round sequencing of 150 bps. Illumina RD2 primer (5′-GTGACTGGAGTTCAGACGTGTGCTCTTCCGATC-3′) was used for second round sequencing of 50 bp, covering the sample barcode and UMI.

The raw sequences were first separated by the sample barcode using a custom Python script. The sequences were then filtered through a custom Python script to separate erroneous sequences under each UMI cluster. Sequences that had a successfully identified sample barcode and were successfully filtered were then processed via the Molecular Identifier Groups-Based Error Correction program (23), which assembles consensus sequences for each UMI cluster and identifies the V, J, and CDR3 sequences via BLAST. Post–Molecular Identifier Groups-Based Error Correction analysis included the removal of contaminated TCR sequences using a custom Python script. A successful TCR had to meet the following criteria: each UMI ≥3 sequence reads, a known V and J gene, and an intact and functional CDR3 amino acid sequence. TCR clonotype was defined as a unique combination of V, J, and CDR3 amino acid sequences. CDR3 length was based on the definition of the International Immunogenetics Information System (31). V-J gene usage, CDR3 length distribution, the percentage of public TCRs, and the overlap of TCR clonotypes between different mice and subsets were calculated using custom Python scripts. For V and J usage analysis, the percentages were calculated based on total distinct TCR sequences and UMI counts. Calculation of overlapping TCR clonotypes relied on the definition of the TCR as V–CDR3–J sequence. We used the DivE (32) R program to estimate the species richness, and the geometric mean of the top models were presented.

Overlap of TCRα and TCRβ sequences between different mice/samples of the same strain of mice or between two strains of mice was analyzed at two levels: 1) overlapping distinct TCR sequences (Eq. 1) 2) and overlapping total sequences based on UMI counts (Eq. 2). The overlapping sequences are presented as percentages for each pair of comparison.

Eq. 1 was used for the calculation of the overlapping percentage of distinct TCR sequences, as follows:

S1 and S2 refer to two individual samples of the same strain of mice or two different strains of mice.

Eq. 2 was used for the calculation of the overlapping percentage of total TCR sequences, as follows:

The sum of UMI counts of TCRs found in both S1 and S2 is calculated as UMI counts of TCRs found in both S1 and S2 in S1 plus UMI counts of TCRs found in both S1 and S2 in S2.

To build an ML classifier to analyze αβ TCR sequences from tTreg and Tconv, we converted the TCR CDR3 sequences into a matrix of length 3 aa (3-mers or trimers) from all three strains of mice used in these studies. We first enumerated all possible 3-mers from all CDR3s (3-mer library) and then embedded each CDR3 into a length L vector that supposed there are L possible 3-mers observed in our data, in which each entry is the number of times that 3-mer appears in a CDR3; all other 3-mers in the library but not found in the analyzed TCR were labeled as 0. To calculate the relative starting location of each 3-mer, we recorded where the 3-mer appeared in CDR3 relative to the 3-mer starting locations that are possible for each CDR3 length (e.g., a length 12 CDR3 has 10 possible 3-mer locations). To determine multiple copies of the same 3-mer in a TCR, we used a Python dictionary to keep track of the number of occurrences of a 3-mer in a particular TCR. These are then combined into an N by L matrix, where N is the number of CDR3s. Furthermore, we vectorized the V/J information and concatenated this with the N × L 3-mer matrix. Suppose there are M possible V genes and K possible J genes; each CDR3’s V gene information is embedded into a length M vector, where the V gene entry that corresponds to the V gene in the TCR is labeled as 1 and the rest as 0 to produce an N × M matrix. The same procedure was done for the length K J gene vector to produce an N × K matrix. To generate the final matrix for ML, we concatenated the three matrices to produce a matrix of size N × (L + M + K) matrix. For two-class classification, we train a random forest binary classifier, which takes in a L vector and predicts the compartment to which a CDR3 belongs. We used 70% of the distinct TCR sequences as a training set and 30% of the sequences as a testing set. Training was performed using the random forest classifier from scikit-learn using 150 trees and default settings for other parameters (33). Model performance was evaluated by calculating the area under the receiver operating characteristic (ROC) using the scikit-learn metrics package. The ROC was plotted using Python’s matplotlib library (https://ieeexplore.ieee.org/document/4160265/). Scripts used for sequence processing and analysis can be made available upon request.

To compare trimer enrichment in the non-V/J portions of the CDR3, we first removed 3 aa from either ends of CDR3 and used the central CDR3 sequences for trimer analysis. Next, to compare the relative abundance of a particular amino acid trimer in tTreg and Tconv, we created a 2 × 2 contingency table for each k-mer and then computed the p value using a χ² test using the Python library SciPy (34). We corrected for multiple comparisons using the Benjamini-Hochberg procedure. The significantly enriched trimers in tTreg were defined as those that met the criteria tTreg/Tconv ratio ≥1.5 and false discovery rate (FDR) ≤ 0.05. The percentages of amino acid usage in the enriched trimers (Supplemental Table I) were calculated by the sum of each amino acid in the enriched trimers multiplied by their respective UMI counts divided by the total number of amino acids based on UMI counts in these enriched trimers. The percentages of amino acid usage in all trimers were calculated by sum of each amino acid in all the trimers multiplied by their respective UMI counts divided by the total number of amino acids based on UMI counts in all trimers unique to either tTreg or Tconv.

The Mann–Whitney U test was used to calculate the significant difference of UMI/TCR ratios between shared and nonshared TCRs with tTreg and Tconv. A p value <0.05 was considered significant.

The authors state that all data generated during this study are included in the article, its supplementary information file, and the Source Data file and are available from the corresponding author upon reasonable request. The raw TCR sequences data were deposited in the National Center for Biotechnology Information BioProject database with accession number 541952 (https://www.ncbi.nlm.nih.gov/bioproject/?term=541952).

To analyze αβ TCR repertoires of tTreg and Tconv, we isolated similar numbers of recently generated tTreg and Tconv from the thymus of Rag-GFP/Foxp3-RFP (25, 26) and Tcrd^CreERZsGreen-Foxp3-RFP mice (27) by cell sorting and applied a high-throughput sequencing method with UMI labeling of TCR mRNA. The Rag-GFP/Foxp3-RFP strain marked newly generated tTreg and Tconv in green fluorescent dye–GFP (25, 26), whereas the Tcrd^CreERZsGreen-Foxp3-RFP strain labeled newly produced tTreg and Tconv in fluorescent dye–ZsGreen after tamoxifen induction (27). These fluorescent markers were used in flow cytometric isolation of tTreg and Tconv to ensure the thymic origin of tTreg and Tconv by excluding the potential contamination of recirculating T cells (Supplemental Fig. 1A, 1B). We analyzed a total of 3.7 × 10⁵ tTreg and 3.0 × 10⁵ Tconv from three Rag-GFP/Foxp3-RFP mice and found that the estimated size of TCRα and TCRβ repertoires, identified as the number of distinct sequences, was comparably diverse between tTreg and Tconv when similar cell numbers were analyzed (Table I). We also isolated tTreg (6.0 × 10⁴ cells) and Tconv (7.2 × 10⁴ cells) from Tcrd^CreERZsGreen-Foxp3-RFP mice by cell sorting and determined their TCRα and TCRβ repertoires. Again, we found that the sizes of estimated TCRα and TCRβ repertoires were comparably diverse in similar numbers of Tconv and tTreg in Tcrd^CreERZsGreen-Foxp3-RFP mice (Table I). However, it should be noted that the total TCR repertoire size of Tconv in a mouse is likely larger than that of tTreg when the actual number of cells in the thymus is adjusted. Consistent with this is the observation that there were higher percentages of overlap in tTreg TCRs (3.4–5.2%) than in Tconv TCRs (1.8–2.1%) between individual mice (Rag-GFP/Foxp3-RFP) or different samples (Tcrd^CreERZsGreen-Foxp3-RFP) (Supplemental Fig. 2A, 2B). Because we analyzed TCRα and TCRβ repertoires separately, it could not be determined from this analysis whether the αβ combinatorial TCR repertoire was also comparable between tTreg and Tconv. V gene usage and CDR3 length distributions of TCRα and TCRβ were also not substantially different between tTreg and Tconv (Supplemental Fig. 3A–D).

Next, we determined the degree to which αβ TCR sequences of tTreg and Tconv of normal non-TCR Tg mice are similar or distinct. To overcome the limited numbers of cells available from each individual mouse, in particular for tTreg populations, we combined TCRα and TCRβ sequences of tTreg and Tconv from three samples of each strain of mice and then compared these pooled sequence sets (Fig. 1). In Rag-GFP/Foxp3-RFP mice, we found that 12% of distinct tTreg TCRα sequences (580 out of total 4906 TCRα sequences from tTreg) were found in Tconv and accounted for 14% of distinct TCRα Tconv sequences. Those shared TCRα sequences were more abundant than nonshared sequences as they accounted for 21 and 25% of total tTreg and Tconv sequences (based on UMI counts), respectively (Fig. 1A). Compared with TCRα, the overlap in TCRβ between tTreg and Tconv was slightly lower, accounting for 11 and 10% of distinct TCRβ in Treg and Tconv, respectively (Fig. 1C). These overlapping TCRβ sequences were also more abundant, accounting for 20 and 26% of total Treg and Tconv (Fig. 1C). TCRα and TCRβ sequences found in both tTreg and Tconv were also observed in Tcrd^CreERZsGreen-Foxp3-RFP mice (Fig. 1B, 1D). Collectively, these findings showed that ∼9–12% TCRα and 2–11% TCRβ sequences of tTreg were found in thymic Tconv, and they accounted for 21–30% of TCRα and 5–20% of TCRβ in total tTreg. With the caveat that these data derived from individual TCRα and TCRβ sequences do not directly measure αβ pairing, these findings indicate that TCR sequence is not the only factor determining tTreg generation in thymus.

To more directly analyze αβ combinatorial TCR expression in tTreg and Tconv, we analyzed TCRα sequences of tTreg and Tconv from TCRα^+/− TCRβ Tg-Foxp3-GFP mice (Table I). Each T cell from these mice expresses a single TCRα (because only one TCRα allele is expressed in these TCRα^+/− heterozygotes) in combination with the Tg TCRβ (the Tg TCRβ accounted for 99.99% of TCRβ sequences) so that the TCRα repertoire reflects overall αβ TCR clonotype expression. Again, we pooled TCRα sequences of three mice and compared tTreg and Tconv. Strikingly, we found that 42% (483 out of total of 1159) of distinct TCRα sequences of tTreg were shared with Tconv, accounting for 71% of total tTreg TCRα sequences (Fig. 2). Shared TCRα sequences represented only 2% of distinct TCRα of Tconv and accounted for 18% of Tconv sequences (Fig. 2). These results showed that approximately half of tTreg expressed αβTCR identical to those expressed by Tconv in this Tg mouse, indicating that the TCR sequence is not the sole determinant of Treg fate for this large proportion tTreg.

To further characterize the αβ TCR sequences shared between tTreg and Tconv, we analyzed their relative abundance by UMI counts of each TCR in two normal mouse strains as well as in a TCRβ Tg. We grouped each of the TCRα and TCRβ sequences into four groups: 1) found only in tTreg or 2) only in Tconv, 3) shared sequences expressed in Treg, and 4) shared sequences expressed in Tconv. We found that shared TCRα sequences were significantly more abundant (two to nine times) than those of distinct TCRα sequences in both tTreg and Tconv in two normal stains of mice as well as in the TCRα^+/− TCRβ Tg-Foxp3-GFP mice (Fig. 3A). This suggests that those shared TCRαβ−expressing Tconv and tTreg may either be derived from more abundant progenitor cells or might have undergone preferential expansion after differentiation. Shared TCRβ sequences were significantly more abundant (three to five times) in Tconv but not in tTreg in two normal strains of mice (Fig. 3B).

Last, if the tTreg and Tconv TCRα sequences of TCRα^+/− TCRβ Tg-Foxp3-GFP are selected based on their sequences, we would expect to see conservation of these TCRα sequences in other strains of mice. To address this, we first pooled sequences from three samples of each strains and then compared tTreg and Tconv TCRα sequences of TCRα^+/− TCRβ Tg-Foxp3-GFP with those from the non-Tg mice (Rag-GFP/Foxp3-RFP and Tcrd^CreERZsGreen-Foxp3-EGFP). Indeed, we found a small overlap in TCRα sequences between TCRα^+/− TCRβ Tg-Foxp3-GFP mice and the two non-Tg strains: 0.14 and 0.08% for distinct and total sequences (based on UMI counts) in tTreg-specific TCRα, 2.7 and 8.3% in Tconv-specific TCRα, and 1.3 and 7.6% in tTreg/Tconv–shared TCRα (Fig. 4).

The identification of αβ TCR sequences found only in tTreg or only in Tconv could reflect true differences in the TCR repertoires of these populations or might be, at least in part, a consequence of the inability to sequence to saturation all TCR in the large Tconv population. We therefore next asked whether the αβ TCR sequences that were not identified as shared between tTreg and Tconv have distinct features that distinguish tTreg and Tconv. We applied an ML algorithm to compare the nonshared αβ TCR (V gene-CDR3 amino acids-J gene) of tTreg and Tconv. We first generated a continuous 3-aa motif (trimer) library found in TCR CDR3 (35) and incorporated V and J gene information. We then used 70% of the combined distinct TCRα and TCRβ sequences of tTreg and Tconv as a training set and 30% as a testing set using a random forest classifier from scikit-learn (33). The model performance on the testing set was calculated using the area under the ROC. We found that nonshared TCRα of tTreg were distinguishable from TCRα of Tconv with ROC = 0.82 (Fig. 5A) and that TCRβ of tTreg were distinguishable from TCRβ of Tconv with ROC = 0.72 (Fig. 5B); thus, both TCRα and TCRβ nonshared sequences of tTreg were distinguished from their counterparts in Tconv.

To further determine the features of tTreg-restricted CDR3 amino acid sequences, we compared the abundance in tTreg and Tconv of specific trimers in the central region of CDR3, which mediates direct contact with Ag–MHC, excluding the N-terminal amino acid and the C-terminal 3 aa of CDR3 (24). We found that a number of trimers were significantly more abundant in tTreg than in Tconv (trimer ratio tTreg/Tconv ≥1.5, FDR < 0.05) in CDR3α (n = 49 found in 2.2% of total tTreg TCRα sequences) and CDR3β (n = 86 found in 1.2% of total tTreg TCRβ sequences) (Supplemental Table I), and the 20 most abundant trimers for CDR3α and CDR3β are presented in Fig. 5C, 5D. Strikingly, 2 aa present in these trimers were highly enriched in the abundant trimers of both CDR3α and CDR3β of tTreg: cysteine (enriched by 6.8- and 3.9-fold compared with the CDR3α and CDR3β of Tconv, respectively) and phenylalanine (enriched by 2.9- and 1.7-fold to the CDR3α and CDR3β of Tconv, respectively) (Fig. 5E), suggesting some common biophysical properties of tTreg TCRs. In addition, lysine was enriched by 2.3-fold in CDR3α, and methionine was enriched by 2.2-fold in CDR3β of tTreg. These findings indicate that TCRα and TCRβ in tTreg express a distribution of amino acids and trimer sequences that differs from those of Tconv.

To determine the degree to which TCRα and TCRβ sequences shared by Treg and Tconv resemble sequences that are found only in tTreg or only in Tconv, we applied the ML algorithm described above. Strikingly, we found that the great majority of TCRα and TCRβ sequences that are shared by Treg and Tconv were classified as Tconv in origin by this algorithm (Fig. 6A, 6B). Therefore, 82.4% of TCRα and 91.9% of TCRβ sequences were classified as Tconv/TCR–based on cutoffs selected from the receiver operator characteristic of the TCRα and TCRβ ML classifiers, respectively. Progenitors that express these shared TCR can thus differentiate to tTreg fate despite the expressions of TCR that are classified as more similar to Tconv.

The factors that determine the selection of Tconv or tTreg fate during thymic development are not completely understood. We designed studies to assess the degree to which TCR sequence determines this lineage choice. We conducted TCRα and TCRβ sequencing of tTreg and Tconv using a UMI-based method. Comparing TCRα and TCRβ sequences between tTreg and Tconv from normal strains of mice, we found that a substantial proportion of TCRα sequences and TCRβ sequences were shared between tTreg and Tconv; in TCRα^+/− TCRβ Tg-Foxp3-GFP mice, shared TCRα sequences or clonotypes were even more abundant. Notably, the TCRα and TCRβ sequences that were not shared between tTreg and Tconv were distinct in tTreg and Tconv as recognized by an ML algorithm and by identification of amino acid trimers more commonly used in CDR3α and CDR3β of tTreg than in Tconv. Finally, ML indicated that the great majority of TCRs that are shared by tTreg and Tconv have features in common with the sequences distinct to Tconv but not with sequences distinct to tTreg. Our TCR sequence analysis identified two populations of tTreg, one in which Treg fate is associated with the unique properties of the TCR and another with TCR properties that are characteristic of Tconv and for which tTreg fate may be therefore influenced by other factors than TCR. As previously described in an instructive model of intraclonal competition, the tTreg fate decision can be influenced by the abundance of tTreg precursors during thymic development (36). Further, the potential presence of precursor tTreg within the Tconv population could affect the analysis of the TCR contribution.

It has been reported that tTreg development can proceed through two progenitor cell pathways, with mature CD25⁺Foxp3⁺ tTreg being generated from either CD25⁻Foxp3^lo or CD25⁺ Foxp3⁻ precursors (37). These two progenitor tTreg produce functionally distinct mature tTreg. tTreg derived from CD25⁺ tTreg progenitors are able to prevent experimental autoimmune encephalitis, whereas tTreg from Foxp3^lo progenitor cells do not, and TCR of CD25⁺ tTreg progenitors have a higher affinity than those Foxp3^lo progenitor cells. Sequence analysis of the Vα2 TCR family in that study revealed both distinct and overlapped sequences between these two progenitor cells. Whether these reported features of the Vα2 TCR family generalize across the entire TCR repertoire and whether the tTreg-specific TCR and the tTreg/Tconv–shared TCR that we have identified in this study have distinct progenitor origins will require further study.

Our study provides a quantitative assessment of the degree of uniqueness and similarity of αβ TCR between tTreg and Tconv. The distinct tTreg and Tconv sequences that we identified are present at a small but consistent percentage among individual mice within the same strain and between different strains, suggesting sequence conservation of tTreg- or Tconv-specific TCR characteristics. Furthermore, these conserved sequences occupy a larger fraction of total TCRα and TCRβ sequences based on UMI counts than on distinct TCR sequences. Several factors could affect the estimated proportion of these shared TCR. Because the number of Tconv that could be used for sequencing was only a fraction of total Tconv in the thymus, the percentage of shared TCR sequences in tTreg may be an underestimate. In contrast, the existence of a potential precursor tTreg in the Tconv population could increase apparent sharing, although it should be noted that the Tconv (Foxp3-RFP⁻) analyzed for Rag-GFP/Foxp3-RFP thymocytes were only ∼1% CD25⁺ and that Tconv for other analyses were >99% Foxp3⁻ CD25⁻. The extensive sharing of αβ TCR between tTreg and Tconv is further confirmed through the analysis of TCRα sequences in a TCRβ Tg mouse, in which we found sharing of TCRα and therefore sharing of identical αβ TCR clonotypes between tTreg and Tconv in 71% of total Treg. This high degree of sharing may result from the substantially reduced size of the TCRα repertoire that can be paired with the single Tg TCRβ. Together, these findings suggest that TCR sequence alone does not determine the fate of tTreg or Tconv for a significant proportion of tTreg.

The UMI marking of TCR mRNA molecules allowed us to calculate the abundance of TCR sequences in tTreg and Tconv more accurately than prior approaches. This strategy led to the observation that those TCRα sequences in normal mice and αβ TCR clonotypes in TCRα^+/− TCRβ Tg-Foxp3-GFP mice that are shared between tTreg and Tconv are significantly more abundant than the nonshared TCRα sequences and αβ TCR clonotypes. It is not clear whether the abundance of those shared TCR occurs prior to or after the fate decision of tTreg and Tconv. A study using fluorescent tracking of thymic mature CD4⁺ and tTreg showed only one cell division postselection (38). This suggests that a high abundance of a given preselection progenitor may simply increase the probability that some progenitors of that clonotype will differentiate to tTreg and others to Tconv, reflecting a role of additional factors that drive Treg differentiation in concert with TCR signaling. However, TCRβ sequences shared by tTreg and Tconv were significantly more abundant than nonshared sequences expressed only in Tconv but not more abundant than nonshared sequences expressed only in tTreg. It remains to be determined whether there is selective expansion occurring after the fate decision in the thymus for Tconv that have shared TCR with tTreg.

ML has become an increasingly powerful tool in biological studies, in particular for those involving large datasets (39, 40). In this study, we applied ML to analyze TCRs by first partitioning each distinct TCR into its component V-CDR3 (multiple continuous tri–amino acids)-J factors and determining whether nonshared TCRα and TCRβ sequences are distinct between tTreg and Tconv. Indeed, a random forest classifier is able to distinguish those TCRα and TCRβ sequences of tTreg from those of Tconv with high accuracy. The identification of specific trimers that occur at a higher frequency in tTreg CDR3α and CDR3β reveals that cysteine and phenylalanine are enriched in both CDR3α and CDR3β of tTreg, suggesting common structural features of at least some tTreg TCR. This ML algorithm can be further improved when more TCR sequences of tTreg and Tconv are generated and can also be modified to test αβ paired TCRs as the potential for single cell sequencing advances. In addition, we observed that lysine is enriched in CDR3α trimers, whereas methionine is enriched in CDR3β trimers of tTreg. Although the roles of these enriched trimers and specific amino acids is currently unknown, cysteine is reported to be enriched in CDR3α and CDR3β in intraepithelial lymphocytes and type A intraepithelial lymphocytes precursors (41). Interestingly, cysteine and phenylalanine are reported to be less frequent in CDR3β of MHC-restricted than in MHC-independent TCR-expressing thymocytes (24), and phenylalanine was found in CDR3β of self-reactive TCR (42). Thus, it is possible that the enriched amino acids and enriched trimers in tTreg CDR3 may be involved in the interaction with self-antigen peptide rather than MHC during tTreg differentiation. Collectively, our findings demonstrate that TCR expressed by tTreg and not by Tconv have distinct CDR3α and CDR3β sequences compared with those of Tconv, supporting the critical role of TCR in the tTreg generation in thymus. Further studies will be needed to characterize the interaction of these TCRs with potentially selecting self-peptides to provide insights into the role of TCR in tTreg generation in the thymus.

Collectively, the results reported in this study have identified features that distinguish the αβ TCR sequences expressed by a significant proportion of tTreg from those expressed by Tconv in the thymus of normal mice and which therefore appear critical to determining differentiation into these lineages. The identification in tTreg of preferentially used trimers and selected amino acids in CDR3α and CDR3β provides molecular features for further understanding of TCR and Ag interaction in tTreg generation. For TCR clonotypes that are shared between tTreg and Tconv, it remains to be elucidated what are the non-TCR factors that drive the same TCR-carrying progenitors into either Treg or Tconv lineage. Candidates for such factors include differential costimulatory signaling, cytokine requirements, and other aspects of APC and thymic environment. The ML classification of tTreg and Tconv subsets by TCR sequence described in this study provides a strategy for dissecting the molecular pathways that mark these lineages at a single cell level. Combining ML and single cell analysis of αβ TCR sequence with transcriptome and other molecular parameters will allow better definition of the selection, function, and activation state of these T cell subsets.

We thank Wanjun Chen, Ethan Shevach, Angela Thornton, and Karen Hathcock for careful reviews and comments on this manuscript, Jeffrey Cifello for independently verifying some calculations used in the results, and Rong-Fong Shen and Wells Wu for help in TCR sequencing.

The authors have no financial conflicts of interest.