All experimental protocols involving animals and their embryos were performed in accordance with the guidelines and regulations of the Korea Institute of Science and Technology. All experimental protocols were approved by the Institutional Animal Care and Use Committee at Korea Institute of Science and Technology (approval number: 2017-038).

Sprague-Dawley rats aged P14, P28 and P45 were used for this study. Rats were anesthetized with halothane (2-bromo-2-chloro-1,1,1-trifluoroethane, Sigma Aldrich), and hippocampi were isolated. The hippocampus, including the dentate gyrus, CA3, CA2, CA1 and subiculum, was equally trisected to divide into a ventral, intermediate and dorsal portion, as shown in Figure 1A. Only the ventral and dorsal parts were collected for the study, and these ventral and dorsal hippocampal tissue lysates were used for total RNA isolation and subsequent RNA-seq and qRT-PCR analyses.

RNA-seq was performed as previously described (Joe et al., 2017). Briefly, total RNA was extracted using an RNeasy kit (Qiagen, Valencia, CA, USA), and the RNA-seq library was then prepared using a TruSeq RNA Sample Prep kit (Illumina, San Diego, CA, USA). The sequencing was performed based on the Illumina NextSeq500 platform to generate 150 bp paired-end reads. The sequenced reads were mapped to the Rat genome (rn5) using TopHat 2, and the raw read counts were obtained using HT-seq (v.0.6.0; Anders et al., 2015). The edgeR (v. 3.12.1; McCarthy et al., 2012) was used to quantify gene expression from the RNA-seq data. Transcript expression levels were calculated as the fragments per kilobase of transcript per million fragments mapped (FPKM). The FPKM value of each gene was scaled via the trimmed mean of M-values (TMM) across all libraries for further analysis. In each analysis, differently expressed genes were selected through edgeR (v. 3.12.1). The criteria for determining a significantly altered gene were a fold change greater than 1.5 and a false discovery rate (FDR) less than 0.05. Heatmaps were constructed using MeV software (Howe et al., 2011). Statistical analyses and graph construction were performed using R (v. 3.3.0) and PYTHON (v. 2.7.6). The mapping statistics for RNA-seq were summarized in Supplementary Table S3.

RNA-seq data have been deposited in the National Center for Biotechnology Information Gene Expression Omnibus (GEO) under the accession number GSE96079; for a link for reviewers: https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?token=qlmjcoacxjonnmx&acc=GSE96079.

For the classification of the up-regulated genes over development, significantly increased genes between P14 and P45 (P45/P14 > 1.5 and FDR < 0.05) were first selected for each dorsal and ventral hippocampus group. Among the genes, the up-regulated genes (P28/P14 > 1 and P45/P28 > 1) were selected by the comparison with the intermediate step, P28. For the classification of the down-regulated genes over development, significantly decreased genes between P14 and P45 (P14/P45 > 1.5, FDR < 0.05) were selected for each of the dorsal and ventral hippocampus group. Among the genes, the down-regulated genes (P14/P28 > 1 and P28/P45 > 1) were selected (Figure 2).

The dorsally or ventrally enriched genes were categorized by applying the same criteria (fold change 1.5, FDR < 0.05) for each time point (P14, P28 and P45), respectively (Figures 3–5). For highly enriched genes dorsally or ventrally (Figure 4D), the genes with a fold change of ventral hippocampus/dorsal hippocampus (VH/DH) ≥4.0 or ≤−4.0 at all developmental stages (P14, P28 and P45) were selected.

The Database for Annotation, Visualization and Integrated Discovery (DAVID) Bioinformatics Resources v6.8¹ (Huang da et al., 2009a,b; Joe et al., 2017) was applied for the GO analysis on the basis of three categories, including biological process, cellular component and molecular function, and also for the Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis (Kanehisa and Goto, 2000; Kanehisa et al., 2014). For the disease analysis, rat genes belonging to the dorsally or ventrally enriched categories were converted to human homologous genes, which were then applied to the KEGG mapper tool², a disease searching tool.

