Improved T-cell Receptor Diversity Estimates Associate with Survival and Response to Anti–PD-1 Therapy

T-cell receptors (TCR) are heterodimers of alpha and beta (TCRαβ; ref. 1) or gamma and delta (TCRγδ; ref. 2) cell surface proteins that confer exquisite sensitivity and specificity of T cells for their cognate peptides presented by MHC molecules on antigen-presenting cells (Fig. 1; ref. 3). T-cell recognition of peptide/MHC complexes through binding interactions with TCR drives the T-cell activation, proliferation, differentiation, cytokine production, and cytolytic activity necessary for protection from pathogens (4). Between 2 × 10⁹ and 6 × 10⁹ unique TCRs may be present in the circulating T-cell population in healthy individuals (5). This vast number of unique clonotypes is produced stochastically by V(D)J recombination in thymic T-cell precursors, followed by positive and negative selection in the thymus to select T cells bearing TCRs that bind MHC but do not bind tightly to self-peptides (4). Thus, T-cell development produces a diverse repertoire of TCRs able to recognize pathogens while minimizing the risk of autoimmunity.

Measuring TCR diversity is important in understanding the fundamental immunology of aging, infectious disease, autoimmune disease, hereditary immunodeficiencies, and cancer. Having originated in the ecology literature, metrics to quantify diversity are functions of species richness (the number of unique individuals in a population) and/or species dominance (the relative abundance of each species; ref. 6). In the context of quantifying TCR repertoires via sequencing methods, “abundance” refers to the total number of sequencing reads, whereas the “species” is a unique TCR, also referred to as a “clonotype.” Quantifying the total number of unique clonotypes, i.e., richness, of an organism would require an impossibly large number of reads (7). Even capturing all TCR in the tumor is a subsampling of the total number of clonotypes that exist in that organism. Diversity, however, informs the degree to which clonotypes are equally represented. Given a sufficient abundance of reads, one could calculate the true diversity of a sample and perhaps glean insight into the diversity of the organism as a whole. Diversity is an important summary statistic for the TCR repertoire because it reflects (i) the potential breadth of antigen specificities, (ii) the degree of antigen-driven clonal expansion present in the underlying T-cell population, and (iii) the degree of access/localization that the T-cell repertoire has to the tissue itself. This clonal expansion is a key consequence of T-cell activation, pathogenesis in autoimmune disease, and protective immunity after vaccination or pathogen exposure. Thus, TCR repertoire diversity may influence and mark specific states of health and disease.

In oncology, multiple groups have reported associations between TCR diversity and clinical outcomes in patients with cancer. High TCR repertoire diversity associates with improved survival in multiple tumor types, response to cetuximab in head and neck cancer, and response to CTLA-4 inhibition in melanoma, hepatocellular carcinoma, and with GVAX (an allogeneic whole cell tumor vaccine) in pancreatic cancer (8–15). In contrast, low TCR repertoire diversity (high clonality) associates with clinical response to PD-1 axis inhibition in melanoma and urothelial carcinoma (16–20). These results highlight the complexity of TCR repertoire biology, along with the promise and challenge of integrating repertoire diversity in biomarker development.

Quantifying TCR species is not as straightforward as counting clearly defined organisms. The sequence region of maximal diversity, complementary determining region 3 (CDR3), responsible for conferring most of the specificity for epitope recognition (21), may include junctional deletion of germline nucleotides and nontemplated nucleotide insertions (Fig. 1; ref. 4). In practice, DNA or RNA is extracted from a sample of interest followed by targeted capture and/or amplification of TCR sequences to prepare next-generation sequencing libraries. Assembling and aligning the sequences thus require dedicated computational tools such as TRUST (22) or MiXCR (23) that must be evaluated. Regardless of the method used, TCR repertoire profiling is heavily biased by sampling depth of the T-cell population. Despite the wide array of methods for generating TCR repertoire profiling data and applications of TCR repertoire diversity in biology and medicine, performance characteristics of algorithms to infer TCR repertoires and their resulting diversity estimates have not been reported on ground-truth datasets. Therefore, optimizing TCR repertoire diversity measurements requires answering three questions: (i) which tools accurately identify and quantify TCR sequences, (ii) which diversity metric is the least sensitive to subsampling, and finally (iii) if no diversity metric is unaffected by subsampling, can a model be developed that provides accurate diversity measurements even with extreme subsampling?

