Identification of T Cell-Restricted Genes, and Signatures for Different T Cell Responses, Using a Comprehensive Collection of Microarray Datasets

Generation of T cell and leukocyte subsets

The microarray experiments used in this study are summarized in Table I⇓. Detailed descriptions of each microarray experiment are provided at 〈http://linkage.garvan.unsw.edu.au/public/microarrays/〉. Generation of Th1 and Th2, T_CM, T_EM, and T_FH cells has been described previously (15). Resting and activated Th0 and Th2 cells derived from peripheral blood naive T cells were generated as follows. Adult CD4⁺ naive T cells were isolated from PBMC using a combination of MACS and FACS sorting. Briefly, CD4⁺ T cells were positively selected with anti-CD4-coated MACS microbeads (Miltenyi Biotec). Cells were then stained with anti-CD45RA Abs, and naive CD45RA⁺ T cells were sorted on a FACSVantage (BD Biosciences). The cells were cultured in RPMI 1640 medium supplemented with 1% Glutamax, 1% sodium pyruvate, 1% nonessential amino acids, 50 μg/ml streptomycin/penicillin (Invitrogen Life Technologies), 5 × 10⁻⁵ M 2-ME, and 5% human serum. Single naive T cells were distributed in 96-well plates by FACSVantage sorting and were stimulated with 10⁵ allogeneic irradiated (40 Gy) PBMC and 1 μg/ml PHA in IL-2-containing medium in the presence of 2 ng/ml rIL-4 and 1 μg/ml neutralizing anti-IL-12 Abs (R&D Systems). Clones were expanded in IL-2-containing medium. A cytokine production assay was performed after ∼2 wk. After 8 wk, selected T cell clones were restimulated under the same Th2 conditions to generate Th2 clones or under opposite Th1 conditions (0.5 ng/ml rIL-12 and 0.5 μg/ml neutralizing anti-IL-4 Abs (R&D Systems)) to generate Th0 clones. Clones were expanded in IL-2-containing medium. A cytokine production assay was performed after 2 wk. RNA was extracted from resting and activated T cell clones after 4 wk using TRIzol (Invitrogen Life Technologies).

For cytokine production assays, T cells were stimulated with 2 × 10⁻⁷ M PMA and 1 μg/ml ionomycin (Sigma-Aldrich) for 4 h. Brefeldin A (10 μg/ml) was added for the last 2 h. Fluorochrome-labeled anti-IFN-γ and anti-IL-4 Abs (BD Pharmingen) were used after fixation and permeabilization, performed using the Cytofix/Cytoperm kit (BD Pharmingen). T cells (10⁷/condition) were left untreated or were activated for 4 h with 1 μg/ml anti-CD3 Ab (clone TR66) and 50 ng/ml phorbol dibutyrate (PdBu; Sigma-Aldrich).

Preparation of cRNA and genechip hybridizations

Total RNA was isolated from cells using the RNeasy Total RNA Isolation Kit (Qiagen) or TRIzol according to the manufacturer’s instructions. cRNA was prepared as previously described (16). Briefly, cDNA was specifically transcribed from 500 ng of mRNA using a poly-T nucleotide primer containing a T7 RNA polymerase promoter (Geneworks). Biotinylated, antisense target cRNA was subsequently synthesized by in vitro transcription using the BioArray High Yield RNA Transcript Labeling kit (Enzo Diagnostics). Fifteen micrograms of biotin-labeled target cRNA was then fragmented and used to prepare a hybridization mixture, which included probe array controls and blocking agents. Hybridization to U133A and B Affymetrix arrays was conducted for 16 h at 45°C and 60 rpm. After hybridization, washing and staining of the hybridized probe array were performed by an automated fluidics station, according to the manufacturer’s protocols. The stained probe array was scanned using the Agilent GeneArray Laser Scanner, and the resultant image was captured as a data image file.

Data analysis