In situ hybridization data were obtained from the coronal or sagittal ABA images of the hippocampal regions, including Ammon’s horn and the dentate gyrus across the long hippocampal axis using mouse ABA Brain Explorer 2 (Lein et al., 2007; Lau et al., 2008). Dorsally enriched genes were marked in cyan, and ventrally enriched genes were marked in red. For all ABA images presented in this study, the intensity was ranged between 124 and >260, and the density was ranged between 0.0248 and >0.1, except for the following genes: Cyp26b1 (Figure 1E) and Cpne2 (Figure 1E) with intensities adjusted between 98 and >260; Tacr3 (Figure 4F, Supplementary Figure S6) with intensity adjusted between 79 and >166 and density adjusted between 0.00414 and >0.1; Kcna1 (Supplementary Figure S6) with intensity adjusted between 159 and >260; Otof (Supplementary Figure S6) with intensity adjusted between 79 and >154 and density adjusted between 0.0152 and >0.0779; Tmie (Supplementary Figure S6) with intensity adjusted between 108 and >154 and density adjusted between 0.0414 and >0.1; and Pfn1, Sbf1, Scn3b and Tuba4a (Supplementary Figure S6) with densities ranging between 0.0655 and >0.1 for clearer visualizations.

The isolated dorsal and ventral hippocampal tissues were homogenized using a Dounce Tissue Grinder (Wheaton), and subsequent total RNA extraction steps were followed as described in the manufacturer’s instructions. Briefly, total RNAs were purified from the isolated tissues using RNAiso Plus (TaKaRa, Japan). The cDNA was reverse transcribed from total RNA using a PrimeScript II 1st strand cDNA synthesis kit (TaKaRa, Japan). The qRT-PCR was performed using Power SYBR Green PCR Master Mix (Thermo Fisher Scientific). The reaction mixture contained 2 μL of cDNA corresponding to 100 ng of total RNA, 300 nM of each gene-specific primer, and 2× Power SYBR Green PCR Master Mix in a total volume of 20 μL. The cycling parameters of the StepOnePlus Real-Time PCR System (Applied Biosystems) were as follows: 95°C for 10 min, followed by 40 cycles of 95°C for 15 s and 60°C for 1 min. Relative expression levels were calculated according to the 2^−ΔΔCT algorithm that was internally normalized to the value of the glyceraldehyde 3-phosphate dehydrogenase (GAPDH) gene, which did not display differential expression among the samples in this study. Experiments were performed with three to five independent samples. The primer sets for the 40 genes used for qRT-PCR validations are listed in Supplementary Table S4.

To obtain a dorsally or ventrally enriched index, the relative expression level (2^−ΔΔCT) of each sample (total of six samples: dorsal and ventral hippocampal samples at the three developmental stages represented on P14, P28 and P45) for each gene was normalized to that of the dorsal hippocampal sample at P14. Then, the normalized expression value in the dorsal hippocampus was subtracted from that in the ventral hippocampus at each developmental age. Genes with index values <0 or >0 were defined as enriched dorsally or ventrally, respectively.

The dorsal and ventral hippocampus have been shown to be functionally segregated (Moser and Moser, 1998; Bannerman et al., 2003; Lein et al., 2007; McHugh et al., 2011), and recent RNA-seq studies have provided comprehensive bioinformatic data supporting the notion of dorsal and ventral heterogeneity in the hippocampus (Zeisel et al., 2015; Cembrowski et al., 2016a,b). However, no RNA-seq–based bioinformatic study of the developing dorsal and ventral hippocampus has been reported. To investigate whether the dorsal and ventral hippocampus exhibit differential transcript profiles during postnatal development, and if so, whether the genes involved are related to biologically functional pathways and/or to neural diseases associated with these brain regions, we performed RNA-seq using 12 samples (in duplicate), which consisted of dorsal and ventral hippocampal tissues at three different postnatal developmental stages, P14, P28 and P45. We first tested whether the samples were appropriately prepared for this study by detecting established dorsal and ventral markers. The gene expression levels of the RNA-seq data were calculated as the fragments per kilobase of transcript per million fragments mapped (FPKM). Gene markers reported in previous studies of the dorsal (Nts, Cyp26b1) and ventral (Nr2f2, Cpne2) hippocampus (Leonardo et al., 2006; O’Reilly et al., 2015; Cembrowski et al., 2016a,b) were readily detected by qRT-PCR (Figure 1B) and RNA-seq (Figures 1E,F) with our experimental preparation and further cross-validated using the in situ hybridization results from the mouse ABA (Figure 1E). In addition, clustering analysis (Figure 1C) and multidimensional scaling (Figure 1D) of RNA-seq data showed sufficient clustering and similarity between the duplicate samples. These data together indicate that our samples are appropriately prepared for the study.

