|
|
||||||||
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||











* Centre for Molecular Medicine and
Arthritis Research Campaign Epidemiology Unit, University of Manchester, Manchester, UK; and
GlaxoSmithKline Research and Development, Stevenage, UK
1 Correspondence: Centre for Molecular Medicine, Stopford Building, Faculty of Medicine, University of Manchester, Oxford Rd., Manchester, M13 9PT, UK. E-mail: david.w.ray{at}manchester.ac.uk
| ABSTRACT |
|---|
|
|
|---|
Key Words: network signaling microarray dexamethasone
| INTRODUCTION |
|---|
|
|
|---|
Evidence indicates that genetic variation within the glucocorticoid receptor gene could explain some of the population variation in glucocorticoid sensitivity (15
16
17)
. We, and others, have associated genetic variation within the glucocorticoid receptor gene, either with differences in response to synthetic glucocorticoids or to physiological variables that may have resulted from altered glucocorticoid sensitivity. However, no clear evidence suggests how this genetic variation may be exerting its effect nor does it appear that this genetic variation is sufficient to explain more than a small percentage of the total population variation seen (17)
.
Microarray technology allows gene profiling to be performed to identify absolute gene expression levels within individuals, comparison between individuals of basal gene expression, and also analysis of the genome response to gradient concentrations of glucocorticoid. This approach has been applied before to analyze the profile of glucocorticoid-regulated genes in human blood-derived T-lymphoblastoid cells, and it appeared that up to 10% of the entire genome was regulated by glucocorticoid exposure (18)
. More recently gene expression profiles predictive of Gc response in asthma have been identified, again in peripheral blood cells (19)
. Peripheral blood is easily accessible, and lymphoblasts derived from peripheral blood afford a tractable model system in which dynamic gene profiling can be performed in response to glucocorticoid. Therefore, we recruited healthy subjects and stratified them into two groups by glucocorticoid sensitivity (responders; those in the lowest decile of post-dexamethasone (dex) cortisol and nonresponders; those in the highest decile of post-dex cortisol) as measured by cortisol response to an overnight low-dose dex suppression test (17
, 20)
. Individuals falling into the upper and lower deciles of response were further analyzed by gene profiling of their peripheral blood T-lymphocyte response to incubation with dex. We sought to determine the transcriptome profile of glucocorticoid responsiveness to identify patterns of genes that were robustly glucocorticoid regulated in all individuals and genes that were regulated by Gc in a manner dependent on the Gc sensitivity of the donors. Finally, we sought to identify a pattern of genes that would correctly predict the glucocorticoid sensitivity of the individuals from whom the blood was taken.
| MATERIALS AND METHODS |
|---|
|
|
|---|
Identification of glucocorticoid regulated genes
Genes were classified as Gc regulated if they showed a linear increase or decrease in expression as the glucocorticoid dose increased. These data were obtained using log10 dose as a linear predictor in an analysis of covariance model.
We selected a probability cutoff of
0.05 to determine significance alongside a slope estimate of
0.2.
Quantitative reverse transcriptase-polymerase chain reaction (RT-PCR) of Gc-regulated genes
A subset of 17 Gc-regulated genes were reanalyzed by using quantitative RT-PCR (Sybr Green, qRT-polymerase chain reaction, using Applied Biosystems, Foster City, CA, USA, methodologies on ABI 96 well Taqman machine). The primer sets used are available on request.
Power analysis for glucocorticoid sensitivity dataset
Power analysis was performed on the cohort. The calculation used the between-subject variability calculated for every gene from a mixed model ANOVA to estimate the number of donors needed in each Gc group to obtain at least 80% power to be able to detect a 1.5-, 2- or 2.5-fold change between donor groups. Percentiles of the distribution of variability were calculated across all the genes and were used to form the power calculations. Additionally a particularly variable gene was chosen and also used to inform the power calculations.
Cluster analysis
A mixed-model ANOVA was used to obtain predicted mean expression for each donor group at each dex concentration for each gene, which had shown glucocorticoid regulation. These means were put into a clustering algorithm to find groups of coexpressed genes. Grouped genetic algorithms developed by Optimal Design, (Brussels, Belgium) as described by Falkenauer, were used to perform Gaussian clustering (21)
.
Discriminant gene analysis
The data were normalized using the RMA algorithm (22)
and log (2-based) transformation. Gene filtering was performed based on P-values from Resolver to remove genes that were not expressed in any experimental condition (Gc groups and doses). Different endpoints, including difference on 100 nM dose, change from baseline on 100 nM, and differences in log-dose-response slope, were considered. The relationship between gene expression and Gc group was explored for each gene using ANOVA, and the results were shown by volcano plots. The genes were ranked based on the P-values from ANOVA, and the best 2000 genes were selected for a discrimination analysis to attempt to classify samples using the expression data.
The partial least-squares based discrimination analysis (PLS-DA) and an approach combining the support vector machine (SMV) (23)
with recursive feature selection (24)
and bootstrap aggregating (Bagging) (25)
were used. Models using different numbers of genes from the top of the ranking list were built to find one with the lowest classification error. To reduce selection bias in estimation of classification error, we used an external cross-validation (CV) method. The importance of using the external validation for SVM was demonstrated (26)
. In this method, data were split into a training set (2/3 of total samples) and a validation set (1/3 of total samples). For the PLS-DA, the genes were ranked using univariate P-values based on the training set and a model of top K genes was built to classify the validation data. These steps were repeated 100 times and the average classification error was calculated. The SVM approach is similar except that a recursive feature elimination approach replaces the P-value based ranking.
Network analysis
The glucocorticoid regulated genes were used to build a literature-derived subnetwork of biological relationships between the genes using a heuristic algorithm and scoring function as described (27)
. The global network used was constructed from biological relationships between genes/proteins using nine data sources: The Ingenuity Pathways Knowledge Base (www.ingenuity.com), Human Protein Reference Database (www.hprd.org), GeneGo (www.genego.com), Jubilant (www.jubilantbiosys.com), NetPro (www.biobase.de), Prime Univ-Tokyo (prime.ontology.ims.u-tokyo.ac.jp), HumanCyc (humancyc.org), TransFac (28)
, www.gene-regulation.com) and Transcription Regulatory Regions Database (29)
.
Transfection analysis
HeLa cells were transfected with either a glucocorticoid reporter gene (TAT3-luc), or the same plasmid but with the glucocorticoid responsive elements deleted (
TAT3-luc), and either empty expression vector or cDNA expression vectors for the full-length p160 coactivator SRC2, all in the pcDNA3 expression vector (30
, 31)
, or an expression vector for BMPRII in pSPORT6, or empty pSPORT6. As the empty pSPORT6 vector gave the same result as empty pcDNA3, for clarity the results only include empty pcDNA3. An additional control was pSPORT6 containing an Fe65 expression cassette, and again results were the same as both empty pSPORT6 and empty pcDNA3. In addition all cells were cotransfected with a cytomegalovirus (CMV)-renilla vector to control for transfection efficiency (32)
. Cells were transfected with FUGENE 6 according to the manufacturers recommendation. Cells were divided posttransfection into wells for treatment with the indicated concentrations of dex. Cells were harvested at 24 h for luciferase and renilla assay.
| RESULTS |
|---|
|
|
|---|
|
|
|
Initial analysis of Dex-regulated genes
Just under 20% (2034/10981) of probe sets showed significant Gc regulation at the 5% level. Genes were identified that showed a strong (estimate of slope [ES]
0.2, P
0.05) dose-dependent regulation by dex; this approach found 98 genes (represented by 131 probe sets). Of these 34 were repressed and 64 induced by exposure to Gc. A number of these genes are known to be Gc regulated, either in PBMC or in a pulmonary epithelial cell line. These include cytokines (IL6 and IL8), enzymatic mediators (thrombospondin), immune regulatory transcription factors (IRF1), and survival genes (IL-7 receptor, and GILZ [TSC22D3]). The full list, functionally grouped, is presented in Table 2
.
|
Independent quantitation of gene expression by qRT-PCR
Seventeen Gc-regulated genes identified by the arrays was reanalyzed by using a quantitative PCR method. In all cases the direction of change was the same as predicted from the array analysis, but in some cases the fold change was greater, as expected from the increased dynamic range of the quantitative PCR method (data not shown).
Network analysis of identified genes
A bioinformatic network analysis was performed using the 98 Gc-regulated genes in order to help interpret the array data from a biological network/pathway perspective. Of these genes, 83 were mapped to a global network of literature-derived biological relationships between genes/proteins. A heuristic algorithm and scoring function approach was used to identify a subnetwork that included as many genes as possible from the query set but with few others. A statistically significant subnetwork was produced comprising 60 nodes, with 54 within the input set, and is shown in Fig. 3
. Six nodes were inferred, the GR itself (NR3C1), histone deacetylase 1 (HDAC1), CCAAT/enhancer binding protein, alpha (C/EBP
), mitogen-activated protein kinase 3 (MAPK3), progesterone receptor (PGR), and epidermal growth factor receptor (EGFR).
|
Analysis of Gc regulated gene responses between Gc-sensitive and resistant cohorts
Identification of the 98 most robustly Gc-regulated genes (represented by 131 Affymetrix probe sets) allowed assessment of differences in response between the Gc-sensitive and -resistant cohorts. It was predicted that the Gc-sensitive group would show an increased amplitude of response to Dex compared to the Gc-resistant group. Initial analysis showed that 21 (represented by 24 probe sets) of the 98 genes had a significant interaction of responder status with dex concentration with P < 0.01, as determined by ANOVA. The individual dose-response curves for each probe set are shown in the supplement (S1). To investigate this phenomenon more comprehensively, the entire 98 gene profile was subjected to cluster analysis. For analysis each probe set was considered as an independent variable, hence multiple probe sets from the same gene are shown where such sets emerged from the data. This approach yielded four meaningful clusters, depicted in Fig. 4
, and the genes listed by cluster in Table 3B
. The first cluster included genes that were up-regulated by Gc but that showed no interaction between dose and responder status. This group included GILZ (TSC22D3), a robustly regulated gene. Two of the 24 apparently interacting genes fell into this cluster (CXCR4 and IRF4). The second cluster, which again included genes up-regulated by Gc exposure, identified genes whose response to Gc differed according to responder status. This group included the majority (20 of the 24) of the probe sets originally identified. The third cluster included genes repressed by Gc, irrespective of responder status, and the fourth cluster genes repressed by Gc, but with nonresponders showing higher basal expression and a parallel dose-response compared to responders. One of the 24 genes fell into this cluster (TRAF1), and the remaining gene (SMARCA2) of the 24 was unclassified. Multiple entries for the same gene in Table 3b
result from each probe set being treated as a separate variable.
|
|
Discriminatory gene profile of the cohort
The ANOVA results showed that the difference between the Gc groups at the 100 nM Dex concentration level had the strongest signal. Therefore, this was the concentration selected for the Gc sensitivity discrimination analysis. The results of ANOVA analysis for 100 nM dose level are summarized in Fig. 5
. This plots the log2 of the fold-difference between the groups for each gene against the minus log10 P-value. Absolute values of 0.58, 1, and 1.32 on the x-axis relate to fold change values of 1.5, 2, and 2.5. A value of 2 or more on the y-axis relates to a P-value of 0.01 or less. Although many genes had very small P-values, most of them showed less than a 2-fold change in expression between the two groups.
|
Discriminant gene identification by predictive model generation
From the array data a predictive model was constructed using two independent methods, PLS-DA and bagging SVM. Using PLS-DA analysis with external CV, the lowest classification error (13%) was achieved by a model with 90 genes. The error rate changed in the range of 13% to 18% for models with 10 genes to 1000 genes. The 20 most frequently selected genes for the discrimination analysis are listed in Table 4
. The most frequently selected gene was selected 65 times out of the 100 runs, indicating no gene had significantly greater contribution to the classification than the other genes.
|
The lowest error rate for Bagging SVM with recursive feature elimination was around 15%. The top 3 most frequently selected genes were the same as those in Table 3
, namely regulator of G protein signaling 14, sialyltransferase 4B, and bone morphogenetic protein receptor, type II. Therefore, the SVM showed consistent results to those of the PLS-DA approach.
The location of the 20 genes shown in Table 4
is annotated on the volcano plot (Fig. 5)
. It is important to note that none of the genes in Table 4
was found to be glucocorticoid regulated.
Principal component analysis of the 20 best discriminant genes
Principal component analysis (PCA) was performed to show the natural separation of the two groups with different Gc sensitivity, based on their transcriptome results. The scores plot (Fig. 6
) shows the separation in the samples, i.e., the distinction between the patient groups at 100 nM dex. The first component (x-axis) accounts for 65% of the variation across the 20 genes. The variation is attributable to the difference between responders to steroid treatment and nonresponders across expression of the 20 genes (listed in Table 4
) at 100 nM.
|
Expression analysis of BMPRII
As BMPRII was found in the top three discriminant genes in two separate analyses, and its known signaling pathway suggested the possibility of functional interaction with the GR, its functional effect was sought. A simple transfection system was used in the Gc-sensitive HeLa cell line. SRC2 was used as a positive control, well known to augment Gc action on transactivation. BMPRII also augmented Gc effects on TAT3-luc expression, an effect dependent on having intact glucocorticoid response elements in the promoter region (Fig. 7
).
|
| DISCUSSION |
|---|
|
|
|---|
Previous work has identified a broad spectrum of genes, 21.5% of those represented on the array, regulated by Gc in ex vivo exposed human peripheral blood mononuclear cells (18)
. The investigators were able to group the regulated transcripts functionally and thus could propose a major role for Gc in up-regulating innate immune-type genes, and down-regulating adaptive immune-type genes. We proposed to use a similar approach to characterize differences in the transcriptome response between individuals, to identify core genes that were robustly and reproducibly regulated, and to determine whether the transcriptome response in vitro would correlate with the in vivo sensitivity to Gc.
Our subjects were stratified by cortisol response to a low-dose, overnight dex suppression test. We, and others, have previously used this approach in healthy volunteers and have shown that the response associates with genetic markers on the GR gene that have been independently found to associate with biological markers of altered Gc sensitivity (15
, 17
, 34
, 35)
. Importantly, no subject in the cohort had been subject to previous treatment with Gc, nor suffered with a Gc remediable condition.
The response to the overnight low-dose dex suppression test appears to be stable over time, supporting a robust underlying biological mechanism (20)
. In our cohort PDC did not correlate with any of the demographic or physiological variables measured. It is possible that we may have missed subtle correlations with, for example body mass index, body composition, or blood pressure, as a result of the cohort size. However, within the cohort was a wide range of dex response, and a very clear separation of the two groups selected for array comparison.
Initial analysis of the array data revealed an unexpected reduction in the number of Gc-regulated genes compared with that predicted from the previous published report (18)
. Some of the genes identified in our study had previously been found using either gene array studies or more traditional approaches to be regulated by Gc (36
, 37)
. Interestingly, some genes, including IFN stimulated gene (20 kD), were robustly stimulated by dex in our study but repressed in the earlier work (18)
. This discrepancy is likely the result of interindividual variability in Gc response, an effect minimized in our study by the number of subjects analyzed and the range of Gc concentrations tested. In addition, the state of cell activation may affect glucocorticoid responses, for example in the earlier study TLR-4 was repressed by dex in CD3-activated cells but stimulated by dex in resting cells (18)
. In our study TLR-4 was not found to be one of the most robustly regulated 98 genes. Importantly, we identified IFN-
and phosphoinositide 3-kinase regulatory subunit (p85) to be Gc regulated, and neither of these two genes were found in the previously published array studies of Gc action (18
, 36
37
38)
.
The table of glucocorticoid regulated genes reveals significant patterns of regulation involving cytokines, cytokine receptors, intracellular signaling mediators, and transcription factors. We built a glucocorticoid gene-regulated network based on literature-derived biological relationships. This analysis resulted in a significant subnetwork containing 60 nodes; 54 were from the query set as depicted in Fig. 3
. This method of network analysis has the very considerable advantage of graphically portraying relationships between the regulated genes and also identifying potential mediator proteins. Of these C/EBP
, MAPK3, and EGFR are the most clinically relevant.
The pattern of Gc response of the 98 genes was examined, and 24 showed differences in response dependent on the Gc sensitivity of the donors. To examine this in an unbiased manner, we applied cluster analysis to the 98-gene dose responses. In this way we could separate the genes into four clusters. Cluster 2 contained 28 genes, and cluster 4 contained 10 genes. Both clusters contained genes with differential expression profiles dependent on Gc sensitivity in vivo. These genes may provide important markers for future clinical studies of Gc response in vivo, a major current deficiency limiting mechanistic analysis of this complex trait.
Initial analysis showed a number of genes with a highly significant P-value for differences in expression between the glucocorticoid-sensitive and glucocorticoid-resistant groups. However, the magnitude of change observed for each of these genes was individually small. The conclusion from this observation is that many genes each make a small contribution to the difference in Gc response between the two groups. Further analysis was undertaken in an attempt to identify a model that could accurately predict the glucocorticoid phenotype of the donor individual. Two complementary approaches were taken; partial least-squares discriminate analysis and also support vector machines (SVM). Importantly, the top 3 most selected genes from the support vector machine were the same as those identified most frequently by the partial least squares discriminant analysis.
We presented the array data by using principle component analysis based on the expression of the 20 most frequently identified genes (Table 4)
. This analysis shows that the two groups can be separated in the first dimension.
Previous microarray analysis of Gc sensitivity has been confined to an asthma cohort in Iceland. This study identified a set of 15 genes whose expression predicted Gc sensitivity with 84% accuracy. In particular basal expression of the RelA component of NFkB was found to be the best predictor of Gc response. We did not find a similar gene set in our study, however, our starting cells were different, being an expanded T lymphoblast population. Also, we excluded subjects with asthma to avoid disease activity effects (19)
.
Analysis of discriminating gene function was informative. A number of these genes are implicated in cell signaling and are therefore good candidates for explaining variation in cell phenotype, and so directly or indirectly glucocorticoid sensitivity. The bone morphogenetic protein (BMP) receptor type II is the ligand-binding partner of the signaling bone morphogenic protein receptor type I. The BMPs act on bone, brain, and in development and are members of the TGFbeta superfamily (39)
. The regulator of G protein signaling 14 (RGS14) gene product is implicated in down-regulation of signaling by heterotrimeric G proteins, and importantly is expressed in brain and spleen (40
, 41)
. The human enhancer of filamentation (HEF1), also known as CRK associated substrate-related protein, or NEDD9 was identified. HEF1 is subject to tyrosine phosphorylation in response to ligation of beta 1 integrin (42
, 43)
.
RAB3 GTPase is the gene associated with Warburg micro syndrome 1, which is characterized by ocular and neurodevelopmental defects and by hypothalamic hypogenitalism (44)
.
The discriminate genes may be markers of differential glucocorticoid sensitivity, but some may also play a causative role in modulating such sensitivity. A good candidate for such an effect is BMP receptor type II (BMPRII). BMPR signaling acts through SMAD proteins, and evidence indicates interaction between Gc and BMPR signaling pathways (45
46
47)
. Definitive evidence for such a role requires functional analysis and so, as proof of principle, we overexpressed BMPRII in HeLa cells. A clear enhancement was found of Gc transactivation of the reporter gene with BMPRII, and also with SRC2, used as a positive control. Deletion of the GR binding sites (
TAT3-luc) abolished the Gc effect, as expected, but also showed that the BMPRII did not activate the minimal promoter. This finding supports the hypothesis that Gc sensitivity modulating genes can be identified by screening using microarray technology in well-characterized human volunteer cohorts.
In summary, we have demonstrated the feasibility of using a gene array approach to profile glucocorticoid regulated genes in individuals with defined in vivo glucocorticoid sensitivity. The proportion of the genome that we have identified to be regulated by glucocorticoids is an order-of-magnitude lower than previously published estimates. It will be important to analyze the gene profile of peripheral blood mononuclear cells from patients with inflammatory diseases with defined responses to exogenous glucocorticoid administration and also individuals with metabolic diseases associated with altered glucocorticoid bioactivity, for example type 2 diabetes, central obesity, and apparently idiopathic osteoporosis.
Received for publication September 12, 2006. Accepted for publication September 22, 2006.
| REFERENCES |
|---|
|
|
|---|
This article has been cited by other articles:
![]() |
R. Newton and N. S. Holden Separating Transrepression and Transactivation: A Distressing Divorce for the Glucocorticoid Receptor? Mol. Pharmacol., October 1, 2007; 72(4): 799 - 809. [Abstract] [Full Text] [PDF] |
||||
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |