Spatially Comprehensive, Genome-Wide, Exon-Level Expression Profiling

We hybridized RNA isolated from the regions illustrated in Figure 1A to Affymetrix GeneChip Human Exon 1.0 ST Arrays ("Exon Arrays") to obtain independent genome-wide “whole transcript coverage” expression data from each of the four microdissected brains (Table S3). Probesets on the Exon Array are classified as Core, Extended, or Full according to the level of annotation of the source sequence(s); all analyses reported here utilized only the highest confidence “core” probesets, which are based on RefSeq and GenBank full-length mRNAs. Probesets are grouped into “transcript clusters” corresponding to all possible isoforms transcribed from a single locus or gene; therefore, for simplicity we refer to transcript clusters as genes throughout these results.

Hybridization of selected samples to both Affymetrix U133 and Exon microarrays revealed comparable differential expression results between the two platforms (R² > 0.5; Figure S4). In addition, we performed principal component analysis (PCA) to assess the consistency of Exon Array data within and between brain regions, in comparison to the variability across subjects or hybridizations (Figure S5). PCA confirmed that brain region, rather than individual differences, contributed the majority of variance to the data. Together, these results supported the validity of the Exon Array for comprehensive expression profiling of the mid-fetal human brain.

Regional Differential Gene Expression and Alternative Splicing

We first analyzed gene expression and splicing differences across the five major embryonic brain divisions (NCTX, HIP, STR, THM, and CBL) represented in our samples. Out of 17,421 core transcripts, just over 76% (13,223) were present in at least one brain region, indicating that the great majority of human genes are significantly expressed in the developing brain (Figure 1B). Moreover, at a stringent false discovery rate (FDR) of 10⁻⁵, 4,369 genes, or about 33% of those present, were differentially expressed (DEX). In addition, 3,755 genes (28%) exhibited differential exon usage suggestive of region-dependent splicing (“differentially alternatively spliced,” DAS), for a combined 44% of expressed genes showing evidence of some form of differential regulation. A total of 2,260 genes (17%) fell into both categories (DEX and DAS). Of particular note are 1,495 genes, or 11% of those present, detected as DAS but not DEX, which represent a cohort of genes whose differential regulation across brain regions would not be detected by older 3’-biased array platforms. For our purposes, AS encompasses alternative promoter usage, polyadenylation site usage, and cassette exon splicing. Numbers of genes present, DEX, and DAS from all analyses are given in Table S4.

To identify candidates for further study, we selected DEX genes (FDR < 10⁻⁵) with a minimum two-fold difference in expression between any two brain regions and performed unsupervised hierarchical clustering on both genes and brain regions. We used the resulting heatmaps to identify groups of genes with the most specific or restricted expression patterns (Figure 2). In addition to large numbers of novel expression patterns, these clusters included transcription factors whose mouse orthologs are critical for the development of region-specific neuronal cell types, including FEZF2, SATB2,SOX5 and TBR1 in cortex, and TITF1 in STR (Sur and Rubenstein, 2005; Chen et al., 2005; Molyneaux et al., 2007; Leone et al., 2008; Britanova et al., 2008; Kwan et al., 2008). This validates our approach to identifying region-specific expression patterns in the developing human brain. Correlations and sizes of these and all subsequent clusters are given in Table S5; complete lists of genes in these clusters are given in Table S6.

Open in a separate window

Figure 2

Unsupervised Hierarchical Clustering and qRT-PCR Validation of Selected Genes Differentially Expressed Between Brain Regions

Genes with FDR<10⁻⁵ and greater than twofold maximum expression difference were clustered using uncentered correlation and average linkage. Selected highly correlated (r>0.75) clusters of genes with the most restricted expression patterns, indicated in the heatmap at left, are shown in detail at right, with representative validated genes labeled in red. Red is higher expression, blue is lower expression. Complete gene lists for these clusters are given in Table S6. Bar graphs show qRT-PCR (black bars) next to Exon Array data (grey bars) for representative validated genes. Bar graphs represent fold-change relative to the average of all samples for qRT-PCR (black bars; mean±SEM) and median normalized log expression level for Exon Array (grey bars).

Consistent with their ontogenetic and phylogenetic closeness and similar cellular composition, NCTX and HIP were the only two brain regions with a positive correlation (r=0.13) across all of the genes clustered, as reflected in the 57 genes in cluster one enriched in both regions (Figure 2 and Table S5 and Table S6). CBL was the most distinct of the brain regions sampled, consistent with previous findings in the adult human brain (Khaitovich et al., 2004; Roth et al., 2006; see also Table S11) and reflecting differences in both its cytoarchitecture and possibly its developmental time course relative to forebrain structures.

