Regional transcriptome analysis of AMPA and GABAA receptor subunit expression generates E/I signatures of the human brain

Theoretical and experimental work has demonstrated that excitatory (E) and inhibitory (I) currents within cortical circuits stabilize to a balanced state. This E/I balance, observed from single neuron to network levels, has a fundamental role in proper brain function and its impairment has been linked to numerous brain disorders. Over recent years, large amount of microarray and RNA-Sequencing datasets have been collected, however few studies have made use of these resources for exploring the balance of global gene expression levels between excitatory AMPA receptors (AMPARs) and inhibitory GABAA receptors. Here, we analyzed the relative relationships between these receptors to generate a basic transcriptional marker of E/I ratio. Using publicly available data from the Allen Brain Institute, we generated whole brain and regional signatures of AMPAR subunit gene expression in healthy human brains as well as the transcriptional E/I (tE/I) ratio. Then we refined the tE/I ratio to cell-type signatures in the mouse brain using data from the Gene Expression Omnibus. Lastly, we applied our workflow to developmental data from the Allen Brain Institute and revealed spatially and temporally controlled changes in the tE/I ratio during the embryonic and early postnatal stages that ultimately lead to the tE/I balance in adults.

. Global and regional transcriptomic analysis of AMPAR subunits (GRIA1-4) in the human brain. (a) Diagram detailing transcriptomic workflow. Publicly available data from the Allen Brain Atlas, Brainspan, and the ADTBI studies were filtered using exploratory factor analysis and analyzed at the global and regional levels for both GRIA1-4 expression and tE/I ratio patterns. (b) Box plots of microarray gene expression from the Allen Brain Atlas (111 substructures per 6 subjects) and corresponding percent contributions at the global level. The median is represented by the line within the box, and the first and third quartiles are represented by the ends of the box. The whiskers extend from each end of the box to the first or third quartile ± 1.5 (interquartile range). Points outside of the whiskers are outliers and color-coded according to the inset. FL frontal lobe, Ins insula, CgG cingulate gyrus, HiF hippocampal formation, PHG parahippocampal gyrus, OL occipital lobe, PL parietal lobe, TL temporal lobe, Amg amigdala, GP globus pallidus, Str striatum, Cl claustrum, Hy hypothalamus, SbT subthalamus, DT dorsal thalamus, Vt ventral thalamus, MES mesencephalon, CbCx cerebellar cortex, CbN cerebellar nuclei, Bpons basal part of the pons, PTg pontine tegmentum, MY myelencephalon. (c) Twoway unsupervised Ward's hierarchical clustering of microarray data from the Allen Brain Atlas for analysis of substructures. Major brain regions were separated by Log2 gene expression, and subunits were clustered according to their regional expression levels. Structures are color coded as in b, and ontogeny is detailed by the inset. For substructure abbreviations, please see Supplementary Data 2.
Scientific RepoRtS | (2020) 10:11352 | https://doi.org/10.1038/s41598-020-68165-1 www.nature.com/scientificreports/ (https ://aging .brain -map.org/), which included a much higher number of individuals (n = 50), but was limited to only four regions: the hippocampus (HIP), temporal cortex (TCx), parietal cortex (PCx), and forebrain white matter (FWM). Together, these complementary analyses created a signature of AMPAR subunit expression in the healthy human brain. To determine whether transcripts for AMPAR subunits and GABA A R subunits were balanced, we then calculated and analyzed the ratio between AMPAR and GABA A R subunit expression at both the whole brain and structural levels. This generated an estimate of the most basic signature of the transcriptional E/I (tE/I) ratio. Finally, to investigate the usefulness of these signatures, we applied our method to single cell data from phenotypically characterized interneurons in mice and to RNA-sequencing data from the Allen Brainspan study that ranged from early embryonic stages to adulthood in humans.