Here, we reported an evaluation of published methods for inferring TCR repertoires from RNA sequencing (RNA-seq) data, a description of bias characteristics in diversity estimation caused by undersampling, and an algorithm derived from statistical learning to provide accurate and precise estimates of TCR repertoire diversity. Applying our method to correct TCR repertoire diversity estimates in data from The Cancer Genome Atlas (TCGA), we found more pronounced correlations between diversity and expected immunogenomics features. Applied to RNA-seq data from patients with melanoma treated with immune checkpoint inhibitors, we found otherwise unobserved survival and response associations with TCR repertoire diversity. The methods and analysis code that generated them were freely available with this publication. We expect the method reported here will be valuable in refining TCR repertoire diversity estimates from tumor RNA-seq data, for use in repertoire profiling studies, and TCR repertoire–based biomarker development going forward.

One thousand TCRαβ repertoires were generated using Synthetic TCR Informatics Generator (STIG; ref. 24). STIG was used to generate 100 million individuals (i.e., virtual cells) from 5,000 clonotypes across a range of diversities. Reads from this population were generated to create paired-end, 48 bp reads to match the TCGA protocol. Base quality scores were provided to STIG to match the same quality scores and error rates as TCGA data. Median insert length and SD were also selected from a range of those commonly found in a subsampling of TCGA samples (153–209 bp and 43–77 bp, respectively). The distribution of TCR clonotypes was set using STIG's χ² distribution with 2 degrees of freedom and a cutoff range from 6 to 8969. The range of modeled diversity falls within the range of 99% of TCGA samples. Matching the distribution of diversity and abundances was performed empirically and required several rounds of repertoire generation, alignment/assembly of repertoire reads, model generation, and diversity calculation. Thymic tissue (THYM) was excluded in the determination of tissue diversity as the diversity of THYM is much higher than others. Read abundances were also matched to the distribution of those found in TCGA data (see “TCGA analysis” methods section). Approximately 10% of STIG reads aligned to the CDR3 region with the remaining unaligned reads representing other regions of the TCR.

TCGA primary tumors were used. See figures for tissue types and counts. Samples were dropped according to the “DO_NOT_USE” column of GDC's TCGA “sample annotations” (see merged_sample_quality_annotations.tsv, https://gdc.cancer.gov/about-data/publications/panimmune).