These data provide a survey of the proportions of the human transcriptome with broad or specific expression in the late mid-fetal brain, as well as identifying a large number of regionally enriched or alternatively spliced genes not previously identified as such. In addition, the detection of many orthologs of important rodent neurodevelopmental genes suggests their human counterparts play evolutionarily conserved roles in establishing regional neuronal identity and connectivity.

Validation of Regional Differential Expression Data

We next sought to validate some of these novel expression patterns by additional methods. Genes were prioritized by expression level, fold change, functional classification, and existing rodent and monkey data. Particular emphasis was placed on previously unknown or uncharacterized genes, and on those that appeared to differ in their regional enrichment from available data in other species. Sixty-five of 68 genes validated by qRT-PCR (bar graphs in Figure 2; Table S6 and Table S9) correlated with Exon Array data (median r=0.95), giving us extremely high confidence in Exon Array gene-level analysis. We further examined over 20 candidate genes by in situ hybridization (ISH) or immunohistochemistry (IHC) in whole sagittal or coronal sections of more than ten additional mid-fetal human brains (Table S6, Figure S6, and data not shown; brain specimens listed in Table S1). These confirmation data provide further confidence in the Exon Array platform, as well as identifying an initial cohort of strong candidate genes for functional analyses.

Genetic Patterning of the Developing Human Neocortex

Next, we used a similar approach to search for areal differences in gene expression and splicing within the NCTX. Previous studies have found few differences between adult human neocortical areas (Khaitovich et al., 2004; Roth et al., 2006). We hypothesized that developing NCTX would show more genetic differences, both due to the possibility of an increased number of genes expressed during development and because of the neuronal differentiation and establishment of connectivity that occur during the late mid-fetal period. Comparing among nine areas of NCTX and using the same FDR threshold as the previous analysis (10⁻⁵), we identified 471 DEX (4.1%) and 496 DAS (4.3%) genes from the total of 11,407 expressed (Table S4). We hypothesized that important gene expression differences among neocortical areas might be smaller in magnitude than those between brain regions, and thus relaxed the FDR to 1% for exploratory cluster analysis, yielding 1,753 DEX (15.4%) and 828 DAS (7.3%) genes (Figure 1C). Clustering revealed a clear division between the four PFC areas and the four non-frontal lobe areas (PAS, TAU, TAS, OCC) (Figure S7). Interestingly, the MS area, which was dissected from the border region of the frontal and parietal lobes, was not highly correlated with either prefrontal or non-frontal areas. Rather, genes enriched in MS were roughly evenly divided into those correlated with PFC and those correlated with more posterior cortical areas (Figure S7), consistent with the mixed frontal/parietal nature of the tissue sample.

In order to investigate both PFC- and non-frontal-enriched genes in greater detail, we performed two additional, targeted intra-NCTX analyses. First, we grouped together the four PFC samples and compared them to the remaining NCTX areas. Hierarchical clustering of DEX genes (Figure 3A and Table S5 and Table S6) revealed, in addition to PFC- and non-frontal-enriched genes, more specific patterns of enrichment including TAU+TAS, PFC+TAS, and OCC. This analysis confirmed previous reports of enrichment of PCDH17 and CNTNAP2 (Figure 3B) in mid-fetal human frontal NCTX (Abrahams et al., 2007) and EPHA3 and EPHA7 (Table S9) in the fetal rhesus macaque monkey occipital and temporal NCTX, respectively (Sestan et al., 2001). However, the vast majority of our results revealed a complexity of expression patterns in the developing NCTX that has not been previously recognized in either humans (see Table S11), non-human primates, or rodents (see Table S12).

Open in a separate window

Figure 3

Unsupervised Hierarchical Clustering and qRT-PCR Validation of Selected Genes Differentially Expressed Between NCTX Areas

(A) With PFC areas grouped, correlated clusters (r>0.8) of gene enrichment included: (1 and 2) PFC; (3) temporal (TAU+TAS) +OCC; (4) temporal lobe (TAS+TAU); (5) PFC+TAS; and (6) OCC.