Results
Global and region-specific expression patterns of AMPAR subunits. We first analyzed the microarray data from the Allen Brain Atlas at the whole brain level. From this data, we selected the most representative probe for each gene (Supplementary Data 1) using exploratory factor analysis ( Fig. 1a; Supplementary Fig. S1 for selection process) and divided the brain into major regions, structures, and substructures following the Allen Brain Atlas' nomenclature (Supplementary Data 2). The initial global analysis after adjusting by age (Supplementary Data 3; Supplementary Fig. S2) revealed high expression levels (Log 2 > 4) for all four subunits (Fig. 1b). GRIA2 was the most expressed AMPAR subunit in the whole brain followed by GRIA4, GRIA1 and GRIA3 (Fig. 1b). However, there was high variability across brain regions in expression levels for all four subunits, particularly for GRIA1 and GRIA4 ( Supplementary Fig. S2). For example, structures such as the hippocampal formation (HiF; CA1-4 and DG) and the internal globus pallidus (GP i ) were outliers at the global expression level for GRIA1. A nested ANOVA analysis of gene expression, where substructures are nested within structures to account for substructure dependencies, confirmed different levels of expression for all 4 subunits across major structures and across some substructures within the same region (Supplementary Data 4). To address this regional variability, unsupervised hierarchical clustering was performed to group together structures with similar expression profiles (Fig. 1c). The strength of these clusters was examined through non-parametric bootstrap resampling to calculate both the approximately unbiased probability (AU) and bootstrap probability (BP) p values ( Supplementary Fig. S3). The structures were first split between structures of telencephalic origin and those of diencephalic, mesencephalic, metencephalic, and myelencephalic origins (AU > 0.95, Fig. 1c). Within the latter group, further divisions based on structure ontogeny were observed. Regions of metencephalic origin exhibited very low levels of GRIA3 expression (Log 2 < 4) and high levels of both GRIA2 and GRIA4 (Log 2 > 9) while regions originating from the diencephalon, mesencephalon, and myelencephalon exhibited lower GRIA2 expression (Log 2 < 9.9; Fig. 1c). There was very strong evidence (AU > 0.95) for clustering of metencephalic structures at nearly all levels. Structures from these three regions were then further split between those with GRIA1 expression levels comparable to metencephalic structures (Log 2 ≈ 8) and those with very low GRIA1 expression levels (Log 2 < 6). However, evidence for clustering of substructures within these structures was not as robust (AU < 0.95). By contrast, the telencephalic cluster exhibited a high degree of homogeneity, particularly within cerebral cortical regions where high levels of GRIA2 and GRIA4 (Log 2 > 9.5) and lower levels of GRIA1 and GRIA3 (Log 2 < 9.5) were observed (Fig. 1c). Some slight variations were introduced with the cerebral nuclei which typically exhibited lower levels of GRIA4 (Log 2 < 9) and a structure-dependent increase in GRIA1 expression (Log 2 > 9). However, this cluster was not as robust as that of the cortical regions with an AU ± S.E.M of only 0.942 ± 0.005 ( Supplementary Fig. S3). The hippocampal formation differed greatly from all other structures with high GRIA1-3 expression (Log 2 > 10) but low GRIA4 expression (Fig. 1c), and the evidence for this cluster was very strong (AU ± SE = 0.999 ± 0.001).
The hippocampal formation contained the substructures with the greatest fold enrichment of GRIA1 (CA4), GRIA2 (DG) and GRIA3 (CA1) expression (Table 1). GRIA1 was particularly enriched in the CA4 with an expression ~ 1,000% larger than the global brain average.

