Genome-Wide Identification of Long Noncoding RNAs in CD8⁺ T Cells

Mice, viruses, and CD8⁺ T cell priming

C57/BL6 (B6) and congenic B6.Ly5.1 female mice were purchased from the Walter and Eliza Hall Institute’s Animal Services Department (Kew, Australia). Animals were generally used at 8–12 wk of age. Experiments were conducted under the auspices and with the approval of the Austin Health Animal Ethics Committee (Heidelberg, Australia) and conformed to the Australian National Health and Medical Research Council’s code of practice for the care and use of animals for scientific purposes.

A/Puerto Rico/8/34 influenza virus was grown in 10-day embryonic chicken eggs. Recombinant vaccinia virus (rVV) expressing full-length OVA protein (rVV-OVA) was generated as described previously (18) and was a gift from Drs. J. Yewdell and J. Bennink (National Institute for Allergy and Infectious Diseases (NIAID), National Institutes of Health, Bethesda, MD). For CD8⁺ T cell priming, mice were injected with either 1 × 10⁷ PFU of influenza virus via an i.p. route or 5 × 10⁶ PFU of rVV-OVA via an i.v. route.

Cell lines and effector CD8⁺ T cell culture

The thymoma cell line EL-4 and dendritic cell-like cell line DC2.4 (a gift from Dr. K. Rock, University of Massachusetts Medical Center, North Worcester, MA) were cultured in RPMI 1640 containing 10% FCS, 50 μM 2-ME, and 2 mM glutamine (RP-10). For growing effector CD8⁺ T cell lines, influenza-primed B6.Ly5.1 mouse spleens were harvested and 3 × 10⁷ splenocytes were cultured in 6-well plates with RP-10 containing 10 U/ml recombinant human IL-2 and stimulated with 3 × 10⁶ autologous splenocytes pulsed with peptide at 10⁻⁹ M. Media were changed every 2 days. Viable T cells were collected through Lymphoprep gradients (Axis-Shield) and enriched by depleting B220⁺ and CD4⁺ cells using M-450 Dynal beads coated with appropriate mAbs. Ten to 14 days after stimulation, peptide-specific effector CD8⁺ T cells generally comprised ∼60–90% of the total cells in bulk culture, as confirmed by tetramer and intracellular cytokine staining (see below).

Abs and peptides

The mAbs used (FITC- or allophycocyanin-labeled anti-IFN-γ, FITC-labeled anti-CD44, FITC-labeled anti-CD45.1, PE-labeled anti-CD62L, YC-Chrome-labeled anti-CD8α, and purified anti-B220 and anti-CD4) were purchased from BD Biosciences. For flow cytometry, mAbs were used at 1/200 dilution in PBS supplemented with 10% FCS. Peptides corresponding to the dominant influenza CD8⁺ T cell determinants NP_366–374 and PA_224–233 were procured and characterized by the Biologic Resource Branch of the NIAID and were gifts from Drs. J. Yewdell and J. Bennink (NIAID). In each case, substances with the predicted mass constituted >95% of the material analyzed. Peptides were dissolved in DMSO at 1 mM and stored at −20°C.

Tetramer analysis and intracellular cytokine staining (ICS)

Influenza-specific CD8⁺ T cells were identified using tetrameric complexes of H-2D^b and NP_366–374 or PA_224–233 peptide (ImmunoID). Cells were stained with PE-conjugated tetramer and Cy-Chrome-labeled anti-CD8α mAb for 30 min at room temperature, washed, immediately assessed by flow cytometry on a FACScalibur apparatus (BD Biosciences), and analyzed using FlowJo software (Tree Star).

For ICS with influenza-specific CD8⁺ effectors, EL4 target cells were used as stimulators and pulsed with peptide at a concentration of 100 nM for 30 min at 37°C, washed thoroughly, and then incubated with 0.5–1 × 10⁵ cultured CD8⁺ effectors at 37°C for 3 h in 200 μl of RP-10 in round-bottom 96-well plates with brefeldin A (Sigma-Aldrich) at 10 μg/ml. For ICS with naive and memory rVV-OVA-specific CD8⁺ T cells, rVV-OVA-infected DC2.4 cells were used as stimulators and infected with saturating amounts of rVV-OVA in PBS and 0.1% BSA for 4 h at 37°C, washed thoroughly, and then incubated with 1.5–2 × 10⁶ splenocytes at 37°C for 4 h in 200 μl of RP-10 in round-bottom 96-well plates with brefeldin A at 10 μg/ml. After incubation, cells were stained with Cy-Chrome-labeled anti-CD8α mAb with or without PE-labeled CD45.1 mAb at 4°C for 30 min, washed, and fixed with 1% paraformaldehyde in PBS at room temperature for 20 min. Fixed cells were then stained with fluorescein-anti-IFN-γ in PBS containing 0.2% saponin (Sigma-Aldrich) and analyzed by flow cytometry.

CD8⁺ T cell subset isolation