(B) We separately analyzed the four PFC areas plus MS, and identified genes enriched within specific frontal lobe areas. Selected clusters include (1) pan-PFC, but not MS; (2) OPFC; (3) orbital and lateral PFC; (4) MPFC; (5) MPFC+MS; (6) VLPFC+MS. Red is higher expression, blue is lower expression. Bar graphs show qRT-PCR confirmation (black bars; mean±SEM) next to array data (grey bars) (axes as in Figure 2) with a median ANOVA p<0.0004 and correlation to array data r=0.9. Complete lists of genes in these clusters are given in Table S7–Table S8.

Areal Differences in the Transcriptome of the Developing Prefrontal Cortex

Next, we hypothesized that small but significant gene expression differences might exist among functionally distinct PFC areas during development. Thus, in our second targeted intra-NCTX analysis, we compared the four PFC areas and the MS sample that partially clustered with them (Figure 3B). This analysis yielded 233 DEX genes (2.1% of those present; FDR 1%; Table S4 and Table S8), representing the first evidence for genetic differences between functionally distinct PFC areas in human or non-human primate developing brain. We found the most specific gene enrichment in OPFC, followed by MPFC; other clusters defined various combinations of PFC areas and MS. Interestingly, VLPFC, which encompasses the prospective Broca’s speech area and its right hemisphere homolog, was more closely correlated to MS than to other prefrontal areas at this developmental stage (r=0.36), suggesting a molecular similarity to the neighboring orofacial motor cortex controlling muscles involved in speech production. Moreover, genes enriched in VLPFC+MS included FOXP2 (Figure 3B, cluster 6), haploinsufficiency of which causes a severe speech and language disorder associated with morphological abnormalities and functional underactivation of Broca's area (Lai et al., 2003; Liégeois et al., 2003). This differential expression of FOXP2 within the human frontal lobe has not been observed in previous studies (Lai et al., 2003; Teramitsu et al., 2004), suggesting that the late mid-fetal stage may be a critical time-point for the role of this gene in establishing speech-related cortical circuitry. In addition, genes that clustered with FOXP2 may be candidates for further investigation in speech and language disorders (e.g., CNTNAP5, BCL6, C6orf206).

Finally, this result prompted us to search more specifically for genes enriched in the perisylvian cortical region that encompasses future speech and language-related areas: VLPFC, MS, PAS, TAS, and TAU. An analysis contrasting these areas with the remaining NCTX samples identified many slightly but significantly enriched genes (Table S10). FOXP2 showed a modest 1.1-fold increase in perisylvian areas, while the greatest enrichment was found for TRPC7, TRPC4, and the unknown locus DKFZp547H025 (1.3- to 1.4-fold enrichment). In addition, NR4A2, which was 1.2-fold enriched in perisylvian areas, has previously been found to be enriched in mid-fetal temporal NCTX (Abrahams et al., 2007). These results suggest that despite being dispersed across the frontal, parietal, and temporal lobes, the combined perisylvian cortical areas may express a developmental genetic signature related to their common involvement in speech and language in humans. Altogether, these patterns of expression suggest genetic programs for the development of PFC and perisylvian areas involved in higher cognitive functions, and thus represent promising candidate genes for evolutionary and functional analyses.

Validation and Detailed Cellular Mapping of Intra-Neocortical Differentially Expressed Genes

Genes identified as DEX within NCTX (Figure 3A) or PFC (Figure 3B) were chosen for further confirmation based on previous association with disorders such as autism, dyslexia, or speech and language impairments (e.g., CNTNAP2, ROBO1, FOXP2; Bakkaloglu et al., 2008; Hannula-Jouppi et al., 2005; Lai et al., 2003) or for apparent divergence from available rodent expression data (e.g., ANKRD32, CPNE8, POPDC3; see Table S12). We validated 75 intra-NCTX DEX genes by qRT-PCR (bar graphs in Figure 3; Table S7–Table S9), with a focus on PFC-enriched candidates, again finding a very high level of correlation between Exon Array and qRT-PCR results (median r > 0.9; median ANOVA p < 0.005). In addition, we confirmed 18 intra-NCTX DEX genes by ISH and/or IHC (Figure 4 and Table S9). These analyses revealed that despite the cellular heterogeneity of NCTX, Exon Array analysis was able to detect areal differences not only in genes exhibiting widespread expression within a neocortical area (Figure 4F) but also in genes whose expression differs in specific cell types, including astrocytes (Figure 4I), subplate neurons and marginal zone cells (Figure 4L), and cortical layer-specific neurons (Figure 4C). Thus, gene expression in various cell types contributes to neocortical areal molecular differences during development.

Open in a separate window

Figure 4

Confirmation and Cellular Mapping of Selected Neocortical Areal Expression Differences