We categorized the RNA-seq–based transcript expression patterns into eight groups (Figure 2A): group 1 (or 2) for genes that increase (or decrease) expression over development in the dorsal hippocampus; group 3 (or 4) for genes that increase (or decrease) expression over development in the ventral hippocampus; group 5 for overlapping genes in groups 1 and 3; group 6 for overlapping genes in groups 2 and 4; group 7 for overlapping genes in groups 2 and 3; group 8 for overlapping genes in groups 1 and 4. To gain a better understanding of the gene expression properties in each group, we used a bioinformatic functional annotation tool, Database for Annotation, Visualization and Integrated Discovery (DAVID), to retrieve annotations from the gene ontology (GO) analysis based on the three categories, including biological process, cellular component, and molecular function (Ashburner et al., 2000) and from the Kyoto encyclopedia of genes and genomes (KEGG) pathway enrichment analysis (Kanehisa and Goto, 2000; Kanehisa et al., 2014) to identify enriched metabolic or signaling pathway gene patterns. Learning-related GO terms, such as learning or memory, long-term synaptic potentiation and learning were selected in groups 1 (DH-up), 3 (VH-up) and 4 (VH-down; Supplementary Figure S1). Genes belonging to groups 2 (DH-down) or 4 (VH-down) were associated with biological processes such as axon guidance, brain development, and/or neurogenesis (Supplementary Figure S1), indicating that transcript levels of genes belonging to those biological processes decrease over development, which is reasonable since axon guidance, brain development, and neurogenesis are predominant processes in early rather than later developmental stages.

Genes belonging to group 5 (DH-up and VH-up; Figures 2A,B) were associated with GO terms, such as oxidation-reduction process, regulation of neuronal synaptic plasticity, and long-term synaptic potentiation (Supplementary Figure S2A). Genes belonging to group 6 (DH-down and VH-down; Figures 2A,C) were associated with GO terms, such as brain development, cell migration, axon guidance, and response to axon injury, while also revealing their association with the KEGG pathways, such as extracellular matrix (ECM)-receptor interaction, axon guidance, microRNAs in cancer, and gap junctions (Supplementary Figure S2B). Only two genes (Mei1, Phkg2) were identified for group 7 (DH-down and VH-up), and no genes were identified for group 8 (DH-up and VH-down; Figure 2A, Supplementary Figure S1E).

We next identified genes with expression levels enriched in the dorsal or ventral hippocampus at the different developmental ages. Genes with fold changes in their FPKM value of ventral to dorsal hippocampus ≥1.5 or ≤−1.5 (FDR < 0.05) were defined as ventrally or dorsally enriched, respectively. A total of 362 and 355 genes at P14, 333 and 381 genes at P28, and 314 and 331 genes at P45 were selected as dorsally and ventrally enriched genes, respectively (Figure 3A). All genes belonging to dorsally or ventrally enriched groups at each developmental age are also listed in Supplementary Table S1.

We then conducted the GO enrichment and KEGG pathway analyses with the dorsally or ventrally enriched genes at individual developmental ages (Figure 3) in an attempt to identify the associated biological processes and pathways (Supplementary Figures S3, S4). Consistent with previous reports on dorsal hippocampus-dependent learning (Moser et al., 1993, 1995; McHugh et al., 2011), GO terms in biological processes, such as learning, long-term synaptic potentiation, and long-term memory at P14 and P28, and visual learning at P28, and the KEGG pathway long-term potentiation (Kouvaros and Papatheodoropoulos, 2016) at P14, P28 and P45 were associated with genes that were enriched in the dorsal hippocampus at individual ages (Supplementary Figure S3). The enrichment of genes belonging to the dorsal hippocampus-enriched KEGG pathway long-term potentiation were further validated by qRT-PCR (Figure 3B). Cholinergic synaptic transmission at P14, P28 and P45, the serotonin receptor signaling pathway at P14, P28 and P45, and morphine and nicotine addiction at P14, P28 and P45 were associated with genes enriched in the ventral hippocampus at each age (Supplementary Figure S4). Cholinergic synaptic transmission (Chrm4, Chrm5, Chrna4, Chrna7, Htr3a, Lypd1) and GABAergic synapse (Adcy8, Gabra3, Gabrb1, Gabrq, Gnai1, Gng13, Hap1, Slc6a1) associated with the ventrally enriched genes were further validated by qRT-PCR (Figure 3C), confirming the ventral hippocampus association with GABAergic and cholinergic synaptic transmission.

