Research ArticleSYSTEMS IMMUNOLOGY

Precursors of human CD4+ cytotoxic T lymphocytes identified by single-cell transcriptome analysis

See allHide authors and affiliations

Science Immunology  19 Jan 2018:
Vol. 3, Issue 19, eaan8664
DOI: 10.1126/sciimmunol.aan8664
  • Fig. 1 Cytotoxicity-related transcripts are enriched in the CD4-TEMRA subset.

    (A) Schematic representation of the CD4+ T cell subsets and their surface markers used for the study. (B) Single-cell RNA-seq analysis showing row-wise z scores of normalized TPM (transcripts per million) counts of cells in each subset (indicated at the top) for each differentially expressed transcript (rows) obtained by pairwise comparison of TEMRA versus TEM, TEMRA versus TCM, and TEM versus TCM [single-cell differential expression (SCDE) and model-based analysis of single-cell transcriptomics (MAST) analysis, Benjamini-Hochberg adjusted P < 0.05 and ≥2-fold change). (C) Violin plots show the single-cell expression pattern of the indicated TEMRA-enriched transcripts (cytotoxicity-related and TFs) in the indicated subsets. The shapes represent the distribution of cells based on their log2(TPM + 1) expression values (y axis). The color scale represents the mean expression. (D) The contour plots show the surface expression of CD45RA, CX3CR1, GPR56, CD244, CD314, and KLRG1 and intracellular expression of GZMB and PRF1 in singlet-gated CD3+CD4+CCR7 cells. The numbers denote the percentage of cells in each quadrant. Bar graphs beneath show the average percentages. Error bars are mean ± SEM from 9 (GZMB), 10 (CX3CR1), and 11 (other indicated proteins) donors. **P < 0.001, ***P < 0.0005, and ****P < 0.0001 from Student’s paired two-tailed t test. (E) GSEA enrichment plots for the indicated gene sets in the transcriptome of TEMRA versus TEM and TCM (see Supplementary Materials and Methods). The top portion of the plot shows the running enrichment score (RES) for the gene set as the analysis walks down the ranked list of genes and reflects the degree to which the gene set is overrepresented at the top or bottom of the ranked list of genes. The middle portion of the plot shows where the members of the gene set (indicated as blue lines) appear in the ranked list of genes. The bottom portion of the plot shows the value of the ranking metric. NES, normalized enrichment score; FDR, false discovery rate. (F) IPA of canonical pathways enriched in the TEMRA subset; P values calculated by Fisher exact test (see Supplementary Materials and Methods). (G) Spearman correlation plot showing the coexpression of the 111 TEMRA-enriched transcripts. Dashed black box shows a cluster of transcripts showing high correlation; the list of these transcripts is shown in the text box on the right, and stars highlight the genes not previously reported to have cytotoxic function.

  • Fig. 2 CD4-TEMRA cells show marked clonal expansion.

    (A) Clonotype network graphs of cells from the TEMRA, TEM, and TCM subset and DENV-specific TEMRA and DENV-specific TEM cells from donor #6. Each circle represents a single cell; the reconstructed TCRα and TCRβ sequences for each cell are depicted as red- and blue-colored bars, respectively, inside each circle. Dark- and light-colored bars represent productive and nonproductive TCRs, respectively. Red connecting lines indicate shared TCRα sequences; blue lines indicate shared TCRβ sequences. Red arrow indicates a precursor cell. Purple dashed line indicates a precursor cell and an effector cell sharing the same unique TCRα and TCRβ chain clonotype. (B) Bar graph shows the percentage of cells with TCRα (left), TCRβ (middle), and TCRαβ (right) chain clonotype frequency ≥2 (clonally expanded cells) in TEMRA, TEM, and TCM subsets. Error bars are mean ± SEM from three donors with the same HLA allele. *P < 0.05 from Student’s paired two-tailed t test. (C) Single-cell RNA-seq analysis shows the row-wise z score of normalized TPM counts for the indicated TEMRA-enriched transcripts (SCDE and MAST analysis comparing TEMRA cells with clonotype frequency ≥2 or =1, Benjamini-Hochberg adjusted P < 0.05 and ≥2-fold change) in clonally expanded cells (unique TCRα and TCRβ chain clonotype frequency ≥2; magenta) and nonexpanded cells (unique TCRα and TCRβ chain clonotype frequency =1; green). The cells are arranged on the basis of the clonotype frequency, and the numbers on the top in black boxes represent clonotype frequency. (D) Violin plots show the single-cell expression pattern of the indicated cytotoxicity-related transcripts in clonally expanded versus nonexpanded TEMRA cells. The shapes represent the distribution of cells based on their log2(TPM + 1) expression values (y axis). The color scale represents the mean expression. (E) Bar graph shows the number of TEMRA cells with a unique TCRα and TCRβ chain clonotype for donor #6. Stars in red color indicate the clonotypes shared by DENV-specific TEMRA. Pie chart shows the percentage of expanded TEMRA cells that share the unique TCRα and TCRβ chain clonotypes of DENV-specific TEMRA cells (pink). (F) Bar graph shows the percentage of the expanded TEMRA cells that share TCRαβ clonotypes of DENV-specific TEMRA cells. Each dot represents a donor. Error bars are mean ± SEM from four donors.

  • Fig. 3 CD4-TEMRA cells are heterogeneous across donors.

    (A) Bar graph shows the percentage of clonally expanded TEMRA cells (cells having a unique TCRα and TCRβ chain clonotype with frequency ≥2) across 12 donors (left) and between donors classified on the basis of their status of previous DENV infection (DENV+ and DENV), geographical location [Sri Lanka, Nicaragua, and San Diego (the Americas)], and CD4-TEMRA proportion (TEMRAhigh, donors #1 to #6; TEMRAlow, donors #7 to #12). Each dot represents data from a single donor. Error bars are mean ± SEM. ***P = 0.001 from Student’s unpaired two-tailed t test; ns, not significant. (B) Correlation between the proportions of the CD4-TEMRA subset in CD4+ T cells and clonally expanded TEMRA cells (cells having a unique TCRα and TCRβ chain clonotype with frequency ≥2). r, Spearman correlation coefficient; P, by Student’s two-tailed paired t test. (C) Pie charts show the percentage of TEMRA cells with the first (clone 1) and second (clone 2) most frequent unique TCRα and TCTβ chain clonotype, as well as other clonally expanded cells (all others) and the rest of the TCRα and TCRβ chain clonotypes present in only one cell (unique clones). (D) Single-cell RNA-seq analysis of TEMRA (blue), TEM (yellow), and TCM (brown) cells shows the row-wise normalized mean TPM counts of 111 TEMRA-enriched transcripts [from Fig. 1 (B and C)] (rows) for each donor (column). The panel above the heat map identifies the cell type and donor number (#1 to #12, ordered on the basis of the frequency of the CD4-TEMRA subset). (E) Violin plots show the single-cell expression pattern of the indicated transcripts in the TEMRA subset across the 12 donors. The shapes represent the distribution of cells based on their log2(TPM + 1) expression values (y axis). The color scale represents the mean expression. (F) Matched data from the same donors showing the proportion of TEMRA cells in the CD4+ T cell subset (blue) and the proportion of CD244+ (red), GPR56+ (green), GZMB+ (pink), and PRF1+ (dark blue) TEMRA cells from 9 to 11 donors. Each dot represents a donor. **P < 0.001, ***P < 0.0005, and ****P < 0.0001 from Student’s paired two-tailed t test.

  • Fig. 4 Four distinct subsets of CD4-CTLs revealed by single-cell RNA-seq analysis.

    (A) t-Distributed stochastic neighbor embedding (tSNE) two-dimensional (2D) plot of single-cell RNA-seq data (n = 917 cells, based on 1000 most variable transcripts) shows four distinct clusters of cells in the TEMRA subset. (B) Bar graph shows the distribution of cells across four clusters in each indicated donor (donor numbers match those in table S1) classified on the basis of geographical location (the Americas and Sri Lanka, left) and the percentage of clonally expanded cells (cells having a unique TCRα and TCRβ chain clonotype with frequency ≥2) in each cluster (right). (C) Single-cell RNA-seq analysis shows row-wise z scores of mean-normalized TPM counts of cells in each cluster for each differentially expressed transcripts (rows) obtained by pairwise comparison of each cluster versus the rest of the clusters (Seurat analysis, Benjamini-Hochberg adjusted P < 0.05) for the four clusters. (D) tSNE 2D plots of single-cell RNA-seq data show the expression of the indicated differentially expressed transcripts (SCDE and MAST analysis, Benjamini-Hochberg adjusted P < 0.05 and ≥2-fold change). B2M is shown as a control housekeeping gene. (E) Violin plots show the single-cell expression pattern of the indicated differentially expressed transcripts across four TEMRA clusters along with TCM. The shapes represent the distribution of cells based on their log2(TPM + 1) expression values (y axis). The color scale represents the mean expression.

  • Fig. 5 CD4-CTL precursor cells in the TEMRA subset.

    (A) tSNE 2D plots of single-cell RNA-seq data show the expression of the indicated differentially expressed transcripts between the four clusters (SCDE and MAST analysis, adjusted P < 0.05 and ≥2-fold change). (B) Violin plots show the single-cell expression pattern of the indicated differentially expressed transcripts across four TEMRA clusters along with TCM. Color scale shows the log2(TPM + 1) expression values (A) and mean expression values (B). The shapes represent the distribution of cells based on their log2(TPM + 1) expression values (y axis) (B). (C) Venn diagram comparing the transcripts enriched in precursor clusters (clusters 3 and 4) versus effector clusters (clusters 1 and 2) and transcripts enriched in TCM versus TEMRA (left); transcripts enriched in precursor clusters (clusters 3 and 4) versus TCM and transcripts enriched in TEMRA versus TCM (right) (gene list obtained from SCDE and MAST analysis, Benjamini-Hochberg adjusted P < 0.05 and ≥2-fold change). Examples of overlapping genes are shown to the right. (D) Percentage of the CD4-TEMRA subset in the CD4+ T cell population across 104 samples from 89 individuals. (E) Contour plots show the surface expression of CD45RA and CCR7 (left) or CD127 (IL-7R) (right) in live and singlet-gated CD3+CD4+ T cells obtained from PBMCs (left). Contour plots show coexpression of CD45RA and CD127 in TEMRA and TCM subsets (right). The gating strategy shows CD127 and CD127high proportion. (F) Proportion of total CD4-TEMRA (left), CD127 CD4-TEMRA (middle), and CD127high CD4-TEMRA cells (right) in CD4+ T cells from 15 donors with two longitudinal samples (V1 and V2). Donor #16, marked in red color, shows marked increase in the proportion of CD4-TEMRA and CD127 CD4-TEMRA (effectors) cells from V1 to V2. (G) Contour plots show the surface expression of CD45RA and CCR7 (top) or CD127 (bottom) in live and singlet-gated CD3+CD4+ cells obtained from PBMCs of donor #16 at two time points (168 days apart) [V1 (left) and V2 (right)]. Bottom panels show the expression of CD127 in the CD4-TEMRA subset for the same donor. ns, not significant; from Student’s two-tailed paired t test. Numbers inside the contour plots show the proportion of the indicated cell type. (H) Bulk RNA-seq analysis of the IL-7Rhigh TEMRA subset (CD4-CTL precursors), IL-7R TEMRA (CD4-CTL effectors), and TCM cells shows the row-wise normalized TPM counts of top 200 (100 up-regulated and down-regulated, based on fold change) differentially expressed transcripts obtained by pairwise comparison of IL-7R TEMRA (CD4-CTL effectors) versus TCM from the five donors (DESeq2 analysis, Benjamini-Hochberg adjusted P < 0.05 and ≥2-fold change). The panel above the heat map identifies the cell type, donor, and visit. (I) Contour plots show the coexpression of CD244 or GPR56 with CD28 in TCM (red), IL-7Rhigh TEMRA (CD4-CTL precursors) (magenta), and IL-7R TEMRA (CD4-CTL effectors) (blue) in singlet-gated CD3+CD4+ T lymphocytes.

  • Fig. 6 CD4-CTL effectors share TCR clonotypes with CD4-CTL precursors.

    (A) Clonotype network graphs of single cells from TEMRA subsets from donors #12 (top) and #3 (bottom). Each circle represents a single cell; the reconstructed TCRα and TCRβ sequences for each cell are depicted as red- and blue-colored bars, respectively, inside each circle. Dark- and light-colored bars represent productive and nonproductive TCRs, respectively. Red connecting lines indicate shared TCRα sequences; blue lines indicate shared TCRβ sequences. Red arrow indicates a precursor cell. Purple dashed line indicates a cluster of cells with intermixed precursor and effector cells sharing the same clonotype. (B) Cell-state hierarchy map constructed for precursor cluster TEMRA (light green), effector cluster TEMRA (orange), TCM (brown), and TEM (yellow) cells from donor #12. The network shows connection (black line) between cells (circle). (C and D) TCR-seq assays in TCM, TEM, IL-7Rhigh TEMRA (CD4-CTL precursors), and the IL-7R TEMRA (CD4-CTL effectors) subset from 14 donors; bar graphs show the Shannon-Wiener diversity index obtained using V(D)J tools (C) and percentage of expanded TCRβ clonotypes (clonotype frequency ≥3) (D); error bars are mean ± SEM from 13 to 14 donors. *P < 0.05, **P < 0.005, and ***P < 0.001 from Student’s paired two-tailed t test; ns, not significant. (E) Pie charts (left) show the distribution of TCRβ clonotypes based on clonal frequency: the most clonally expanded (top expanded clone; blue) and the next most clonally expanded (second most expanded clone; orange), as well as the rest of the expanded (≥3) (all other expanded clones; yellow) and nonexpanded clonotypes (frequency ≤2; gray). The pie charts (middle) show the distribution of shared (overlapping) expanded effector TCRβ clonotypes (frequency ≥3) within the other indicated subsets. The graphs (right) show the percentage overlap of the most expanded TCRβ clonotype from effector subset with the other indicated subsets. (F) Pie charts (left) show the distribution of TCRβ clonotypes based on clonal frequency: the most clonally expanded (top expanded clone; blue) and the next most clonally expanded (second most expanded clone; orange), as well as the rest of the expanded (≥3) (all other expanded clones; yellow) and nonexpanded clonotypes (frequency ≤2; gray) in V2 samples from the longitudinal samples. The pie charts (middle) show the distribution of shared (overlapping) expanded V2 effector clonotypes (frequency ≥3) within other subsets from V1 (V2 samples). The graphs (right) show the percentage overlap of the most expanded TCRβ clonotype from the effector subset at V2 with other subsets at V1.

  • Fig. 7 Variable number of CD4-CTL precursors across donors.

    (A) tSNE 2D plot of single-cell RNA-seq data obtained using the Seurat software package for donors #2 (left) and #1 (right) from the 3′ transcriptome of single cells (RNA-seq carried out using the 10× genomics platform) for the most variable transcripts (1799, donor #2; 1569, donor #1). (B) tSNE 2D plots (left) and violin plots (right) of single-cell RNA-seq data show the expression [log2(UPM + 1) expression; UPM, UMI per million] of the indicated differentially expressed transcripts (SCDE and MAST analysis, Benjamini-Hochberg adjusted P < 0.05 and ≥2-fold change) for the indicated donors. Precursor cluster is marked with a black dashed line. The shapes in violin plots represent the distribution of cells based on their log2(UPM + 1) expression values (y axis). The color scale represents the expression of log2(UPM + 1) (tSNE 2D plots) and mean expression (violin plots). (C) Violin plots show the single-cell expression pattern of the indicated differentially expressed transcripts between IL-7Rhigh TEMRA (CD4-CTL precursors) and IL-7R TEMRA (CD4-CTL effectors) subsets. The shapes represent the distribution of cells based on their log2(TPM + 1) expression values (y axis). Color scale shows the mean expression values.

Supplementary Materials

  • immunology.sciencemag.org/cgi/content/full/3/19/eaan8664/DC1

    Materials and Methods

    Fig. S1. Most of the expanded TCR clonotypes are DENV-specific.

    Fig. S2. TEMRA cells cluster into four distinct clusters.

    Fig. S3. Batch analysis of TEMRA cells.

    Fig. S4. TEMRA clusters express distinct set of transcripts.

    Fig. S5. Clusters 1 and 2 are enriched for transcripts involved in cytotoxicity-related pathways.

    Fig. S6. Precursor cells in TEMRA subset.

    Fig. S7. CD4-CTL precursors share clonotypes with CD4-CTL effectors.

    Fig. S8. Smart-seq2 is more sensitive compared with droplet-based approach (10× genomics).

    Table S1. Summary of study donors.

    Table S2. List of differentially expressed genes between memory subsets [for data in Fig. 1 (B and C)].

    Table S3. Summary of TCRα and TCRβ chains recovered from full-length single-cell transcriptomes of all single cells in the 15 donors.

    Table S4. TCRα and TCRβ chain reconstruction from full-length single-cell transcriptome of every single cell sequenced from the 15 donors.

    Table S5. List of differentially expressed genes between clusters (for data in Fig. 4C and fig. S4).

    Table S6. List of differentially expressed genes between TCM and CD4-CTL effectors (for data in Fig. 5H).

    Table S7. TCRα and TCRβ chain sequences derived from TCR-seq analysis of 14 donors [for data in Fig. 6 (E and F) and fig. S7].

    Table S8. TCRα and TCRβ chain clonotype sharing derived from TCR-seq analysis of longitudinally collected samples from five donors (for data in Fig. 6F and fig. S7).

    Table S9. List of differentially expressed genes between the two clusters in donors #1 and #2 [for data in Fig. 7 (A and B)].

    Table S10. List of antibodies used in the study.

    Table S11. Primers used for TCR-seq.

    Table S12. List of gene sets used for GSEA (for data in Fig. 1E).

    References (6473)

  • Supplementary Materials

    Supplementary Material for:

    Precursors of human CD4+ cytotoxic T lymphocytes identified by single-cell transcriptome analysis

    Veena S. Patil, Ariel Madrigal, Benjamin J. Schmiedel, James Clarke, Patrick O?Rourke, Aruna D. de Silva, Eva Harris, Bjoern Peters, Gregory Seumois, Daniela Weiskopf, Alessandro Sette, Pandurangan Vijayanand*

    *Corresponding author. Email: vijay{at}lji.org

    Published 19 January 2018, Sci. Immunol. 3, eaan8664 (2017)
    DOI: 10.1126/sciimmunol.aan8664

    This PDF file includes:

    • Materials and Methods
    • Fig. S1. Most of the expanded TCR clonotypes are DENV-specific.
    • Fig. S2. TEMRA cells cluster into four distinct clusters.
    • Fig. S3. Batch analysis of TEMRA cells.
    • Fig. S4. TEMRA clusters express distinct set of transcripts.
    • Fig. S5. Clusters 1 and 2 are enriched for transcripts involved in cytotoxicity-related pathways.
    • Fig. S6. Precursor cells in TEMRA subset.
    • Fig. S7. CD4-CTL precursors share clonotypes with CD4-CTL effectors.
    • Fig. S8. Smart-seq2 is more sensitive compared with droplet-based approach (10× genomics).
    • Table S1. Summary of study donors.
    • Table S3. Summary of TCRα and TCRβ chains recovered from full-length single-cell transcriptomes of all single cells in the 15 donors.
    • Table S10. List of antibodies used in the study.
    • Table S11. Primers used for TCR-seq.
    • References (64—73)

    Download PDF

    Other Supplementary Material for this manuscript includes the following:

    • Table S2 (Microsoft Excel format). List of differentially expressed genes between memory subsets for data in Fig. 1 (B and C).
    • Table S4 (Microsoft Excel format). TCRα and TCRβ chain reconstruction from full-length single-cell transcriptome of every single cell sequenced from the 15 donors.
    • Table S5 (Microsoft Excel format). List of differentially expressed genes between clusters (for data in Fig. 4C and fig. S4).
    • Table S6 (Microsoft Excel format). List of differentially expressed genes between TCM and CD4-CTL effectors (for data in Fig. 5H).
    • Table S7 (Microsoft Excel format). TCRα and TCRβ chain sequences derived from TCR-seq analysis of 14 donors for data in Fig. 6 (E and F) and fig. S7.
    • Table S8 (Microsoft Excel format). TCRα and TCRβ chain clonotype sharing derived from TCR-seq analysis of longitudinally collected samples from five donors (for data in Fig. 6F and fig. S7).
    • Table S9 (Microsoft Excel format). List of differentially expressed genes between the two clusters in donors #1 and #2 for data in Fig. 7 (A and B).
    • Table S12 (Microsoft Excel format). List of gene sets used for GSEA (for data in Fig. 1E).
    • Source data (Microsoft Excel format)

    Files in this Data Supplement:

Navigate This Article