Confirmation of array results by qRT-PCR, and the detection of areal- and cell type-restricted expression patterns by ISH and IHC.

(A–C) NPY enrichment in non-frontal areas is confirmed by qRT-PCR, with highest expression in temporal (TAS, TAU) and occipital lobes (OCC). (B) ISH on a whole sagittal section of 24 wg human brain confirms NPY enrichment in OCC. NPY enrichment in temporal cortex is not visible in this very medial tissue section. (C) Higher magnification reveals specific enrichment in the middle of the occipital cortical plate (CP), and high expression in scattered cells throughout the subplate (SP) (red arrowheads).

(D–F) CBLN2 enrichment in OPFC and lateral PFC is confirmed by qRT-PCR (D) and ISH on a whole sagittal section of 24 wg brain (E). (F) Higher magnification reveals that CBLN2 is enriched throughout the prefrontal CP and SP, but absent from the marginal zone (MZ).

(G–I) CNTNAP2 is selectively enriched in OPFC and lateral PFC areas. (H–I) IHC reveals specific diffuse enrichment in orbitofrontal SP and high expression in scattered MZ cells. Triple-immunofluorescent staining (I, middle panel) reveals colocalization of CNTNAP2 with astrocytic marker GFAP but not neuronal marker NeuN, suggesting differential expression of CNTNAP2 in SP astrocytes.

(J–L) FOXP2 is differentially expressed within the frontal cortex, and enriched in perisylvian cortex. IHC in coronal 22 wg brain sections (K–L) suggests that these differences are accounted for by a combination of higher cellular expression levels in the CP, particularly in OPFC (lower insets), and greater numbers of FOXP2-immunopositive SP cells, especially in VLPFC and PAS. Interestingly, strongly FOXP2-positive cells were present in the MZ in VLPFC and OPFC, but were completely absent from the MZ in other areas (upper insets). Bar graphs are mean±SEM. Scale bars in C, F, I, L: 500 µm.

Global Left-Right Symmetry of Late Mid-Fetal Neocortical Gene Expression

The human brain exhibits structural and functional left-right differences, a prominent example of which is the interhemispheric asymmetry of perisylvian NCTX underlying the functional lateralization in hand preference and speech and language processing (Galaburda et al., 1978). Subtle neocortical structural asymmetries first become evident during the late mid-fetal period analyzed in this study (Chi et al., 1977). Furthermore, a recent study by Sun et al. (2005) identified left-right differences in expression of LMO4 and other genes in the perisylvian regions of the human neocortex during the early mid-fetal period. To investigate whether such molecular correlates of cortical structural asymmetry persist into the late mid-fetal period, we compared gene expression and alternative exon usage between the left and right cortical hemispheres. Our Exon Array analysis did not detect any population-level hemispheric bias of gene expression or AS in individual neocortical areas or other analyzed brain regions. For example, we analyzed perisylvian areas (VLPFC, PAS, TAS, TAU) involved in speech and language processing, but were unable to detect significant differences either in perisylvian cortex as a whole (Figure 5B) or in individual areas (Figure S8). To further validate these findings, we performed qRT-PCR for over 120 genes displaying nonsignificant trends towards neocortical interhemispheric asymmetry on the Exon Arrays (data not shown). Only four of these genes showed a significant hemispheric bias (p<0.05; THBS1, SUMO2, CXorf1, and FAM5C), and all four of these results were driven by large but inconsistent asymmetry in specific neocortical areas, possibly reflecting inter-individual differences in these genes’ spatial or temporal regulation. For example, THBS1, a member of a gene family implicated in synaptogenesis (Christopherson et al., 2005) and human neocortical evolution (Cáceres et al., 2007), was ~2-fold enriched in the right VLPFC in two of the four brains tested. Overall, our results reveal a population-level global interhemispheric symmetry of gene expression in the late mid-fetal NCTX.

Open in a separate window

Figure 5

Population-level Global Left-Right Symmetry of Gene Expression in the Late Mid-Fetal Prefrontal and Perisylvian Neocortex

Volcano plots depicting the range of fold-differences and uncorrected p-values in representative perisylvian (A–B) and prefrontal (C–D) areas of NCTX. A large number of genes show significant differences in expression between NCTX areas (A and C). In contrast, no genes showed significant differences between left and right hemispheres within areas (B and D). Dashed lines indicate a corrected p-value of 0.01 (step-up FDR).

Spatial Regulation of Alternative Splicing