We next selected the genes with fold changes in their FPKM value of ventral to dorsal hippocampus >1.5 (ventrally enriched; FDR < 0.05) or <−1.5 (dorsally enriched; FDR < 0.05) at all developmental stages (P14, P28 and P45). A total of 120 and 144 genes were selected as dorsally and ventrally enriched genes, respectively, at all ages (Figure 4A; see also Supplementary Table S2). Among these genes, those with a fold change >4.0 or <−4.0 at all developmental stages were selected as highly enriched genes in the ventral or dorsal hippocampus, respectively (Figure 4D). From these groups of genes, 7 (cyan font; Figure 4D) and 11 (red font; Figure 4D) genes were potentially identified as dorsal and ventral markers, respectively. Among the 7 cyan and 11 red font genes in Figure 4D, to further ensure the correct identification of markers, genes that were cross-validated with ABA in situ hybridization data were selected as novel dorsal (Figure 4E) and ventral (Figure 4F) markers for all developmental stages. A total of four genes as dorsal markers and six genes as ventral markers were discovered (Figures 4E,F). Further bioinformatic analysis revealed that biological processes, such as long-term memory (Adcy1, Camk4, Rgs14), the ephrin receptor signaling pathway (Epha5, Epha7, Ntrk1), and circadian rhythm (Adcy1, Ntrk1, Kcnh7, Tph2; Molero-Chamizo, 2013) are associated with the dorsal hippocampus (Supplementary Figure S5), and the genes associated with these processes were validated by qRT-PCR (Figure 4B). In addition, the ventral hippocampus was shown to be significantly associated with biological processes, such as neuropeptide signaling pathways (Grp, Galr1, Npy2r, Nppa, Oprl1, Penk, Rxfp3, Sstr1) and behavioral fear responses (Htr2c, Lypd1, Npy2r, Penk), and the genes belonging to these processes were validated by qRT-PCR (Figure 4C).

Functional differentiation between the dorsal and ventral regions of the hippocampus has been experimentally and bioinformatically demonstrated in previous studies, strongly suggesting that hippocampal regions may have differential roles in distinct diseases. Indeed, schizophrenia, anxiety and depression have been associated with altered function of the ventral hippocampus (Tseng et al., 2009; Fanselow and Dong, 2010; Brooks et al., 2011; Small et al., 2011; O’Reilly et al., 2015).

Therefore, we further analyzed the genes dorsally or ventrally enriched in the postnatally developing hippocampus using the KEGG mapper tool, a disease searching tool (Kanehisa, 2013). Most diseases selected by the KEGG mapper tool, except hypogonadotropic hypogonadism, were distinctly segregated between the dorsal or ventral hippocampus (Figure 5A). Consistent with previous reports (Tseng et al., 2009; Brooks et al., 2011; O’Reilly et al., 2015), schizophrenia was associated with the genes ventrally enriched in the hippocampus (Figure 5A). We further identified that Charcot-Marie-Tooth disease, episodic ataxias spinocerebellar ataxia and amyotrophic lateral sclerosis are associated with gene profiles in the dorsal hippocampus, whereas early infantile epileptic encephalopathy, autosomal recessive deafness, Usher syndrome and Joubert syndrome are associated with those in the ventral hippocampus (Figure 5A). In addition to these neural diseases, type II diabetes mellitus and atrial fibrillation were associated with genes in the dorsal hippocampus while Ehlers-Danlos syndrome, primary ciliary dyskinesia and glycogen storage diseases were associated with genes in the ventral hippocampus (Figure 5A). We next cross-validated the expression patterns of the disease-related genes in Figure 5A with the ABA in situ hybridization data (Supplementary Figure S6). Based on the genes that were validated using the ABA data and our information on the associated diseases (Figure 5A), we generated a model of disease-related molecular domains specific to the dorsal and ventral hippocampus (Figure 5B, Supplementary Figure S7).