To purify naive CD8⁺ T cells, T cells were initially isolated by magnetic separation from the spleens of naive B6 mice using Thy1.2 MACS beads (Miltenyi Biotec) as per the manufacturer’s instructions, labeled with anti-CD8α and anti-CD62L mAb at 4°C for 30 min, and washed. Naive CD8⁺ T cells (CD8⁺CD62L^high) were then isolated via FACS on a high-speed MoFlo sorter (DakoCytomation). To purify memory CD8⁺ T cells, T cells were isolated by magnetic separation from the spleens of rVV-OVA-primed B6 mice, labeled with anti-CD8α and anti-CD44 mAb, and washed. Memory CD8⁺ T cells(CD8⁺ CD44^high) were then isolated via FACS. To purify activated effector CD8⁺ T cells, cells in bulk culture were initially activated by coculturing for 3 h at 37°C with EL-4 target cells prepulsed with the relevant peptide at 10⁻⁷ M. The cells were then labeled with anti-CD8α and anti-CD45.1 mAbs and washed. Activated effector CD8⁺ T cells (CD8⁺CD45.1⁺) were then isolated via FACS. To purify nonactivated, effector CD8⁺ T cells, cells in bulk culture were initially cocultured for 3 h at 37°C with unpulsed EL-4 target cells. Cells were then labeled with anti-CD8α and anti-CD45.1 mAbs and washed. Nonactivated effector CD8⁺ T cells (CD8⁺CD45.1⁺) were then isolated via FACS. After sorting of each CD8⁺ T cell subset, cell purity was checked by flow cytometry and found to be ∼91–98% of the desired population.

RNA isolation and labeling

RNA was purified from CD8⁺ T cell subsets using an RNeasy mini kit (Qiagen), and the quality of purified total RNA samples was assessed with an RNA 6000 Nano assay kit using the Agilent 2100 Bioanalyzer (Agilent Technologies) according to the manufacturer’s instructions. For microarray experiments, RNA was amplified and labeled using the Amino Allyl MessageAmp II kit (Ambion) following the instructions provided by the manufacturer. Amplified aminoallyl RNA was labeled with either Cy3 or Cy5 monoreactive dyes (GE Healthcare) according to the MessageAmp II protocol (Ambion). Satisfactory dye incorporation was confirmed by spectrophotometric analysis, and the quality and quantity of amplified RNA samples were assessed using the Agilent 2100 Bioanalyzer.

Custom microarray expression analysis