Proportional contributions of AMPAR subunits and inter-individual variability.
To determine the contribution of all available transcripts for AMPAR subunits per brain region, we calculated the average percentage of expression for each subunit as was previously done for GABA A Rs subunits 16 ; that is to calculate the sum of probe intensities across AMPAR subunits and then determine to proportional contribution of the intensity of each probe to this sum. This proportional contribution is an estimate that represents the available Table 1. Human brain substructures with greatest anatomical enrichment. The expression (log2) represents the average per brain region. The fold enrichment is the substructure's value divided by the global average. HiF Hippocampal Formation, VT ventral thalamus. www.nature.com/scientificreports/ pool of AMPAR subunits mRNA in each brain region, structure, or substructure, and normalizes distinct levels of expression between different brain areas (Fig. 2). As before, the cerebral cortex was highly homogenous with a greater contribution from GRIA2. GRIA1 contributed more in the hippocampus and cerebral nuclei, and noncortical regions were more variable and dominated by GRIA4 (Fig. 2). Further analysis revealed opposing expression patterns between some subunits. For example, the hippocampus showed a gradual decrease in the proportional contribution of GRIA2 from the DG to CA4 that was paralleled by a gradual increase in the contribution of GRIA1 (Fig. 2). Similarly, GRIA4 and GRIA2 followed opposing patterns of expression in the myelencephalon and thalamus compared to the cortex (Fig. 2). It is important to note that microarray data is highly influenced by technical factors, and intensity levels are typically compared across experimental units for the same probe, not across probes within the same experimental unit as we are doing in this study; therefore we compared the proportional contribution of GRIA subunits in the temporal cortex using microarray (n = 6 subjects with 12 temporal lobe substructures per subject) with that using RNA-Seq of the ADTBI study (n = 50 subjects), which is a better measure of mRNA expression. Proportional contributions in both datasets showed similar relationships across GRIA1, GRIA2 and GRIA3 subunits providing support for our microarray analysis of these subunits; the proportion of GRIA4 seems to be overrepresented in microarray compared to RNA-Seq ( Supplementary Fig. S4), a probable consequence of choosing the probe with higher correlations across all GRIA4 probes in the whole brain, that happened to also have the highest intensity, highlighting some of the technical limitations aforementioned ( Supplementary Fig. S1).
In addition to differences in level of gene expression of isolated subunits, we also explored the variability in the collective organization of AMPAR subunits between individuals. To do so, we followed the previous approach used for GABA A Rs 16 . First, inter-individual variability, within and across substructures of the microarray study, www.nature.com/scientificreports/ was quantified by calculating the Euclidian distances (d) between the expression levels of GRIA subunits in each region and correlating these results. For example, the correlation coefficient (R) between the superior frontal gyrus across individuals (di) was higher (di SFG-SFG = 0.94 ± 0.05; Mean ± SD) than in the dentate gyrus (di DG-DG 0.81 ± 0.18), indicating that the dentate gyrus is more variable across individuals than cortical regions (Supplementary Fig. S5). To extend the analysis and compare brain regions at the level of substructures the frontal operculum (fro), which is the frontal-most structure displayed in the Allen Atlas, was chosen to construct a reference (consensus) against which all other 110 substructures were compared (dc Fro ). Each substructure from each subject was then measured against this reference ( Fig. 3a and Supplementary Fig. S6). Cerebral cortical structures were mostly homogenous with correlation coefficients close to one, indicating a high similarity to other cerebral structures and to the frontal operculum ( Fig. 3a and Supplementary Data 6 for statistical analyses). Hippocampal structures, as expected, differed greatly from other cerebral regions (r < 0.6) and were more variable between individuals (Fig. 3a). However, the microarray study included only six subjects, so we turned to the ADTBI study (n = 50 subjects) to examine inter-individual variability more accurately. For this analysis, the parietal cortex was chosen as the reference structure against which the other three structures were compared. The parietal and temporal cortices were identical with virtually no global variability between all 50 subjects (Fig. 3c).
Surprisingly, unlike what was previously observed for GABA A R subunits 16 , the forebrain white matter was very similar to the cortices with low inter-individual variability (Fig. 3c). The hippocampus, in agreement with the microarray data, was very different from the other structures and highly variable between subjects (Fig. 3c). This difference arose primarily out of a higher proportional contribution of GRIA1 in the hippocampus as opposed to the other three regions.