A major advantage of the Exon Array platform is the ability to probe individual exons within a transcript and thus test for tissue-specific expression of alternative isoforms. Analysis of probeset-level expression data identifies genes with a significant interaction between exon expression and tissue as candidates for DAS (Figure S9). From the 28% of genes significantly DAS across brain regions (Figure 1B), we selected a number of promising candidates, some with known alternative isoforms and some previously unknown to be alternatively spliced, and used exon- and isoform-specific qRT-PCR for validation (Figure 6). One of these, NTRK2, encodes multiple isoforms, including the full-length NTRK2a and the truncated NTRK2b (also known as TrkB-T1), which through distinct signaling pathways promote cortical neurogenesis or astrogliogenesis, respectively, in response to BDNF in mice (Cheng et al., 2007). We found that the truncated isoform, NTRK2b, is drastically and specifically downregulated in NCTX (Figure 6A), suggesting that the transition from neurogenesis to astrogliogenesis is delayed in the NCTX compared to other brain regions.

Open in a separate window

Figure 6

Validation of Selected Region-Specific Alternative Splicing Patterns

Validation of differential AS across brain regions (A–D) or neocortical areas (E–F) by exon- or isoform-specific qRT-PCR.

(A) Confirmation of a specific reduction of truncated NTRK2b in NCTX, whereas full-length NTRK2a is more evenly expressed throughout the brain.

(B) The LIMK2b splice variant, which lacks the N-terminal LIM protein-binding domain, is predominantly enriched in THM and CBL.

(C–D) CPVL and SKAP2 are each predicted to encode a previously unknown short isoform as a result of alternative promoter usage. Exon 2 of CPVL is drastically reduced in NCTX, while exons 5–6 are expressed in a complementary pattern (C). A more than 10-fold enrichment in STR of SKAP2 exon 13 is consistent with expression of a novel variant composed of exons 8–3 (D).

(E) ROBO1 encodes two known isoforms, the full-length ROBO1a and the alternative short isoform ROBO1b. In the developing human NCTX, ROBO1a is highly enriched in temporal lobe, while ROBO1b is slightly enriched in PFC.

(F) ANKRD32 is a previously uncharacterized gene that appears to encode two splice isoforms, ANKRD32a, which is evenly expressed across the NCTX, and ANKRD32b, which is significantly enriched in PFC. Bar graphs are mean±SEM. P-values represent the interaction between brain region or NCTX area and splice isoform.

Another DAS gene, LIMK2, is one of two LIM kinases that regulate actin dynamics and may be involved in neurite morphogenesis (Endo et al., 2007). LIMK2 has three known isoforms, at least two of which, LIMK2a and LIMK2b, are regulated by distinct promoters, encode distinct proteins, and display expression differences between human brain and other tissues (Nomoto et al., 1999). We found that within mid-fetal brain, LIMK2b is specifically enriched in THM and CBL (Figure 6B). Although functional differences between these two splice variants have yet to be characterized, region-specific splicing within the developing brain suggests that they may play distinct roles in neuronal morphogenesis. Interestingly, Nomoto et al. (1999) found evidence of a role for RORA in transcriptional regulation of LIMK2b, but not LIMK2a. Our data also supported such a relationship, with RORA and LIMK2b showing parallel enrichment in THM and CBL (Figure S10), suggesting the potential for Exon Array data to help predict transcriptional regulatory relationships at the level of individual splice variants.

We also identified DAS candidates with previously unknown alternative splice variants. For example, Exon Array data indicated differential expression of novel variants of CPVL and SKAP2, two genes with unknown functions in brain development (Figure S9). Our qRT-PCR data confirmed a specific reduction of CPVL exon 2 in NCTX (Figure 6C), and a more than 10-fold enrichment of SKAP2 exon 13 in STR (Figure 6D). Together, Exon Array and qRT-PCR data predict novel short isoforms of CPVL and SKAP2 consisting of exons 4–13 and 8–13, respectively, of their full-length transcripts (Figure 6 and Figure S9).

Finally, we identified and confirmed several intra-NCTX DAS candidates. These included ROBO1 (Figure 6E), a key axon guidance gene whose mouse ortholog is required for the formation of major cortical axonal projections (Andrews et al., 2006). It encodes two main isoforms, ROBO1a and ROBO1b (Dutt1), recently shown to be differentially expressed in the embryonic mouse brain (Nural et al., 2007). We found that human ROBO1a is enriched in temporal lobe, while ROBO1b is enriched in PFC (Figure 6E), suggesting that ROBO1 AS may be involved in patterning human intracortical connectivity.