Design of the custom microarray has been described elsewhere (10), and it is available from the ArrayExpress data warehouse (European Molecular Biology Laboratory-European Bioinformatics Institute (EMBL-EBI; ArrayExpress accession no. A-MEXP-1070; www.ebi.ac.uk/microarray- as/ae/). In summary, the microarrays contained 22,038 65-mer oligonucleotide probes from the Compugen mouse oligo library and 2,118 70-mer oligonucleotide probes that were designed to specifically target ncRNAs, including known mouse pre-miRNAs from miRBase (19), longer mouse ncRNAs from RNAdb (20), and “high confidence” ncRNAs identified from the FANTOM3 project (21). Probes were printed on PowerMatrix slides (Full Moon BioSystems) at the Special Research Centre Microarray Facility of the University of Queensland (Brisbane, Australia), and the quality of the print run was verified by hybridizing random 10-mer oligonucleotides to the first and last slides of the run.

Naive and memory CD8⁺ T cells were compared by competitively hybridizing 5 μg of Cy3- or Cy5-labeled RNA from each subset to individual microarrays. Activated and nonactivated effector CD8⁺ T cells were similarly compared. Two biological and three technical replicates (including dye swaps) were performed for each comparison. Blocking, hybridization, and washing were performed according to the manufacturer’s instructions (Full Moon BioSystems). Slides were scanned at 5-μm resolution using a DNA microarray scanner (Agilent Technologies). Feature extraction was performed with ImaGene software (BioDiscovery) using manual grid adjustment, autospot finding, and segmentation. Data was exported from ImaGene and analyzed using the LIMMA (Linear Models for Microarray Data) software package via the R Project for Statistical Computing (www.r-project.org). Data were background corrected, normalized both within and between arrays (22), and differential expression analysis was performed by fitting a linear model of the data to the experimental design matrix and then calculating Bayesian statistics (B statistics; posterior log odds) adjusted for multiple testing using Benjamini-Hochberg analysis (23). Based on previous analyses (10), an A value of >7 was selected as a conservative indicator of probe expression, and a B statistic of >3 and a fold change of >1.5 were selected as the threshold for significant differential expression. Raw and processed microarray data is available at the ArrayExpress data warehouse (EMBL-EBI; ArrayExpress Accession nos. E-TABM-491 and E-TABM-492; www.ebi.ac.uk/microarray-as/ae/).

Microarray probe classification

Our microarray probe classification has been reported previously (10, 24) but has been recently updated and is described here in detail. Although the Compugen mouse oligo library probe set was predominantly designed to recognize protein-coding transcripts, several thousand probes targeted miscellaneous transcripts whose coding status was not well characterized at the time this commercial probe set was first produced. A computational pipeline was therefore designed to reannotate the entire probe set. Sequences for all probes were mapped to the February 2006 (National Center for Biotechnology Information (NCBI) build 36) assembly of the mouse genome using BLAT (BLAST-like alignment tool) (25), and probes that mapped uniquely to the genome were analyzed further. Probes that targeted protein-coding transcripts were initially identified based upon genomic overlap with all known protein-coding mRNAs in mouse and other mammals as defined by NCBI reference sequence (RefSeq) (26), the Mammalian Gene Collection (27) and the University of California Santa Cruz (UCSC) Known Genes (28). Next, we applied the CRITICA (coding region identification tool invoking comparative analysis) algorithm to assess the protein-coding capacity of all remaining transcripts (29). Probes corresponding to transcripts that CRITICA predicted to be protein coding were annotated as such. Probes that did not target protein-coding genes were then classified as either “noncoding” if they targeted a transcript of >200 nt designated by CRITICA as a ncRNA or “unknown” if they targeted a genomic region without any associated transcript. Notably, many individual probes were found to jointly recognize both protein-coding and noncoding transcripts. This was hardly surprising given the highly interlaced nature of the mammalian transcriptome (30, 31), but given the potential ambiguity of using such probes, they were excluded from all subsequent analyses.

Genomics Institute of the Novartis Research Foundation (GNF) SymAtlas analysis

Publicly available GNF Gene Expression Atlas data (symatlas.gnf.org) (32) were downloaded for mouse and human tissues. GNF probes were classified using the same computational pipeline described above. Absent/present/marginal calls were obtained from the above GNF web site and used to identify probes with detectable expression. Hierarchical clustering of expressed probes was based upon the Pearson correlation and the single linkage method using Cluster software (33) and visualized via Java TreeView (34).

Gene Ontology (GO) term analysis

The GO project (www.geneontology.org/) provides a controlled vocabulary for functionally annotating genes. By studying GO terms associated with a large set of genes, one obtains useful information about the types of genes represented. Lists of differentially expressed protein-coding genes identified from the custom ncRNA microarray experiments were uploaded into GOstat (gostat.wehi.edu.au/) (35). Statistically over-represented GO terms in the biological process and molecular function categories were obtained by applying a p value cutoff < 0.05, correcting for multiple testing with the Benjamini false discovery rate, and using the Mouse Genome Informatics (MGI) database (www.informatics.jax.org) as the reference set.

Quantitative real-time PCR (qRT-PCR)

cDNA preparation and quantitative RT-PCR analysis were performed as previously described (10). The following primers were used to quantitate Ptpre and Ptpreas isoform levels: Ptpre long, 5′-TCCCAGCCGCTGCTCA-3′ (forward) and 5′-AGAACTTCTTGGGGCCTGAT-3′ (reverse); Ptpre short, 5′-AAGAACTTTTCCCGGCTCAC (forward) and 5′-AGAACTTCTTGGGGCCTGAT-3′ (reverse); AK017962, AK020025, and AK08864, 5′-AGCCTCCAGCACACTCATTT-3′ (forward) and 5′-GAGAAGCAGGCTTTTGATGG-3′ (reverse); and AK148045, 5′-CCTTGTGAAGCCCTTGAGAG-3′ (forward) and 5′-GTGAACGAACATCACCATCG-3′). In all quantitative RT-PCR experiments, a minimum of three replicates was performed and RNA levels were quantified relative to Tubulin-D1 RNA expression using the following primer: 5′-GCTTGGCATGTCGATAGTGA-3′ (forward) and 5′-GAACAGAGTCCCTGGGAACA-3′ (reverse).

ncRNA genomic context analysis

Determination of the genomic context of ncRNAs (relative to protein-coding genes) has been described previously (9, 36) but was updated for this study. In summary, cis-antisense ncRNAs were defined where the corresponding targeted transcript was mapped to the opposite strand of a RefSeq-annotated exon (this includes 5′-untranslated region (UTR), coding, and 3′-UTR); intronic ncRNAs were defined where the corresponding probes mapped within the intron of a protein-coding gene; bidirectional ncRNAs were defined where the corresponding probes targeted transcripts that were oriented head-head to a protein-coding gene at a distance of <1000 bp; and promoter-associated ncRNAs were defined where the targeted transcripts overlap the RefSeq-annotated transcription start site.

ncRNA promoter analysis

Promoter regions were classified as described previously (37). Briefly, transcripts with a 500-bp interval within −0.5 to +2 kb of the transcription start site with a GC fraction ≥ 0.55 and an observed-to-expected ratio ≥ 0.6 were classified as high CpG promoters. Promoters where all of the 500-bp intervals within −0.5 to +2 kb of the transcription start site have a CpG observed-to-expected ratio ≤ 0.4 were classified as low CpG promoters. CpG observed-to-expected ratios were calculated as described previously (38). Promoter methylation status was obtained from publicly available DNA methylation maps previously generated for mouse CD8⁺ T cells (39). These maps were constructed using high-throughput, reduced representation bisulfite sequencing and provided >5 million high-quality sequencing reads covering a majority of CpG islands within the genome.

Evolutionary conservation and secondary structure analysis of ncRNAs

