Expression quantitative trait loci and genetic regulatory network analysis reveals that Gabra2 is involved in stress responses in the mouse

Abstract

Previous studies have revealed that the subunit α 2 (Gabra2) of the γ-aminobutyric acid receptor plays a critical role in the stress response. However, little is known about the gentetic regulatory network for Gabra2 and the stress response. We combined gene expression microarray analysis and quantitative trait loci (QTL) mapping to characterize the genetic regulatory network for Gabra2 expression in the hippocampus of BXD recombinant inbred (RI) mice. Our analysis found that the expression level of Gabra2 exhibited much variation in the hippocampus across the BXD RI strains and between the parental strains, C57BL/6J, and DBA/2J. Expression QTL (eQTL) mapping showed three microarray probe sets of Gabra2 to have highly significant linkage likelihood ratio statistic (LRS) scores. Gene co-regulatory network analysis showed that 10 genes, including Gria3, Chka, Drd3, Homer1, Grik2, Odz4, Prkag2, Grm5, Gabrb1, and Nlgn1 are directly or indirectly associated with stress responses. Eleven genes were implicated as Gabra2 downstream genes through mapping joint modulation. The genetical genomics approach demonstrates the importance and the potential power of the eQTL studies in identifying genetic regulatory networks that contribute to complex traits, such as stress responses.

Keywords::

• Anxiety
• expression QTL
• Gabra2
• genetical genomics
• genetic regulatory network
• hippocampus

Introduction

γ-Aminobutyric acid (GABA) is the primary inhibitory neurotransmitter in the central nervous system (CNS). Nineteen GABA_A receptor subunits have been identified in mammals by cloning, including six α, three β, three γ, one δ, one ε, one π, one θ and three ρ (Barnard et al. 1998; Bonnert et al. 1999). GABA_A receptors are modulated by various clinically important compounds including ethanol, benzodiazepines, and barbiturates. Previous studies have shown that the GABA_A receptor α 2 (Gabra2) gene is related to the risk for alcohol dependence (Edenberg et al. 2004; Soyka et al. 2008). Mice with a Gabra2 knock-in point mutation are insensitive to anxiolytic effects of benzodiazepines (Low et al. 2000). More recently, Enoch et al. (2006) reported that a Gabra2 haplotype alters risk for alcoholism via the mediating action of anxiety, which generally indicates stress.

The hippocampus is a crucial brain structure involved in spatial memory formation (Kandel 2001). The effects of stress on memory have been well studied in various species (Bremner and Narayan 1998). The hippocampus is enriched with receptors for corticosterone which is the principal glucocorticoid released by the rodent adrenal cortex in response to stress (Freund and Buzsaki 1996). In addition, stress has been shown to affect synaptic plasticity, dendritic morphology, neurogenesis, and neurotoxicity within the hippocampus (Bremner et al. 1993; Kim and Wu 1996; Gould et al. 1998). The Gabra2 gene has been shown to be highly expressed in the hippocampus and hypothalamus (Laurie et al. 1992).

The stress response is a complex trait that is influenced by both genetic and environmental factors (Dumas et al. 2000). In order to study complex traits, several mouse genetic models have been developed, such as the BXD recombination inbred (RI) panel. The parental strains of BXD mice, C57BL/6J (B6), and DBA/2J (D2) are strikingly different in their stress responses (Tarricone et al. 1995). Quantitative trait locus (QTL) analysis has identified numerous chromosomal regions that significantly influence the stress response, such as Ssrq1–Ssrq11 (Thifault et al. 2008). Although molecular techniques have been used to study the stress response (Steckler 2001), the genetic regulatory network remains unclear.

The genetical genomics approach, combining global expression data and linkage analysis, has been used to reconstruct genetic networks by integrating gene expression and linkage analysis. Since Jansen and Nap (2001) proposed the genetical genomics concept, it has been successfully applied in several organisms (Brem et al. 2002; Schadt et al. 2003; Yvert et al. 2003; Bystrykh et al. 2005; Chesler et al. 2005; West et al. 2007). Identification of an expression QTL (eQTL) for a gene suggests that the gene transcript levels are regulated by genetic factors. Recent advances in the genetical genomics approach are promising in the identification of causal genes influencing complex traits, and ultimately the pathways and gene networks underlying these traits (Mehrabian et al. 2005). In this study, we explore the molecular genetic mechanisms that may contribute to the stress response through the construction of a genetic regulatory network and the identification of downstream genes for Gabra2.

