Introduction
Dysregulation of the striatum is linked to multiple neuropsychiatric diseases, including Huntington’s (HD), Parkinson’s, X-linked dystonia-parkinsonism, addiction, autism, and schizophrenia. The dorsal striatum comprises the caudate and putamen in humans but consists of a single nucleus in the mouse. This nucleus is a key component of cortical and subcortical circuits regulating movement, reward, and aspects of cognition, including speech and language. 85 to 95% of striatal neurons are medium spiny (MSNs), the striatal projection neuron. They are morphologically homogeneous, but phenotypically heterogeneous, and adult subtypes are distinguished by unique transcriptomes (Anderson et al., 2020; Gokce et al., 2016; Ho et al., 2018; Märtin et al., 2019; Muñoz-Manchado et al., 2018; Ortiz et al., 2020; Stanley et al., 2020; Zeisel et al., 2018). MSNs are equally distributed between direct neurons (dMSNs), which express the dopamine D1 receptor (D1R) and project directly to the substantia nigra (SN) or to the internal segment of the globus pallidus, and indirect neurons (iMSNs), which express the dopamine D2 (D2R) and adenosine 2A receptors (A2aR) and project to the external segment of the globus pallidus or to the subthalamic nucleus (Keeler et al., 2014).
The striatum is comprised of two main compartments named striosomes (patch) and matrix. Striosomes occupy 10–15% of the volume and are dispersed in a continuum throughout the 80–85% occupied by the matrix. Importantly, an imbalance between striatal compartments likely contributes to movement disorders (Cazorla et al., 2015; Crittenden and Graybiel, 2011; Keeler et al., 2014), and compartmentation also appears to be required for non-motor functions, for example, speech and language and discrimination learning (Campbell et al., 2009; Deriziotis and Fisher, 2017; Feuk et al., 2006; Konopka et al., 2009; Lai et al., 2001; MacDermot et al., 2005). In X-linked dystonia-parkinsonism, degeneration is initiated in the striosomes (Beste et al., 2018; Goto et al., 2005). Compartmentation is altered in autism (Kuo and Liu, 2020), and HD features early, preferential loss of striosome neurons (Victor et al., 2014). Moreover, dopaminergic signaling has opposing effects on D1R activation in a compartment-specific manner, which is important for task dependent behaviors (Prager et al., 2020). The distinction between the striosome and matrix in adults is based on differences in gene expression, the origins of afferents from cortical regions (e.g. sensorimotor, limbic, and associative), and to some extent, by the destination of their efferents (Brimblecombe and Cragg, 2017; Crittenden and Graybiel, 2011; Fujiyama et al., 2019). They both contain direct and indirect neurons, but the striosome dMSN content regionally varies (Cirnaru et al., 2019; Miyamoto et al., 2018). Currently, no robust differentiation protocol is available for the generation of striosome cells for in vitro studies and replacement therapies (Arber et al., 2015; Golas, 2018; Kemp et al., 2016; Richner et al., 2015; Telezhkin et al., 2016; Victor et al., 2014), and the molecular mechanisms governing compartmentation are incompletely defined.
Striosome and matrix compartments differ in their developmental timelines, as there are two waves of striatal neurogenesis (Matsushima and Graybiel, 2020). Striosome MSNs complete mitosis by embryonic day 13 (E13) and project toward substantia nigra pars compacta by E18 (Fishell and van der Kooy, 1987; Fishell and van der Kooy, 1989). Matrix MSNs are born in the second wave and are postmitotic by E18–E20 and do not project toward substantia nigra until the first postnatal week (Fishell and van der Kooy, 1987; Fishell and van der Kooy, 1989). This period is essential for the specification of striatal compartmentation, and most of the markers that distinguish striosome from matrix do not reach their final adult distribution until at least the second postnatal week. For example, dopamine- and cAMP-regulated neuronal phosphoprotein (DARPP-32), encoded by
Figure 1.
Schematic of the developmental cascade for the striosome and matrix compartments and our experimental approach.
Upper panel shows the cascade of the two compartments in the developing mouse. Lower panels show our experimental approach. We FACs sorted EGFP+ cells from Nr4a1 EGFP mice at PND3, performed RNAseq, identified top transcription factors in the striosome, validated them in vivo (mice) and in vitro (primary mouse cultures and human iPSCs).
Multiple transcription factors (TFs) are required for induction and differentiation into MSN subtypes (Arlotta et al., 2008; Fjodorova et al., 2019; Golas, 2018; Ivkovic and Ehrlich, 1999; Long et al., 2009; Marin et al., 2000; Martín-Ibáñez et al., 2012; Martín-Ibáñez et al., 2017; Precious et al., 2016; Victor et al., 2014; Wang et al., 2011b; Zhang et al., 2016), and some are useful for direct conversion of human fibroblasts into MSNs (Hedreen and Folstein, 1995; Lawhorn et al., 2008). Importantly, a widely used protocol yielded only a matrix phenotype (Victor et al., 2014), and several others yielded largely calbindin-positive neurons, with no mention of a striosome phenotype (Adil et al., 2018; Arber et al., 2015; Golas, 2018; Kemp et al., 2016; Telezhkin et al., 2016). Only a handful of TFs and signaling systems contributing to striosome/matrix compartmentation have been identified, including
We used unbiased transcriptomic and epigenetic assays, RNA-seq and ATAC-seq, to identify TFs and open chromatin regions (OCRs) associated with striosome maturation on postnatal day 3, specifically after they are first formed but not yet mature (Figure 1). Genes expressed at this point may also include terminal differentiation effectors as defined by Hobert, 2016, which initiate and maintain the adult identity of neurons. We used the GENSAT
Results
Transcriptional analysis of
The TF
Figure 2.
Identification of striatal striosome and matrix compartment-specific transcription factors using RNA-seq analysis in FACs EGFP+ and EGFP– cell populations from PND3
(A) Coronal section of PND3
Figure 2—figure supplement 1.
Validation of
(A) Spontaneous EGFP fluorescence (green) co-localizes with immunolabeling of the striosome marker Oprm1 (magenta) in coronal sections of
Table 1.
List of striosome and matrix markers.
| Gene symbol | Gene name | PND3 | Adult |
|---|---|---|---|
|
| Calbindin | Striosome | Matrix |
|
| Ephrin A5 | Striosome | |
|
| Eph Receptor A4 | Matrix | Matrix |
|
| Eph Receptor A7 | Striosome | |
|
| Forkhead Box P1 | Striosome | pan-MSN |
|
| Forkhead Box P2 | Striosome | Striosome |
|
| Nuclear Receptor Subfamily 4 Group Member 1 | Striosome | Striosome |
|
| Mu opioid receptor | Striosome | Striosome |
|
| Protein phosphatase 1 regulatory subunit 1B | Striosome | pan-MSN |
|
| CalDAG GEFII | Striosome | Striosome |
|
| CalDAG GEFI | Matrix | Matrix |
Next, we carried out RNA-seq, transcriptomic analysis on both the EGFP+ and EGFP– cells (Figure 2—source data 1). A total of 9124 genes were differentially expressed between the two groups, with 4714 enriched in the EGFP+ cells (positive log2FC) and 4410 enriched in the EGFP– cells (negative log2FC) (p<0.01, Figure 2—source data 1), but this entire list includes many genes with small differences between compartments. Twenty-three differentially expressed genes were selected with an arbitrary threshold of log2FC greater than seven and mean normalized counts greater than 40, which along with
The
To further highlight compartmentalized differentially expressed TFs, we used the mouse TF database (Figure 2—source data 1). The top 60 TFs enriched in each population according to their log2FC [i.e., EGFP+ population relative to the EGFP– (upper panel), and in the EGFP– population (lower panel), without regard to mean normalized counts] are shown in Figure 2F. TFs enriched in the EGFP+ group include
The ongoing maturation of the two compartments was further highlighted by the localization of other mRNAs recently associated with striosome, matrix, and exopatch (Ortiz et al., 2020; Saunders et al., 2018; Smith et al., 2016). The fact that some of these are not yet compartmentalized and/or not present in the data base confirms that both neuronal maturation and compartmentation are incomplete on PND3. Thus, adult and/or PND9 striosome markers
As the
Analysis of striosome-enriched transcription networks define
To identify potential master regulators driving the transcription program of the striosome, we conducted TF and co-expressor enrichment analysis with Enrichr (Chen et al., 2013; Kuleshov et al., 2016). We input the list of either differentially expressed genes or differentially expressed transcription factors to obtain the upstream regulators and co-expressors for both gene sets (p<0.05). The intersect of the enrichment results and the differentially expressed TFs were further used as input to generate striosome- and matrix-specific TF co-expression networks, using GeneMANIA (Warde-Farley et al., 2010) as a prediction tool (Figure 3A,B). This allowed the generation of network reconstruction and expansion for striosome and matrix differentially expressed TFs. From the network analysis, we observed that striosome-enriched
Figure 3.
Transcription networks and enrichment maps of striosome and matrix.
(A,B) TF networks enriched in either striosome or matrix show interactions and co-expression among the regulatory genes. An arbitrary threshold of absolute log2 fold-change > five was set to designate hub genes. (C,D) Co-expression networks derived from striosome and matrix RNA-seq data show targets of
In addition to defining the transcriptional network for the striosome and matrix, we analyzed the differentially expressed genes for overrepresented biological annotations or pathways, including GO term and KEGG pathway enrichment analysis (The Gene Ontology Consortium, 2017; Figure 2—source data 1). Striosome enriched genes were associated with general development, for example, integrins, angiogenesis or cell motility, and nodes include ‘regulation of pre-synapse organization’. Matrix-enriched genes were associated with neural system development and neurotransmission (Figure 3E). Both
ATAC-seq analysis defines compartment-specific OCRs in the striatum
To map chromatin accessibility of striosome and matrix cells, we used FACS, followed by preparation of eight ATAC-seq libraries from 12 animals pooled in four independent samples for the EGFP+ and EGFP– populations. Overall, we obtained 322 million (average of 40.3 million) uniquely mapped reads after removing duplicates and those aligning to the mitochondrial genome (Figure 4—source data 1). To quantitatively analyze differences between striosome and matrix cells, we generated a consensus set of 69,220 OCRs by taking the union of peaks called in the individual cells (Materials and methods). We next quantified the number of reads that overlapped each OCR. Uniform Manifold Approximation and Projection for Dimension Reduction (UMAP)-based clustering using the normalized read counts clearly separated striosome from matrix samples (Figure 4A). Comparison of striosome and matrix peaks resulted in 44% of OCRs that were significant after multiple testing corrections (30,799 differentially modified OCRs at FDR ≤ 0.05) (Figure 4—source data 2). Among these, 16,963 were striosome, and 13,836 were matrix.
Figure 4.
Open chromatin region (ATAC-seq) analysis identifies OCRs differentially located in striosome or matrix.
(A) Clustering of the individual samples using UMAP. (B) Distribution of the location within the gene (see Materials and methods, Annotation of OCRs) of all and differentially distributed OCRs. (C) Correlation of log2 fold-changes in RNA-seq and ATAC-seq analyses. OCRs within TSS were considered for this analysis, and only genes with adjusted p-value <0.01 from RNA-seq analysis were included. (D)
Figure 4—figure supplement 1.
Analysis with ANANSE integrating RNAseq and ATACseq identifies key TFs required for striosome MSN differentiation.
(A) ANANSE analysis defines key transcription factors in striosomes. ANANSE is a network-based method that uses properties of enhancers and their gene regulatory networks to predict key TFs. The influence score vs. log 2Fc demonstrates
We examined the location of OCRs with respect to the distance from transcription start sites (TSS) and genic annotations. As expected, OCRs are in the vicinity of TSSs (Figure 4B) but are also enriched for non-promoter regulatory elements, suggesting a more important role for long-range regulation of gene expression in developing striosome and matrix. We tested the concordance of striosome and matrix-specific genes and OCRs, based on the RNA-seq and ATAC- seq analyses, respectively, and found a moderately high correlation (Pearson’s R = 0.34; p-value < 2.2x10-16) (Figure 4C). To determine if the ATAC-seq identified regulatory regions relevant to
Figure 5.
(A) Representative confocal micrographs of Foxf2 (green), DARPP-32 (magenta) and DAPI (blue) immunolabeling on 16-µm-thick coronal sections from wild-type PND3 showing the localization of Foxf2 within the DARPP-32-immunopositive striosomes. Scale bars correspond to 20 µm. (B) RNAscope analysis in 16-µm-thick coronal sections from wild-type PND3 mice shows
Figure 5—figure supplement 1.
(A) Low-power magnification of RNAscope assay visualization of
Figure 5—figure supplement 2.
Foxf2 overexpression in medium spiny neurons promote their maturation in vitro.
(A) Representative DARPP-32 immunostaining in DIV9 WT primary striatal neurons transduced at DIV5 with ADV-GFP (control) or ADV-
Conversely, endogenously regulated, early overexpression of
We also assayed the effect of Foxf2 on MSN maturation in vitro. Similar to early postnatal striatum in vivo, adenoviral-mediated overexpression of human FOXF2 promoted an increase in the percentage of DARPP-32-immunopositive cells in primary mouse neuronal cultures derived from embryonic striatum (Figure 5—figure supplement 2A) There was also an increase in striosome (
Next, we evaluated the role of
Figure 6.
Olig2 is required for striosome compartmentation and a downstream intergenic
(A) Confocal micrographs of Olig2 (green), DARPP-32 (magenta) and DAPI (blue) immunolabeling of WT PND3 16-µm-thick coronal sections showing localization of Olig2 within the DARPP-32-positive striosomes (unfilled arrowhead). An examples of a single Olig2-immunopositive and DARPP-32-negative cell is indicated by a filled overhead. Scale bars = 20 µm (B) RNAscope visualization of
Figure 6—figure supplement 1.
Impact of neonatal AAV2 driven Olig2 overexpression on the striatum of PND21 mice and analysis of Olig2 OCR-driven expression in adult founders.
(A) Low- and high-power magnification of Olig2 (green), DARPP-32 (magenta) and DAPI (blue) immunolabeling of WT PND3 coronal sections (16 μm) using a second antibody against Olig2 (gift of Dr. Wichterle) (Wichterle et al., 2002) showing the localization of Olig2 within the DARPP-32-positive striosomes. Scale bars = 200 μm and 20 μm, respectively. (B) Map of the p688/mCherry plasmid used for the cloning of the
Figure 6—figure supplement 2.
Olig2 overexpression in medium spiny neurons promote their maturation in vitro.
(A) Representative DARPP-32 immunostaining in DIV9 WT primary striatal neurons transduced at DIV5 with ADV-GFP (control) or ADV-
To determine if the ATAC-seq identified functional, compartment-specific OCRs distinct from promoters, we sought to validate a putative enhancer peak in a gene preferentially expressed in the striosomes. We compared the
To determine if Olig2 is required for striosome compartmentation, we utilized a mouse in which a Cre recombinase cDNA inserted into the
Next, we sought to determine if expression of OLIG2 in mouse primary MSNs promotes maturation toward a striosome phenotype (Figure 6—figure supplement 2A,B). Adenovirus-mediated expression of OLIG2 does not increase the number of DARPP-32-immunopositive cells, but importantly, it increases levels of the mRNAs for striosome markers
Transcription network analysis defines the STAT pathway as a key regulator of striosome development
We used the arithmetic mean of the rank of each of the TFs at each node to quantitatively assess the importance of individual TFs in the striosome network (Figure 3A). The ranking was based on the number of edges expanding from each node and the number of overlapping TFs shared between nodes that interact with
Figure 7.
STAT1 overexpression in MSNs in vitro promotes maturation and increases levels of
(A) Network analysis indicates that Stat1 may be a ‘master’ regulator in striosomes as it modulates the levels of the greatest number of TFs enriched in that compartment. (B,C) GFP and DARPP-32 immunolabeling in DIV9 WT primary striatal neurons either non-transduced (NT) or transduced for 96 hr with ADV-STAT1-V5 or ADV-GFP showing that STAT1 overexpression increases the number of DARPP-32-immunopositive cells. Scale bars = 50 µm. n=four images from four individual cultures, t-test ***p<0.001. (D) RT- qPCR assay shows increase of mRNA for striosome markers
TFs drive MSN fate in iPSCs-derived neural stem cells (NSCs)
MSN subtypes are critically needed for disease modeling, and protocols and TFs that specifically drive a striosome phenotype are just beginning to emerge (Cirnaru et al., 2019). To determine if the FOXF2, OLIG2, and STAT pathways are key for human development of MSNs and represent useful TFs for disease modeling of MSN subtypes, we expressed the TFs alone or in combination in NSCs derived from HD patient-induced pluripotent stem cells (Figure 8; An et al., 2012; Naphade et al., 2017; Ring et al., 2015; Voisin et al., 2020; Zhang et al., 2010). STAT1 overexpression, with or without FOXF2, increases the number of OPRM1-immunopositive cells (Figure 8A), with OLIG2 producing a trend in the same direction. STAT1 and OLIG2 together also increase the expression of PPP1R1B (Figure 8A,B). Alone, only OLIG2 increases
Figure 8.
FOXF2, OLIG2, and STAT1, alone and in combination, promote MSN differentiation in NSCs from human HD induced pluripotent stem cells.
(A) HD72-NSCs transduced for 4 days with ADV-
Figure 8—figure supplement 1.
Expression of FOXF2, OLIG2 and STAT in HD NSCs resulted in an increase expression of pan neuronal marker Nestin in HD NSCs.
HD NSCs transduced for 4 days with ADV-
Figure 8—figure supplement 2.
FOXF2, OLIG2, and STAT1, alone and in combination, promote MSN differentiation in NSCs from human isogenic control C116-induced pluripotent stem cells.
C116-NSCs transduced for 4 days with ADV-
OLIG2 alone also increases the matrix marker
Discussion
In the first postnatal week, striosome neurons are already compartmentalized, and their maturation is relatively advanced as compared to the matrix, based on their description as ‘dopamine islands’ and their expression of
The use of these new data bases to construct a transcriptional regulatory network also led to the identification of
The transcription regulatory networks generated for each cell compartment have defined other potential TFs and pathways regulating striosome and matrix cell fate. The co-expression enrichment analysis highlighted a set of TFs that coexist with
Consistent with our hypothesis, the EGFP– cells were indeed enriched in markers of immature MSNs. Most notably, these included multiple members of the
In summary, we present two data bases derived from PND3 striatum that compare gene expression and OCRs in developing striosome and matrix, demonstrating clear gene expression and epigenetic distinctions. Using these data sets, we describe a novel pathway via which Stat1 appears to regulate a transcriptional hierarchy that includes Foxf2 and Olig2, critical for striosome compartmentation and phenotypic maturation. Further, we demonstrate that these TFs can be utilized to better model striatal striosome MSNs in hiPSC systems and to develop genetic tools to direct expression to neuronal subsets for their eventual in vivo manipulation and study.
Materials and methods
Animals
Animal procedures were conducted in accordance with the NIH Guidelines for the Care and Use of Experimental Animals and were approved by the Institutional Animal Care and Use Committee of our institutions (LA09-00272, 16–0847 PRYR1). The
Enzymatic dissociation and FACS purification of striatal neurons
PND3
Tissue preparation and immunofluorescence
PND3 mice were rapidly euthanized by decapitation, and brains were removed, washed in ice-cold PBS, and post-fixed for 24 hr at 4°C in 4% PFA. The brains were then incubated in 30% sucrose/1X PBS for 24 hr at 4°C and cryopreserved in OCT embedding medium (4583, Tissue-Tek Sakura). Serial coronal sections (16 µm) were cut on a Leica cryostat, collected on Superfrost Plus microscope slides (Fisher Scientific) and frozen at −20°C. Immunofluorescence was performed as described (Cirnaru et al., 2019). Sections were incubated with mouse anti-DARPP-32 (1:250, sc-271111, Santa Cruz Biotechnology), sheep anti-Foxf2 (1:2000, AF6988, R and D), rabbit anti-Olig2 (1:500, ab136253, Abcam), rabbit anti-Olig2 (1:8000, gift from Dr. Wichterle) (Wichterle et al., 2002), mouse anti-Olig2-Alexa488 (1:1000, MABN50A4, Millipore) rabbit anti-Irx1 (1:2000, PA5-36256, Thermo Fisher Scientific) or rabbit anti-tyrosine hydroxylase (1:1000, OPA1-04050, Thermo Fisher Scientific) antibodies. The respective secondary antibodies included: anti-mouse Alexa 488 (1:400, A-11008, Thermo Fisher Scientific), anti-mouse Alexa 594 (1:400, A-11005, Thermo Fisher Scientific), anti-rabbit Alexa 488 (1:400, A-11034, Thermo Fisher Scientific), anti-rabbit Alexa 594 (1:400, A-11012, Thermo Fisher Scientific), or anti-sheep Alexa 594 (1:400, A-11016, Thermo Fisher Scientific). Sections were sealed with Vectashield hard-set mounting medium (H-1400, Vector Laboratories). Images were acquired using an Olympus BX61 microscope or a Confocal Zeiss LSM 510.
Fluorescent in situ hybridization (FISH) using RNAscope technology
Postnatal day 3 WT brains were fixed in freshly prepared, ice-cold 4% PFA for 24 hr at 4°C, followed by equilibration in 10% sucrose gradient, then 20% and finally 30%, each time allowing the tissue to sink to the bottom of the container. The tissue was embedded in OCT compound (4583, Tissue-Tek Sakura) and stored at −80°C until sectioned. 16 µm-thick sections were cut using a Leica cryostat, collected onto Superfrost Plus slides maintained at −20°C during the sectioning. The slides were then stored at −80°C. RNAscopeProbe murine Mm-DARPP-32-C1 Ppp1r1b (NM_144828.1, bp590-1674, 405901), Mm-Olig2-C3 (NM_016967.2, bp865-2384, 447091-C3)and Mm-Foxf2-C2 (NM_010225.2, bp846-2316, 473391-C2) were purchased from Advanced Cell Diagnostics probe catalog (ACD). For signal detection, we used Opal 520 and Opal 690 TSA plus fluorophores (Akoya Biosciences). The RNAscope Multiplex Fluorescent Reagent Kit v2 (ACD, 323100) used here provides the target retrieval solution, hydrogen peroxide, protease III, amplification reagents (Amp1-3), HRP reagents, DAPI, TSA buffer and wash buffer. We used a modified version of the manufacturer’s protocol for sample preparation, probe hybridization, and signal detection. Briefly, the fresh frozen sections on slides were retrieved from −80°C and briefly immersed in 1X PBS to wash off the O.C.T and then baked at 60°C for 30 min. Slides were then post-fixed in fresh 4% PFA for 1 hr at room temperature (RT). After fixation, the sections were dehydrated in a series of ethanol solutions (5 min each in 50%, 70%, and two changes of 100% ethanol) at RT and left to dry for 5 min at RT. Sections were treated with hydrogen peroxide for 10 min at RT and washed twice with distilled water. Subsequently, target retrieval was performed by boiling the slides for 5 min in 1X Target Retrieval Reagent (ACD), washed in distilled water, immersed in 100% ethanol, and air-dried for 5 min at RT. A hydrophobic barrier was created around the section using an ImmEdge Pen (ACD, 310018) and completely dried at RT before proceeding to the next step. Sections were then treated with protease III for 5 min at 40°C in the pre-warmed ACD HybEZ II Hybridization System (ACD, 321721) inside the HybEZ Humidity Control Tray (ACD, 310012) and washed twice with distilled water. The
Primary neuronal cultures
E16.5 embryos were obtained from wild-type (WT) Swiss Webster timed bred females purchased from Charles River. E16.5 striatum was removed from by microdissection in cold Leibovitz’s medium (L-15) (Gibco-Invitrogen,11415064) and primary medium spiny neuronal cultures were prepared as described (Cirnaru et al., 2019). Briefly, the tissue was incubated in Ca2+/Mg2+-free Hanks’ balanced salt solution (Sigma, 55021C) for 10 min at 37°C. The incubation mixture was replaced with 0.1 mg/ml papain in Hibernate E/Ca2+ (BrainBits), incubated for 8 min, and rinsed in Dulbecco’s minimum essential medium (Gibco-Invitrogen, 21013024) with 20% fetal bovine serum (Gibco-Invitrogen, 10438026) and twice in Leibovitz’s medium (L-15). The tissue was then suspended in Dulbecco’s minimum essential medium with 10% fetal bovine serum, glucose (6 mg/mL) (Sigma, G7021), glutamine (1.4 mM) (Gibco- Invitrogen, 25030081) and penicillin/streptomycin (100 U/mL) (Gibco-Invitrogen, 15140122). Cells were triturated through a glass bore pipette and plated onto either Lab Tek eight-well slides (1.25 x105 cells/well) for immunocytochemistry or 24-well plates (1x106 cells/well) for RT-PCR analysis, each coated with polymerized polyornithine (0.1 mg/mL in 15 mM borate buffer, pH 8.4) and air-dried. One hr later, the medium was replaced with Neurobasal (Gibco-Invitrogen, 21103049), supplemented with B27, (Gibco-Invitrogen, 17504044), GLUTAMAX (Gibco- Invitrogen, 35050061) and penicillin/streptomycin. The medium was changed every 2 days, and the cells were assayed on day in vitro (DIV) 9.
Neuronal adenovirus (ADV) transduction
ADV-CMV-
Mouse neuronal immunocytochemistry
Cells were fixed in 4% paraformaldehyde in 0.1 M phosphate buffer, pH 7.4, and immunolabeled with mouse anti-DARPP-32 (1:250, sc-271111, Santa Cruz Biotechnology), rabbit anti-Olig2 (1:500, ab136253, Abcam) or rabbit anti-STAT1 (1:400, 14994S, Cell Signaling Technology) followed by anti-mouse Alexa 488 (1:400, A-11008, Thermo Fisher Scientific) and anti-rabbit Alexa 594 (1:400, A-11012, Thermo Fisher Scientific). To identify the total number of cells the nuclei were stained with DAPI (4', 6-diamino-2-phenylindole dihydrochloride) (1:10,000, Millipore-Sigma). Images were acquired using Olympus BX61 microscope and analyzed using Fiji software (ImageJ).
Real-time qPCR of cultured mouse neurons
RNA from DIV9 primary MSNs was extracted with the miRNeasy micro kit (Qiagen), according to the manufacturer’s instructions. RNAs (500 ng) were reversed-transcribed using the High Capacity RNA-to-cDNA Kit (Applied Biosystems). Real-time qPCR was performed in a Step-One Plus system (Applied Biosystems) using All-in-One qPCR Mix (GeneCopoeia). Quantitative PCR consisted of 40 cycles, 15 s at 95°C and 30 s at 60°C each, followed by dissociation curve analysis. The ΔCt was calculated by subtracting the Ct for the endogenous control gene GAPDH from the Ct of the gene of interest. Mouse primer sequences are listed in Table 2. Relative quantification was performed using the ΔΔCt method and expressed as a fold-change relative to control by calculating 2-ΔΔCt.
Table 2.
RT-qPCR murine primer sequences.
| Gene | Primer forward 5'−3' | Primer reverse 5'−3' |
|---|---|---|
| GAPDH | AACGACCCCTTCATTGACCT | TGGAAGATGGTGATGGGCTT |
|
| GAAGAAGAAGACAGCCAGGC | TAGTGTTGCTCTGCCTTCCA |
|
| CCCTCTATTCTATCGTGTGTGT | AGAAGAGAGGATCCAGTTGCA |
|
| AAGCAGCTTGCCTTTGCTAAG | GGATTGAATGTATGTGTGGCTGA |
|
| GGACCTACCAAGAACTGGAAC | GATCCCAGTAAACCCGTCTG |
|
| ACTCTCAAACTAGCCGCTGCA | TCAGCGTCGAAATGAAGCC |
|
| CTTGGACCAGAACCAGGATG | GTGGCAGTTCACACCACAAG |
|
| TGCTCCGCTTTGCACACACAGG | TAAGTTCTCAATAATGGACCAGCAC |
|
| TCGTGGTCATTCTCATTG | TCTCTTCATCTGCTTCTTG |
|
| CGATAGAACCAAGATAATACT | TAGAATCAGAGGACTCAG |
Human-induced pluripotent stem-cell-derived NSC culture
All work on human iPSCs was approved by the Buck Institute institutional stem cell and the human ethics committee (Approval S1002). The isogenic C116 and HD iPSCs were previously published (An et al., 2012; Naphade et al., 2017; Ring et al., 2015; Voisin et al., 2020; Zhang et al., 2010), were sequenced, karyotyped and were negative for mycoplasma contamination. Human-induced pluripotent stem cells (iPSCs) were differentiated into prepatterned activin A-treated neural stem cells (NSCs) by the following protocol. Briefly, iPSC colonies were detached using 1 mg/ml collagenase (Type IV, Thermo Fisher Scientific, 17104019) in Gibco KnockOut DMEM/F-12 medium (Thermo Fisher Scientific, 12660012), and the resulting cell clumps were transferred to a 0.1% agarose (Sigma- Aldrich, A9414) coated low-attachment petri dish in embryonic stem (ES) culture medium Gibco KnockOut DMEM/F12 supplemented with 20% Gibco KnockOut Serum Replacement (Thermo Fisher Scientific, 10828028), 2.5 mM L-glutamine (Thermo Fisher Scientific, 25030081), 1 X Non- Essential Amino Acids (NEAAs) (Thermo Fisher Scientific, 11140050), 15 mM HEPES (Thermo Fisher Scientific, 15630106), 0.1 mM β-mercaptoethanol (Thermo Fisher Scientific, 31350010), 100 U/ml penicillin-streptomycin (Thermo Fisher Scientific, 15140122). Every 2 days, 25% of ES medium was replaced by embryoid body (EB) differentiation medium [DMEM (Corning, 10–013- CV) supplemented with 20% FBS (Thermo Fisher Scientific, 16000036), 1 X NEAA, 2 mM L- glutamine, 100 U/ml penicillin-streptomycin]. At day 8, 100% of the culture medium was EB medium. At day 10, the embryoid bodies were attached to dishes coated with poly-L-ornithine (1:1000 in PBS; Sigma-Aldrich, P3655) and laminin (1:100 in KnockOut DMEM/F-12; Sigma- Aldrich, L2020), and cultured in neural induction medium [DMEM/F12 supplemented with 1 X N2 (Thermo Fisher Scientific, 17502001), 100 U/ml penicillin-streptomycin] and 25 ng/ml βFGF (Peprotech, 100-18B) and 25 ng/ml Activin A (Peprotech, 120–14P). Media change was performed every 2 days. Rosettes were harvested after 7–10 days, plated on poly-L-ornithine- and laminin-coated dishes, and cultured in Neural Proliferation Medium [NPM; Neurobasal medium (Thermo Fisher Scientific, 21103049), B27-supplement 1 X (Thermo Fisher Scientific, 17504001), GlutaMAX 1 X (Thermo Fisher Scientific, 35050061), 10 ng/ml leukemia inhibitory factor (Peprotech, 300–05), 100 U/ml penicillin-streptomycin] supplemented with 25 ng/ml β-FGF and 25 ng/ml activin A. The resulting NSCs were passaged and maintained in this same medium. These prepatterned activin A-treated NSCs were validated by immunofluorescence analysis and labeled positively for putative NSC markers, namely Nestin, SOX1, SOX2, and PAX6.
Adenovirus transduction of human NSCs
For the ADV transduction experiments human iPSC- derived NSCs were plated at 700,000 cells per well of a six-well plate in 2 ml of NPM supplemented with 25 ng/ml βFGF (Peprotech, 100-18B) and 25 ng/ml activin A (Peprotech, 120–14P). Two day after plating, cells were transduced at MOI of 10 with ADV-CMV-
Cell immunofluorescence of human NSCs
Cells were fixed using 4% paraformaldehyde in 0.1 M phosphate-buffered saline (PBS), pH 7.4 (Corning, 21–040-CV) for 30 min. After three washes in PBS, cells were permeabilized and blocked for 1 hr at RT using 0.1% Triton X-100 (Thermo Fisher Scientific, 28313) and 4% donkey serum in PBS. Primary antibodies were added in the presence of blocking buffer overnight at 4°C. Secondary antibodies (1:500) were added after three PBS washes in blocking buffer at RT for 1 hr. The following primary antibodies were used for the immunofluorescence studies: rabbit anti-DARPP-32 (Santa Cruz, sc-11365, 1: 100) and rabbit anti-Opioid Receptor-Mu (Millipore, AB5511, 1:500). The secondary antibodies were donkey anti- rabbit IgG conjugated with Alexa-488 (Invitrogen, A12379) or Alexa-647 (Invitrogen, A22287). Images were acquired using a Biotek Cytation five microscope and were prepared using Fiji software (ImageJ).
Quantitative real-time PCR of human NSCs
For qRT-PCR analysis of prepatterned activin A- treated human NSCs, total RNA was isolated using the ISOLATE II RNA Mini Kit (Bioline, BIO- 52072). cDNA was prepared from 300 ng of RNA in a total reaction volume of 20 µl using the Sensi-FAST cDNA synthesis kit (Bioline, BIO-65053). RT-PCR reactions were set up in a 384-well format using 2X SensiFAST Probe No-ROX Kit (Bioline, BIO-86005) and 1 µl of cDNA per reaction in a total volume of 10 µl. RT-PCR was performed on the Roche LightCycler 480 instrument. Quantitative PCR consisted of 40 cycles, 15 s at 95°C and 30 s at 60°C each, followed by dissociation curve analysis. The ΔCt was calculated by subtracting the Ct for the endogenous control gene β-actin from the Ct of the gene of interest. Human primer sequences are listed in Table 3. Relative quantification was performed using the ΔCt method and is expressed as a fold- change relative to control by calculating 2- ΔCt.
Table 3.
qRT-PCR human primer sequences.
| Gene | Primer forward 5'−3' | Primer reverse 5'−3' | Universal probe |
|---|---|---|---|
|
| CACACCACCTTCGCTGAAA | GAAGCTCCCCCAGCTCAT | 82 |
|
| AGAAACAGCAGGAGCTGTGG | ACCGAGACTTTTCGGGTTC | 30 |
|
| CCCAGAGGGAGCTCATCAC | TTTGACACTGGCCACAGGT | 45 |
|
| CACAGCCTCACAGTTTTTCG | CCTTTCCTTCCAGGTAACCA | 36 |
|
| ATGGTGATGGGATATTTTCTCC | GCCATGCCCTTTGTCTTAAC | 46 |
|
| GCAACTACGTGGGCGACT | CGAGGCACGTACTTGTGAGA | 78 |
|
| GAGCCAAAGATCTGCTCCAT | GGTCCGATCCTTACTCTCCTC | 71 |
|
| TGAGCCACAGCTCCATCTC | CCGTCACTAGTTCCGTGAGAC | 75 |
|
| AGCAGCCACTCAGGCAAC | ACGAAAATAGGGCGAAATAGAA | 51 |
|
| CCAACCGCGAGAAGATGA | CCAGAGGCGTACAGGGATAG | 64 |
Generation of RNA-seq libraries
For RNA-seq, EGFP-positive and -negative cells were sorted into low-binding tubes containing Arcturus PicoPure Extraction buffer. RNA was isolated in accordance with the PicoPure RNA Isolation kit manufacturer’s instructions, which included a DNase treatment step. Samples were eluted in RNase-free water and stored at −80°C until preparation of RNA-sequencing libraries using the Takara Clontech Laboratories SMARTer Stranded Total RNA-Seq Pico Kit, according to the manufacturer’s instructions. After construction of the RNA-seq libraries, libraries were analyzed on an Agilent High Sensitivity D1000 TapeStation, and quantification of the libraries was performed using the KAPA Library Quantification Kit.
RNA-seq
RNA-seq was carried out at New York University using Illumina HiSeq 4000 Paired-End 150 Cycle Lane from purified RNA PicoPure RNA Isolation Kit KIT0204 Arcturus (ThermoScientific) and SMARTer Stranded Total RNA-Seq Kit - Pico Input Mammalian (250 pg–10 ng RNA) (635005, Clontech). The low-quality base (quality score lower than 20), as well as the adapters of the raw reads from the sequencing experiments, was removed using Trim Galore! (https://www.bioinformatics.babraham.ac.uk/projects/trim_galore/). External and internal rRNA contamination was filtered through SortMeRNA 2.1b (Kopylova et al., 2012). Then the filtered raw reads were then mapped to the Genome Reference Consortium Mouse Build 38 striosome release 6 (GRCm38.p6) assembly by GENCODE using STAR 2.7.2b (Dobin et al., 2013). The counts of reads mapped to known genes were summarized by featureCounts, using GENECODE release M22 annotation (GSE143276). The data discussed in this publication have been deposited in NCBI's Gene Expression Omnibus GEO Series accession number, GSE143276. Next, R Bioconductor package DESeq2 (Love et al., 2014) was used to normalize raw read counts logarithmically and perform differential expression analysis. Differentially expressed genes were based on an arbitrary cutoff of adjusted p-value less than 0.01 (Figure 2—source data 1). We found the regulation pattern of a gene with the same EGFP status is more likely to be the same across replicates, whereas it is more likely to be different when the EGFP status is different between replicates with the exception of
Terminology enrichment analysis and pathway enrichment analysis
Enrichment analysis was performed on gene clusters in specific databases to determine if a specific biological annotation could be considered as significantly represented under the experiment result. Both terminology enrichment analysis and pathway enrichment analysis were conducted by clusterProfiler (Yu et al., 2012), a Bioconductor package. In our analysis, biological process (BP), molecular function (MF), and cellular component (CC) terms in gene ontology (GO) (The Gene Ontology Consortium, 2017), as well as pathway annotations derived from Kyoto Encyclopedia of Genes and Genomes (KEGG), were chosen to identify predominant biological processes of the differentially expressed gene clusters and differentially expressed transcription factor clusters involved in the development of the striosome neurons. We conducted both analyses on the differentially expressed gene clusters with the arbitrary cutoff of adjusted p-value less than 0.01 and the absolute Log2 fold-change greater than 0, 1, and 2, respectively, and we conducted both analyses on the differentially expressed transcription factor clusters with the arbitrary cutoff of adjusted p-value less than 0.01 and the absolute Log2 fold- change greater than 0 and 1, respectively (Figure 2—source data 1).
TF enrichment analysis and co-expressor enrichment analysis
Both TF enrichment analysis and co-expressor enrichment analysis were conducted by Enrichr (Kuleshov et al., 2016), a comprehensive online tool for doing enrichment analysis with a variety of biologically meaningful gene set libraries (Figure 2—source data 1). In our analysis, ChEA (Lachmann et al., 2010) and ENCODE (Frankish et al., 2019) databases were chosen to identify the significant upstream TFs regulating genes and other TFs differentially expressed in striosome cells and matrix cells, respectively, and the ARCHS4 database was chosen to identify the significant co-expressors of those differentially expressed genes and transcription factors. An arbitrary cutoff of adjusted p- value less than 0.01 and the absolute log2 fold-changes greater than 0 and 1 were chosen (Figure 2—source data 1).
GeneMANIA gene regulatory network analysis
GeneMANIA (Warde-Farley et al., 2010) is an online tool using published and computational predicted functional interaction data among proteins and genes to extend and annotate the submitted gene list by their interactive biological pathways and visualize its inferred interaction network accordingly. We used GeneMANIA to conduct the interaction network inference analysis to TFs enriched in either the striosome or matrix compartments with the arbitrary cutoff of adjusted p-value less than 0.01 and the absolute log2 fold-change greater than 1.
Gene regulatory network inference through data curation
A gene regulatory network links TFs to their target genes and represents a map of transcriptional regulation. We used all the TFs and their target gene data curated by ORegAnno (Lesurf et al., 2016) to build the network. To simplify the network, we only chose the compartmental differentially expressed TFs that are high on the hierarchy. In other words, only the differentially expressed TFs that served as a regulator of other differentially expressed TFs were chosen as the candidates of a gene regulatory network.
Gene set enrichment analysis
GO (Subramanian et al., 2005) was performed using ranked list of differential gene expression with parameters set to 2000 gene-set permutations and gene-set size between 15 and 200. The gene-sets included for the GSEA analyses were obtained from Gene Ontology (GO) database (GOBP_AllPathways), updated September 01, 2019 (http://download.baderlab.org/EM_ Genesets/). An enrichment map (version 3.2.1 of Enrichment Map software Merico et al., 2010) was generated using Cytoscape 3.7.2 using significantly enriched gene-sets with an FDR <0.05. Similarity between gene-sets was filtered by Jaccard plus overlap combined coefficient (0.375). The resulting enrichment map was further annotated using the AutoAnnotate Cytoscape App.
Data processing
The preprocessing of ATAC-seq data involved the following steps:
Alignment
Sequencing reads were provided by the sequencing center demuxed and with adaptors trimmed. Reads from each sample were aligned on GRCh38-mm10 reference genome using the STAR aligner (Dobin et al., 2013) (v2.5.0) with the following parameters:
alignIntronMax 1
outFilterMismatchNmax 100
alignEndsType EndToEnd
outFilterScoreMinOverLread 0.3
outFilterMatchNminOverLread 0.3
This produced a coordinate-sorted BAM file of mapped paired-end reads for each sample. We excluded reads that: (1) mapped to more than one locus using SAMtools (Li et al., 2009), (2) were duplicated using PICARD (v2.2.4; http://broadinstitute.github.io/picard; Picard, 2016), and (3) mapped to the mitochondrial genome.
Quality-control (QC) metrics
The following quality-control metrics were calculated for each sample: (1) total number of initial reads, (2) number of uniquely mapped reads, (3) fraction of reads that were uniquely mapped and additional metrics from the STAR aligner, (4) Picard duplication and insert metrics, and (5) normalized strand cross- correlation coefficient (NSC) and relative strand cross-correlation coefficient (RSC), which are metrics that use cross-correlation of stranded read density profiles to measure enrichment independently of peak calling. Figure 4—source data 1 describes the main QC metrics. The bigWig tracks for each sample were manually inspected. None of the libraries failed QC and visual inspection and eight libraries were subjected to further analysis.
Selection and further processing of samples meeting quality control
We subsequently subsampled samples to a uniform depth of 10 million paired-end reads and merged the BAM-files of samples from the same cell type. We called peaks using the Model-based Analysis of ChIP-Seq (MACS) (Zhang et al., 2008) v2.1 (https://github.com/macs3-project/MACS; Liu, 2015). It models the shift size of tags and models local biases in sequencability and mapability through a dynamic Poisson background model. We used the following parameters:
—keep-dup all
—shift −100
—extsize 200
—nomodel
We created a joint set of peaks requiring each peak to be called in at least one of the merged BAM-files. That is, if a peak was identified in just one or more samples it was included in the consensus set of peaks. If two or more peaks partially overlapped, the consensus peak was the union of bases covered by the partially overlapping peaks. After removing peaks overlapping the blacklisted genomic regions, 69,229 peaks remained. We subsequently quantified read counts of all the individual non-merged samples within these peaks, again, using the feature counts function in RSubread (Liao et al., 2014) (v.1.15.0). We counted fragments (defined from paired-end reads), instead of individual reads, that overlapped with the final consensus set of peaks. This resulted in a sample by peak matrix of read counts, obtained using the following parameters:
allowMultiOverlap = F,
isPairedEnd = T,
strandSpecific = 0,
requireBothEndsMapped = F,
minFragLength = 0,
maxFragLength = 2000,
checkFragLength = T,
countMultiMappingReads = F,
countChimericFragments = F
Differential analysis of chromatin accessibility
To identify genomic regions with significant regional differences in chromatin structure among the two cell types, we performed a statistical analysis of chromatin accessibility. Here, chromatin accessibility was assessed by how many ATAC-seq reads overlap a given OCR: the higher the read count, the more open the chromatin is at a given OCR. For this, we performed the following steps:
Read counts
As a starting point, we used the sample-by-OCR read count matrix described in the previous section (eight samples by 69,229 OCRs). From here, we subsequently removed nine OCRs using a filtering of 0.5 CPM in at least 50% of the samples, resulting in our final sample-by-OCR read count matrix (eight samples by 69,220 OCRs). Next, we normalized the read counts using the trimmed mean of M-values (TMM) method (Robinson et al., 2010).
Covariate exploration
To explore factors affecting the observed read counts, we examined several biological and technical sample-level variables. For these covariates (e.g. number of peaks called in the sample, chrM metrics, RSC and NSC, and Picard insert metrics), we normalized to the median of the cell. We next assessed the correlation of all the covariates with the chromatin accessibility values in the normalized read count matrix to determine which of these variables should be candidates for inclusion as covariates in the differential analysis. We did this using a principal component analysis of the normalized read count matrix and by examining which variables were significantly correlated with the high-variance components (explaining > 1% of the variance) of the data. We did not identify any significant association even when we used a lenient false discovery rate threshold of 0.2.
Differential analysis
We used the edgeR package (Robinson et al., 2010) to model the normalized read counts by negative binomial (NB) distributions. The estimateDisp function was used to estimate an abundance-dependent trend for the NB dispersions (McCarthy et al., 2012). To normalize for compositional biases, the effective library size for each sample was estimated using the TMM approach as described above. For each open chromatin region, we applied the following model for the effect on chromatin accessibility of each variable on the right-hand side:
chromatin accessibility ~ cell type.
Then, for each OCR, the cell type coefficient was statistically tested for being non- vanishing. A quasi-likelihood (QL) F-test was conducted for each OCR using the glmQLFTest function (Lund et al., 2012) from the edgeR package, with robust estimation of the prior degrees of freedom. p-values were then adjusted for multiple hypothesis testing using false discovery rate (FDR) estimation, and the differentially accessible regions of chromatin were determined as those with an estimated FDR below, or at, 5%.
Annotation of OCRs
We used the gene annotations form the org.Mm.eg.db (version 3.8.2) package for all analyses in this paper. We assigned the closest gene and the genomic context of an ATAC-seq OCR using ChIPSeeker (Dixon et al., 2015). The genomic context was defined as promoter (+/- 3 Kb of any TSS), 5’-UTR, 3’-UTR, exon, intron, distal intergenic and downstream. TF binding motif analysis of ATAC-seq data was performed using HOMER suit function findMotifsGenome.pl tool. Differential footprinting analysis was performed using TOBIAS (https://doi.org/10.1101/869560) using TF motifs from Jaspar database (doi: 10.1093/nar/gkx1126).
© 2021, Cirnaru et al. This work is published under https://creativecommons.org/licenses/by/4.0/ (the “License”). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.