In addition, many intra-NCTX DAS genes have not been previously charaterized. For example, we found that PFC enrichment of ANKRD32 (Figure 3A, cluster 2) reflects differential expression of a novel short isoform, ANKRD32b, whereas full-length ANKRD32a is expressed at lower levels uniformly across the NCTX (Figure 6F). Altogether, our evidence for regional differential splicing and predictions of novel isoforms provide insights into specific functional mechanisms of human brain development and suggest a wealth of biological hypotheses for future work. More generally, our Exon Array AS data illuminate a new level of complexity in the transcriptome of the developing human brain.

Regulatory Evolution of Genes Expressed in the Developing Human Brain

Recently, the completed genome sequences of humans, chimpanzees, and other species have been leveraged to identify highly conserved noncoding sequences (CNSs) that often act as transcriptional cis-regulatory elements (Pennacchio et al., 2006). A subset of these elements appears to have undergone human-specific accelerated evolution (haCNSs) (Prabhakar et al., 2006). We tested whether DEX genes are disproportionately associated with these haCNSs. We first identified the list of 6,921 genes nearest to CNSs, 10% of which are near haCNSs. We found that out of 1,045 intra-NCTX DEX genes near CNSs, 181 (17%) were near haCNSs (Table S13a), representing a highly significant enrichment of accelerated evolution in regulatory regions (p=7.15×10⁻¹⁶; Table 1). In contrast, in a control set of 985 genes equally highly expressed in NCTX but with no differential expression across areas, only 8% were found near haCNSs, a significant decrease relative to the genome-wide rate (p=0.0056). Thus, for NCTX-expressed genes, those with differential or restricted expression were twice as likely to be associated with accelerated human evolution of cis-regulatory elements. Consistent with this effect, a disproportionate 20% of genes most significantly enriched in perisylvian areas involved in speech and language processing were associated with haCNSs, including FOXP2, PCDH9, LPHN2, and SORCS3 (Table S10). In order to determine whether this effect was specific to NCTX or was more generally related to fine spatial control of gene expression during human brain development, we tested the top regionally DEX genes, as well as another control list of highly expressed, non-DEX genes. Again, DEX genes were more than twice as likely to be associated with haCNSs as non-DEX genes (16% versus 7%; Table 1 and Table S13b). Finally, because many genes are also associated with chimpanzee-accelerated CNSs (caCNSs; Prabhakar et al., 2006), we tested whether this trend was preserved among genes near haCNSs but not caCNSs. We again found a significant, though lesser, enrichment of haCNSs near DEX genes (Table 1, right-hand columns). Together, these results suggest that accelerated evolution of putative cis-regulatory elements is a feature of a subset of genes with highly specific expression in the developing human brain.

Biological Themes of Differential Gene Expression

To characterize the major biological themes present in regional DEX genes, we assigned Gene Ontology, protein domain, pathway, and other annotations to the five clusters illustrated in Figure 2, and tested for significant enrichment of functionally related terms (Table S14). Not surprisingly, most of the enriched terms involved neural development or function. In the cortex (HIP+NCTX; cluster 1), transcriptional regulation-related annotations were both the most significant (top four terms, all p<0.004) and most numerous (79% of enriched genes). These included transcription factors known to be important in the development of cortical projection neurons such as FEZF2, SATB2, and SOX5, as well as functionally uncharacterized factors such as TSHZ3 and ZNF238. Another over-represented class of genes enriched in cortex were those involved in axon guidance (e.g., SEMA3A, SEMA3C, EPHB6, UNC5D), which may reflect the transition of afferent axons from the subplate to the deep cortical plate at this stage (Kostovic, 1990). In contrast, in STR (cluster 3), the most enriched annotation clusters involved synaptic transmission or specific receptor signaling pathways (e.g., dopamine receptor and GPCR signaling), consistent with the more developed and neurochemically diverse striatal circuits present at this stage (Sajin et al., 1992).

We next annotated the intra-NCTX DEX gene clusters highlighted in Figure 3 and identified similar themes, including signal transduction, transmembrane proteins, and neurotransmitter receptor activity (Table S14). Notably, there was again a highly significant enrichment of axon guidance-related genes (e.g., NTN4, ROBO1, EPHA4, PLXND1, NTNG1), consistent with the notion that differences in neuronal connectivity are critical for the definition of distinct functional areas. However, there was little overlap between the axon guidance genes uniformly enriched in NCTX and those DEX between NCTX areas, suggesting that while the former may be responsible for the general establishment of cortical afferent and efferent projections, the latter may be involved in more specific targeting of neocortical area-specific circuitry.