ncRNA conservation and secondary structure composition were determined by intersecting their chromosomal positions with those of genome-wide PhastCons elements (40) and RNAz predictions (41), respectively. The PhastCons program uses genome alignments across multiple vertebrate species to identify conserved genomic regions based upon a phylogenetic hidden Markov model (40). The RNAz program predicts functional ncRNAs by identifying genomic regions that contain thermodynamically stable and conserved RNA secondary structures (41). The significance of these predictions is represented by the class probability P, and conserved RNA secondary structure predictions were considered significant at a threshold of P > 0.5 as previously reported (42). We also used RNAfold software to predict and visualize the secondary structure of selected ncRNAs (43).

Genome-wide expression analyses reveal presence of hundreds of long ncRNAs in CD8⁺ T cells

To investigate the extent to which long ncRNAs are expressed in CD8⁺ T cells, we analyzed three independent microarray expression data sets. The first two data sets, from human and mouse, were derived from the GNF Symatlas (symatlas.gnf.org/) (32). The third was based on a custom-designed microarray (10), and the results from each analysis are presented below.

The GNF Symatlas is a large-scale, publicly available gene atlas that features expression data for 44,775 human and 36,182 mouse transcripts using whole genome Affymetrix arrays. In total, the atlas profiles 79 human and 61 mouse tissues, including a variety of immune-related organs and cell types (e.g., CD8⁺ T cells). Although the probe set used by GNF was designed primarily to target the protein-coding transcriptome, we had informally observed that some probes appeared to target long ncRNAs. We therefore designed a computational pipeline to systematically classify the GNF probe set (see Materials and Methods) and found that 1,234 human and 5,622 mouse probes uniquely recognized long ncRNAs. Of these, 268 (22%) and 1,211 (22%) were expressed in human and mouse CD8⁺ T cells, respectively (supplemental Tables S1A and S1B).⁵

We recently designed custom ncRNA microarrays whose probes uniquely profile 4,329 high-confidence long ncRNAs (and 14,316 protein-coding transcripts) from the mouse genome (10). To complement the GNF data, we used these arrays to generate a third data set that examined the expression of CD8⁺ T cells at different stages of differentiation and activation. After naive CD8⁺ T cells encounter Ag, a minority of them differentiate into Ag-specific memory cells over subsequent weeks. This transition is typically characterized by down-regulation of L-selectin (CD62L) and up-regulation of Pgp-1 (CD44) (44). We used FACS to isolate highly purified populations of naive (CD62L^high) and memory (CD44^high) CD8⁺ T cells from the spleens of unprimed and vaccinia-infected B6 mice, respectively (supplemental Fig. S1A). Consistent with the expected functional differences between these two subsets, only memory CD8⁺ T cells were able to produce IFN-γ upon exposure to viral Ag (supplemental Fig S1B). To obtain effector CD8⁺ T cells, we harvested splenocytes from mice infected >30 days earlier with influenza virus and exposed the cells to influenza-specific peptides in vitro. After 10–14 days, effector cells were activated via coculture with peptide-pulsed target cells and isolated to high purity by FACS (supplemental Fig. S1A); nonactivated control effector cells were similarly purified by FACS following coculture with unpulsed target cells. As expected, only activated CD8⁺ T cells produced IFN-γ in an ICS assay (supplemental Fig. S1B) and, based upon tetramer staining, down-regulated TCR surface expression (supplemental Fig. S1C). Following isolation of each of these CD8⁺ T cell subsets, RNA was competitively hybridized to custom ncRNA microarrays. Experiments were designed to directly compare either naive and memory cells or activated and nonactivated effector cells. Of 4,329 noncoding probes examined, 1,106 (26%) were expressed in at least one of these experiments (supplemental Table S1C). Interestingly, the majority of these were common to both experiments (740/1106; 67%), whereas the remainder were expressed only in naive/memory cells (97/1106; 9%) or effectors (269/1106; 24%).

Differential expression of long ncRNAs in CD8⁺ T cells

We expected that some ncRNAs expressed in CD8⁺ T cells would fulfill general cellular housekeeping functions and were thus likely to be expressed ubiquitously; alternatively, other ncRNAs might have specific roles in T cell immunity and be expected to display lymphoid-restricted and/or stage-specific expression. To investigate these possibilities, we performed hierarchical clustering of the mouse GNF Symatlas ncRNA data set. The majority of noncoding probes (60%; 3354 of 5622) were expressed in at least one of 61 tissues examined. Of these, 33% (1110 of 3354) were present in the majority of tissues (supplemental Table S2A) whereas the remainder showed more limited expression patterns (including testis-specific, oocyte-specific, CNS-specific, and thyroid-specific clusters) (Fig. 1⇓A). Of particular interest was the presence of 96 lymphoid-specific ncRNAs whose expression was generally either shared between or restricted to particular lymphoid organs (spleen, thymus, bone marrow, or lymph node) or lymphocyte subsets (CD8⁺, CD4⁺, or B cells) (Fig. 1⇓B; supplemental Table S2B). In total, 29 of these lymphoid-specific ncRNAs were expressed in CD8⁺ T cells. Hierarchical clustering of the human GNF ncRNA data set produced similar clusters of ubiquitous and lymphoid-specific ncRNAs (supplemental Fig. S2).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 1.

