Expression of Dopamine-Related Genes in Four Human Brain Regions

A better understanding of dopaminergic gene expression will inform future treatment options for many different neurologic and psychiatric conditions. Here, we utilized the National Institutes of Health’s Genotype-Tissue Expression project (GTEx) dataset to investigate genotype by expression associations in seven dopamine pathway genes (ANKK1, DBH, DRD1, DRD2, DRD3, DRD5, and SLC6A3) in and across four human brain tissues (prefrontal cortex, nucleus accumbens, substantia nigra, and hippocampus). We found that age alters expression of DRD1 in the nucleus accumbens and prefrontal cortex, DRD3 in the nucleus accumbens, and DRD5 in the hippocampus and prefrontal cortex. Sex was associated with expression of DRD5 in substantia nigra and hippocampus, and SLC6A3 in substantia nigra. We found that three linkage disequilibrium blocks of SNPs, all located in DRD2, were associated with alterations in expression across all four tissues. These demographic characteristic associations and these variants should be further investigated for use in screening, diagnosis, and future treatment of neurological and psychiatric conditions.


Introduction
Multiple neurological and psychiatric diseases result from alterations in the production, degradation, or improper signaling of dopamine in the brain. In order to develop better therapeutic targets for these diseases, we must first better understand the organization and connectivity of dopaminergic circuitry in the brain. Within the mesocortical and mesolimbic pathways, the intricate balance of dopaminergic gene expression controls much of the reward and motivation machinery in the brain. Such a balance in the nigrostrial pathway influences movement and procedural learning, while the hippocampus is involved in long-term memory, spatial processing, and navigation. An imbalance of dopamine expression in any of these areas can create any number of heterogenous behavioral outcomes, including many common neurological and psychiatric conditions. Such an imbalance is often due to genotypic variations that alter expression. Although much work has been done in this area using animal models, it is not yet clearly understood in humans how dopaminergic expression (a) changes across the regions of such circuitry and (b) varies by given genotype.
One way to gain a better understanding of such information is through the use of the National Institutes of Health's Genotype-Tissue Expression project (GTEx) database. The GTEx dataset includes genotype, expression, and basic demographic data for nearly 1000 donors of up to 53 different tissues. Within the GTEx database, a web-based portal can be used to check sample and tissue distributions, to browse and search data by various parameters, and to query expression Quantitative Trait Loci (eQTL). While the portal information is useful for general inquiry, raw data access (genotype, expression, and a limited set of demographic and phenotypic data) can be granted for more detailed projects by application through the National Center for Biotechnology Information's Database of Genotypes and Phenotypes (dbGAP; https://dbgap.ncbi.nlm.nih.gov/). Here, we applied and obtained such access to the controlled access dataset to investigate eQTL in seven dopamine pathway genes (including ankyrin repeat and kinase domain containing 1  , and dopamine transporter 1/solute carrier family 6 member 3 [SLC6A3]) across four human brain tissues, selected for their importance in dopaminergic reward, movement, and memory processes: prefrontal cortex, nucleus accumbens, substantia nigra, and hippocampus. A wide variety of psychiatric and neurological disease processes are associated with these genes and brain regions. In particular, single nucleotide polymorphisms (SNPs) in ANKK1 and DRD2 have been associated with addictive behaviors [1,2] and schizophrenia [3,4], while variants in DBH, DRD1, and DRD3 have been associated with Alzheimer's Disease, Parkinson's Disease, and cognitive performance in a number of different clinical populations [5][6][7]. Additionally, variants in DRD5 and SLC6A3 have been associated with attention deficit hyperactivity disorder (ADHD) [8]. The purpose of this work was to investigate genotype by expression changes across these tissues, with the hope that this exploratory inquiry will build a foundation for future work in finding disease risk variants and potential targets for therapeutic intervention. To our knowledge, this work is the first of its kind in summing the pathway effects of genetic variation for expression in dopaminergic genes for human brain tissue samples.

Materials and Methods
After receiving access to the controlled dataset and obtaining approval of this protocol under the designation of non-human subject research by the Institutional Review Board of the University of Tennessee Health Science Center (#19-06937-NHSR), our team downloaded expression data (gene read counts or Transcript Per Kilobase Million [TPMs]), sample attributes, and subject demographic and phenotype data from the GTEx portal. The results reported here are based on GTEx data phs000424.GTEx.v8.p2.c1.GRU with serial reformatting, genotype data GTEx_Analysis_2017-06-05_v8_WholeGenomeSeq_838Indiv_Analysis_Freeze.vcf, and expression data phASER_GTEx_v8_matrix.gw_phased.txt from phe000037.v1.GTEx_v8_RNAseq. expression-data-matrixfmt.c1. We are interested in and limited our inquiries to seven dopaminergic genes (ANKK1, DBH, DRD1, DRD2, DRD3, DRD5, and SLC6A3) across the following four tissues as described in the database: brain-nucleus accumbens (basal ganglia), brain-substantia nigra, brain-hippocampus, and brain-frontal cortex (BA9; representing prefrontal cortex), with the purpose to identify SNPs associated with corresponding gene expression across the four brain tissues. The genotype variants within 10 kilobases up-and down-stream to the seven genes were extracted and coded as the number of alternative alleles (0, 1, and 2). The total RNAseq expression values of the seven genes in dopamine pathway were also extracted for each of these four tissues. Any expression value of 0 was excluded from any association analysis, as this could mean no expression, or the expression level was under the detection threshold level. The association of each gene expression in a brain tissue was tested using rank-based Kruskal-Wallis test with sex and race. The association of expression with number of alternative alleles of variants annotated to the corresponding gene in dopamine pathway was measured by Spearman correlation, a rank-based conservative approach. The magnitude of the rank-based correlation reflects the degree with which the SNP effects expression between brain tissues. False discovery rate was controlled by Benjamini-Hochberg method, with q value cutoff of 0.2 in each of the 4 brain tissues [9]. As only the variants annotated to the seven genes in dopamine pathway were interrogated, the selection criteria were less stringent than would be required by a whole genome approach. All analyses were performed in R-3.6.1 [10].

Single Tissue Genotype x Expression Analyses
A total of 6602 variants were annotated to the selected seven genes in the dopamine pathway and within 10 kilobases up-and down-stream of the gene. The variants with minor allele frequency (MAF) less than 0.01 or call rate less than 0.9 were excluded from analysis. In total, 1367 variants were interrogated in the association analysis with corresponding gene expression. Table 1 demonstrates the distribution of significant associations of variants by tissue for each gene at false discovery rate (FDR) of 0.2. Information on the results for each individual variant is provided in Supplementary Table S1. The overlap of significant variants for expression between tissues is demonstrated in the Venn diagram of Figure 1. A total of 6602 variants were annotated to the selected seven genes in the dopamine pathway and within 10 kilobases up-and down-stream of the gene. The variants with minor allele frequency (MAF) less than 0.01 or call rate less than 0.9 were excluded from analysis. In total, 1367 variants were interrogated in the association analysis with corresponding gene expression. Table 1 demonstrates the distribution of significant associations of variants by tissue for each gene at false discovery rate (FDR) of 0.2. Information on the results for each individual variant is provided in Supplementary  Table S1. The overlap of significant variants for expression between tissues is demonstrated in the Venn diagram of Figure 1.

Single Tissue Expression Analyses by Demographic and Phenotypic Characteristics
For each individual tissue type, we next investigated general demographic and phenotype characteristics (such as age, sex, race, and body mass index [BMI]) for associations with expression to uncover patterns of variation in gene expression across tissues. The association of each gene with at least 30% non-zero expression in a brain tissue with sex and race was tested using Kruskal-Wallis test. The association with age and BMI was measured by Spearman correlation. A limitation of this dataset is that few subjects had full genetic and RNAseq data for all four tissues, and so sample size and demographic associations with genotype x expression varied between 51 and 177 individuals, for a total of 249 unique persons in the dataset. Demographic characteristics for this group of subjects can be found in Supplementary Table S2, but in brief this group was 71.89% male and 90.73% white, had a median age of 61 years old (range: 55-66 years), and had a median BMI of 27.44 (range: 24.99-31.32). We found that age was inversely associated with expression within nucleus accumbens for DRD1 (n = 145; r = −0.1977; p = 0.017) and DRD3 (n = 136; r = −0.1879; p = 0.029), in hippocampus for DRD5 (n = 89; r = −0.2648; p = 0.012), and in prefrontal cortex for DRD1 (n = 123; r = −0.2212; p = 0.014) and DRD5 (n = 103; r = −0.3179; p = 0.001; these associations and demographic information can be found in Supplementary Table S3). Increased expression was associated with female sex within substantia nigra for DRD5 (n = 51; p = 0.018) and SLC6A3 (n = 100; p = 0.044), but with male sex in hippocampus for DRD5 (n = 89; p = 0.05; Supplementary Table S3). Although race was found to be significantly associated with expression for DRD3 in nucleus accumbens (p = 0.0002), DRD2 in substantia nigra (p = 0.047), and DBH in prefrontal cortex (n = 161; p = 0.018), the majority (>90%) of each sample was Caucasian, limiting the conclusions that could be drawn from this information. There were particularly limited numbers of samples collected from Asian subjects. Within substantia nigra, expression was associated with BMI for DRD2 (n = 87; p = 0.05) but there were no other associations in other tissues or genes with this characteristic.

Four Tissue Genotype x Expression Analyses
Next, we are interested in the SNPs associated with gene expression across all four brain tissues. Out of 1367 variants interrogated, there are 28 variants in three linkage disequilibrium (LD) blocks that were associated with alterations in expression of DRD2 across all four tissues (see the quad overlapping portion of Venn diagram in Figure 1, LD structure for the three haplotype blocks is presented in Supplementary Figure S1) at FDR of 0.2 cutoff. The first LD block contains intron variant rs7103679 in DRD2. The second LD block contains two variants: rs6277 and rs10891549, with correlation of 0.99, and average correlation of −0.42 with the third LD block, which contains 25 variants with smallest correlation of 0.93. Detail on these associations in each tissue is provided in Table 2, with shading indicating respective LD blocks.
The variants in LD block 1 and 3 are all intronic or fall outside of the coding region of DRD2. For rs10891549, the reference allele incorporates a thymine and the alternate allele incorporates a cytosine at chr11:113407725 (GRCh38.p12). The minor allele frequency (MAF) from 1000 genomes estimated 0.26 for the alternate allele, although the variant is less common in African and Asian populations with MAF < 0.13 (information from dbSNP). In rs6277, a synoynmous variant, the incorporation of the amino acid proline does not change whether the reference allele guanine or alternate allele adenine is present at chr11:113412737 (GRCh38.p12). Similarly to rs10891549, MAF from 1000 genomes is estimated at~0.24 for the alternate allele, while the variant is less common in African and Asian populations (MAF < 0.07). For full detail on all DRD2 SNPs in each tissue, please see Supplementary  Table S1. As there was no significant association between total RNAseq expressions and sex or race in each tissue, the association tests of expression with variants was not adjusted for demographics or population structure for the target pathway. Detail on the effect of the genotype x expression associations for each variant and tissue in our three LD blocks is detailed in Table 2, whereas a negative correlation indicates that the number of alternate alleles is associated with lower expression of DRD2 for that tissue (i.e., the number of C alleles in rs7103679 is negatively associated with DRD2 expression in all four tissues, with the association in hippocampus the strongest correlation (r = −0.339, p = 0.0002) and the weakest correlation occurring in prefrontal cortex (r = −0.252, p = 0.0048) for these four tissues).

Discussion
Our work is some of the first of its kind to investigate dopaminergic genotype by expression changes in human brain tissue samples. Our results identified demographic and phenotypic characteristics associated with alterations in gene expression, and also identified individual variants within each gene and within each tissue that are associated with these expression changes. Such information may provide additional evidence for the association of various disease states with particular demographic and phenotypic characteristics, and may provide a foundation for further inquiry into the circuitry of gene expression alterations in each disease state. Furthermore, our results demonstrate that there are three haplotype blocks, for a total of 28 DRD2 SNPs that are associated with expression across our four tissues. In regard to demographic and phenotypic characteristics, we found associations of gene expression with age, sex, and BMI. Each is discussed in turn.

Age
We noted that age inversely alters expression of DRD1 in the nucleus accumbens and prefrontal cortex, DRD3 in the nucleus accumbens, and of DRD5 in the hippocampus and prefrontal cortex. Clinically, such findings hint at the potential for age-related changes in expression to be associated with the development or worsening of neurological and psychiatric conditions; older age is a common risk factor for such changes. However, very little work has been done in human brain tissues to determine the mechanism by which such associations occur. We do know that upregulation of DRD1, DRD3, and DRD5 expression in blood, along with alterations in the expression of other dopaminergic signaling genes, has been associated with lower resilience to stressful events and other psychological parameters, including personality, depression, anxiety, and intelligence, in age, sex, and race-matched individuals [11]. We also know that increasing binding of D1 receptors has been associated with younger age in suicide [12], and DRD5 variants have been associated with age of onset in pediatric attention-deficit/hyperactivity disorder (ADHD) patients [13]. Such information further supports our findings.
One study investigated expression in DRD1, DRD2, and DRD5 (among other dopaminergically relevant genes) in dorsolateral prefrontal cortex samples (Brodmann's area 46) taken from 69 patients, grouped into seven different age groups (considering categorically neonates, infants, toddlers, young adults, and adults) [14]. The sample from this study was significantly younger than our sample at hand (our sample age range was from 55-66 years at time of death; theirs was from 0.1-50 years old), but found that DRD1 expression was correlated positively with age, with adults having the highest expression, while DRD2 and DRD5 expression levels were inversely correlated with age [14]. While this information is in agreement with our findings on DRD5 expression, it is in disagreement with our findings for DRD1, and DRD2 expression was not found to be associated with age for our variants of interest in our sample. These discrepancies could be from analyzing our data specifically in relationship to genotypic variation, in treating age as a continuous rather than a categorical variable, or it could be due to the different age ranges of the subjects at time of death. Differences in sampling locations within the prefrontal cortex (Brodmann's area 9 in our sample) could also contribute to the discrepancies, as these areas may exhibit differential expression patterns.
In humans, age-related expression associations by tissue are likely to be highly complex and may be due to a complement of environmental factors that occur over the lifespan. More controlled animal models have found some mixed results. Rat models of aging did not show this variation in DRD1 expression by age, though other genes that contribute to dopamine activity (such as NURR1, GCH1, NTRK2, and GFRA1) showed a relationship to aging [15]. Mouse models of aging demonstrate a relationship with DRD2 expression but not with other dopamine receptors [16]. Future work should include epigenetic analysis to determine the mechanism for expression alteration in aging, work that has already shown promise in relationship to the dopamine transporter (DAT) in rat models of aging [17].

Sex
Our results demonstrate that sex was associated with expression of DRD5, with higher expression in females in substantia nigra and higher expression in males in the hippocampus. Higher expression of SLC6A3 was found in females in substantia nigra. This result is interesting, as sex has long been known to be a risk factor for several different neurological and psychiatric disease processes. In particular, males are at higher risk for Parkinson's disease (PD; many of the symptoms of this disease are associated with reduced dopamine activity in substantia nigra), ADHD (hippocampal volume is reduced in persons with this disorder), and autism (also with hippocampal effects related to long-term memory), all of which may be linked to lower expression of dopaminergically active genes across our tissues of interest [18]. Furthermore, brain morphology and cognitive function is known to be sexually dimorphic in schizophrenia, a condition which is highly associated with dopaminergic genetic variation and expression [19]. We are only just beginning to understand the impact of such biological influences of sex. However, expression of sex chromosome genes, such as sex-determining region on the Y chromosome (SRY), may interact with expression of dopaminergic genes to create such differences [20], though this idea should be investigated further.

BMI
Although obesity is often associated with dopaminergically active genes (such as in binge eating disorder [21]), in this inquiry, BMI was only associated with expression of DRD2 in our sample within substantia nigra. There were no other associations with this phenotype in other brain tissues or genes studied in this work. Still, our work adds to the growing body of literature on the involvement of DRD2 in weight and obesity [21][22][23], and is in agreement with PET imaging of obese subjects that demonstrated increased D2 availability in substantia nigra [24]. Work in animal models has demonstrated altered DRD2 gene expression in animals that were resistant to calorie-restriction weight loss [25]. Further inquiry is needed to uncover the potential of DRD2 targets for weight loss management.

Genotype x Expression Associations
Apart from the demographic and phenotypic inquiry, we found that there were three haplotype blocks of SNPs within the DRD2 gene that were associated with expression across all four brain tissues. In brief, the DRD2 gene codes for the dopamine receptor type 2, a G-protein-coupled receptor that has a very strong affinity for binding dopamine and functions to inhibit adenylyl cyclase activity. The gene is located on Chromosome 11q23.2, with the genomic coordinates 11:113, 409, 594-113, 475, 397. Variants in this gene are associated with schizophrenia [3,4], and addictive behaviors, including alcoholism and nicotine use [1,2], substance abuse disorder [26,27], and weight/binge eating behavior [21][22][23].
The first SNP, rs7103679, is intronic and is not reported to have clinical significance in the literature to date. However, variants in the second and third blocks have been associated with various neurological and psychiatric diseases in the literature. We will discuss these associations by haplotype block in turn.

Haplotype Block 2: rs10891549 and rs6277
The first variant within this block is rs10891549. This is a silent variation with no amino acid consequence, but still the variant has been associated with alcohol consumption in the literature [28]. No other results have been indexed in PubMed, but this SNP is in high LD with a much more well-studied variant, rs6277.
The second variant within DRD2 found to be associated with expression across the four brain tissues of interest is rs6277. This variant (also known as C957T) is more widely reported in the literature. Although the structure of the amino acid chain is not altered, this polymorphism does still alter mRNA expression overall, resulting in altered protein folding, which reduces translation and receptor synthesis [29]. Alone, the variant has been shown to be associated with higher fiber tract integrity for connections from basal ganglia to frontal cortex [30], a fact that may result from altered expression of the alternate allele, and a fact that supports the current findings. However, the variant is known to be in linkage disequilibrium not only with rs10891549, but also with other well-studied variants Taq1A (rs1800497) and -141C ins/del (rs1799732), located in the promoter region of this gene. Indeed, when the haplotype block CCGCCGTT (rs6277-rs1076560-rs2283265-rs2734833-rs2075652-rs1079596-rs4436578-rs11214607) was considered, phenotypic associations of the block were found to be associated with heroin dependence behavior [26]. As a clinically actionable target, a metanalysis has supported the hypothesis that rs6277 potentially alters response to risperidone treatment of schizophrenia [31], possibly doing so through altered expression, although a definitive conclusion could not be found at this time.

Haplotype Block 3: 25 Variants in DRD2
The 25 variants that make up this haplotype block are shown in Table 2. Several of these have been reported in the literature, either alone or in portions of the block demonstrated here. Overall, variants in this block are associated with alterations in the D2 short/long isoform expression ratio, schizophrenia and schizotypal symptomatology [32], altered fronto-striatal networks [33], responsiveness to neuroactive medications, and substance abuse disorders. To highlight a few associations, rs2283265 has demonstrated a differential effect by female sex for cocaine addiction (along with a 48 bp VNTR), but in mixed-sex groups the variant has been shown to confer risk of cocaine overdose leading to death (OR~3 [34]). Sex-based effects have also been found for ADHD symptomotolgoy for rs1079727 and rs1124491 [35].
Additionally, within this same block, the rs1076560 T allele has been associated with responsiveness to anti-psychotic medications [36] and reduction in symptomatology in pharmacologic treatment of Parkinson's Disease [37], but was most recently demonstrated to be associated with the personality characteristics of neuroticism and anxiety in polysubstance dependent males [38]. However, other studies have not found an association for this variant and substance use disorder [39,40]. In conjunction with rs2242593, the two variants together have been shown to result in higher dopamine binding in the ventral striatum for healthy adults; our work here indicates that the phenomenon may continue in other regions of the brain [41].

Limitations
Despite these interesting results, we had several limitations in the conduction of this first-of-its-kind inquiry. First, due to the nature of this large-scale but preexisting database of human tissue samples, there were varying sample sizes between tissues and genes. More problematically for clinical implications of this work, detailed phenotypic data were not available for diseases and processes that would be of high interest. Data on subject health histories, underlying diseases/co-morbidities, and medication histories are not available. Considering the prior work in the genetics of our variants of interest, we feel that any information that could have been provided on schizophrenia, Parkinson's disease, ADHD, drug use, and connectivity (structural integrity) of these regions would greatly enhance the work, but is unfortunately not available through this resource. Educational level would also be of relevant interest for this work but was similarly not included in the dataset. Furthermore, this dataset is racially quite homogenous (>90% Caucasian). Although we had significant results for expression variation by race, such lack of variability is likely to have skewed the result. Racial/ancestral marker associations should be investigated further in other projects.

Conclusions
Still, despite these limitations, we found that age alters expression of DRD1 in the nucleus accumbens and prefrontal cortex, DRD3 in the nucleus accumbens, and DRD5 in the hippocampus and prefrontal cortex. Sex was associated with expression of DRD5 in substantia nigra and hippocampus, and SLC6A3 in substantia nigra. We found that three blocks of SNPs, all located in DRD2, were associated with alterations in expression across all four tissues of interest. These demographic characteristics and these variants should be further investigated for use in screening, diagnosis, and future treatment of neurological and psychiatric conditions. Supplementary Materials: The following are available online at http://www.mdpi.com/2076-3425/10/8/567/s1, Table S1: Genotype x expression variants by tissue (detailed information). Table S2: Demographic and phenotypic characteristics for participants (N = 249). Table S3: Association of gene expression with demographic and phenotypic characteristics, by tissue. Figure S1. Linkage disequilibrium structure for 3 haplotype blocks of DRD2.
Author Contributions: Conceptualization, A.G.S. and X.C.; methodology, A.G.S. and X.C..; formal analysis, X.C.; writing-original draft preparation, A.G.S. and X.C.; writing-review and editing, A.G.S. and X.C.; visualization, X.C.; funding acquisition, A.G.S. and X.C. All authors have read and agreed to the published version of the manuscript.