Absolute levels of expression of genes were determined and scaled to 150 using algorithms in MicroArray Analysis suite 5.0 (MAS5) software (Affymetrix). The signal value represents the level of expression of a transcript. The signal log ratio is the change in expression level of a transcript expressed as the log₂ ratio (a signal log₂ ratio of 1 is equal to a fold change of 2). Annotations were extracted using the NetAffx analysis tool (〈www.affymetrix.com〉). All the microarray data presented in this study are publicly available from 〈http://linkage.garvan.unsw.edu.au/public/microarrays/〉.

In-depth analyses and clustering of data were conducted using GeneSpring software (Silicon Genetics). After data transformation (to convert any negative value to 0.01), normalization was performed using a per chip 50th percentile normalization and per gene median normalization method. Genes that were consistently absent or below the noise level were excluded from analysis. To identify genes with statistically significant differences between T cell subsets or cell types, one-way ANOVA with a p value cutoff of 0.05 and the Benjamini and Hochberg false discovery rate as multiple testing correction were performed. The Student-Newman-Keuls post-hoc test was used to identify the specific groups in which significant differential expression occurred. Hierarchical clustering was performed on both genes and individual experiments, with Pearson correlation as a measure of similarity to group genes and samples with similar expression patterns. Data points were arranged in a hierarchy and were displayed in a phylogenetic tree of clusters of genes in a hierarchically ordered relationship. Branch lengths represent the degree of similarity between sets. Gene expression profiles that were similar across the experimental samples were clustered together. Principal component analysis using GeneSpring software was preformed using all genes that were expressed in at least two samples, to identify the components responsible for the greatest variability.

Data were imported as a Microsoft Excel file into Spotfire for graphical representation of gene expression patterns in T_CM and T_EM cells. Spotfire software was used to map genechip analysis results where fluorescence intensity detected on the genechips is represented by a color scale. Genes that showed a change of 2-fold or greater were considered differentially expressed.

Real-time PCR validation of gene expression

Total RNA was isolated from purified cells as described. RT was performed on 200 ng of total RNA per reaction using the Reverse-IT RTase Blend Kit (ABgene) according to the manufacturer’s instructions. After cDNA synthesis, semiquantitative real-time RT-PCR was performed with FastStart DNA Master SYBR Green I Reagent by LightCycler (Roche), as previously described (17). Primers were designed using Primer3 software (18). mRNA expression was normalized to GAPDH, and the fold change was calculated relative to Th1 gene expression.

Comprehensive leukocyte subset profiling allows the identity of T cell-restricted genes

T cells differ from other leukocyte subsets in their phenotypic and functional properties, which is reflected by different gene expression patterns. We have produced a large dataset of genechip expression profiles for all the major leukocyte subsets, including various subsets of B cells, plasma cells, NK cells, eosinophils, neutrophils, basophils, dendritic cells, macrophages, mast cells, and all major effector and memory T cell subsets (Table I⇑).