Hierarchical clustering of long ncRNAs expressed in mouse tissues. Affymetrix MAS5 call data for individual mouse GNF ncRNA probes were obtained for 61 different tissues examined in the GNF Symatlas. A, Three thousand, three hundred fifty-seven probes were expressed in at least one tissue and hierarchically clustered based upon their tissue distribution. Individual rows represent a single ncRNA probe, individual columns represent different tissue types (supplemental Table S2A), and individual pixels are colored to reflect expression (white, present; gray, marginal; black, absent). Many ncRNA probes show widespread expression, whereas others form tissue-specific clusters. B, Lymphoid-specific ncRNAs were clustered into two main blocks representing either shared (bottom, orange) or restricted (top, blue) lymphoid expression in the indicated tissues.

To identify ncRNAs that might be important during CD8⁺ T cell differentiation and activation, we examined the custom ncRNA array data set for ncRNAs whose expression was dynamically regulated. One hundred six of 1106 expressed noncoding probes (10%) were differentially expressed under the experimental conditions tested. Of these, 21 were significantly altered during naive to memory cell differentiation, 81 changed during effector cell activation, and four were differentially expressed during both transitions (supplemental Table S2C). Importantly, analysis of protein-coding probes on the custom array confirmed that our microarray experiments had successfully identified biologically relevant transcripts: many mRNAs with known roles in CD8⁺ T cell biology were dynamically regulated (supplemental Table S2D) in conjunction with long ncRNAs, and GO terms such as T cell differentiation and T cell activation were significantly enriched among the differentially expressed mRNAs (p < 0.05).

Expressed long ncRNAs share loci with protein-coding genes

Long ncRNAs have previously been shown to originate from complex loci that contain interlaced networks of noncoding and protein-coding transcripts (31, 36). In these instances, the exact nature of the genomic relationship between long ncRNAs and their overlapping neighbors often has important functional consequences, especially in relation to local epigenetic regulation (5, 6, 8, 45, 46, 47). Therefore, we characterized long ncRNAs as cis-antisense, intronic, bidirectional, and promoter-associated ncRNAs based upon their orientation to local protein-coding genes (see Materials and Methods, Table I, and supplemental Tables S3–S5). In doing so, we discovered that many ncRNAs associated with protein-coding genes, many of which fulfill previously characterized roles in T cell biology. A number of compelling examples are detailed below.

cis-Antisense ncRNAs

cis-Antisense transcripts are prevalent throughout the mammalian genome (48, 49) and are believed to regulate expression of neighboring mRNAs via a range of mechanisms, including chromatin modification, alternative splicing, endogenous siRNA generation, translational interference, and small RNA promoter targeting (47, 50, 51, 52, 53, 54, 55, 56). Within our data sets, we found that 16–22% of expressed ncRNAs was organized antisense to protein-coding exons (Table I⇓ and supplemental Tables S3A–S5A). For example, a ncRNA (AK029296, Lef1as RNA) was identified antisense to the lymphoid enhancer binding factor 1 Lef1 (Fig. 2⇓A). Lef1 is a downstream effector molecule in the Wnt signaling pathway and acts as a transcription factor to critically regulate T cell differentiation. Lef1 exists in two isoforms; the longer isoform acts as a transcriptional enhancer while the shorter one, missing a critical 5′ domain, has dominant negative function (57). In naive CD8⁺ T cells, the shorter inhibitory isoform predominates but is down-regulated following Ag encounter (58). The antisense ncRNA, which we named Lef1as, overlaps exclusively with the longer isoform and is predominately expressed in naive T cells (Fig. 2⇓A), consistent with a possible role in suppression of the longer isoform. This organization is also apparent in the human Lef1 locus, which indicates that this cis-antisense pairing is conserved and therefore likely to be functional.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 2.

Expressed ncRNAs share loci with CD8⁺ T cell-related, protein-coding genes. Long ncRNAs expressed in CD8⁺ T cells were categorized based upon their orientation to neighboring protein-coding genes. Examples of cis-antisense (A and B), intronic (C), bidirectionally oriented (D), or promoter-associated (E) ncRNAs are shown. Each panel illustrates the organization of the ncRNA (red, transcription initiation direction indicated by arrow) to associated protein coding gene (blue). For clarity, only selected mRNAs from the GenBank track are shown. A, Lef1as (AK029296) and Lef1 genomic locus. B, Top panel and inset, Ptpreas (short isoforms: AK017962, AK020025, and AK088640; long isoform: AK148045) and Ptpre genomic locus. Lower panel, Expression profiling of the indicated transcripts by quantitative RT-PCR (see Materials and Methods for details). C, Il2ra intronic ncRNAs (M21979, M121980, and M21891) are found adjacent to the exon boundaries of Il2ra (see detailed inset). This intronic organization is also conserved in humans as shown. D, The Snhg1 ncRNA is organized bidirectional to the Slc3a2 gene with which it shares a CpG island. Multiple C/D box snoRNAs (Snord22–31) are nested within introns of Snhg1. E, A promoter-associated ncRNA (AK009498) overlaps the promoter and 5′-UTR of the Flip gene in the sense orientation.