Materials and methods

BXD RI mice

The BXD RI panel, originally bred by Taylor in the late 1970s, was recently expanded and now includes 80 RI strains (Williams et al. 2001b; Peirce et al. 2004). In both the original and extended BXD RI panels, the progenitor strains, B6 and DBA/2J (D2), were mated, followed by more than 20 generations of sibling-matings. This resulted in a panel of inbred strains with fixed genotypes at each locus, with the parental B6 and D2 alleles segregated among the strains. For all of the described experiments, mice were weaned at 25 days of age and housed in same-sex cages, 2–5 mice per cage until the time of testing. Mice had free access to standard laboratory chow and water. Mice were maintained on a 12:12 h light/dark cycle (lights on from 6 a.m. to 6 p.m.). Room temperature ranged from 20 to 24° C. The study was approved by the University Review Boards and was performed in accordance with animal model research guidelines.

RNA isolation

Naïve mice of 67 BXD RI strains and their parental strains (B6 and D2) were analyzed and almost all strains were represented by matched male and female samples. The mice were killed at 60 days of age by cervical dislocation without anaesthesia. The brain was quickly removed and stored in RNAlater (Ambion, Inc., Austin, TX, USA) over night at 4°C, then kept at − 80°C until use. The brain was dissected into discrete brain structures including the hippocampus, neocortex, striatum, ventral midbrain, and hypothalamus. Brain dissections have been described elsewhere in detail (Lu et al. 2001; Williams et al. 2001a; Peirce et al. 2003). The extraction of RNA from hippocampal tissues was carried out with the RNA STAT-60 reagent (Tel-Test, Inc., Friedenswood, TX, USA). Tissue samples were washed three times in phosphate-buffered saline (GibCo BRL, Grand Island, NY, USA) to remove blood contamination. Tissue (100 mg) was homogenized in 1 ml of RNA STAT-60 reagent. Following homogenization, the homogenate was stored for 5 min at room temperature. Next, 0.2 ml of chloroform per 1 ml of RNA STAT-60 was added and the mix was vigorously shaken for 15 s and centrifuged at 12,000 g for 15 min at 4°C. After centrifugation the aqueous phase was transferred to a fresh tube and mixed with 0.5 ml of isopropanol per 1 ml of RNA STAT-60 as used for the homogenization. The precipitate was washed twice in 75% ethanol, air-dried and re-diluted in Nuclease-free water (Ambion, Inc.). The purity of the extracted RNA was determined with a NanoDrop spectrophotometer (NanoDrop Technologies, Inc., Rockland, NC, USA). The 260/280 nm ratio of the samples was within the range of 1.9–2.1. RNA integrity was assessed on an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA), and the RNA integrity number values were required to be >8.0.

Microarray analysis

RNA samples were amplified using the TotalPrep RNA amplification kit. The TotalPrep RNA amplification kit uses T7 oligo(dT) primers to produce higher yields of first strand of cDNA followed by a second strand synthesis. cRNA was generated from the double-stranded cDNA by in vitro transcription by T7 RNA polymerase. Then, total cRNA sample was hybridized to the Affymetrix microarray (MOE430V2). Post-hybridization staining and washing were performed according to the manufacturer's protocols. Finally, the arrays were scanned and images were quantified. Normalization was performed by the robust multichip average (RMA) method. All normalizations were done using default analysis parameters. We generated probe set data using RMA (Gautier et al. 2004), obtained the base 2 log of each probe set and standardized using Z scores. We doubled the Z scores and added eight to produce a set of Z scores with a mean of 8, a variance of 4 and a standard deviation (SD) of 2. The advantage of this modified Z score is that a twofold difference in expression level corresponds approximately to a 1-unit difference. Expression levels below six are usually close to background noise levels (Chesler et al. 2005).

eQTL mapping