Furthermore, upon recommendation of TCGA personnel: TCGA-06-0156-01A-02R-1849-01 was used for TCGA-06-0156-01, TCGA-06-0211-01B-01R-1849-01 was used for TCGA-06-0211-01, TCGA-21-1076-01A-01R-0692-07 was used for TCGA-21-1076-01, TCGA-23-1023-01A-02R-1564-13 was used for TCGA-23-1023. Finally, for the remaining duplicated samples, the last sample (ordered alphabetically) was used. The TCGA data types used were batch-corrected RNA-seq expression data (EBPlusPlusAdjustPANCAN_IlluminaHiSeq_RNASeqV2.geneExp.tsv, https://gdc.cancer.gov/about-data/publications/pancanatlas) and TCR diversity metrics (see “Diversity of subsampled populations” methods section) were calculated from MiXCR outputs (see “MiXCR and TRUST comparisons” methods section).

This initial subsampling was performed on the ground-truth counts generated by STIG prior to read generation. Chao1 (25) and Shannon entropy (26) were calculated using the R Package vegan. Evenness was calculated as Shannon entropy/Log (Richness; ref. 27). Evenness produces a comparable result to the “productive clonality” as reported by Adaptive Biotechnologies. Their “clonality” is 1 - evenness, although they use log base 2 in both their entropy and normalization (personal communication). d25, d50, and d75 (dXX) indices were calculated as described by VDJTools (28) GitHub repository: “The estimate equals to 1 - n/N, where n is the minimum number of clonotypes accounting for at least XX% of the total reads and code N is the total number of clonotypes in the sample.”

MiXCR, TRUST, and TRUST 3.0 Dockerfiles and example commands can be found here: https://github.com/Benjamin-Vincent-Lab. As TRUST does not align sequences, we used Bowtie2 (29) to align reads prior to using TRUST for these analyses. Morisita–Horn was calculated using the vegdist function in the vegan R package and shown as a similarity index (1 - Morisita–Horn dissimilarity index). Mismatch Index was calculated as the median number of mappings to a mismatch repertoire over the number of mappings to the correct repertoire. With this index, a number greater than one indicates the tool matched a CDR3 in the wrong repertoire more commonly than the correct repertoire. Basic local alignment (BLAST) searches were performed using BLAST+ (30) v2.2.30 with the STIG ground-truth repertoire sequences being used as the database for each sample.

Elastic net model generation (31) was performed using Monte Carlo cross-validation. One thousand simulated TCR repertoires were split into discovery (n = 667 repertoires) and validation (n = 333 repertoires) sets. The discovery set was repeatedly split into 3 folds with 2 folds for training and 1 fold for testing. Elastic net regression models were fit with repertoire features (abundance, Chao1, d25 index, evenness, inverse Simpson, richness, and Shannon entropy) used as predictor variables and ground-truth Shannon entropy as the response variable. Cross-validation was performed on the discovery set 100 times across 10 α values, with the selected model having the lowest root mean squared error. Samples were balanced by the number of events for both the discovery-validation and test-training set splits. For the calculation of modeled entropy, analyzed samples were required to have an abundance of 8 or higher for both α and β TCR chains.

TCGA data (32) and Riaz and colleagues' SKCM immune checkpoint data (33) can be accessed through their original publications. The TCR datasets for RNA and amplicon sequencing can be generated using STIG: https://github.com/Benjamin-Vincent-Lab/stig, see “TCR repertoire generation” and “Code availability” method sections for details.

Our tool for TCR repertoire simulation, STIG (24), was written in the python programming language, and our method for correcting diversity estimates was written in the R statistical programming language. Both are freely available, along with code necessary to reproduce our analysis, at https://github.com/Benjamin-Vincent-Lab. Code used to calculate modeled entropy may be obtained using the following R package: devtools::install_github(“Benjamin-Vincent-Lab/binfotron”, ref = “0.2–14”). Modeled entropy may be calculated using the functions: predict_true_entropy_from_counts or predict_true_entropy_from_diversity. Help documentation is provided in the package. Ground-truth TCR repertoires may be generated using the STIG software provided here: https://github.com/Benjamin-Vincent-Lab/stig.

Analyses were performed and visualized using R (34) v3.5.2 through RStudio (35) with the R packages ggplot2 v3.1.1 (36), ggsurvplot v0.4.3, ggrepel v0.8.2, ggseqlogo v0.1 (37), glmnet v2.0-18 (38), parallel, stats v3.5.2, survival v2.43-3 (39), survminer v0.4.3, tidyr v0.8.3, vegan v2.5-3, and viridis v0.5.1. FDR corrections were performed via the Benjamini–Hochberg method (40) using the R-package stats::p.adjust. P value thresholds were 0.05, 0.01, and 0.001 (*, **, and ***, respectively). “Significance” was calculated by -log₁₀ (P value).

TCR receptor repertoire diversity is typically calculated by assembling and counting the sequences that encode their CDR3 receptors (Fig. 1). It is becoming more and more common to use bulk RNA-seq data to assess this diversity (Fig. 1B), as RNA-seq data are more available than amplicon-based data for many applications. However, RNA-seq datasets often have far lower CDR3 counts, and thus more subsampling, than TCR amplicon datasets. We first sought to understand the extent to which this subsampling affected mathematical indices of diversity. Starting with a hypothetical population of 50 individuals, we observed the variance in diversity estimation using the metrics of Shannon entropy and evenness (1 – “productive clonality”), a commonly used metric for normalizing Shannon entropy. The variance increased as the subsampling became more severe for both methods (Fig. 2A). To evaluate if subsampling only affected diversity metrics at lower population sizes, we simulated 1,000 repertoires at varying degrees of subsampling from 2 to 100 million individuals and assessed the variation in diversity estimation using common diversity metrics. The ratios of ground-truth diversity index values (i.e., the diversity measurement taken at 100 million individuals) to subsampled values for each index are plotted as a function of abundance (Fig. 2B). The error in diversity was not limited to Shannon entropy and evenness, nor was it confined to those samples that were severely undersampled. As abundance increased, calculated diversity indices did eventually converge to accurate estimates of true diversity. However, this happened at abundances well above those of the relatively low coverage of TCR loci in RNA-seq experiments, including those found in TCGA data (32). Even evenness, a metric commonly used to correct for subsampling, failed to provide a consistent metric for diversity. Conversely, for counts over about 2⁷, evenness was further from the ground-truth value than using an unnormalized metric like Shannon entropy. Although it has been reported that evenness tends to overestimate diversity at lower abundances (41), we were surprised to see that this overestimation was much more pronounced in samples with lower initial diversity (Fig. 2C). Not only would this systematic overestimation of evenness obfuscate correlations with diversity, but in datasets where samples with higher true diversity also have higher abundance, this normalization could flip the direction of a correlation.

We next investigated the rigor of commonly used diversity metrics. Prior to comparing these measurements, we first needed to assess TCR repertoire inference methods for correct inference of individual TCR sequences. MiXCR (23) and TRUST (22) are the two primary tools used for assembling TCR clonotype sequences from paired-end reads. Because the only published comparison of these tools was done by the MiXCR group themselves (42), an independent confirmation was desired. TCR repertoires and reads were generated using Synthetic TCR Informatics Generator (STIG; ref. 24), which allowed us to know the true sequences and distributions of the TCR repertoires for each sample. Our findings indicated MiXCR outperformed TRUST methods with respect to sensitivity, specificity, and precision for both chain identity and CDR3 sequence (Supplementary Fig. S1). These findings were consistent regardless if either exact match or 95% sequence matching by BLAST was used (Q95; Supplementary Fig. S2). In addition, we assessed the performance of MiXCR amplicon and RNA-seq methods on STIG-simulated repertoires of amplicon reads (Supplementary Fig. S3). Amplicon-based sequencing produces higher quality and abundance reads by targeting the area of sequencing to one specific region of the genome, such as the V(D)J segments of TCR (43). TRUST methods could not be tested on these amplicon repertoires as the tool either failed or took over a week to run. Both MiXCR methods produced nearly identical results in all respects. As expected, due the higher quality input data, the performance of MiXCR on amplicon-simulated reads in determining V/J region and CDR3 sequence far exceeded that of the RNA-seq simulated dataset. We next determined the capacity of these TCR methods to support accurate estimation of true diversity (i.e., the known diversity of the STIG samples at the full 100 million individuals). Regardless of the TCR chain or diversity index used, no measurements supported accurate estimation of true index values from RNA-seq–simulated reads (Supplementary Fig. S4A–S4D). An ideal diversity index should be stable across a range of abundance and entropy values. Surprisingly, even at the higher abundances of amplicon-simulated reads, evenness still failed to track well with either the ground-truth evenness or Shannon entropy (Supplementary Fig. S4E and S4F). Only Shannon entropy was consistent across a wide range of entropy values and abundances for these amplicon-simulated datasets (Supplementary Fig. S4G and S4H, respectively).

Given the stability of Shannon entropy as a diversity index across a wide range of parameters, we sought to derive a model to predict the true TCR repertoire entropy when abundances are low (i.e., in conditions of undersampling such as tumor RNA-seq experiments). To accomplish this, we utilized elastic net regression (31) to model the ground-truth Shannon entropy. This statistical learning was chosen as it (i) includes parameters for penalizing and removing features that would needlessly add model complexity and variance, and (ii) can perform well in cases where predictor variables are co-correlated. Although direct calculation of Shannon entropy showed modest correlation with ground-truth values, albeit with high variance (Fig. 3A), elastic net models showed tighter correlations with less variance, whether models were trained and tested on the same chain (Fig. 3B), trained on one type of chain and tested on the other (Fig. 3C), or trained and tested in a chain agnostic fashion (Fig. 3D). All of these models favored evenness and the d25 index in their feature selection and beta coefficient weight (Supplementary Fig. S5). When abundance was sufficiently high, entropy was accurately estimated by direct calculation, and model-based correction was not needed (Fig. 2B; Supplementary Fig. S4E). To determine the exact conditions under which direct calculation is insufficient and model-based correction improves entropy estimation, we analyzed subsampled repertoires for accuracy of direct calculation of entropy based on the measured value of that entropy calculation (Fig. 3E). Here, we determined the abundance needed at a measured Shannon entropy to be within 5% of the ground-truth Shannon entropy. If a sample's measured Shannon entropy and abundance put the sample over this threshold (blue line, Fig. 3E), the measured Shannon entropy was taken prima facie, and no correction was made for our model entropy calculation (Supplementary Fig. S6).

Validating the model-based corrected diversity measurements on biological samples is a difficult task given that there is no ground-truth reference in this context. There are, however, multiple features believed to correlate with diversity that can be tested for association with diversity estimates. TCR diversity is generally accepted to inversely correlate with age, as the pool of TCRs available for selection is gradually lost over time (44). Similarly, diversity may inversely correlate with tumor mutations, as the number of mutations is associated with the number of neoantigens that can be targeted by T cells, driving successfully binding T cells to clonally expand and attack their targets (4). Therefore, we hypothesized that a more accurate diversity metric would have a stronger inverse correlation with age, SNV neoantigens, and InDel neoantigens. Analyzing data from TCGA project (primary tumors with TCR abundance over 8 reads and neoantigen counts available), the corrected diversity for TCRα and TCRβ chains showed coefficients closer to -1 with age and neoantigens than other diversity metrics with low variance and higher significance than any other metric (Fig. 4A; Supplementary Data S1). Conversely, we would expect diversity to be higher in the thymus than any other tissue type because thymus tumor purity is less than 100% and the thymus is the source of T-cell generation (4). Again, the corrected diversity metric outperformed the other diversity metrics in this regard with low variance, a coefficient close to 1 and higher significance than any other metric (Fig. 4A; Supplementary Data S1).

A robust estimate for diversity should also be more internally consistent, providing similar predictions across a range of TCR abundances. To test the internal consistency of diversity metrics, we randomly divided all TCGA tissues into three quantiles based on sample abundance and tested for association between overall survival and TCRβ diversity metrics using Cox proportional hazards (CoxPH) regression (45). A fourth group was made by sampling equally across the three abundance ranges. Precision of the predictions was estimated as the difference between the log₁₀ of the upper and log₁₀ of the lower confidence intervals for the HRs of each tissue (Fig. 4B and C). An internally consistent metric should provide a narrow width that is consistent across abundance ranges. Although confidence interval widths were not significantly different across abundance ranges for evenness and modeled entropy, Shannon entropy predictions were significantly different between abundance groups (Kruskal–Wallis, P = 0.0049). Of note, the variance of the survival predictions based on evenness was so large that the confidence interval widths needed to be placed on an axis 20x larger than that of Shannon entropy and modeled entropy (Fig. 4B; see Supplementary Fig. S7A for same axis comparison). Similarly, the variability in evenness resulted in large coefficients but low significance in correlating with other biological correlates of diversity (Fig. 4A). Together, these results demonstrated that the model-based diversity estimates showed improved correlations with true biological diversity and better internal consistency.

With model-based corrected diversity providing better correlations with biological diversity than other metrics, we then asked: Does this improved metric of diversity improve prediction of survival of patients with cancer? We note that even a perfect measurement of diversity may not predict survival in cancer. That said, we hypothesized that a decrease in diversity should prognosticate better survival as it would indicate the immune system is mounting an effective response against the tumor (16–18). To our surprise, in many TCGA tissue types and in all TCGA samples as a whole, higher corrected diversity associated with increased survival. In addition, we tested the association of pretreatment tumor-associated TCR repertoire diversity with survival in patients with melanoma treated with nivolumab, a monoclonal antibody that inhibits programmed cell death receptor 1 (PD-1). We restricted our analysis to the Riaz and colleagues' dataset (33) as other publicly available RNA-seq datasets associated with immune checkpoint inhibition trials were generated using formalin-fixed paraffin-embedded (FFPE)–derived RNA. Repertoire inference from RNA-seq data is untrustworthy when degraded RNA derived from FFPE tissue is used as template for sequencing library preparation (46). Increased TCR repertoire diversity was associated with improved survival in patients with melanoma treated with nivolumab (Fig. 5A; Supplementary Data S2). The magnitude of the effect was at least that of survival association with diversity in TCGA melanoma samples. Only the significance of TCRα for all TCGA samples was revealed if evenness was used as the diversity normalization (Supplementary Fig. S7B). A similar result was seen looking at the Kaplan–Meier plots for TCGA-treated samples (Fig. 5B; Supplementary Data S2) and the Riaz and colleagues' immune checkpoint–treated samples (Fig. 5C). This effect was largest when both alpha and beta chains had high modeled entropy. When comparing TCR diversity metrics according to patient response, responders (mResist score of partial or complete response) showed significantly higher diversity when looking at modeled entropy (Fig. 5D; TCRα P = 0.0083, TCRβ P = 0.0168, two-sided Mann–Whitney test). Evenness showed no significance (TCRα P = 0.3702, TCRβ P = 0.6273). Nonprogressors (mResist score of partial response, complete response, or stable disease) showed a similar although less significant difference as responders, with nonprogressors showing higher TCRα-modeled entropy than progressors (Fig. 5E; P = 0.0176). The difference when comparing evenness was not significant (P = 0.0985). As mentioned previously, there was a statistically significant association of TCRα-modeled entropy with overall survival in the Riaz and colleagues' immune checkpoint–treated samples. This difference was not observed in TCGA SKCM dataset, in spite of it having a larger sample size (283 TCGA SKCM samples vs. 33 for Riaz and colleagues; Fig. 5A). To explore this further, we bootstrapped TCGA and Riaz SKCM samples and quantified the percentage of samplings that indicated high model predictions of TCRα entropy predict survival (Fig. 5F). In addition, we calculated the median HRs for these samplings (Fig. 5G). Both of these results indicated that fewer immune checkpoint inhibitor–treated samples were needed to obtain the same predictive power as untreated samples.

Taken together, these results indicate our modeled diversity metric improved discrimination of responses in TCGA dataset and in patients with melanoma treated with PD-1 inhibition.

TCR repertoire diversity reflects the degree of clonal expansion in a T-cell population and thus is a valuable biological readout as well as a potential biomarker of health or disease. Here, we showed that accurate measurement of TCR repertoire diversity required correction for undersampling that may occur as a result of biological limitations, such as testing only a small portion of a blood or tissue-resident T-cell population of interest, or technical limitations such as low sequencing read coverage for a given sample. Optimal estimation of true population diversity also required use of the most accurate algorithm for TCR repertoire inference. Based on our results, we recommend that MiXCR be used due to improved sensitivity, precision, and sequence length independence regardless of sequence identity threshold used to call individual TCR clonotypes.

Despite the good performance of MiXCR in reconstructing TCR repertoires from RNA-seq reads, existing methods of calculating diversity from these repertoire counts were inaccurate under typical RNA-seq conditions. We developed a method for correcting diversity estimates that can be applied to any dataset using software supplied with this article. Of note, models built from TCRα data performed well when applied to TCRβ datasets and vice versa. Even “chain agnostic” models trained on subsampled population counts performed well on TCR data, suggesting that our method may be generally extensible to entropy estimation in the context of undersampling beyond TCR repertoire data as long as the expected diversity of the sampled population is within the range provided to the model. In datasets where sample entropy estimates fall outside the range of entropies upon which the model was trained (0.98 to 8.11), we recommend training a new entropy correction model on the desired true Shannon entropy range.

Although our method of modeling the true Shannon entropy from RNA-seq data should effectively correct entropy measurements with different abundances, it should only be used to correct entropy within a batch of samples. Sequencing methods that are better or worse at detecting low abundance sequences in low-quality reads could produce very different values. For example, amplicon data could be very sensitive to detecting low abundance sequences and would therefore produce much higher diversity measurements. Given the differences in sensitivity across sample preparations and sequencing methods, it would be very difficult to make a generalizable tool that corrects across sequencing methods. For these high-quality high-abundance sequencing methods, we recommend using Shannon entropy as this metric approaches the ground-truth Shannon entropy when the number of reads is in the tens of thousands range or higher. The use of evenness and “productive clonality” (1-evenness) is strongly discouraged.

We applied our method of improving the accuracy of entropy estimation to TCR repertoires derived from RNA-seq datasets in TCGA and a study of pretreatment samples from patients with melanoma undergoing treatment with the PD-1 inhibitor nivolumab. In the Riaz and colleagues' melanoma study, we found an association between diversity and survival where the magnitude of effect size was at least as great as the survival associations in TCGA (i.e., independent of immunotherapy treatment as well as in the small number of patients with melanoma in TCGA that were treated with immunotherapy after sample collection). Given our analysis of power by bootstrapping, we believe the overlapping effect sizes were a function of far greater numbers in the nontreated TCGA dataset; however, we cannot rule out that TCR repertoire diversity was prognostic rather than predictive in the case of nivolumab treatment. Future work on additional datasets will be required to answer this question.

Our use of modeled entropy indicated higher TCR diversity was associated with better survival. This observation was in contrast to our expectation. We previously hypothesized that better outcomes would correlate with reduced diversity, indicating a clonally restricted T-cell population primed to respond to tumor antigens. There is evidence that lower pretreatment diversity is associated with response to PD-1 axis inhibition (17–20); however, other studies have found associations between higher TCR diversity and more favorable clinical outcomes (8–15). It is possible that the association of low diversity and better outcomes is unique to mechanisms of action that reverse T-cell exhaustion, related to necessity of a clonally expanded T-cell population, reactive against tumor antigens but exhausted, to be present in the tumor prior to treatment. In other therapeutic contexts (i.e., those not driven by the mechanism of reversing T-cell exhaustion), the necessity for a diverse TCR repertoire to recognize tumor antigens may be more important. In addition, there is emerging evidence that immune checkpoint inhibition efficacy is driven by clonal expansion that happens after initiation of therapy (15, 47, 48). Thus, we expect that more accurate measurement of TCR repertoire diversity, as well as measurements of repertoire diversity while on therapy, will provide better biomarkers of response. We also expect that the magnitude and direction of associations between TCR repertoire diversity and clinical outcomes will vary by therapeutic mechanism of action.

If, as our results suggest, clonal expansion occurs after administering immune checkpoint blockade, then monitoring the changes in diversity upon receiving ICP therapy could be an early indicator of how well a patient is responding to immunotherapy. Our method of accurately estimating diversity even with low counts may make such diagnostics readily affordable and aid in numerous fields that must rely upon assessing diversity from subsampled populations.

D.S. Bortone reports personal fees from GeneCentric Therapeutics Inc. (for assay creation and biomarker development) outside the submitted work. B.G. Vincent reports other from GeneCentric Therapeutics (equity) outside the submitted work. No disclosures were reported by the other authors.

D.S. Bortone: Conceptualization, data curation, software, formal analysis, validation, investigation, visualization, methodology, writing–review and editing. M.G. Woodcock: Software, writing–review and editing. J.S. Parker: Conceptualization, supervision, writing–review and editing. B.G. Vincent: Conceptualization, resources, supervision, funding acquisition, methodology, writing–original draft, project administration, writing–review and editing.

The authors thank members of the Vincent Lab in the UNC Lineberger Comprehensive Cancer Center for helpful discussion of the article.

Funding was provided by the UNC University Cancer Research Fund (to B.G. Vincent and J.S. Parker), the UNC-Duke Immunotherapy Training Program (T32-CA211056 to M.G. Woodcock), the Susan G. Komen Foundation (Career Catalyst Research Award to B.G. Vincent), and the V Foundation for Cancer Research (to B.G. Vincent).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.