Global and region-specific patterns of transcriptional excitation-inhibition ratio markers.
To combine these results with the expression patterns in GABA A Rs subunits we previously uncovered, we used the microarray gene expression data as a marker of the tE/I ratio by dividing the sum of probe intensities across AMPAR subunits by the sum of probe intensities across GABA A R subunits (Supplementary Data 7). This is the most basic and non-subunit specific estimation of the transcriptional foundations of the E/I balance based on gene expression for the principal gates of excitation and inhibition. Because age correction may change GABA A Rs and AMPARs differently and lead to distinct ratios in each subject, we compared the tE/I ratio www.nature.com/scientificreports/ using gene expression levels with and without age correction. Both ratios were linearly correlated (R 2 = 0.789; p < 0.0001) and no significant difference in each substructure was observed between both ratios, although more dispersion was observed using non-corrected data (p > 0.29 in all cases by Welch's testing means equal, allowing SD not equal, Fig. 4 and Supplementary Fig. S7). For simplicity, we decided to use the tE/I based on non-adjusted gene expression values. Expression patterns across cerebral cortical regions were relatively homogenous with tE/I ratios close to one (Fig. 4a). Similarly, cerebellar regions exhibited highly homogenous tE/I ratios though these were lower than those found in the cerebral cortex (Fig. 4a). Hippocampal structures exhibited a higher tE/I ratio overall compared to the rest of the brain, apart from the reticular nucleus (Supplementary Data 7). Within the hippocampus, the tE/I ratio gradually decreased from DG to CA2 and then increased from CA3 to CA4 where AMPAR subunit expression was nearly double that of GABA A Rs subunit expression (Fig. 4a). Cerebral nuclei and thalamic regions were more variable with more caudal regions tipping towards greater inhibition (Fig. 4a). The structure with the greatest difference was the reticular nucleus in the ventral thalamus (tE/ I mean = 3.8, Fig. 4a). These observations were further supported by the ADTBI RNA-Seq dataset. The tE/I ratios using non-adjusted by age gene expression in the parietal and temporal cortices were remarkably homogeneous and both were significantly different from the hippocampus, which was both higher and more spread ( Fig. 4b; Supplementary Data 8 for statistical analysis). No differences in the tE/I ratio in the temporal and parietal cortices, hippocampus or forebrain white matter by sex were observed (p > 0.11 double tailed t test; Supplementary Fig. S8). Even though the tE/I ratio in the forebrain white matter was more variable across individuals, the tE/I ratios from these structures were remarkably constant overall, indicating that AMPARs and GABA A Rs were highly correlated within brain structures, although with different slopes (Fig. 4c). Such strong correlations between AMPARs and GABA A Rs suggest that transcriptional correlations start at the cellular level. Indeed, the www.nature.com/scientificreports/ electrophysiological E/I (eE/I) balance has been demonstrated in excitatory pyramidal neurons as well as inhibitory interneurons [31][32][33] . However, analyzing the tE/I at the single cell level by single-cell RNA-Seq (sc-RNA-Seq) in humans is challenging due to the inherent low coverage of the method (approx. 10-20% of the transcriptome) and the potential mRNA degradation by postmortem interval and agonal factors 34 . Therefore, we used a high resolution sc-RNA-Seq dataset from six phenotypically identified cortical interneurons from mice 35 to test for intra-and inter-variability of the tE/I ratio at the cell-type level (Fig. 5). For this analysis, we used the dataset from Paul et al. 35  These cells demonstrated strong linear correlations between AMPARs and GABA A Rs with tE/I ratios averaging between 0.66 to 1.8 ( Fig. 5b; Table 2). Only MNCs were significantly different with a tenfold higher tE/I ratio (P < 0.0001; ANOVA followed by Tukey-Kramer HSD), and CHC1 cells were significantly different from CCKC cells (P = 0.022). No other cellular types were significantly different from the others (P > 0.05). te/i ratio expression patterns change during development. To demonstrate the flexibility and capabilities of our approach, we applied our workflow to developmental data gathered by Allen Brainspan focusing on the cortex since it has more data than the other structures studied in that dataset (Supplementary Data 9). The www.nature.com/scientificreports/ tE/I ratio across cortical regions was initially spread across a large range, but all structures reached a relatively stable tE/I ratio by one postnatal year (Fig. 6a). Interestingly, the tE/I ratio in the cortex using the Brainspan dataset was smaller than the one in the ADBTI study, although both used RNA-Seq methods (tE/I = 0.31 ± 0.0061 (mean ± S.E.M) in Brainspan vs 0.5842 ± 0.0045 in ADTBI; double tailed t test, p < 0.001). During pre-natal development, some structures, namely the dorsolateral and inferolateral temporal cortices, auditory primary cortex, and posteroventral parietal cortex exhibited high tE/I ratios that rapidly declined approaching birth (Fig. 6a).
Other structures, such as the primary somatosensory and anterior cingulate cortices rose slightly during early development and peaked at around 19 post-conception weeks (pcw) before declining towards the steady state observed in adulthood (Fig. 6a). Overall, the cortical tE/I ratio stabilizes after birth. Closer inspection of the data revealed that the developmental reduction of tE/I resulted from a gradual increase in GABA A R subunit expression levels while AMPAR subunit levels remained stable ( Supplementary Fig. 9). In addition to the tE/I ratio, we also analyzed the patterns of expression of the sodium-potassium-chloride cotransporter (NKCC1, SLC12A2) and the neuron-specific chloride-potassium symporter 5 (KCC2, SLC12A5) which have been demonstrated to change with development and determine whether GABA A Rs are inhibitory or excitatory 36 . KCC2 expression increases, while NKCC1 expression decreases during development. To do so, we compared the ratio of SLC12A5 to SLC12A2 with the tE/I ratio and found that during prenatal development, the higher tE/I ratio corresponded to a low SLC12A5/SLC12A2 ratio while postnatal development exhibited no relationship between the two ratios ( Fig. 6c). KCC2 expression is the main driver of variability in the KCC2/NKCC1 ratio. Interestingly, KCC2 gene expression was linearly correlated with the global expression of GABA A R and AMPAR subunits.

Discussion
Our clustering analysis indicates that the co-expression patterns of AMPARs subunits is region specific according to the ontogenic origin of the brain structures. It is quite remarkable that with the expression of only 4 genes, major brain regions could be separated by unsupervised clustering analysis. This suggests that AMPARs may act as strong gene markers in synaptic components previously described by bioinformatics analysis [37][38][39][40] . Our microarray analysis fit the clustering scheme created by Gold et al. in the rodent model using in situ hybridization 41 , though with some differences. For instance, while the higher contribution of GRIA1 in the amygdala and hippocampus relative to other cerebral structures was mirrored in our analysis, it was to a much lesser extent, only 20-30% as compared to the 50% reported. Further, the cerebellum exhibited a high GRIA4 profile that is inconsistent with any of the reported clusters. On the other hand, cortical regions aligned fully with the pattern of ~ 50% GRIA2 contribution. While data from mice showed GRIA3 and GRIA4 enrichment relative to GRIA1 in the thalamus, our analysis found only GRIA4 enrichment. Both differences may be related to the high anatomical specificity of structures in hybridization studies compared to the microarray data that we used for our analysis, or it could be due to divergent characteristic between species 42 . One potential drawback of our microarray analysis is the low number of samples which reduces the statistical validity of our findings. However, our main aim for this microarray data was to generate a transcriptional signature of structural clusters based on the relative expressions of the different GRIA subunits. For this purpose, we utilized bootstrapping to measure the strength of our clusters. As a result, we discovered very strong evidence for upper-level divisions between cortical regions of telencephalic origin and regions of myelencephalic, metencephalic, and diencephalic origin, as well as the existence of distinct structural clusters, such as the hippocampal formation. While these results are promising, additional samples would boost the power of this study.
Our observed region-specificity in AMPAR subunit expression correlates with the differing complexity in pathways and cytoarchitecture and with embryonic origin. Indeed, previous studies have found that a large selection of genes exhibit region-specific expression largely dependent on their embryonic structure of origin 38,43 . For brain structures with recurrent cytoarchitecture, such as the cerebral and cerebellar cortices 44 , AMPAR subunit expression patterns are highly stable well into adulthood. This shared pattern may result from a shared embryonic origin such as the pallial regions for cortical glutamatergic cells 45 and the rhombic lip for the glutamatergic cells in the cerebellar cortex 46 . In contrast, the more heterogenous patterns observed in deep cerebral nuclei, thalamic nuclei, and mesencephalic structures are congruent with more heterogenous cytoarchitecture of these regions and may underly the diverse modulatory functions of these nuclei 47 . Recent studies have begun to demonstrate the importance of using functional rather than neuroanatomical principles for the analysis of synaptic proteins 48 ; www.nature.com/scientificreports/ www.nature.com/scientificreports/ similar analytical approaches with gene expression may help to explain the heterogenous expression of AMPA receptor subunits in the thalamic nuclei and mesencephalic structures in future studies. Interestingly, the hippocampal formation differed greatly from the cerebral cortices, exhibited the highest expression levels of all AMPAR subunits compared to all the other brain structures, and showed the largest variability across individuals when using the correlation of Euclidian distances between subunits. The heightened expression may result from the unique neurogenesis in the hippocampus 49 , the rapid regulation of synaptic AMPA receptors necessary for its role in learning [50][51][52] and the high neuronal density in this area. Previous studies have shown the necessity of GRIA1 subunit in long-term potentiation and contextual learning at the CA3-CA1 synapse 52,53 and the rapid subunit specific turnover of AMPARs in response to alterations in neuronal activity 54,55 . Moreover, the subgranular zone of the dentate gyrus is one of the only regions in the human brain where adult neurogenesis occurs 56 . As a part of this process, neural stem cells enter hippocampal circuits after modulation by various environmental and cell-intrinsic factors 57,58 . This intense plasticity suggests that hippocampal AMPAR subunit expression varies greatly depending on the subject and environment. Our analyses of both the microarray and RNA-seq studies support this hypothesis as inter-subject variability in AMPAR subunits expression is particularly high in the hippocampal formation while AMPAR subunits expression in the cerebral cortex is highly stereotypical across individuals.
A major interest in our analysis was to determine whether the mRNA for the receptors involved in establishing the E/I ratio were balanced across brain structures in human. Because many subunits of GABA A Rs and AMPARs follow complementary or opposing patterns of expression, we used the most conservative approach which includes all potential receptor configurations by adding transcripts for all subunits when estimating the tE/I ratio. Our analysis found that the total amount of AMPARs and GABA A Rs is highly correlated in a regionspecific manner. Therefore, the tE/I ratio remains highly stable across substructures within major regions such as the cerebral and cerebellar cortices but varies gradually within the hippocampus and amygdala, suggesting different balances within these regions. Our results using the ADTBI study did not find sex differences for the tE/I ratio, suggesting that in healthy individuals, sex differences for particular GABA A Rs subunits like higher cortical α1 in men 59 , may be balanced by compensatory mechanisms regulating the expression of AMPARs.
Our study unfortunately cannot distinguish between synaptic and non-synaptic receptors; therefore, our inferences are limited to estimates of total transcription levels per brain structure in the microarray and ADTBI analyses, or to the cell-type level in the single cell analysis. Nevertheless, the electrophysiological activity of synaptic and non-synaptic receptors GABA A Rs strongly correlates within single neurons 60 , and AMPARs show high lateral mobility between non-synaptic and synaptic regions 61 suggesting that whole transcription levels may at least partially reflect the activity-dependent need for AMPAR and GABA A R subunit transcription. Recent studies have also shown that electrophysiological measures of global eE/I by recording cortical synaptic AMPARs and GABAARs are highly correlated 62 . Future studies integrating transcriptomic with electrophysiological measures of the E/I will help to understand deviations of this balance in non-physiological states.
Our analysis also shows that the tE/I ratio is established near birth. While AMPAR subunits expression is relatively stable during the embryonic stages, GABA A R subunit expression gradually increases until stabilizing near birth reducing the tE/I ratio to a stable value. This indicates that GABA A Rs are highly important in setting the tE/I ratio. In fact, development of GABA A Rs from immature GABA A α2-and GABA A α3-enriched to mature GABA A α1-enriched receptors has been observed to play a major role in synaptic plasticity in the cat visual cortex during the early post-natal phases 63 and a fundamental role in the establishment of circuits and neuronal-firing properties within the mouse cortex 64 . On the other hand, although the amount of AMPARs remains quite stable during the development, the subunit expression and receptor assembly profile change markedly leading to a shift in receptor properties. As described by Henley and Wilkinson 65 , subunit adjustment is directly related to function; for instance, GRIA1 is developmentally restricted while expression of GRIA2, which produces Ca 2+ -impermeable receptors, increases after birth guiding the activation of silent glutamatergic synapses and the consolidation of the synaptic neural network. Our findings using the Brainspan dataset are in line with these descriptions (Supplementary Fig. S10), suggesting that the relationship between GRIA1 and GRIA2 are important in stabilization and regulation during brain development.
In conclusion, we show that AMPARs, as well as GABA A Rs, follow stereotypical patterns of expression in healthy controls; therefore, the tE/I ratio is remarkably stable across brain regions with recurrent cytoarchitecture. Additionally, this balance is established early after birth and remains stable through adulthood. Our analysis also shows that the tE/I ratio is determined at the cellular level, highly constant within neurons of the same type, and very similar across transcriptionally different cell types, with some exceptions. The stability of the tE/I ratio suggests that many redundant mechanisms control the concerted abundance of AMPAR and GABA A R subunit transcripts, and disruption of this balance may have severe impacts on proper neuron and brain function.

Methods
Microarray and RNA-sequencing databases. Three publicly available databases from the Allen Institute were used. For global brain analysis, normalized microarray transformed data (Log2) from the Allen Brain Atlas was downloaded from the webpage (https ://human .brain -map.org). This database was generated with brain samples obtained from six subjects (five males and one female) between the ages of 24 and 57 years of age with no known neuropsychiatric or neuropathological history. Brain tissue was collected after obtaining informed consent from decedent's next-of-kin. Institutional Review Board (IEB) review and approval was obtained for collection of tissue and non-identifying case information at the tissue banks and repositories that provided the tissue for the project. Case qualification and donor profiles can be seen in the Allen website: https ://help.brain -map.org/displ ay/human brain /Docum entat ion. Only brain regions where data for all AMPAR and GABA A R subunits were measured were used for the analysis (111 brain structures). For gene expression measurements www.nature.com/scientificreports/ the Allen Institute used a custom design (by Beckman Coulter Genomics) Agilent 8 × 60 K array that includes the 4 × 44 K Agilent Whole Human Genome probe set supplemented with an additional 16,000 probes. The upper limit for detection was 2 and was used as part of the initial quality control as described in the white papers for microarray survey found at https ://help.brain -map.org/displ ay/human brain /Docum entat ion. Normalization methods are also documented extensively in the white papers; briefly the data was normalized within each dissection batch through multivariate local regression fitting, each brain via quantile-quantile mapping of averages, and across all brains by aligning control samples. For analysis of variability across individuals, normalized RNA-Seq Fragments Per Kilobase Million (FPKM) data from the Aging, Dementia, and Traumatic Brain Injury study (https ://aging .brain -map.org/downl oad/index ) covered four regions from 56 healthy control subjects (35 males and 21 females) between the ages of 78 and 99 years of age. Detailed documentation of the ADTBI study, including tissue collection, tissue processing, and RNA-Sequencing and quantification, can be found at: https ://help.brain -map.org/displ ay/aging /Docum entat ion. Briefly, RNA-Sequencing was done through Illumina TruSeq (random hexamer first strand cDNA synthesis with rRNA depletion and fragmentation) on Illumina HighSeq 2,500 using v4 chemistry, producing a minimum of 30 M 50 bp paired-end clusters per sample. FPKM gene quantifications were normalized via TbT normalization and corrected for RNA quality and batch effects.
Normalized developmental RNA-seq Reads Per Kilobase Million (RPKM) data from the BrainSpan atlas (Gencode v10 summarized to genes; https ://www.brain span.org/stati c/downl oad.html) covered up to 26 regions from 42 subjects ranging from eight weeks post conception to 40 years of age. Detailed documentation of the developmental data, including Institutional Review Board approval for the tissue collection, donor and sample metadata, and transcriptome profiling methods, can be found at: https ://help.brain -map.org/displ ay/devhu manbr ain/Docum entat ion. Briefly, RNA-Sequencing was done with poly(A) selection and normalized by the addition of spike-ins.
For single cell analysis we used the processed high resolution sc-RNA-seq unique Transcripts Per Million (uTPM) dataset downloaded from GEO expression omnibus accession number: GSE92522 35 . This dataset comprises single cell transcriptomes of anatomically and physiologically characterized GABAergic neurons. Paul et al. developed combinatorial recombinase driver lines to identify and isolate 6 types of neurons. Detailed methods, including the approval for the study, can be found in the author's manuscript 35 . Sequencing was performed with poly(A) linear amplification followed by Illumina TruSeq and normalized to a single cell's total unique counts across all genes; approximately 9 k genes per cell were detected 35 . Our study focused on the analysis of the databases described above. All the methods to generate the original data were performed in accordance with the relevant guidelines and regulations of the respective authors' institutions as described in the documentation accessed by the provided links.
Analysis of AMpAR subunits. All analysis was done as previously stated for GABA A Rs 16 . The normalized transformed data (Log 2 ) microarray data from the Allen Brain Atlas used 17 probes to measure the expression of 4 AMPARs genes, 4 probes for each GRIA1-3 and 5 probes for GRIA 4. To avoid redundant clustering due to collinearity between probes for the same gene, only the most representative probe for each gene as determined via principal axis Exploratory Factor Analysis with no rotation in JMP 14Pro were used. While the effects of sex and ethnicity were negligible on any of the datasets analyzed, corrections were made for age effects via linear regression in JMP 14. As before, all expression patterns were represented by the Mean (M) ± Standard Deviation (SD) unless otherwise stated. The proportional contributions of each subunit were expressed as percentages of untransformed (non-Log2) probe intensity level of each subunit to the sum of probe intensities across all subunits. Unsupervised hierarchical clustering was done using Ward's minimum variance method in JMP 14. Bootstrap probability (BP) and approximately unbiased (AU) p-values were calculated through 10,000 bootstraps with pvclust in RStudio using ward.D2 and Euclidean distances as the clustering method and distance measures, respectively 66 . Pearson product-moment correlations were used to measure both inter-structure and inter-subject variance. Euclidean distances of each structure or subject (d i ) were correlated against a chosen consensus (d c ). For subjects, d c was calculated from the average of all subjects; for structures, the frontal operculum (fro) and parietal cortex (PCx) were used. The Euclidean distances were calculated in RStudio as follows: where d i is the Euclidean distance between two subunits, given by x and y, per subject or structure, represented by the integer i, and d c is the consensus Euclidean distance, for all n subjects or structures.
Analysis of the excitation-inhibition ratio. The tE/I ratio was calculated by dividing the sum of all AMPAR subunit expression by the sum of all GABA A R subunit expression. Pearson product-moment correlations were then run on these ratios to measure inter-structure variance. Data from the microarray (Log2, uncorrected and corrected for age), ADTBI (normalized FPKM, uncorrected and corrected for age), single-cell (uTPM), and BrainSpan (normalized RPKM) studies were all analyzed separately. In addition to the tE/I ratio, the expression of NKCC1 (SLC12A2), NKCC2, (SLC12A1) and KCC2 (SLC12A5) were also analyzed in a fashion similar to the AMPARs and GABA A Rs. (1)