Abstract
The transcription factor FOXP3 is essential for the formation and function of regulatory T cells (Tregs), and Tregs are essential for maintaining immune homeostasis and tolerance. This is demonstrated by a lethal autoimmune defect in mice lacking Foxp3 and in immunodysregulation polyendocrinopathy enteropathy X-linked syndrome patients. However, little is known about the molecular basis of human FOXP3 function or the relationship between direct and indirect targets of FOXP3 in human Tregs. To investigate this, we have performed a comprehensive genome-wide analysis for human FOXP3 target genes from cord blood Tregs using chromatin immunoprecipitation array profiling and expression profiling. We have identified 5579 human FOXP3 target genes and derived a core Treg gene signature conserved across species using mouse chromatin immunoprecipitation data sets. A total of 739 of the 5579 FOXP3 target genes were differentially regulated in Tregs compared with Th cells, thus allowing the identification of a number of pathways and biological functions overrepresented in Tregs. We have identified gene families including cell surface molecules and microRNAs that are differentially expressed in FOXP3+ Tregs. In particular, we have identified a novel role for peptidase inhibitor 16, which is expressed on the cell surface of >80% of resting human CD25+FOXP3+ Tregs, suggesting that in conjunction with CD25 peptidase inhibitor 16 may be a surrogate surface marker for Tregs with potential clinical application.
Regulatory T cells (Tregs) are essential “master regulators” of tolerance (1), allowing effective responses to infection and non-self Ags while preventing autoimmunity and immune pathology (2). Although there are several subsets of cells with regulatory capacity (reviewed in Ref. 3), it is clear that the CD4+CD25+FOXP3+ natural Tregs (nTregs) is a key participant in the orchestration of tolerance, because without these cells an early onset aggressive autoimmune condition arises in both mouse and man (4–6). Furthermore, a number of autoimmune diseases involve a perturbation in Treg numbers or function (reviewed in Ref. 7), and control of the local response to acute tissue damage appears to require Treg function (8). The ability to monitor and manipulate Tregs may improve diagnosis or therapeutic intervention for autoimmune diseases. To achieve this, sensitive and specific biomarkers are required, but currently, however, there are limited markers available. This genome-wide study of human Tregs takes an unbiased approach for identification of novel markers and gene regulation networks.
The forkhead (FKH) transcription factor FOXP3 is expressed constitutively in natural CD4+ Tregs and has been shown to be essential for nTreg function (9). FOXP3 is a member of the P subfamily of the FKH family of transcription factors sharing a number of common features including a C2H2 zinc finger, a leucine zipper, and a FKH DNA-binding domain (10). These domains are involved in DNA binding, nuclear transport (11), homomeric and heteromeric complex formation (12), and transcriptional repressor activity (12, 13). Mechanisms of action of FOXP3 in Tregs have been identified in mouse, but less is known about its role in the development and function of human Tregs (reviewed in Ref. 12). In particular, the existence of two major human-specific FOXP3 isoforms (14) with different regulatory capacities (15, 16) and the transient induction of FOXP3 following activation of CD4+CD25− non-Tregs in humans (17) suggest a degree of species-specific function by, and regulation of, FOXP3.
In mouse and human Tregs, FOXP3 not only represses but also can induce gene expression, although the mechanism by which FOXP3 can negatively regulate transcription of some genes, such as IL2, but positively regulate others, such as CD25 and CTLA4, is still unclear (reviewed in Ref. 10). FOXP3 interacts with multiple transcription factors known to be involved in activation, differentiation, and response of CD4+ T cells to TCR stimulation, including NFAT, NF-κB, Runx1, retinoic acid receptor-related orphan receptors (RORs) (RORα and RORγT), IFN regulatory factor 4, STAT3, and Jun (15, 16, 18–24). FOXP3 may act as a transcriptional corepressor, because FKH domain-independent inhibition of the trans-activation activity of NF-κB, CREB, and RORα has been reported (16, 21). FOXP3 is also likely to modulate gene expression through epigenetic mechanisms, such as chromatin remodeling and histone deacetylation (25, 26). It is evident that the composition and activity of the FOXP3 regulatory complex is dynamic and is influenced by external stimuli acting partly through the posttranslational modification of FOXP3 (27, 28). These observations demonstrate the complex nature of the interaction of FOXP3 with its target genes.
In this article, we identify molecular targets of human FOXP3 in nTregs using a combination of chromatin immunoprecipitation array profiling (ChIP-on-chip) and expression profiling. Defining the direct binding sites for FOXP3 reveals target genes that may be regulated by FOXP3 and allows a detailed analysis of FOXP3-bound promoter elements, providing important information relating to FOXP3 transcription factor function. We have also identified direct targets of FOXP3 that are differentially expressed in human nTregs, and analysis of the FOXP3 target genes and differentially regulated genes has uncovered a number of transcription factors, microRNAs, signaling molecules, and surface molecules that have potentially important roles in Treg formation or function. In the search for positive cell surface markers for human Tregs, we report for the first time that peptidase inhibitor 16 (PI16) is expressed on the surface of a significant proportion of human Tregs, pointing to the potential utility of this molecule as a new biomarker for nTregs.
Materials and Methods
Isolation, in vitro expansion, and characterization of cord blood T cell populations
Cord blood was obtained with approval from the donor and the Children’s, Youth and Women’s Health Service Research Ethics Committee. Mononuclear cells were isolated from cord blood postpartum as previously described (29). Briefly, cord blood CD4+CD25+ (Treg) and CD4+CD25− (Th) cells were isolated from purified mononuclear cells using a Regulatory CD4+CD25+ T Cell Kit (Dynabeads; Invitrogen, Carlsbad, CA). Ex vivo expansion of isolated T cell populations (1 × 106 cells per well in a 24-well plate) were performed in X-Vivo 15 media supplemented with 5% human AB serum (Lonza, Walkersville, MD), 20 mM HEPES (pH 7.4), 2 mM l29).
ChIP, ChIP PCR, and whole-genome array
Expanded cord blood Tregs were cultured overnight in X-Vivo 15 media with 100 U/ml recombinant human IL-2 prior to a 2 h restimulation with 1 μM ionomycin (Sigma-Aldrich, St. Louis, MO) before cross-linking. Isolated T cells were subjected to cross-linking for 10 min in 1% formaldehyde solution (50 mM HEPES KOH [pH 7.5], 100 mM NaCl, 1 mM EDTA, 0.5 mM EGTA, and 11% formaldehyde). Formaldehyde was quenched by the addition of glycine to a final concentration of 112 mM. Cell lysis, ChIP, and DNA isolation steps for human FOXP3 ChIP-on-chip experiments were carried out essentially as described in Ref. 30, using 5 × 107 cells per immunoprecipitation with either a rabbit anti-Foxp3 IgG (Novus Biochem, Littleton, CO) or ChIP-grade control rabbit IgG Sera (Abcam, Cambridge, U.K.). Chromatin shearing conditions consisted of 12 repeats of 30 s pulses with an amplitude of 90% with 90 s between pulses with a S-4000 ultrasonic processor (Misonix, Farmingdale, NY) with a microtip. Amplification of immunoprecipitated and input chromatin (10 ng) for labeling and hybridization to Affymetrix Human Tiling 2.0R arrays was by whole-genome amplification (Sigma-Aldrich GenomePlex Whole Genome Amplification Kit) (31) adapted for Affymetrix labeling systems by the inclusion of 0.11 mM dUTP in the amplification reaction. Amplified material was purified using a PCR clean-up kit (Qiagen, Valencia, CA). Data from two independent ChIP-on-chip experiments were analyzed as replicates using Model-based Analysis of Tiling-array (32). Labeling and hybridization of amplified material to Affymetrix human tiling arrays was carried out at the Biomolecular Resource Facility (John Curtin School of Medical Research, Australian National University). Data files for this experiment are lodged on Gene Expression Omnibus (accession no. GSE20995; www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE20995).
ChIP-on-chip data analysis
Data from two independent FOXP3 ChIP whole-genome array sets were analyzed as replicates using Model-based Analysis of Tiling-array (32). A false discovery rate (FDR) of 0.5% was used to identify significantly enriched regions in FOXP3-immunoprecipitated material relative to input chromatin. Gene accessions were assigned to an individual ChIP region if the peak of the enriched region was within 20 kb upstream of a transcription start site (TSS) and 20 kb downstream of the transcription end site. Annotation was performed using gene accessions from both University of California Santa Cruz and National Center for Biotechnology Information, encompassing RefSeqs, mRNAs, and miRs, and where more than one accession could be associated with a ChIP peak, all possible annotations were recorded. Species conservation was analyzed using VISTA software (http://genome.lbl.gov/vista/index.shtml) (33).
RNA preparation and expression array
Differential expression analysis was carried out using Affymetrix Human Exon 1.0 ST arrays. Total RNA was isolated from expanded Th cells and nTregs that were rested for ~60 h prior to treatment for 2 h with either ionomycin or vehicle (DMSO). Total RNA, including small m.w. RNA, was isolated using QIAshredder and a miRNeasy Mini Kit (Qiagen). RNA quality was assayed using an Agilent Systems Bioanalyzer (Santa Clara, CA). Labeling and hybridization were carried out according to the manufacturer’s protocols at the Biomolecular Resource Facility (John Curtin School of Medical Research, Australian National University). Data files for this experiment are lodged on Gene Expression Omnibus (accession no. GSE20934; www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE20934).
RNA expression array data normalization
Microarray data were processed using RMA (34) and version 11.0 of an Entrez Gene-centric cdf (35). All of the processing was conducted using the statistical software R (R Development Core Team) under the aroma.affymetrix framework (36). Log fold change in a four-way comparison between resting and stimulated treatments of nTregs and Th cells was estimated within individual donors after Loess normalization. Final expression analysis was performed using Limma (37), and raw p values for each term were adjusted globally to provide an estimate of the FDR (38).
Differential gene expression group assignment and analysis
Differentially expressed genes were divided into seven mutually exclusive groups (T1–T6 and A1), based upon their behavior across all four comparisons. For genes with a significant adjusted p value in only one comparison, the accepted FDR was raised to p < 0.1 in all of the other comparisons to enable more accurate group assignment. Clustering within each subgroup was performed using the AGNES algorithm (M. Maelcher, unpublished observations), based on log fold-change and whether there was an associated FOXP3 binding site. Overrepresentation of ChIP hits within individual groups was assessed by using Fisher’s exact tests (all p values <10−7).
Pathway mapping
FOXP3 ChIP hits and the sets of differentially expressed genes comprising groups T1–T6 were assessed for overrepresentation of genes from canonical pathways and known biological functions either individually or together using Fisher’s exact tests within the Ingenuity Pathways Analysis software package (Ingenuity, Redwood, CA).
Transcription factor binding motif analysis
Transcription factor binding site analysis was carried out on 400-bp sequences centered on the ChIP peak associated with either the top 1000 FOXP3-bound regions (ranked on p value) or genes differentially regulated in Tregs. Scans for overrepresented motifs were performed on these regions using WEEDER and MatInspector (Genomatix, Munich, Germany) (39–42). Motifs identified by WEEDER were matched against known transcription factor binding motifs in JASPAR (43). For the 739 genes that were ChIP hits and members of groups T1–T6, 400-bp genomic sequences centered at the ChIP peak were scanned for the presence of transcription factor binding motifs using the default settings of MatInspector. Of these sequences, the 681 identified as containing one or more binding motifs belonging to the FKH family were further examined for co-occurrence of other transcription factor family motifs. A reference set of sequences was generated by randomly selecting 1000 genes from the Entrez Gene database that were not FOXP3 ChIP hits and selecting 400-bp sequences around the TSS using the observed distribution for the differentially regulated ChIP hits, as described above. A further set of 701 sequences containing a FKH motif was extracted as a subset of this reference set. Transcription factor families present in at least 30% of the differentially expressed ChIP genes with an FDR-adjusted p value < 0.05 were considered significantly overrepresented.
Quantitative real-time PCR
Validation of differentially regulated target genes identified by expression array analysis was performed on RNA from five donors. Random primed cDNA was prepared from total RNA using a QuantiTect Reverse Transcription Kit (Qiagen). Relative gene expression analysis was carried out using SYBR Green quantitative real-time PCR (qPCR) on a Rotor-Gene 6000 PCR machine (Corbet Research, Quiagen, Valencia, CA). Gene-specific primers pairs for qPCR were selected from PrimerBank (44). A primer set specific for RPL13a was used as an internal control (primer sequences listed in Supplemental Table X). Results were analyzed using a Rotor-Gene 6000, QGene software (45), and R. Validation of FOXP3 binding to genomic regions was carried out by ChIP qPCR. For ChIP, PCR primers to amplify target regions were either designed using Primer 3 Plus (www.bioinformatics.nl/cgi-bin/primer3plus/primer3plus.cgi) or purchased from SABiosciences, Frederick, MD as validated ChIP PCR primers (ChampionChIP PCR assays). Reactions were performed using a RT2 SYBR Green/ROX qPCR Master Mix (SABiosciences). The relative enrichment of target regions in FOXP3-immunoprecipitated material relative to input chromatin analysis was carried out using a ChampionChIP data analysis spreadsheet (SABiosciences) using the 2−ΔΔCT method. MicroRNA expression level comparisons were performed on total RNA (miRNeasy kits; Qiagen) isolated from freshly purified Tregs or Th cells by real-time PCR using specific TaqMan miR assays (Applied Biosystems, Carlsbad, CA). Data were normalized to a housekeeper (RNU6B) and non-FOXP3-regulated miRS (miR-16/miR-24). RNA from three independent donors was used to confirm the differential expression.
Validation of cell surface molecule expression
For surface molecule analysis, freshly isolated adult peripheral blood CD4+CD25+ Tregs or CD4+CD25− Th cells were either assayed immediately or stimulated overnight in the presence of CD3/CD28 beads (bead-to-cell ratio, 1:1) and 100 U/ml IL-2 prior to four-color flow cytometry using Abs against CD4, CD25, FOXP3, and test Ag. For PI16 detection, a mouse polyclonal Ab was used (Abnova, Walnut, CA; catalog no. H00221476-B01P)
Results
Isolation and expansion of human cord blood Tregs
We chose to use ex vivo expanded cord blood nTregs to dissect the FOXP3-driven events within human Tregs because large numbers can be grown for ChIP-on-chip studies and these cells display strong FOXP3 expression and a robust suppressive ability in vitro (29, 46, 47). We routinely obtained 100- to 200-fold expansion of cord blood CD4+CD25+ Tregs when cultured in vitro using anti-CD3/CD28 beads. These cells maintained a Treg phenotype, with >90% of the expanded cells staining CD4+, CD25hi, and FOXP3+ (Fig. 1A). These cells were also CD127dim (data not shown). The Tregs retained regulatory function because they were able to robustly suppress the proliferation of CD4+CD25− cells in vitro in an unmatched donor mixed leukocyte suppression assay (Fig. 1B). Thus, having demonstrated that we could expand large numbers of functional Tregs ex vivo, these cells were used for the ChIP-on-chip and expression profiling experiments.
Expanded cord blood nTregs express FOXP3 and have suppressive activity. A, Flow cytometric analysis of surface CD4 and CD25 and intracellular FOXP3 expression. Isotype control of CD4+ and CD25+ population (dotted line), FOXP3 expression in CD25− cells (dashed line), and FOXP3 expression in CD4+ and CD25+ cells (solid line). B, Functional activity of expanded cells in a MLR. Expanded cord blood nTregs were cultured for 5 d with unmatched T cells, APCs, and soluble anti-CD3 prior to a 16 h pulse of [3H]thymidine. Robust suppression of proliferation of CD4+CD25− responder cells by expanded nTregs was observed from a responder-to-nTreg ratio of 1:1 to 10:1.
Identification of human FOXP3 binding sites using ChIP-on-chip assays
The human FOXP3 ChIP protocol was established using the IL7Ralpha (CD127) and IL2 genes (Supplemental Fig. 1), because both promoters have been previously reported to contain binding sites for human FOXP3 (23, 48). FOXP3-bound DNA isolated by ChIP was amplified using a modification of the whole-genome amplification method (31), and enrichment of the same target genes postamplification was confirmed (data not shown) prior to labeling and hybridization to Affymetrix Human Tiling 2.0R arrays. With a FDR of 0.5%, a total of 8308 unique genomic regions were determined as significantly enriched in the FOXP3 ChIP material (Supplemental Table I). As a first step in validating the human FOXP3 ChIP-on-chip data set, a search for the presence of FKH DNA-binding motifs was performed on the top 1000 FOXP3-bound genomic regions (ranked on p value) using MatInspector (Genomatix). FOXP3-bound regions were significantly enriched for FKH motifs when compared against the Genomatix promoter database, with 922 of 1000 sequences (92%; p = 2.19 × 10−33) containing at least one FKH motif. Using the modeling algorithm WEEDER, we identified a motif (consensus [(G/A)(T/A)AAA(C/A)AA]) that was significantly overrepresented in these human FOXP3 ChIP regions, and comparison with transcription factor position weight matrices in the JASPAR database indicated that this represented a FKH binding site [(A/G)(C/T)(A/C)AA(C/T)A] (49). Our human FOXP3 target sequence (Fig. 2) is similar to the published mouse Foxp3 binding motif [(A/G)(T/C)AAACA] (50). This is consistent with localization of human FOXP3 occurring primarily through binding to a FKH-like motif.
The human FOXP3 binding motif in the top 1000 human FOXP3 ChIP hits derived from sequence information using WEEDER.
Location and functional annotation of FOXP3 binding regions
We mapped FOXP3-bound regions relative to known and predicted genes annotated in University of California Santa Cruz and National Center for Biotechnology Information to identify potential FOXP3 target genes. With a window beginning 20 kb upstream of annotated TSSs and ending 20 kb downstream of the annotated transcript end site of a gene, 7012 of the 8308 unique FOXP3-bound regions were associated with a gene. The remaining 1296 regions (13%) fell outside this window and were labeled intergenic. All of the possible gene annotations for a given ChIP peak were recorded, resulting in 5579 genes associated with at least one FOXP3-bound region (Supplemental Table I). Within this data set, FOXP3 binding sites were found in close proximity to 63 loci encoding miRS, including binding sites near loci encoding mir-155, mir-21, mir-22, and mir-142 that were are also mouse Foxp3 targets (50, 51) (Supplemental Table II). Examples of FOXP3-bound regions annotated to the human genome are shown in Fig. 3A. FOXP3 appears to have a preference for engaging in short-range promoter interactions to regulate the transcription of target genes because we observed a clear peak of FOXP3-bound regions starting 2 kb upstream and extending 3 kb downstream of the TSS (Fig. 3B).
Annotation, validation, and distribution of FOXP3 target genes identified by ChIP-on-chip. A, Annotation of three genes enriched for FOXP3 binding by ChIP showing the mouse/human conserved regions of the locus, the scale in kilobases of the region, the chromosome number and location, the individual FOXP3 binding regions, and the gene identification and the gene structure. B, Distribution of FOXP3 binding sites annotated to genes in relation to the distance from the TSS. C, The fold enrichment of target regions in FOXP3-immunoprecipitated material, normalized to control IgG immunoprecipitations, relative to input chromatin was determined for selected FOXP3 target genes by qPCR. Displayed is the mean fold enrichment ± SD; n = 3 for each region.
We selected 17 ChIP regions and four control regions for confirmation by ChIP PCR (Fig. 3C). Included in this subset are species-conserved FOXP3 targets and FOXP3 targets unique to our human data set. Greater than 2-fold enrichment relative to input chromatin (normalized for nonspecific chromatin immunoprecipitation) was detected at the 17 sites tested by qPCR, whereas no enrichment (≤1-fold) was observed for four different genomic regions predicted not to bind FOXP3 (intergenic 12p13, AFM intron 1, PDE3B intron 10, and human IGX1A negative control, SA Biosciences).
Analysis of the potential FOXP3 target genes using Ingenuity Pathway Analysis software demonstrated a substantial enrichment for genes involved in the proliferation, activation, and control of cell death in T lymphocytes in normal immune responses and in disease states, such as inflammatory and autoimmune diseases. Molecules involved in the response to TCR signaling and immunoregulatory cytokines/growth factors were significantly overrepresented in the human FOXP3 target gene data set (Supplemental Table III).
Differential gene expression in resting and activated nTregs
To correlate FOXP3 binding with the transcriptional regulation of target genes, expression profiling was carried out on expanded cord blood nTregs and Th cells. Expression profiling of donor-matched Treg and Th cell populations in both a stimulated and a resting state (Fig. 4A) identified three classes of genes: those with an nTreg intrinsic expression pattern, those that responded to stimulation in a cell type-specific manner, and those that display common responses to T cell activation. Genes displaying an nTreg-specific response were then clustered into six groups based on patterns of differential expression by cell type and/or stimulation (Table I). In total, 1851 genes displayed significant differential expression in a Treg compared with a Th cell (T1–T6), with a further 746 genes showing a similar response to activation in both cell types (A1) (Fig. 4B; described in Table I; gene identifications in Supplemental Table IV). A total of 739 genes in these groups were also FOXP3 ChIP targets (gray bar in Fig. 4B). These FOXP3 ChIP hits were statistically overrepresented in all T groups (p < 1 × 10−7), but we did not observe any statistically significant difference in the distribution of FOXP3 ChIP hits between the upregulated and downregulated genes.
A, Schematic representation of the cells and comparisons used in the expression profiling experiment. B, Heat map representation of significant differentially expressed genes in all four expression array comparisons. The normalized log2 fold changes for repressed (blue) and induced (red) transcripts across each comparison are shown. FOXP3 target genes identified by ChIP-on-chip analysis are indicated in gray on the right-hand side. C, RT-PCR validation of differentially expressed genes in expanded cord blood nTregs versus Th cells. Displayed is the mean log fold change in the relative level of the genes in expanded nTregs compared with Th cells in stimulated (open bar) and resting (filled bar) cells (n = 6 donors). The log fold change determined from analysis of the exon arrays for the nTreg versus Th cell comparisons in resting (hatched bar) and stimulated (stripped bar) cells is also shown. Dotted lines represent fold change of ±2. D, RT-PCR validation of potential FOXP3 target miRs. The relative expression of the miR in matched freshly isolated cord blood nTregs and Th cells was determined on n = 3 donors in triplicate using TaqMan miR assays. Statistical significance (p > 0.05) using a paired t test was achieved for all but those marked with an asterisk.
Of note, groups T1 and T2, which represent the steady-state human nTreg signature, also contained a higher proportion of FOXP3 ChIP hits when compared with those of T3–T6 (50% versus 35% p = 3.064 × 10−9), which is consistent with a role of FOXP3 in stabilizing and maintaining the Treg phenotype (52). Also, within the group of genes that were differentially expressed upon stimulation of the cells (T3–T6), there were clear differences, with a number of genes responding to stimulation in a Th cell but remaining unchanged in a Treg (T5), whereas a further set of genes responded to stimulation in an nTreg but remained unchanged in a Th cell (T6). qPCR performed on a group of candidate genes identified in our array confirmed the differential expression and clustering analysis (Fig. 4C). Similarly, TaqMan microRNA (miR) assays performed on 20 of the miRs predicted to be FOXP3 targets in nTregs demonstrated that five were differentially expressed in nTregs (nine miRs shown in Fig. 4D). To identify potential targets regulated by these miRs, a search for miR binding sites in the 3′ UTRs of downregulated FOXP3 target genes was performed using the microRNA target prediction programs in miRGen (www.diana.pcbi.upenn.edu/miRGen.html). Using stringent criteria in which a miR binding site must be found by three separate prediction programs (miRanda, PicTar4, and TargetScan), we identified a number of high confidence gene targets for these validated microRNAs (Supplemental Table V). Many of these predicted miR target genes (21 of 33) were also FOXP3 ChIP hits, identifying a group of genes that are potentially subjected to multiple FOXP3-dependent mechanisms to ensure their downregulation in nTregs.
Transcription factor motifs in differentially expressed FOXP3 target genes
The molecular mechanism by which FOXP3-dependent regulation of gene expression occurs is not well characterized, but a growing list of transcription factors has been shown to interact with FOXP3. To identify other potential transcription factors that may coregulate gene expression with FOXP3, we analyzed sequences flanking the FOXP3 ChIP peak for the 739 genes that were also differentially expressed. Position weight matrices from Genomatix corresponding to the binding sites for transcription factor families AP-1, BCL6, HAML (Runx family), HNF1, MYBL, NFAT, OVOL, PAX2, PRDF (PRDM1/Blimp1 family), and STAT transcription factor families were found to be significantly overrepresented in FOXP3-bound regions, whereas the CAAT motif only appears significantly overrepresented in groups T1 and T2 (FDR-adjusted p = 0.0462) (Fig. 5A). In the STAT family (Fig. 5B), the STAT5 consensus binding motif was the most significantly enriched, suggesting a preferential interaction of FOXP3 with STAT5 within nTregs. These results have therefore uncovered several new transcription factors that potentially directly interact with, or coregulate gene expression by, FOXP3.
Overrepresentation of transcription factor binding sites in FOXP3 target gene promoters. A, Analysis of differentially expressed genes in the Treg groups that are also ChIP hits for the corepresentation of transcription factor family binding sites. Red bar, T1 + T2; blue bar, T3–T6; white bar, control promoters. B, Analysis of the STAT family of transcription factors demonstrating that STAT5 binding sites are significantly overrepresented in FOXP3 target genes.
Comparison of human and mouse FOXP3/Foxp3 ChIP data
We next compared our ChIP data set with the previously identified murine targets (50, 51). Because substantial differences existed between the two mouse data sets, the mouse data sets were each compared separately with our human data (Table II). The mouse data were then combined for clarity (Fig. 6). Comparison of the human FOXP3 targets with mouse Foxp3 targets indicated that 23% (888 genes) of our human FOXP3 ChIP targets were also Foxp3 targets in at least one mouse data set (Supplemental Table VI, Table II) (50, 51), with 107 genes common to all three data sets. These overlaps were statistically significant (p < 10−5), and a number of genes previously implicated in Treg function were identified in these overlapping groups. On the basis of our gene expression profiling, we found a significantly higher proportion of differentially expressed genes (20–24%) among the conserved FOXP3 targets (Table II), suggesting that these genes may represent a core group of FOXP3-dependent genes involved in the maintenance of the regulatory phenotype in both species. Pathway and biological function analysis of these 888 species-conserved FOXP3 target genes using Ingenuity Pathway Analysis software revealed that members of critical T cell signal transduction pathways including TCR, CD28, and CTLA4 signaling and the intracellular MAPK signaling pathway were overrepresented among these conserved targets. These pathways were also found to be significantly overrepresented among the human-specific ChIP genes (Supplemental Tables VII, VIII), suggesting that although the overall regulation of essential pathways and biological processes may be conserved in human and mouse Tregs, the genes within a pathway used to achieve this may vary between species.
Overlap between human and mouse genes associated with a FOXP3-bound region. The numbers represent target genes common (center) and unique to each species (outer). All of the comparisons were statistically significant (p < 0.005).
Analysis of cell surface molecules in nTregs
When we applied functional annotation to differentially expressed genes in the data set, we discovered a highly significant (p = 3.36 × 10−16) enrichment for surface molecules in groups T1 and T2, with 42.6% of genes in these groups being cell surface related versus 23.8% in groups T3–T6. This included known CD molecules (Supplemental Table IX) as well as less well characterized surface molecules. However, FACS-based analysis failed to identify a CD molecule uniquely expressed on >50% of the FOXP3+ cells (data not shown). We therefore further mined our expression array data set to identify additional surface proteins not assigned a CD designation that were selectively upregulated on the nTregs, and from this we focused on PI16. The differential expression of PI16 in resting and stimulated cord blood nTregs versus Th cells was confirmed by real-time PCR (Fig. 4C).
We next examined PI16 expression in adult peripheral blood CD4+ T cell subsets. PI16 was readily detectable on the surfaces of resting and stimulated nTregs (Fig. 7A, 7B). In particular, on resting Tregs, a significant portion (~85%) of the CD4+CD25brightFOXP3+ cells were PI16+, and the majority (~70%) of the PI16+ cells coexpress FOXP3, with only ~15% FOXP3+PI16− (Fig. 7Aii). In contrast, although ~50% of the resting CD25− cells are PI16+, the majority of these (80%) are FOXP3− (Fig. 7Aiii). Upon isolation of CD4+CD25+ cells and overnight stimulation, ~80% of the PI16+ cells coexpress PI16 and FOXP3, although ~40% of the activated FOXP3+ cells are PI16− (Fig. 7Bii). Importantly, stimulated CD4+CD25− do not express significant amounts of PI16 or FOXP3, suggesting that PI16 expression is not linked to activation but rather to the nTreg phenotype. It is interesting to note that although the proportion of FOXP3+ PI16+ cells is similar between resting and stimulated overnight Tregs (7B, Fig. 7Bii) the proportion of FOXP3+PI16− cells increases upon stimulation, suggesting that stimulation alters surface expression of PI16. This is not due to loss of PI16 transcription, as confirmed by the array and quantitative RT-PCR data in resting and activated Tregs. There are also FOXP3+PI16− cells at a lower frequency, which may represent a distinct Treg subset. In our chip analysis, PI16 was not a direct target of FOXP3 in human or mouse, suggesting that it is regulated by genes downstream of FOXP3. Together, these data indicate that PI16 warrants further investigation as a FOXP3 surrogate marker.
PI16 and FOXP3 expression in nTregs and Th cells. A, Validation of PI16 as a surrogate for FOXP3 in resting adult Tregs. When CD4+ cells were gated on the top 1% CD25+ (i), >85% of the FOXP3+ cells are PI16+ (ii). When gated on CD25−, few of the PI16+ cells are FOXP3+ (iii). B, CD25+ and CD25− CD4 cells were isolated according to the gates shown (i) and then stimulated overnight with anti-CD3/CD28 beads and stained for PI16 and FOXP3 expression. Sorted stimulated CD25+ cells stained robustly with the anti-PI16 Ab, and this population contains >60% of the FOXP3+ cells (ii). The majority of the activated CD25− cells do not express PI16 or FOXP3 (iii). Representative data from n = 3 donors.
Discussion
Using expanded cord blood nTregs to generate FOXP3 ChIP-on-chip and expression profiling data, we show, for the first time, human FOXP3-bound genes and their differential expression profile in human nTregs. We detected 8308 potential FOXP3 binding sites in the human genome, which were associated with 5579 genes. Of these potential FOXP3-regulated genes, 13% (739) displayed differential expression between donor-matched expanded nTregs and Th cells. We selected an early time point for the expression analysis to identify genes directly downstream of this altered signaling. Expression profiles at later time points may uncover additional FOXP3 target genes whose expression level changes more slowly following stimulation. Nonetheless, our finding that a subset of FOXP3-bound regions can be linked to the regulation of the neighboring genes is consistent with other genome-wide ChIP studies where a similar level of transcription factor binding and gene expression changes in nearby genes has been reported (reviewed in Ref. 53).
The comparison of FOXP3-bound genes from this study in humans with published studies in mice (50, 51) provided a means of identifying genes with a conserved function in Tregs. Although >50% of the FOXP3 targets are bound in a species-specific manner, 23% of FOXP3-bound genes identified in our human study were present in one of the published mouse Treg studies. However, upon pooling of the Foxp3-bound genes identified in the two mouse studies, 45% of FOXP3 targets identified in mouse Tregs were also bound by FOXP3 in human Tregs. A caveat to these comparisons is the small overlap of the Foxp3 target genes between the two mouse studies (7% overlap; data not shown). However, this degree of conservation of FOXP3 target genes in human and mouse cells is consistent with other large scale transcription factor binding and expression studies, where a substantial divergence in transcription factor binding events in human and mouse cells has been reported (54–57). These species-conserved FOXP3 targets contained proportionately more differentially regulated genes (20% of conserved targets versus 16.2% in all targets, p = 0.0104), suggesting that the curation of conserved targets is a valid approach for identifying proteins associated with the key role of FOXP3 within the Treg. The proportion of differentially regulated genes increased to ~25% (p = 1.66 × 10−4) when only the Foxp3-bound genes from primary mouse Tregs (51) were compared with our human FOXP3 targets. In contrast, when the group of mouse Foxp3-bound genes derived using forced expression in hybridoma cells was compared (50), there was no significant enrichment for differentially regulated genes (18%, p = 0.224). This may be a consequence of the differences in intracellular signaling pathways and transcription factor levels in hydridomas compared with those of primary mouse T cells or a result of overexpression of FoxP3 and highlights the difficulties in comparing transcription factor binding in different cell systems. In support of this, expression profiling performed by Marson et al. (50) indicated that Foxp3 overexpression almost overwhelmingly downregulated gene expression in an activation-dependent manner in hybridoma cell lines, whereas 30–50% of the Foxp3-bound genes were upregulated in mouse (50) and human (this study) primary Tregs.
Our ChIP-on-chip assay also identified FOXP3 binding regions in proximity to loci encoding 63 microRNAs (miRs), which is consistent with FOXP3 partly exerting its regulatory program by modifying the expression of these posttranscriptional rheostats of protein expression. In mice, a critical role of miRs in nTreg biology has been demonstrated by the T cell-targeted deletion of either of the two key enzymes involved in the generation of mature miRs [DICER (58–60) or DROSHA (61)]. Deletion of either of these enzymes led to a spontaneous lethal inflammatory disease in mice due to a severe reduction in the stability and function of nTregs. MicroRNA signatures for mouse and human nTregs have been reported (59, 62) and several miRs have been shown to be important in nTreg function, including miR155 in mouse (63, 64) and miR21 in human cells (62). We have confirmed that both of these miRs are located near FOXP3 binding sites by ChIP PCR and that both are differentially upregulated in human nTregs, suggesting that these are conserved targets for FOXP3 regulation. In addition, we have identified three further miRs potentially regulated by Foxp3 (miR-146a, miR-101, and miR-7) that are upregulated in human nTregs. Thus, in addition to the four miRs also identified as potential Foxp3 targets in the mouse (miR-7, miR-21, miR-22, and miR-155) (50, 51) and 10 microRNAs recently identified as differentially expressed in human Tregs (62), we find that the FOXP3 targets miR-101 and miR-146a are also significantly upregulated. These results suggest that other miRs are likely to contribute to the Treg gene program and that further work is required to establish their identity and role. Knockdown of FOXP3 in human Tregs to assess the role of individual microRNAs underway.
Analysis of the human FOXP3-bound regions that were associated with differentially regulated genes identified a number of transcription factor motifs that appeared to co-occur with FKH motifs, including AP-1, Runx, NFAT, and STAT binding motifs. The co-occurrence of these motifs is consistent with known roles for these transcription factors in nTregs (18, 65–69) and, in the case of AP-1, NFAT, and Runx family members, with their known physical interaction with FOXP3 (18, 20, 23). Several of the other transcription factor family motifs that are overrepresented in these sequences, such as BCL6 and PRDF, also contain members known to have a role in T cells, such as BCL6, PRDM1/Blimp-1, and c-myb (70–73). The finding of potential PRDM1 motifs within FOXP3-bound regions is of particular interest given its relatively high level of expression in natural Tregs and its requirement for T cell homeostasis, nTreg function, and self-tolerance in the mouse (74, 75), suggesting that FOXP3 and PRDM1 may co-operate to regulate gene expression in human Tregs. Furthermore, the PRDM1 gene is a conserved FOXP3/Foxp3 target in humans and mice (Ref. 51 and this study), suggesting a regulatory loop involving these two transcription factors.
A major role for FOXP3 is to modulate the response of a T lymphocyte to external stimuli, such as TCR and cytokine signaling, by altering the level and activity of downstream transcription factors (76). Several of the transcription factors that are able to bind to sites enriched in FOXP3-bound regions are themselves potential targets of FOXP3 and are differentially expressed. For example, NFATc2, NFATc3, STAT1, and STAT4 are all potential ChIP targets of FOXP3 and are downregulated in Tregs. A direct role for FOXP3 downregulating the level of NFAT proteins has recently been supported by the finding that forced FOXP3 expression in human cell lines inhibits the activation-induced expression of NFAT2 (NFATc/NFATc1) in T cells by competing with NFAT1 (NFATc2) for binding to a site within the NFAT2 promoter (77). Among the binding sites identified as belonging to the STAT family of transcription factors, there was a strong overrepresentation of STAT5 binding sites at FOXP3-bound regions. This enrichment was not observed for the binding sites of other STAT family members. Thus, the apparent enrichment of STAT5 sites within FOXP3-bound regions, coupled with a decrease in the expression of other STAT family members, such as STAT1 and STAT4, suggests an essential requirement for and increased sensitivity to the IL-2/STAT5 signaling pathway for Treg development and homeostasis (66, 67, 69).
Clustering of our gene expression data showed clear differences in the response to stimulation, with over one third of the differentially expressed genes showing cell type-specific responses to stimulation (groups T5 and T6; 35% = 897/2597). On top of these activation-influenced signatures, another group of genes was revealed that was differentially expressed in both resting and stimulated Tregs, suggesting that these genes represent a steady-state human nTreg expression profile (T1 + T2; 21% = 552/2597). Together these data may be helpful for the further development of a gene signature to investigate the relative contributions of nTregs and other CD4 subsets in disease pathology (78, 79).
We detected a significant enrichment of differentially regulated cell surface molecules in the steady-state nTreg signature (T1 + T2), consistent with a role for surface molecules in Treg function, including the homing response (80, 81) and suppressive function (82). In particular, two molecules, LRRC32/GARP and PI16, were significantly upregulated in human nTregs. LRRC32/GARP has recently been reported to be specifically upregulated in human nTregs upon stimulation where it has been proposed to play a role in stabilizing FOXP3 expression and the localization of latent TGF-β on the surface of nTregs upon TCR stimulation (83, 84). Significantly, we have found that PI16, unlike LRRC32/GARP, is differentially expressed on resting Tregs. Hence, we propose it as a potential biomarker for human nTregs, where in conjunction with CD25 and LRCC32/GARP it may be feasible to identify both resting and stimulated Treg subsets. Little is known about the function of PI16, which was initially identified as a serum binding partner of prostate secretory protein 94 (85) and more recently has been shown to have growth inhibitory properties (86, 87). The Treg expression profile and the potential growth suppressive activity of PI16 indicate that the contribution of this surface protein to Treg function clearly warrants further investigation.
In conclusion, we have generated a gene profile of human Tregs and identified putative direct targets of FOXP3, which allows us to model direct and indirect FOXP3 target genes. Pathway mapping and transcription factor over representation analysis provide novel information on the mechanisms by which human nTregs function, and we can now undertake functional studies of the cell surface, secreted, and other molecules, which may be critical for Treg function. This may hence provide new leads for diagnosis or therapeutic intervention.
Acknowledgments
We thank Sylvie Nobbs and Dr. Randall Grose for expert assistance with flow cytometry and Dr. Doreen Krumbiegel for critical reading of this manuscript.
Disclosures The authors have no financial conflicts of interest.
Footnotes
The microarray data presented in this article have been submitted to Gene Expression Omnibus under accession numbers GSE20995 and GSE20934.
The online version of this article contains supplemental material.
Abbreviations used in this paper:
- ChIP-on-chip
- chromatin immunoprecipitation array profiling
- DE
- differentially expressed
- FDR
- false discovery rate
- FKH
- forkhead
- miR
- microRNA
- mRNA
- messenger RNA
- nTreg
- natural regulatory T cell
- PI16
- peptidase inhibitor 16
- qPCR
- quantitative real-time PCR
- ROR
- retinoic acid receptor-related orphan receptor
- Treg
- regulatory T cell
- TSS
- transcription start site.
- Received January 12, 2010.
- Accepted May 3, 2010.
- Copyright © 2010 by The American Association of Immunologists, Inc.
References
In this issue