Spatial and Functional Organization of the Human Mid-Fetal Brain Transcriptome

These analyses reveal novel developmental gene expression patterns corresponding to known anatomical and functional subdivisions of the human brain. While some spatial expression specificities may reflect underlying structural differences, others may reflect transient cellular developmental events, temporal neurogenic and maturational gradients, or spatially restricted signaling centers.

Our analysis has identified more DEX genes in the fetal human brain than have previously been reported in studies of adult human or embryonic mouse brain (Funatsu et al., 2004; Khaitovich et al., 2004; Roth et al., 2006; Kudo et al., 2007; Mühlfriedel et al., 2007). Multiple factors most likely contribute to this difference, including the increased sensitivity and genomic coverage of the Exon Array and other methodological differences. In addition, due to both the evolutionary differences of neocortical development between rodent and human (Preuss, 1995) as well as differences in developmental time-points surveyed, the larger number of DEX genes found in the present study may represent in large part a greater molecular diversity of human cortical areas and cell types. Furthermore, although methodological differences once again contribute, our finding of roughly two orders of magnitude more gene expression differences compared to the adult human NCTX suggests that prenatal differences in gene expression are more robust and complex than those present in the adult human brain.

Global Population-level Interhemispheric Genetic Symmetry of the Late Mid-fetal Neocortex

Structural left-right asymmetry is a prominent feature of the human NCTX and first appears during the late mid-fetal stage (Chi et al., 1977; Galaburda et al., 1977). However, we have found global population-level symmetry of gene expression during this period. Using the SAGE technique, Sun et al. (2005) reported significant asymmetric expression of the transcription factor LMO4 at 12 and 14 wg, but this difference was reduced at 16 to 17 wg and not detectable by 19 wg. Together, these data indicate that significant interhemispheric asymmetry of gene expression is likely a transient feature of the embryonic and early fetal NCTX, thus preceding by several weeks the structural asymmetry visible at late mid-gestation. If gene expression asymmetries are present in the late mid-fetal NCTX, they are likely limited to small differences in a few genes, and in a limited set of cell types.

Spatial Regulation of Alternative Splicing

Although the importance of AS in nervous system development has by now been well established, it has not previously been studied on a genome-wide scale in the developing human brain, nor has the prevalence of differential splicing between brain regions been comprehensively addressed. Our study has uncovered spatial patterns of enrichment of both known and novel splice variants. These data include, to our knowledge, the first evidence for intra-NCTX differential splicing, including enrichment of specific variants in the developing human PFC. One example is the axon guidance gene ROBO1, which is required for formation of major cortical connections in mouse (Andrews et al., 2006). Interestingly, a translocation breakpoint in the ROBO1 gene that has been associated with developmental dyslexia, a disorder linked to alterations in cortical circuits, results in loss of ROBO1a transcription while leaving ROBO1b unaffected (Hannula-Jouppi et al., 2005). In addition, recent work implicates differential splicing of mouse Robo3 in the midline crossing of spinal cord axons (Chen et al., 2008). Intriguingly, cortical midline (callosal) projections exhibit a rostrocaudal developmental gradient and prominent areal differences in the fetal rhesus macaque monkey (Dehay et al., 1988). Thus, although a role for AS in cortical axon guidance has not previously been identified, our finding of ROBO1 intra-NCTX DAS, together with several other lines of evidence, suggests that differential areal expression of axon guidance gene splice variants is likely an important mechanism of cortical circuit formation.

Our splicing analysis also identified 19 members of the protocadherin family of cell adhesion molecules (data not shown), consistent with previous reports of extensive splicing in this gene group (Wu and Maniatis, 1999). In contrast, neither DSCAM, famous for extensive AS in Drosophila (Schmucker et al., 2000), nor the related gene DSCAML1 appeared to be DAS in our analysis. We expect that the public availability of these Exon Array data will enable discovery and functional analysis of many more AS patterns and variants, generating a more complete picture of the transcriptional and posttranscriptional complexity of the developing human brain transcriptome.

Transcriptional Landscape of the Developing Prefrontal Cortex