We also identified an antisense ncRNA that lies opposite to the protein tyrosine phosphatase receptor type E (Ptpre) gene in a highly complex locus (Fig. 2⇑B). Ptpre exists as two isoforms generated by alternative promoter usage; the longer form is membrane bound and the shorter one is cytoplasmic and selectively inhibits the Jak-Stat pathway critical for cytokine signaling in T cells (59). The antisense ncRNA, named here Ptpreas, itself exists as several different isoforms that vary in the extent to which they overlap with the different isoforms of Ptpre. Although all isoforms overlap with a common Ptpre coding exon, only the three shorter isoforms (AK017962, AK020025, and AK088640) overlap with the transcription start site and first exon of the shorter cytoplasmic Ptpre transcript (Fig. 2⇑B, see inset), because the sequence complementary to this region is spliced from the longer ncRNA isoform (AK148045). We hypothesized that the different Ptpreas ncRNAs might function to regulate the balance of the two Ptpre variants in CD8⁺ T cells. However, specific profiling of the different Ptpreas and Ptpre isoforms by quantitative RT-PCR failed to establish any direct correlation between them (Fig. 2⇑B, lower panel).

Intronic ncRNAs

At complex loci, the most common organizational arrangement we observed involved expressed CD8⁺ T cell ncRNAs hosted within the introns of protein-coding genes (31–46%; Table I⇓ and supplemental Tables S3B–S5B). In some cases, long ncRNAs were completely nested within individual introns. For example, we identified a number of ncRNAs (M21979, M21980, and M21981) hosted within the introns of the IL 2 receptor IL2Rα (Fig. 2⇑C). These ncRNAs directly flank the seventh and ninth exon of the Il2ra gene and were previously considered simple artifacts of unprocessed Il2ra mRNA (60). However, we observed strong and up-regulated expression of one of these ncRNAs (M21981) upon effector cell activation, although the ncRNA was not expressed during the differentiation of naive to memory cells. In contrast, the Il2ra gene itself is strongly up-regulated in both transitions as previously reported (61, 62), suggesting that the fate of the mRNA and intronic ncRNA might be independently regulated. Furthermore, homologous ncRNAs are apparent within the human Il2ra gene, suggesting this unusual organization is conserved within the mammalian lineage (Fig. 2⇑C).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 3.