QTL mapping (single marker regression and interval mapping) was performed using web-based complex trait analysis (www.genenetwork.org) which uses QTL reaper software (Chesler et al. 2004). Single marker regression across all chromosomes was performed where a hypothetical QTL was evaluated at the location of 8222 informative markers (http://www.genenetwork.org/genotypes/BXD.geno). Interval mapping evaluates potential QTL at regular intervals and estimates the significance at each location with a graphical representation of the likelihood ratio statistic (LRS). A permutation test establishes genome-wide significance criteria of 5% for the trait.

Genetic network construction

We computed Pearson product–moment correlations of strain means for each pair of probe sets on the array. Pearson product–moment correlations can be computed using GeneNetwork (www.genenetwork.org). The network structure was drawn using the software Cytoscape 2.1 (http://www.cytoscape.org).

Results

Expression variation and QTL mapping for Gabra2

DNA microarrays were used to measure the abundance of mRNA that is complementary to the DNA in the array substrate. Recent studies have shown that gene expression can be influenced by the single nucleotide polymorphisms (SNPs) on the probe affecting the hybridization efficiency which could cause artifact cis-QTLs, and by the position of probe on the target sequence (Li et al. 2005). There are seven probe sets in the Affymetrix MOE430V2 array for Gabra2 that were used to analyze its mRNA expression in the hippocampus (Table I). Four of them (1420299_at, 1420300_at, 1446805_at, and 1449807_x_at) were excluded because of low expression and targeting at non-coding regions of Gabra2. We further examined the location of SNPs in the remaining three probe sets and found that no SNPs in the probe sets were detected. Their expression in BXD RI strains averaged 11.025, and ranged from a low of 10.744 to a high of 11.503 (Table I). The expression level of these three probe sets was twofold higher in D2 than that in B6 strains (Figure 1; P < 0.01, t-test).

Table I.  Summary of probe sets for Gabra2 in the Affymetrix MOE430V2 array.

Probe set name	Expression level^*	Location of probe set^†
1420299_at	5.423 ± 0.095	3′ Inter-genic region
1420300_at	5.091 ± 0.259	3′ Inter-genic region
1421738_at	10.828 ± 0.451	Last exon and proximal 3′ UTR
1443865_at	10.744 ± 0.536	Far 3′ UTR
1446805_at	9.324 ± 0.455	Intron 3
1449807_x_at	5.783 ± 0.146	3′ Inter-genic region
1455444_at	11.503 ± 0.527	Far 3′ UTR

^*

The expression level is the mean for the hippocampus 67 BXD mouse strains, shown as mean ± SD

^†

Location of probe set indicates the location of the sequences that the probe set maps to in the Gabra2 gene.

Figure 1  The normalized relative expression level of three probe sets for Gabra2 mRNA between parental mouse strains, B6, and D2. The expression levels of all three probe sets are higher in DBA/2J than in B6 mice (P < 0.01).

Interval mapping was performed to define the chromosome regions associated with expression differences for these three probe sets. The analysis used BXD RI panel and 8222 markers as genotypes through GeneNetwork (www.genenetwork.org). We detected the same interval between the marker rs13478294 and rs3722245 for all three probe sets (1421738_at, 1443865_at, and 1455444_at). These probe sets showed highly significant linkage (LRS) scores of 33.6 for 1455444_at, 32.1 for 1443865_at, and 26.6 for 1421738_at, respectively. All three probe set LRS scores exceeded the significance threshold level (P = 0.01). DBA/2J (D2) alleles at this interval had an additive effect of 0.3. These three transcripts were associated with a QTL mapping to within ± 3 Mb of the gene itself, and so are indicative of a cis-acting QTL (Figure 2).

Figure 2  LRS scores for Gabra2 mRNA eQTL in the BXD mouse hippocampus. The solid blue line indicates LRS across the genome. A positive additive coefficient (green line) indicates that DBA/2J alleles increase trait values, whereas a negative additive coefficient (red line) indicates that B6 alleles increase trait values. The y-axis denotes the level of linkage and the x-axis denotes genetic markers along the different chromosomes (autosomes 1–19 and X).

Genes co-regulated with Gabra2

One probe set, 1421738_at, was selected for correlation and network analysis, because this probe set not only had the highest LRS (33.6) in eQTL mapping results but also is located on the last exon and proximal 3′ UTR. Using Pearson correlations between 1421738_at and whole mouse genome transcripts based on the gene expression abundance from 67 BXD strains, the top 500 genes (Supplemental Table I) whose expression levels were highly correlated with 1421738_at were selected (r>0.53; P < 0.001). To reduce this list of 500 genes still further while retaining a high chance of biological plausibility, these genes were further subsetted based on their literature correlation with Gabra2. The average literature correlation of the top 29 genes was 0.49, suggesting that an association between those 29 genes has been reported in a large number of studies. Subsequently, we treated the expression level of correlated genes as the perturbation to construct the gene co-regulatory network for these 29 genes (Figure 3).

Figure 3  Gene co-regulatory network for Gabra2 was constructed using genome-scale gene expression in hippocampus of the BXD RI mice. Genes with Pearson correlation coefficient >0.4 are shown linked to others.

To examine whether these 29 genes are related with stress responses, we used the Chilibot literature searching database (http://www.chilibot.net/index.html). Ten genes were identified when we searched the literature correlation with stress, including Gria3, Chka, Drd3, Homer1, Grik2, Odz4, Prkag2, Grm5, Gabrb1 and Nlgn1. Another GABA_A receptor subunit, Gabrb1, is particularly noteworthy because it not only directly interacts with Gabra2 in the regulatory network, but also is in a gene cluster with Gabra2 and Gabrg1 on chromosome 5. We could not detect any significant or suggestive cis-QTL for Gabrb1 and Gabrg1. This suggests that Gabrb1 may be regulated by Gabra2 at this locus. Glutamatergic mechanisms are thought to be involved in stress-induced alterations of brain function, especially in the hippocampus (Schwendt and Jezova 2000). This network includes two ionotropic glutamate receptor genes (Gria3 and Grik2) and one metabotropic glutamate receptor (Grm5). In addition, Homer1 is a critical element in glutamatergic signaling and the mRNA or protein expression can be regulated by stress during development (Ary et al. 2007). Dopamine D3 receptor (Drd3) plays an important role in mediating stress induced reinstatement (Xi et al. 2004). Previous study has demonstrated that Odz4 is induced by endoplasmic reticulum (ER) stress response in the C/EBP homologous protein pathway (Wang et al. 1998). Choline kinase α (Chka) gene induction seems to be associated with certain cell stress or cell defence responses (Janardhan et al. 2006). Neuroligin1 (Nlgn1) which acts as a postsynaptic cell surface molecule and binding partner for Neurexins, can trigger presynaptic synaptogenesis under stress (Scheiffele et al. 2000). Study in humans showed that a heterozygous mutation in PRKAG2 can prevent activation of the heterotrimer by metabolic stress in intact cells (Akman et al. 2007). Collectively, the genes that are presented in this co-regulated network are directly or indirectly associated with stress responses.

Downstream genes affected by Gabra2

We mapped the top 500 correlated genes (Supplemental Table I) using the ClusterMap tool at www.GeneNetwork.org. Forty-seven transcripts/genes with trans-regulatory QTL bands reside in the locus of chromosome 5 where Gabra2 is located (Figure 4). Gabra2 is a cis-regulated gene at this locus, making it a candidate genetic modifier for trans-regulatory transcripts.

Figure 4  Cluster QTL map on chromosome 5 for the 47 co-correlates of Gabra2. Notice the strong hotspot on chromosome 5 indicating most of these genes share a common trans-QTL. Colors are used to encode the significance of the LRS (brightness), as well as the additive effect of polarity of alleles at each marker (red or blue color). The red box denotes the locus of Gabra2 on chromosome 5.

In order to exclude correlated transcripts due to linkage disequilibrium, we next examined the rank Spearman correlation between Gabra2 and those 47 trans-regulatory genes using 28 different tissue expression data sets from GeneNetwork web site (www.genenetwork.org). Eleven genes, Elf2, Herc1, Slc13a3, Slc4a7, Rod1, Pcdhga6, Dlg2, Serpinb9d, Art3, Fus, and Svs5 were statistically significantly correlated with Gabra2 (Figure 5; P < 0.01). This suggests that 11 genes are regulated by Gabra2 as a consequence of joint modulation and high tissue expression correlation with Gabra2. Two of them, Elf2 and Fus, are associated with stress responses. The Elf2 valosin–containing protein (VCP) pathway might be important for cell survival and proliferation under cytokine stress (Zhang et al. 2007). Fus mutants are hypersensitive to oxidative stress and rapidly lose viability in the presence of hydrogen peroxide (Macvanin et al. 2004). These two genes are most likely to be potential downstream genes of Gabra2. This result also indicated that Gabra2 may play a pivotal role in stress responses through these 11 downstream genes.

Figure 5  Genes statistically significantly correlated with Gabra2 using expression data sets for 28 different tissues from www.genenetwork.org web site.

Discussion

The genetical genomics approach has been demonstrated as a powerful tool for identification of candidate genes and construction of regulatory genetic networks. In this study, we performed eQTL analysis for Gabra2 using BXD RI mouse strains and found it is a cis-regulated gene. We further demonstrated that the evidently co-regulated genes are known to be involved in different aspects of stress responses.

Several regions of the CNS that are most important in emotion include the thalamus, hypothalamus, the limbic system and cerebral cortex. The hippocampus is one part of the limbic system and Papez's circuit, which are important for mood modification. Moreover, the hippocampus is a brain area primarily affected by stress (Sala et al. 2004). The SymAtlas database (http://symatlas.gnf.org/SymAtlas/symform) showed that Gabra2 (probe set: 207014_at) is highly expressed in brain specific regions, particularly in hippocampus, amygdala, hypothalamus, and cerebral cortex (Figure 6), which are closely involved in stress responses. Many genes in the co-regulated network (Figure 3), such as Gria3, Grm5, and Gabrb1, are also highly expressed in hippocampus.

Figure 6  The expression level of Gabra2 (probe set: 207014_at) among 61 tissues from the SymAtlas database (http://symatlas.gnf.org/SymAtlas/symform). The y-axis denotes the relative expression level and the x-axis denotes different tissues. Gabra2 is highly expressed in the hippocampus, amygdala, hypothalamus, and cerebral cortex.

Amino acid changes within the coding sequence that affect the level of protein produced may lead to a change in gene expression. Such local variation affecting the protein product would act in trans mode to alter the expression of both alleles (Doss et al. 2005). Alternatively, polymorphisms in the promoter region that affect chromatin structure or transcription factor binding sites, or polymorphisms in the 3′ untranslated region (3′ UTR) that affect mRNA stability, would be expected to act in cis mode, altering the abundance of the transcript in an allele-specific manner in a diploid (Ronald et al. 2005). In this study, Gabra2 was identified as a cis-regulated gene (Figure 2). Thus, we attempted to find polymorphism in the promoter region and 3′ UTR in Gabra2 by searching known SNPs and Indels. However, we were unable to detect single nucleotide and insertion or deletion polymorphism in the 3′ UTR and upstream 5000 bp of Gabra2. By contrast, three SNPs (mCV23399336, mCV24189554, and mCV23399622) were detected in intron 6 and intron 7. Therefore, expression regulation of Gabra2 is likely to be driven by the SNPs in the introns (Tsukada et al. 2006).

We identified a group of genes closely related to stress responses through transcript and literature correlation analysis in BXD RI mouse hippocampus. A genetic regulatory network for Gabra2 and its potential role in stress could be initially constructed with this group of genes, indicating a component of the molecular architecture of the stress response. Ten of these genes co-regulated with Gabra2, including Gria3, Chka, Drd3, Homer1, Grik2, Odz4, Prkag2, Grm5, Gabrb1, and Nlgn1, have known links to stress responses (Figure 5). Although, from literature searches or correlations, the remaining genes in the co-regulated network are not directly associated with the stress response, they represent a biologically plausible novel extension of the network.

Acknowledgements

This research was founded by National Natural Science Foundation of China (Nos 30771200, 30770666, and 30700517), key project of Jiangsu Natural Science Foundation (No. BK2007703) and Jiangsu Natural Science Foundation (No. BK2007065) and National Institute on Alcohol Abuse and Alcoholism (NIAAA) of USA (No. U01-AA014425). We thank David C. Airey for helpful comments on manuscript drafts.

Declaration of interest: The authors report no conflicts of interest. The authors alone are responsible for the content and writing of the paper.