We have identified more than 200 genes with putative expression differences within the mid-fetal human frontal lobe, many of which appear to be absent from or uniformly expressed in the developing mouse cortex (Table S12). These expression patterns may reflect species-specific differences in functionally specialized prefrontal areas, such as the complex social and emotional processing of the OPFC, and may also suggest new hypotheses regarding the genetic mechanisms controlling arealization of the human PFC. Many of these genes encode members of the same or related families of proteins, e.g., cerebellins; contactin associated proteins; and cadherins, suggesting a particular relevance of specific pathways or functions. In addition, several of the genes identified have been previously implicated in disorders that are thought to involve alterations of human PFC circuitry. These include CNTNAP2, which is not only related to language delay in autism, but is a target of the language-related transcriptional repressor FOXP2, and has a more general role as a susceptibility factor for specific language impairment (Vernes et al., 2008). Our data elaborate on these results, showing specific co-enrichment of CNTNAP2 and FOXP2 in OPFC and VLPFC. Thus, our study has uncovered complex spatial patterns of gene expression and AS that may reflect the underlying developmental, cellular, and species-specific differences between distinct PFC areas.

Implications for Clinical Research

Our data reveal previously unknown spatial expression patterns for many human disease-relevant genes. Furthermore, our data can help evaluate the results of genome-wide association or linkage studies by narrowing the focus to those genes that are specifically expressed or restricted to a relevant brain circuit during development. Finally, we contribute a rare resource in the form of whole-genome genotyping and expression data from the same individuals, enabling correlateion of copy number variation to expression levels across different regions of the developing human brain.

RNA Isolation, Processing, and Microarray Hybridization

Total RNA was extracted using TRIzol (Invitrogen), followed by treatment with RNeasy Mini Kit (Qiagen) to exclude smaller RNAs. The quality of total RNA was evaluated by 2100 Bioanalyzer (Agilent) and RNA 6000 Nano Kit (Agilent) before being processed with the Affymetrix GeneChip Whole Transcript Sense Target Labeling Assay and hybridized to the Affymetrix Exon 1.0 ST Arrays following the recommended Affymetrix protocols. Hybridized arrays were scanned on an Affymetrix GeneChip Scanner 3000 and visually inspected for hybridization artifacts.

Exon Array Data Analysis

Exon Array data were preprocessed using standard RMA normalization, DABG, and probeset summarization methods in either Partek Genomics Suite (Partek, Inc.) or the Excel Array Analysis software (XRAY; Biotique Systems, Inc.). Principal component analysis, left-right hemisphere ANOVAs, and t-tests were performed in Partek using gene summary values for all core transcript clusters. Global DEX and DAS ANOVAs were performed in XRAY using default parameters. All ANOVAs included brain specimen and date of hybridization as co-factors, to eliminate batch effects and variations due to individual genetic differences. All p-values were corrected for multiple comparisons using the FDR step-down method. Unsupervised hierarchical clustering was performed in Cluster 3.0 (bonsai.ims.u-tokyo.ac.jp/~mdehoon/software/cluster), and heatmaps were generated using Java Treeview (jtreeview.sourceforge.net). For further details, see the Supplemental Experimental Procedures.

Analysis of Patterns of Cis-Regulatory Evolution and Gene Expression

Genomic coordinates of all ~110k CNSs identified by Prabhakar et al. (2006) were mapped to the human genome (hg18; NCBI build 36.1) on the UCSC Genome Browser (genome.ucsc.edu) and cross-referenced with the RefSeq Genes track using Galaxy (main.g2.bx.psu.edu) to identify the nearest human RefSeq gene to each CNS. After removing duplicates, this yielded 6921 genes, of which 694 were near CNSs reported as showing evidence of accelerated evolution in the human lineage (haCNSs). We then intersected these gene lists with lists of DEX and non-DEX genes, calculated the proportion of haCNSs for each condition, and assigned a p-value according to the hypergeometric distribution.

Gene Ontology Annotation Analysis

Annotation analysis was performed using the web-based DAVID software (david.abcc.ncifcrf.gov; Dennis et al., 2003). Intra-NCTX DEX gene clusters from Figure 3 were grouped together to control for the length of annotated gene lists and allow direct comparison with the annotation of Figure 2 clusters. See Supplemental Experimental Procedures for details.

We thank Bradford Poulos for assistance with tissue acquisition, Aiping Lin for help with microarray platform comparisons, James Noonan for advice on analysis of haCNSs, Fuying Gao for statistical analysis, Donna Karolchik and Andy Pohl for assistance in creating tracks for the Genome Browser, and many colleagues for their help and comments. This work was supported by NIH grants NS054273, MH081896 (to N.S.), NS056276 (to M.W.S.), NS051869 (to S.M.M.), MH060233 (to D.H.G.), Kavli Foundation, NARSAD, and the James S. McDonnell Foundation Scholar Award (to N.S.).