Long ncRNAs expressed in CD8⁺ T cells overlap small RNAs. Long ncRNAs that were expressed in CD8⁺ T cells and overlap small RNAs were identified. Each panel illustrates the organization of the ncRNA (red; transcription initiation direction indicated by arrow) to annotated small RNA (green) or protein coding gene (blue). Blue histogram indicates conservation within vertebrates according to MultiZ alignment (mm8; UCSC Genome Browser). For clarity, only selected mRNAs from GenBank track are shown. A and B, We identified ncRNAs that overlap previously annotated miRNAs such as mmu-mir-22 (AK008813) (A) and mir-142–5p and mir-142–3p (AK020764) (B). C, We also identified examples of long ncRNAs (AK010403) that overlap novel small RNAs (see detailed inset illustrating red, purple and blue small RNAs from mouse embryonic stem cells according to Barbiarz et al. (82 ). D, These small RNAs correspond to stems, as colored, within RNA secondary structure as predicted by RNAz and folded according to CONTRAfold (97 ).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 4.

Evolutionary conservation of long ncRNAs expressed in CD8⁺ T cells. Long ncRNAs that were expressed in CD8⁺ T cells and show evidence of either primary sequence or secondary structure conservation were identified. NcRNAs (red; transcription initiation direction indicated by arrow) and protein coding gene (blue) are shown above histogram (blue), indicating conservation within vertebrates according to MultiZ alignment (mm8; UCSC Genome Browser). A, Tug1 and an alternative isoform (DQ517434) that exhibit strong conservation are widely overlapped with PhastCons elements (black). B, An ncRNA (AK154685) overlaps the entire Ccl14 gene and two inverted LINE elements (green) that are strongly predicted to form an extended RNA duplex (see inset). Additional predicted RNA secondary structures are shown (pink).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 5.

Genomic context of expressed ncRNAs^a

^a The percentage of expressed transcripts from the total targeted by the array in each class is indicated in parentheses. For detailed information on each expressed transcript, see supplemental Tables S3–S5.

Bidirectional ncRNAs

More than 1000 bidirectional transcript pairs have been previously identified within both mouse and human genomes that encompass ∼10% of known genes (36, 63). Unlike cis-antisense transcripts, bidirectional pairs do not share any genomic overlap but instead are separated by <1000 bp, are divergently transcribed, and typically share a common CpG island promoter. Among expressed CD8⁺ T cell ncRNAs, <10% formed bidirectional pairs (Tables I and supplemental Tables S3C–S5C). As an example, the long ncRNA Snhg1, previously identified as a host gene encoding multiple small nucleolar RNAs (snoRNAs) within its introns (64), is bidirectional to the amino acid transporter, Slc3a2 (Fig. 2⇑D). Also known as Ly-10 or CD98, Slc3a2 was originally identified in lymphocytes and subsequently considered a marker of T cell activation (66, 66). We observed an up-regulation of Slc3a2 upon CD8⁺ effector cell activation in mouse as expected, as well as an up-regulation of Snhg1. Interestingly, the bidirectional orientation of these two transcripts, arising from a common CpG island promoter, is conserved in both human and mouse, which suggests that the observed coupling is functional.

NcRNAs overlapping promoters

Transcription of long ncRNAs across promoters of protein-coding genes that lie downstream on the same genomic strand has been shown to regulate mRNA expression both positively and negatively (13, 67). For example, cascading ncRNA transcription across a yeast promoter was recently shown to convert chromatin to an open structure and thereby increase accessibility to transcription cofactors and RNA polymerase II binding (68). ncRNA transcription has also been shown to inform the promoter choice of the TCR (13). Indeed, in our profiling experiments we identified a long ncRNA (AK030929) expressed from the Vα-Jα region at the 3′-end of the TCRα locus, consistent with this previously reported role of noncoding transcription regulating Jα promoter activity.

We therefore searched systematically for long ncRNAs that overlapped promoters and transcription start sites of protein-coding genes and found that such ncRNAs represented <5% of expressed ncRNAs (Table I⇑ and supplemental Tables S3D–S5D). Many of these ncRNAs overlapped the 5′-end of protein-coding genes that function in a T cell development. For example, a long ncRNA (AK009498) overlapped the transcription start site of Flip (Fig. 2⇑E), a gene essential for the development and survival of mature T lymphocytes (69, 70). Interestingly, we observed that Flip is up-regulated during the transition of T cells to an activated or memory phenotype, whereas AK009498 is in contrast down-regulated during both of these transitions in keeping with a potential negative regulatory function.

Expressed long ncRNAs frequently overlap small RNAs

Small ncRNAs such as miRNAs and snoRNAs can be processed from longer primary ncRNAs (51, 72, 73, 74) as we noted earlier with Snhg1 (Fig. 2⇑D). To explore the possibility that other long ncRNAs expressed in mouse CD8⁺ T cells might also act as primary transcripts for small RNAs, we systematically searched for genomic overlap between expressed long ncRNAs and known small RNAs.

Mature miRNAs bind to the 3′-UTRs of mRNAs and, in doing so, regulate the stability and fate of thousands of protein-coding transcripts (reviewed in Ref. 75). We found that 18 long ncRNAs expressed in mouse CD8⁺ T cells overlap annotated miRNAs (supplemental Table S6A), and a further 21 overlap with snoRNAs (supplemental Table S6B). For example, a long ncRNA (AK008813) up-regulated during effector cell activation contains mmu-mir-22 within its second exon (Fig. 3⇑A). This direct overlap suggests that the long ncRNA is likely to be processed into and function via this smaller RNA, and indeed mir-22 has recently been detected at a low copy number in CD8⁺ effector cells (76). Moreover, the most strongly predicted targets of mmu-mir-22, based on the TargetScan algorithm (77), include multiple genes that function in CD8⁺ effector cells (supplemental Table S6C). One of the most strongly up-regulated ncRNAs during effector cell activation (AK053349) partially overlaps a well-characterized cluster of miRNAs (mir-17/18/19a/20a-1/19b/92-1). Interestingly, this cluster is often amplified in human lymphomas (78, 79) and its forced overexpression in mouse lymphocytes results in lymphoproliferative disease and autoimmunity, most likely via enhanced lymphocyte proliferation and survival following activation (80). Lastly, two miRNAs, mir-142-5p and mir-142-3p, recently shown to be among the most highly expressed miRNAs in naive, memory, and effector CD8⁺ T cell populations (76), are hosted within the first intron of a long ncRNA (AK020764) that is also strongly expressed in CD8⁺ T cells (Fig. 3⇑B). When mir-142-5p is ectopically expressed in hematopoietic progenitor cells, the proportion of cells differentiating down the T lymphoid lineage increases (81), further supporting a role for mir-142, and by extension AK020764, in T cell biology.

Recently, it has been demonstrated that endogenous siRNAs exist within a variety of mouse tissues, including oocytes and embryonic stem cells (51, 82, 83). These reports also revealed a diversity of mammalian small RNA-generating pathways and an abundance of other novel small RNAs, the functions of which are still largely unknown. To estimate the extent to which our long ncRNAs might be processed into siRNAs and other novel small RNAs, we examined the genomic locations of these recently identified small RNAs and compared them to those of the long ncRNAs expressed in CD8⁺ T cells. Although the small RNAs we examined were derived from embryonic stem cells and might not include the full range of small RNAs found in T cells, we nevertheless found that ∼20% of long ncRNAs expressed in mouse CD8⁺ T cells on the custom arrays directly overlapped these small RNAs (supplemental Table S6D). In some cases, these small RNAs overlapped stable secondary structures within the long ncRNAs. For instance, a long ncRNA (AK010403), abundantly expressed in T cells, was identified directly downstream of Nbl4a, a component of the β-catenin/Tcf pathway (84). This ncRNA is predicted to fold into a highly stable structure whose extended dsRNA regions are likely to be processed into numerous small RNAs (Fig. 3⇑, C and D). Another prominent example of a long ncRNA overlapping small RNAs was MALAT1 (supplemental Fig. S3), a long ncRNA previously associated with nuclear S35 splicing domains (85). We found that MALAT1 is associated with abundant small RNAs previously identified from mouse oocytes and embryonic stem cells (51, 82). Interestingly, these small RNAs are not significantly altered after Dicer knockdown, suggesting they are not produced by the canonical Dicer-dependent dsRNA cleavage pathway (82). Consistent with this observation, MALAT1 was recently shown to be cleaved by RNase P to generate a small tRNA-like cytoplasmic ncRNA (86), and our findings suggest that additional small RNAs might be generated from non-Dicer RNases.

Promoter analysis and DNA methylation status of expressed ncRNAs

Mammalian gene promoters can be classified into distinct categories based on their CpG dinucleotide content. The majority of protein-coding genes (∼64%) have high CpG content (HCG) (37) and typically play “housekeeping roles” within the cell (87). In contrast, low CpG (LCG) content is seen with ∼17% of genes (37) whose functions are generally tissue specific (87). In this way, promoter class can provide clues to gene function. We therefore classified the promoters of long ncRNAs expressed in mouse CD8⁺ T cells. Of those identified on the custom arrays, surprisingly only 22% had HCG content, whereas 32% had LCG content. Similar results were obtained for CD8⁺ T cell ncRNAs identified from the mouse GNF atlas (23% HCG and 36% LCG). Taken at face value, these results indicate that a greater proportion of long ncRNAs are likely to exert tissue-specific effects, as suggested by previous studies (9). However, our results should be treated with some caution, because the 5′-ends and promoter regions of long ncRNAs are currently not as well annotated as those of protein-coding genes.

DNA methylation helps determine local chromatin architecture and gene expression, with unmethylated CpG dinucleotides indicative of an open chromatin conformation and regional transcription. To determine the methylation status of ncRNA promoters, we compared previously published, genome-wide DNA methylation maps for mouse CD8⁺ T cells with the 5′ start sites of long ncRNAs expressed on the custom arrays (39). For this analysis, we focused on 208 ncRNAs with HCG promoters, because their annotation (and therefore methylation status) was likely to be the most reliable. Of these ncRNAs, 29 (14%) had HCG promoters that were specifically unmethylated in CD8⁺ T cells (and methylated in other tissues) (supplemental Table S7), providing further evidence that the expression of these transcripts is regulated in CD8⁺ T cells.

Expressed ncRNAs exhibit evolutionary conservation of primary sequence and predicted secondary structure

Across evolutionary time, the purifying selection of functional genomic elements results in the presence of sequences that exhibit high levels of conservation across multiple species, providing a useful indicator of function. We therefore identified elements (known as PhastCons elements) within expressed CD8⁺ T cell ncRNAs that are strongly conserved across species (40). We found that the majority of human and mouse ncRNAs expressed in CD8⁺ T cells contained PhastCons elements: 180 of 268 human (67%), 964 of 1211 mouse GNF (80%), and 700 of 1106 mouse custom array (63%) probes targeted ncRNAs with ≥20 PhastCons bases (supplemental Table S1A—S1C). One illustrative example was the Tug1 ncRNA, whose DQ517434 isoform was comprised of 32% PhastCons bases (Fig. 4⇑A). Although Tug1 was originally identified as having an essential role in retinal differentiation, our observation of Tug1 expression in both mouse and human CD8⁺ T cells is consistent with previous reports of splenic Tug1 expression and induction by taurine (88), the predominant free amino acid in lymphocytes that attenuates activation-induced T cell death (74, 89).

Functional ncRNAs often contain conserved secondary structures that may be processed into small RNAs (see above) or serve specific functions, such as binding protein partners. Stereotyped selection signatures can facilitate the identification of these functional RNA secondary structures. We used the RNAz prediction algorithm to identify conserved RNA secondary structures in the mouse and human genome (42) (see Materials and Methods). We found that a significant number of long ncRNAs contained predicted secondary structures. Thirty-nine of 268 human (15%), 417 of 1211 mouse GNF (34%), and 278 of 1106 mouse custom array (25%) probes targeted ncRNAs with ≥20 RNAz conserved bases (supplemental Table S1A–S1C). In addition to identifying a number of previously annotated secondary structures such as miRNA stem loops (e.g., miR-296) and snoRNAs (e.g., snoRNA-HACA39/HACA60), our analysis also unearthed a number of high-confidence novel structures (supplemental Table S8). For example, the ncRNA (AK154685) is a noncoding isoform of Ccl4, a chemokine important in T cell adhesion and migration (90). Up-regulated in both memory and activated T cells, this ncRNA overlaps the entire Ccl4 pre-mRNA and then extends ∼750 bp downstream (Fig. 4⇑B) to encompass two inverted L1 repeats that are strongly predicted to fold into a long, stabilized, double-stranded RNA hairpin (Fig. 4⇑B, inset).