The first question we asked was whether this comprehensive dataset of genechips for various leukocyte types might allow us to identify T cell-specific transcripts. We applied one-way ANOVA to identify genes that most reliably discriminate between T cells and other leukocytes (Fig. 1⇓A). We also used the publicly available SymAtlas database (〈http://symatlas.gnf.org/SymAtlas/〉) (19) to assess the expression of gene transcripts in many nonlymphoid tissues. We identified a large number of genes that were preferentially expressed by T cells compared with non-T cell leukocytes (>1200 probe sets on the Affymetrix U133 A+B arrays; p < 0.0001). The top ∼100 genes most intimately associated with T cells are listed (Fig. 1⇓A), and ∼400 additional genes selectively expressed by T cells (p < 0.000001) are provided as supplementary information (supplementary Table I). As expected, some T cell-restricted genes were expressed outside of the hemopoietic lineage, whereas others were highly specific to T cells (Fig. 1⇓A). Interestingly, we found that many genes were shared between T and NK cells, indicating that the origin and functional properties of these two leukocyte types are closely related. Importantly, many of the genes identified through this approach (Fig. 1⇓A) were well-characterized, T cell-specific genes, such as CD3 chains, CD28, ICOS, CD40L, TCR, Lck, Zap70, and GATA3, which confirmed the validity of our approach. This gave greater credence and reliability to the numerous additional genes that have not yet been necessarily associated with T cell-selective roles. We selected several genes and confirmed their T cell-specific expression using real-time PCR (Fig. 1⇓B). These genes included IL-32 (previously known as NK4), a recently identified inflammatory cytokine (20); BCL11b, a transcriptional repressor that controls T cell development (21, 22); EOMES, a paralogue of T-bet important during embryogenesis (23), but also able to regulate CD8⁺ effector function (24); and UBASH3A (also known as suppressor of TCR signaling-2), which is a regulator of signaling pathway involved in T cell activation. We also confirmed T cell-specific expression of two GTPase members of the immunity-associated protein (IMAP) (or immunity-associated nucleotide (IAN)) family GTPase, IMAP family member (GIMAP)5 (hIAN5 and IAN4L) and GIMAP7 (hINA7). Interestingly, we also noted preferential expression of GIMAP6 (hIAN2) by T cells (supplementary Table I). This family is relatively poorly characterized to date; however, members of this family seem to be highly expressed in the spleen and lymph nodes, and there is increasing evidence of their importance in T cell development and maturation (25, 26, 27). Fig. 1⇓A shows that in addition to genes specifically expressed by T cells, there were also many genes that were specifically absent from T cells compared with other leukocyte subsets. Not surprisingly, these included several genes associated with Ag processing and presentation, endocytosis, complement receptors, and numerous other genes. It is beyond the scope of this study to discuss all the novel genes identified as T cell specific; however, the important point of this work is that bioinformatic strategies and large gene expression datasets provide a novel way of identifying T cell-biased molecules.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 1.

Identification of T cell-selective gene expression using comprehensive gene microarray profiling. A, cRNA was generated from numerous leukocyte types (Table I⇑), and microarray gene expression profiles were generated using Affymetrix U133A+B genechips and analyzed using GeneSpring software. One-way ANOVA was used to identify the genes that were significantly different among the numerous T cell subsets and non-T cell leukocyte subsets. Approximately 100 genes most significantly (lowest p value) preferentially expressed in T cells are listed. SymAtlas (19 ) was used to assess the expression of these genes in nonhemopoietic tissues, and a summary is provided in the right column. B, Differential expression of genes identified from microarray analysis was validated using semiquantitative real-time PCR. RNA was isolated from purified populations of different human cell types, and relative mRNA expression of selected genes was determined. Lane 1, Th1 cells; lane 2, Th2 cells; lane 3, blood αβ T cells; lanes 4 and 5, eosinophils; lanes 6 and 7, neutrophils; lane 8, mast cells; lane 9, monocytes; lane 10, macrophages; lane 11, immature DC. EOMES, eomesodermin homologue; BCL11B, B cell CLL/lymphoma 11B; UBASH3A, ubiquitin associated and SH3 containing. GIMAP5, GTPase IMAP family member 5; GIMAP7, GTPase IMAP family member 7. C, Genes involved in TCR complex signaling and T cell costimulation were selected from the genes preferentially expressed by T cells and mapped using the software application GenMapp2. Genes that were preferentially expressed by T cells compared with all other leukocytes are shown in yellow boxes. Downward blue arrows signify that the gene was down-regulated after activation; upward red arrows indicate up-regulation. An asterisk next to the gene name indicates that multiple probes showed variable expression patterns. Gene expression was normalized at ∼1 for all experiments, and the scale shown in this figure is the same for the following figures. The color scale indicates normalized expression levels, a value of 1 indicates normal expression, 0 indicates low expression, and values >1 indicate high expression.

One of the processes unique to T cells is TCR signaling. We found >150 TCR signaling or costimulation related genes that were significantly preferentially expressed in T cells (supplementary Table II). Selected genes involved in the development of the immunological synapse and downstream signaling in T cells are shown in Fig. 1⇑C. The majority of genes involved in T cell signaling cascades were restricted to T cells, although some were expressed by other leukocytes. Interestingly, many of the T cell signaling molecules were expressed at a higher level in resting, rather than activated, T cells, which was noted in previous studies (28, 29) and could be linked to down-modulation of the TCR signaling apparatus after stimulation (30, 31).

Gene expression signatures of T cell subsets

After identification of all the genes that distinguish T cells from other leukocytes, we focused on the differences between the T cell subsets. Our dataset of genechip profiles contains a diverse array of T cell subsets, including resting and activated Th0 and Th2 cells; naive, central, and effector memory CD4⁺ T cells; CD8⁺ T cells; and αβ and γδ T cells (Table I⇑). A one-way ANOVA with Student-Newman-Keuls post-hoc testing was used to identify the genes that distinguish particular subsets of T cells from all others. NK cells were also included in this study, because we noted previously that T cells share many genes with NK cells. Unsupervised clustering of these genes (Fig. 2⇓A) revealed easily identifiable and distinct gene expression signatures that distinguished each subset. For instance, the T_FH-specific signature consisted of well-known markers of this subset, such as CXCR5 and CXCL13, and many genes not previously associated with this subset (2, 3). Of note, gene expression signatures were not only characterized by highly preferentially expressed genes, but also by the genes whose expression was particularly low or absent in a particular subset, compared with others. This was especially true among closely related subsets, where few genes were expressed exclusively by a single subset. For instance, we identified few genes that were specific to T_EM or T_CM cells, and these memory T cell populations appeared to be highly related based on the limited differences in gene expression between the two. We have listed select genes from signature clusters, including genes up- and down-regulated in T_FH cells (Fig. 2⇓A, i and vi), down-regulated in γδ T cells (Fig. 2⇓A, ii), and down-regulated in T_CM (Fig. 2⇓A, v) and T_EM (Fig. 2⇓A, vii) and others. Genes contributing to Th1 and Th2 signatures have been discussed previously (15). Importantly, gene expression profiling of the numerous T cell subsets allowed us to distinguish clear and discriminating signatures based on the presence and the absence of certain transcripts.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 2.

Distinct gene expression signatures characterize effector/memory T cell subsets. Gene expression profiles of different T cell subsets were analyzed using one-way ANOVA with Benjamini multiple testing correction and Student-Newman-Keuls post-hoc tests to identify genes with significant changes in expression among the subsets. A, The union of all genes that distinguish the various cell types was created by combining genes that differentiate a single subset from all others. Unsupervised clustering using Pearson similarity measure was used to arrange the genes (y-axis) and microarrays (x-axis). B, A selection of differentially expressed genes (highly significant) that reliably distinguish T cell subsets. One or two genes that most reliably (lowest p value) differentiate between one particular T cell subset and all the others were selected and combined into the predictor list. The classification value of this list was validated when the correct T cell subset was predicted using the k-nearest neighbors method (p < 0.05 for most subsets).

A signature for γδ T cells

A perplexing question in T cell biology is the relevance and functions of the γδ T cell subset (32). Certain transcriptional profiling approaches such as serial analysis of gene expression (33, 34) have identified preferentially expressed genes in the mouse; however, to date, sheep and cattle are the only species in which γδ T cells are known to express a distinctive marker, T19 (35). We sought to identify genes specific to human γδ T cells that might underlie a function unique to this enigmatic subset. Surprisingly, we found few genes expressed specifically by γδ T cells. A number of genes highly expressed by γδ T cells were also expressed by NK cells, in particular, a number of KIR genes (Fig. 2⇑A, iv). The expression of these genes, however, was not confined to NK and γδ T cells, because some αβ effector T cells also expressed transcripts for these molecules. The gene expression signature of γδ T cells could be distinguished more easily by low or absent expression of a number of αβ T cell transcripts, rather than by the expression of any γδ-specific genes (Fig. 2⇑A, ii). Even the γ-chain of the TCR was expressed in some αβ T cells. Among the genes down-regulated in γδ T cells were several genes involved in pre-mRNA splicing processes, including transformer 2β, suppressor of white apricot homologue 2 (splicing factors, arginine/serine rich 10 and 16, respectively), and pre-mRNA splicing factor 16 (the conserved motif Asp-Glu-Ala-Asp/H box polypeptide 38). The down-regulation of these splicing factors, compared with that in αβ T cells, may be due to γδ T cells not requiring extensive alternative splicing of the TCR. It should be noted that to date, only blood γδ T cells have been assessed, not those from epithelial or lymphoid tissues.

T cell subset predictor genes

As an extension of our signature gene analysis for the different T cell subsets, we next sought to identify small numbers of predictor genes that could reliably distinguish between all the different T cell subsets in this study. We selected one or two genes with the lowest p value from each subset that most reliably distinguished that particular T cell subset from all other T cells, creating a list of 13 genes (Fig. 2⇑B). We found that using this small subset of genes it was possible to correctly identify the T cell type used in each microarray experiment using a k-nearest neighbors method. As few as 13 genes were sufficient to distinguish between 11 different T cell subsets and NK cells, demonstrating the remarkable power of this approach for cell type classification. Interestingly, the genes that made up the optimal predictor set were not necessarily well-known markers and included intracellular signaling molecule (such as phosphatases (DUSP4)) and obscure uncharacterized genes (such as BC006146 and AI766311, which do not show homology to any known gene families). These results demonstrate the utility of gene expression profiling for T cell type classification, similar to the way in which predictor genes have been used to distinguish cancer subtypes (36, 37, 38, 39). An important qualification, however, is that a much larger number of replicates will be required for a more reliable classification of unknown subsets of T cells, and probably the predictor set of genes will be refined with increased sample numbers for each subset and additional subsets, such as regulatory T cells. Nevertheless, our analysis demonstrates the potential of microarray dataset analysis for identification of unexpected molecules as biomarkers for T cell subset classification.

T cell activation induces an extensive transcriptional program

T cell subsets examined in this study represent many of the major paradigms of T cell biology, for instance, TCR γδ vs αβ T cells, Th1 vs Th2, and T cell activation and memory cell development. Principal component analysis was undertaken to identify the most distinctive expression patterns in all our T cell subsets. Principal component analysis involves a mathematical procedure that transforms a number of possibly correlated variables (in this case, gene expression data points) into a smaller number of uncorrelated variables, called principal components. The first three principal components, which accounted for more than half the variability in gene expression, clearly separated T cells based on the isolation method and activation state. Thus, T cell activation and culture were the two major influences on gene expression patterns in the T cell subsets in this study (Fig. 3⇓). This was not totally surprising, because previous studies have documented the profound changes in T cell transcription upon activation (28).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 3.

Principal component analysis (PCA) identifies key parameters affecting gene expression patterns in T cell subsets. PCA was performed on all T cell arrays using GeneSpring software. Only probe sets that were expressed in at least two T cell arrays were used in the analysis. This analysis identified the predominant patterns or parameters responsible for most of the variability in the T cell gene expression profiles. The first three principal components, which accounted for over half the variability in gene expression, separated T cells based on the isolation method and activation state.

To identify transcriptional changes induced by T cell activation, we analyzed gene expression in Th0 and Th2 cells before and after activation in vitro. TCR engagement (anti-CD3 and PdBu treatment) was accompanied by a vast change in gene expression in Th0 and Th2 cells. The majority of genes affected by activation were regulated in a similar manner in both Th0 and Th2 cells (Fig. 4⇓A). Activation induced an extensive change in gene expression, with >7000 genes significantly different between resting and activated T cells. A detailed description of many (or even some) of the individual genes is not possible in this study; however, global analyses that reflects changes in biological processed can be more informative, particularly when very large numbers of gene transcripts are regulated. We used gene ontology in the GeneSpring software program to identify the biological processes that were up-regulated in either resting or activated T cells (Fig. 4⇓B). The p value for the random overlap allowed us to identify the processes that were over-represented in our gene list compared with what would be expected by chance. As expected, activated T cells up-regulated a number of processes associated with higher metabolic and effector activities, including genes involved in transcription and translation, and also cell-cell signaling and cell adhesion. Not surprisingly, genes involved in immune effector function were significantly regulated; for instance, a large number of cytokines and chemokines showed major changes in expression. Indeed, the transcript levels for almost all cytokines and chemokines were increased after activation (not shown). Interestingly, IL-10R was down-regulated after activation, suggesting that T cells are less susceptible to suppression after TCR engagement.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 4.

Gene expression profiles of resting and activated T cells. Th2 clones generated from naive CD4⁺ T cells were cultured under Th1- and Th2-polarizing conditions to generate Th0 and Th2 clones, respectively, and RNA was extracted from resting and activated T cell clones after 4 wk. T cells (10⁷/condition) were left untreated (NA) or were activated (ACT) for 4 h with 1 μg/ml anti-CD3 Ab and 50 ng/ml PdBu. Gene expression of resting and activated T cells was analyzed using U133A+B Affymetrix microarrays (three different donors). A, Genes that showed significant changes in expression between resting and activated Th0 and Th2 cells were identified using one-way ANOVA. The numbers of genes regulated by activation in Th0, but not Th2; in Th2, but not Th0; and in both are shown by Venn diagram. B, The Gene Ontology biological process was created using GeneSpring software, and processes that showed a significant overlap with genes that increased or decreased after T cell activation are listed.

A number of processes involving proliferation, cell cycle, and apoptosis were regulated in both resting and activated cells. Resting T cells expressed higher levels of transcripts for genes involved in cytoskeleton organization, including dynein; myosin VA and IXB; capping protein; tubulin-specific chaperones c, e, and d; and others. Cytoskeleton rearrangements are an integral part of the T cell activation process (40) and are early events in the formation of the immunological synapse. Down-regulation of these transcripts in T cells stimulated for 4 h may serve as a feedback mechanism.

Gene expression in memory T cell subsets: T_CM vs T_EM

Gene signatures for T cell subsets are not necessarily informative for revealing all the differences between two particular subsets when a proportion of the genes is shared by numerous other T cell types. One interesting subdivision of the memory pool is the T_CM vs T_EM distinction defined by the presence or the absence of CCR7 (9). We compared T_CM and T_EM cells directly and identified differentially expressed genes using Affymetrix software. Surprisingly few genes were reproducibly differentially expressed between these two memory T cell populations (Fig. 5⇓), suggesting that the phenotypes of these two subsets are quite similar. Predictably, T_CM cells expressed higher levels of CCR7 and L-selectin, which are defining molecules for this subset and facilitate secondary lymphoid organ homing. T_EM cells expressed a number of molecules more closely associated with effector function, including chemokines CCL4, CCL5, and XCL1 and CTL-related molecules such as granulysin and granzymes A and B (41). Because the CD45RO⁺CCR7⁻ T_EM subset might be quite heterogeneous (10), there may well be circulating Th1- and Th2-like effector cells and possibly also some NK T and γδ T cells in this population. Nevertheless, our genechip results suggest that T_EM are somewhat related to and resemble effector T cells, although the distinction between T_EM and T_CM is also relatively small. Whether CD45RO⁺CCR7⁻ T_EM cells in tissues display a phenotype even closer to that of effector T cells is a likely possibility.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 5.

Genes differentially expressed by T_CM and T_EM cells. Genes that differed between T_EM and T_CM by >2-fold (as determined using MAS5.0 software) in two replicate experiments are shown. The color scale indicates the level of expression. The signal log ratio (SLR) is the log of fold change.