The cell bodies of midbrain dopaminergic (mDA) neurons are localized in the substantia nigra and ventral tegmental area, where these neurons play critical roles in the central regulation of motor and motivational behaviors (Elsworth and Roth, 1997; Wallen and Perlmann, 2003). Progressive degeneration of dopaminergic (DA) cells leads to the cardinal symptoms of Parkinson’s disease (Shastry, 2001; Wallen and Perlmann, 2003). Similar to the development of other neuron types in the central nervous system (CNS), maturation of mDA neurons may require neuron migration, cell positioning, and axon wiring, all processes for establishing neural networks in different neuronal systems.

A previous report described a series of events that are important for mesencephalic DA cell differentiation in mice (Andersson et al., 2013; Blakely et al., 2011; Deng et al., 2011; Di Salvio et al., 2010; Inestrosa and Arenas, 2010; Panman et al., 2014; Prakash and Wurst, 2006; Smits et al., 2006; Theofilopoulos et al., 2013; Van den Heuvel and Pasterkamp, 2008; Wallen et al., 1999; Yang et al., 2013; Zhang et al., 2013). Initially, the mDA progenitor appears in the ventricular-most cell layers within the ventral mesencephalon at embryonic day (E)10.5, undergoes neurogenesis, and gives rise to mDA neuroblasts. The mDA precursor cells then migrate from the ventricular zone to a medial location in the mantle zone and exit the cell cycle (become post-mitotic) between E11.5 and E13.5. At this stage, expression of nuclear receptor-related 1 (Nurr1) protein, a transcription factor, is strongly upregulated.

This transcription factor is a member of a family of nuclear receptors that are critical for the development and survival of DA neurons. During the development of DA neurons, Nurr1 expression is induced in early postmitotic DA progenitors and then maintained during differentiation and in the adult brain (Perlmann and Wallen-Mackenzie, 2004). Of interest, in mice lacking the Nurr1 gene, DA neurons fail to differentiate. In addition, in contrast to wild-type animals, in the knockout mice, DA progenitor cells show no lateral migration in the midbrain, fail to innervate the striatal target area, and become apoptotic (Castillo et al., 1998; Saucedo-Cardenas et al., 1998; Zetterstrom et al., 1997). In adult DA neurons, Nurr1 also maintains fiber integrity, and Nurr1 ablation results in a progressive pathology associated with reduced striatal DA, impaired motor behaviors, and dystrophic axons and dendrites (Kadkhodaei et al., 2013). The early development of mDA neurons has received much attention in recent years, particularly with regard to morphogenesis, progenitor specification, mDA differentiation, and neuritogenesis (Andersson et al., 2013; Blakely et al., 2011; Deng et al., 2011; Di Salvio et al., 2010; Inestrosa and Arenas, 2010; Prakash and Wurst, 2006; Smits et al., 2006; Theofilopoulos et al., 2013; Van den Heuvel and Pasterkamp, 2008). Little is known, however, about the downstream targets of Nurr1 that are involved in this pathway.

Few reports are available that characterize the gene expression profile of Nurr1 overexpression (Jacobs et al., 2009a; Sousa et al., 2007), and no reports focused on human development are available. The dramatic differences between mice and humans in brain development—including size, events at specific stages, and structural proportions—likely reflect meaningful differences in developmental gene expression (Bohland et al., 2010; Hawrylycz et al., 2012; La Manno et al., 2016; Lein et al., 2007; Myers et al., 2015; Ng et al., 2009). In the absence of human developmental data, gaps in understanding of gene expression profiling are likely, given the limitations of mouse-derived data. Furthermore, the two available reports describe Nurr1 as regulating many genes in the murine mesencephalon MN9D cell line (Jacobs et al., 2009b) and in meso-diencephalic dopamine neurons (Sousa et al., 2007). Neither group, however, confirmed whether the affected genes are direct targets of Nurr1 regulation. For this reason, using human neural stem cells (hNSCs), we sought to address this gap and identify direct targets of NURR1. Here, we show that NURR1 modulates sets of genes implicated in cell migration, synapse wiring, and postmitotic events.

Cell cultures

Immortalized hNSC lines (HB1.F3, HB1.F5, and HB1.A4) (Kim, 2004) and the human embryonic kidney cell line HEK293 were maintained and passaged on uncoated culture dishes in Dulbecco’s modified Eagle medium (Gibco – ThermoFisher Scientific, USA) with 10% fetal bovine serum (Hyclone), and 10 µg/ml penicillin-streptomycin (Gibco, USA). All cells transduced with the plasmids pLPCX or pLPC-NURR1 were maintained in the same culture condition. Cells from the mouse amphotropic retrovirus packaging cell line PA317 were cultured in RPMI 1640 medium with 10% fetal bovine serum. All cells were incubated at 37°C with 5% CO₂.

Retrovirus-mediated gene transfer

We used an amphotropic replication-incompetent retroviral vector to infect target cells, as previously described (Kim et al., 2013). A vector encoding the human NURR1 gene was generated using pLPCX to infect the PA317 retrovirus packaging cells. Cells then were infected with a retrovirus encoding human NURR1 (pLPC-NURR1) or with the control vector (pLPCX).

Reverse transcription polymerase chain reaction (RT-PCR) and Northern blot analysis

Total RNA was extracted from all cells using the RNeasy Mini Kit (Qiagen, Germany). Total RNA (2 μg) was subjected to reverse transcription with Superscript II reverse transcriptase (Invitrogen, USA) according to the manufacturer’s instructions. Amplification reactions were performed with 1/10 volume of the reverse-transcribed product in a final volume of 25 ml using recombinant Taq DNA polymerase (Invitrogen), as described previously (Kim et al., 2013).

For Northern blotting, 20 mg of total RNA was separated on a 1% formaldehyde-agarose gel and transferred to polyvinylidene fluoride membranes (Millipore, USA). Blots were hybridized with a 1.8-kb BamHI–MluI fragment of the NURR1 gene labeled with [α-³²P]-dCTP (3000 Ci/mmol; Amersham, USA). X-ray film was exposed for 2 days at –80°C and developed.

DNA microarray analysis

The GeneChip Human Genome U133 Plus 2.0 Array (Affymetrix, USA) was used to analyze differential gene expression profiles modulated by NURR1 protein expression. Total RNA was extracted according to the manufacturer’s protocol (Affymetrix). For extraction of RNA from hNSCs overexpressing NURR1, we harvested the cells within one week of puromycin selection after transduction to avoid shutdown of NURR1 expression. Hybridization to the U133A DNA microarray, washing, and scanning were performed according to the manufacturer’s protocol, and expression patterns were compared between samples.

To normalize signals, we divided each probe by the average value of the chip to avoid differences between different chips and experiments. NCBI UniGene cluster nomenclature was used to describe uncharacterized sequences. All of the raw microarray data are available from Gene Expression Omnibus (GEO) (series record No. GSE58475). For hierarchical clustering of the list of differentially expressed genes using Hamming distance correlation, we used software made publicly available by the lab of Michael Eisen at the University of California, Berkeley.

Quantitative real-time PCR

Quantitative real-time RT-PCR assays were carried out using gene-specific double-labeled fluorescent probes and sets of specific primers in an ABI PRISM 7700 Sequence detection system (PE Applied Biosystems, USA). The primers and probes were obtained for 12 selected genes and the glyceraldehyde phosphate dehydrogenase (GAPDH) gene control from Assays-on-Demand system, which are ready-made and validated assays for the human genome (Applied Biosystems). To compensate for differences in cell number and/or RNA recovery, the copy number of hNSC mRNA was determined relative to human GAPDH mRNA simultaneously using the comparative ΔCt method.

RNA interference knockdown

Four sets of human NURR1 small interfering (si) RNA duplexes (19 nucleotides) were purchased from Dharmacon RNA Technologies (USA), and used as a cocktail. As a negative control, the commercially available control siRNA (Non-specific Control Duplex IX; Dharmacon RNA Technologies, USA) was used. Transfection of the siRNAs was performed according to the manufacturer’s instructions. Briefly, cells were transfected with 100 nM NURR1 siRNA by Lipofectamine Plus (Invitrogen). Forty-eight hours after transfection, cells were harvested for real-time PCR analysis.

Cloning of the 5'-flanking region of human NMU and G0S2 genes and mutagenesis

We cloned the 5'-flanking region of the human NMU (neuromedin U) and G0S2 (lymphocyte G0/G1 switch) genes using a PCR-based method of genomic DNA from SH-SY5Y cells. The following sets of primers were designed: forward, NMUp (5'--3') and reverse NMUp (5'-gggggtaccACATTCCAAGCAGCCTGGTTCA-3'); and forward, G0S2p (5'-gggaagcttCTCGGCGGCTGCGGGAGGCT-3') and reverse, G0S2p (5'gggggtacccatttctgccagttatcagaggt-3'). To facilitate subcloning, the enzyme site (indicated with lowercase letters) was added to the 5'-end of each primer. The PCR conditions were 25 cycles at 94°C for 1 min, at 60°C for 45 s, and at 72°C for 3 min. PCR products were directly cloned into a pGEM-T Easy vector (Promega, USA) and confirmed by sequence analysis. To generate luciferase reporter constructs, each promoter fragment (human NMU, 2102 bp; human G0S2, 2082 bp) was subcloned upstream of the firefly luciferase reporter gene of pGL3-Basic vector (Promega) employing KpnI/HindIII and KpnI/XhoI sites.

The QuickChange Mutagenesis kit (Stratagene, USA) was used to perform site-directed mutagenesis of the NurRE element within the human NMU promoter. The human NMU-luc construct was used as a template. Oligonucleotides with 34 to 36 nucleotides containing the desired point mutations were created as follows, with the mutated residues underlined: mNurRE, 5'-GTTCCCTGTTTTTCAAAAACAGGTCAAATA-3. The site-directed mutations were confirmed by sequence analysis.

Electrophoretic mobility shift assay (EMSA)

Nuclear extracts from HB1.F3 cells were prepared as described earlier (Kim et al., 2013). Sense and antisense oligonucleotides were annealed and then end-labeled with [γ-³²P]ATP (Amersham) and T4 polynucleotide kinase. Labeled probes were purified on 19% non-denaturing polyacrylamide gels. The DNA–protein binding reaction was performed in a final volume of 20 μl reaction buffer containing 10 mM Tris (pH 7.6), 50 mM NaCl, 0.5 mM dithiothreitol, 0.5 mM EDTA, 1 mM MgCl₂, 5% glycerol, and 250 μg of poly(dI–dC) per milliliter. Nuclear extract (20 μg of protein) was added to the reaction buffer in the absence or presence of unlabeled competitor DNA and pre-incubated for 15 min on ice. Radioisotope-labeled probes (50,000 cpm) were added, and the mixture was incubated for a further 30 min at room temperature. To resolve DNA–protein complexes, electrophoresis was performed on a 5% non-denaturing polyacrylamide gel. Gels were fixed, dried, and visualized by autoradiography. The oligonucleotides employed were as follows (only sense strands presented):

NMU-NurRE, 5'-GTTCCTCACCTTTCAAAGGGAGGTCAAATA-3';

NMU-mtNurRE, 5'-GTTCCCTGTTTTTCAAAAACAGGTCAAATA-3'; and

G0S2-NBRE2, 5'-CATCACTGACCTTTGCAATT-3'.

Chromatin immunoprecipitation (ChIP)

Soluble chromatin from HB.F3 cells was prepared as described previously (Kim et al., 2013) and immunoprecipitated with antibody against NURR1(TransCruz, sc-991X). Specific primer pairs were designed to amplify the promoter region of human NMU (forward, 5'-GGCTTAATTCCTAGCTGACGGGT-3'; and reverse, 5'-TCTAACTGGTCTGGAGCGCTGGT-3') and human G0S2 (forward, 5'-AGTGGGACCTTCGCGTGCAC-3'; and reverse, 5'-CCAGCTAATCTTGGGGAGGGCTC-3') from human genomic DNA. PCR conditions were optimized to allow for semiquantitative measurement, and PCR products were visualized on 1% agarose-Tris-acetate EDTA gels

Transfection and luciferase reporter assays

Transfections were performed using Lipofectamine PLUS reagent (Invitrogen) according to the manufacturer’s instructions. A total of 3 mg of DNA was used per transfection, as previously described (Kim et al., 2013), and briefly as follows: 1 to 2 mg of pLPCX or pLPCX-NURR1, 1 mg of reporter vector, and 0.5 μg of the internal control plasmid, pSV-β-galactosidase (Promega). All transfections were carried out three times in triplicate. Promoter activities were determined using the Single-Luciferase Reporter Assay System (Promega). Luciferase activities were normalized based on β-galactosidase activity in each well.

Cell cultures

Immortalized hNSC lines (HB1.F3, HB1.F5, and HB1.A4) (Kim, 2004) and the human embryonic kidney cell line HEK293 were maintained and passaged on uncoated culture dishes in Dulbecco’s modified Eagle medium (Gibco – ThermoFisher Scientific, USA) with 10% fetal bovine serum (Hyclone), and 10 µg/ml penicillin-streptomycin (Gibco, USA). All cells transduced with the plasmids pLPCX or pLPC-NURR1 were maintained in the same culture condition. Cells from the mouse amphotropic retrovirus packaging cell line PA317 were cultured in RPMI 1640 medium with 10% fetal bovine serum. All cells were incubated at 37°C with 5% CO₂.

Retrovirus-mediated gene transfer

We used an amphotropic replication-incompetent retroviral vector to infect target cells, as previously described (Kim et al., 2013). A vector encoding the human NURR1 gene was generated using pLPCX to infect the PA317 retrovirus packaging cells. Cells then were infected with a retrovirus encoding human NURR1 (pLPC-NURR1) or with the control vector (pLPCX).

Reverse transcription polymerase chain reaction (RT-PCR) and Northern blot analysis

Total RNA was extracted from all cells using the RNeasy Mini Kit (Qiagen, Germany). Total RNA (2 μg) was subjected to reverse transcription with Superscript II reverse transcriptase (Invitrogen, USA) according to the manufacturer’s instructions. Amplification reactions were performed with 1/10 volume of the reverse-transcribed product in a final volume of 25 ml using recombinant Taq DNA polymerase (Invitrogen), as described previously (Kim et al., 2013).

For Northern blotting, 20 mg of total RNA was separated on a 1% formaldehyde-agarose gel and transferred to polyvinylidene fluoride membranes (Millipore, USA). Blots were hybridized with a 1.8-kb BamHI–MluI fragment of the NURR1 gene labeled with [α-³²P]-dCTP (3000 Ci/mmol; Amersham, USA). X-ray film was exposed for 2 days at –80°C and developed.

DNA microarray analysis

The GeneChip Human Genome U133 Plus 2.0 Array (Affymetrix, USA) was used to analyze differential gene expression profiles modulated by NURR1 protein expression. Total RNA was extracted according to the manufacturer’s protocol (Affymetrix). For extraction of RNA from hNSCs overexpressing NURR1, we harvested the cells within one week of puromycin selection after transduction to avoid shutdown of NURR1 expression. Hybridization to the U133A DNA microarray, washing, and scanning were performed according to the manufacturer’s protocol, and expression patterns were compared between samples.

To normalize signals, we divided each probe by the average value of the chip to avoid differences between different chips and experiments. NCBI UniGene cluster nomenclature was used to describe uncharacterized sequences. All of the raw microarray data are available from Gene Expression Omnibus (GEO) (series record No. GSE58475). For hierarchical clustering of the list of differentially expressed genes using Hamming distance correlation, we used software made publicly available by the lab of Michael Eisen at the University of California, Berkeley.

Quantitative real-time PCR

Quantitative real-time RT-PCR assays were carried out using gene-specific double-labeled fluorescent probes and sets of specific primers in an ABI PRISM 7700 Sequence detection system (PE Applied Biosystems, USA). The primers and probes were obtained for 12 selected genes and the glyceraldehyde phosphate dehydrogenase (GAPDH) gene control from Assays-on-Demand system, which are ready-made and validated assays for the human genome (Applied Biosystems). To compensate for differences in cell number and/or RNA recovery, the copy number of hNSC mRNA was determined relative to human GAPDH mRNA simultaneously using the comparative ΔCt method.

RNA interference knockdown

Four sets of human NURR1 small interfering (si) RNA duplexes (19 nucleotides) were purchased from Dharmacon RNA Technologies (USA), and used as a cocktail. As a negative control, the commercially available control siRNA (Non-specific Control Duplex IX; Dharmacon RNA Technologies, USA) was used. Transfection of the siRNAs was performed according to the manufacturer’s instructions. Briefly, cells were transfected with 100 nM NURR1 siRNA by Lipofectamine Plus (Invitrogen). Forty-eight hours after transfection, cells were harvested for real-time PCR analysis.

Cloning of the 5'-flanking region of human NMU and G0S2 genes and mutagenesis

We cloned the 5'-flanking region of the human NMU (neuromedin U) and G0S2 (lymphocyte G0/G1 switch) genes using a PCR-based method of genomic DNA from SH-SY5Y cells. The following sets of primers were designed: forward, NMUp (5'--3') and reverse NMUp (5'-gggggtaccACATTCCAAGCAGCCTGGTTCA-3'); and forward, G0S2p (5'-gggaagcttCTCGGCGGCTGCGGGAGGCT-3') and reverse, G0S2p (5'gggggtacccatttctgccagttatcagaggt-3'). To facilitate subcloning, the enzyme site (indicated with lowercase letters) was added to the 5'-end of each primer. The PCR conditions were 25 cycles at 94°C for 1 min, at 60°C for 45 s, and at 72°C for 3 min. PCR products were directly cloned into a pGEM-T Easy vector (Promega, USA) and confirmed by sequence analysis. To generate luciferase reporter constructs, each promoter fragment (human NMU, 2102 bp; human G0S2, 2082 bp) was subcloned upstream of the firefly luciferase reporter gene of pGL3-Basic vector (Promega) employing KpnI/HindIII and KpnI/XhoI sites.

The QuickChange Mutagenesis kit (Stratagene, USA) was used to perform site-directed mutagenesis of the NurRE element within the human NMU promoter. The human NMU-luc construct was used as a template. Oligonucleotides with 34 to 36 nucleotides containing the desired point mutations were created as follows, with the mutated residues underlined: mNurRE, 5'-GTTCCCTGTTTTTCAAAAACAGGTCAAATA-3. The site-directed mutations were confirmed by sequence analysis.

Electrophoretic mobility shift assay (EMSA)

Nuclear extracts from HB1.F3 cells were prepared as described earlier (Kim et al., 2013). Sense and antisense oligonucleotides were annealed and then end-labeled with [γ-³²P]ATP (Amersham) and T4 polynucleotide kinase. Labeled probes were purified on 19% non-denaturing polyacrylamide gels. The DNA–protein binding reaction was performed in a final volume of 20 μl reaction buffer containing 10 mM Tris (pH 7.6), 50 mM NaCl, 0.5 mM dithiothreitol, 0.5 mM EDTA, 1 mM MgCl₂, 5% glycerol, and 250 μg of poly(dI–dC) per milliliter. Nuclear extract (20 μg of protein) was added to the reaction buffer in the absence or presence of unlabeled competitor DNA and pre-incubated for 15 min on ice. Radioisotope-labeled probes (50,000 cpm) were added, and the mixture was incubated for a further 30 min at room temperature. To resolve DNA–protein complexes, electrophoresis was performed on a 5% non-denaturing polyacrylamide gel. Gels were fixed, dried, and visualized by autoradiography. The oligonucleotides employed were as follows (only sense strands presented):

NMU-NurRE, 5'-GTTCCTCACCTTTCAAAGGGAGGTCAAATA-3';

NMU-mtNurRE, 5'-GTTCCCTGTTTTTCAAAAACAGGTCAAATA-3'; and

G0S2-NBRE2, 5'-CATCACTGACCTTTGCAATT-3'.

Chromatin immunoprecipitation (ChIP)

Soluble chromatin from HB.F3 cells was prepared as described previously (Kim et al., 2013) and immunoprecipitated with antibody against NURR1(TransCruz, sc-991X). Specific primer pairs were designed to amplify the promoter region of human NMU (forward, 5'-GGCTTAATTCCTAGCTGACGGGT-3'; and reverse, 5'-TCTAACTGGTCTGGAGCGCTGGT-3') and human G0S2 (forward, 5'-AGTGGGACCTTCGCGTGCAC-3'; and reverse, 5'-CCAGCTAATCTTGGGGAGGGCTC-3') from human genomic DNA. PCR conditions were optimized to allow for semiquantitative measurement, and PCR products were visualized on 1% agarose-Tris-acetate EDTA gels

Transfection and luciferase reporter assays

Transfections were performed using Lipofectamine PLUS reagent (Invitrogen) according to the manufacturer’s instructions. A total of 3 mg of DNA was used per transfection, as previously described (Kim et al., 2013), and briefly as follows: 1 to 2 mg of pLPCX or pLPCX-NURR1, 1 mg of reporter vector, and 0.5 μg of the internal control plasmid, pSV-β-galactosidase (Promega). All transfections were carried out three times in triplicate. Promoter activities were determined using the Single-Luciferase Reporter Assay System (Promega). Luciferase activities were normalized based on β-galactosidase activity in each well.

Profiling NURR1-regulated gene expression by microarray in hNSCs

Loss of Nurr1 resulted in ventral mDA agenesis (Zetterstrom et al., 1997). Despite the many approaches that have been used to identify Nurr1 target genes in rodents (Jacobs et al., 2009a; Sousa et al., 2007), the full range of Nurr1 target genes is unlikely to have been identified. In particular, NURR1 targets have not been systematically studied in DA neurogenesis through human embryo stages, despite the requirement for Nurr1 during this time for specification of ventral mDA neurons. Thus, to identify additional NURR1 targets in human systems, we performed microarray expression profiling of three hNSC lines overexpressing NURR1 (Fig. 1A).

These cell lines—HB1.F3, HB1.A4, and HB1.F5—represent separately isolated clones from primary human mesencephalon tissues immortalized with the v-myc gene (Kim, 2004; Kim et al., 2002) and express several markers, including nestin (Fig. 1B). In addition, microarray data showed that these cell lines highly express genes such as ALDH2, TNC, Nestin, and FOXA2 that are associated with DA progenitors, involved in DA differentiation, or characterize mature phenotypes (Supplementary Table S1). Based on these data, we determined to use the human NSC lines to investigate NURR1-related gene expression that is relevant for DA development. To trigger NURR1 overexpression in the hNSC lines, we transduced them with pseudovirus harboring plasmids overexpressing NURR1 cDNA (pLPCX or pLPC-NURR1). To avoid detecting secondary gene expression changes resulting from downstream effects of a direct NURR1 target, we treated cells with puromycin for only one week to select transduced cells and checked the surviving cells for NURR1 expression (Supplementary Fig. S1A). We used Affymetrix U133 Plus 2.0 oligonucleotide microarrays to analyze mRNA gene expression in RNA extracted from these cells (Supplementary Fig. S1).

Software made available by the lab of Michael Eisen at the University of California, Berkeley, was used to filter the microarray data for each cell line. Of the > 40,000 gene queries, we found 2,437 that were upregulated by NURR1 overexpression (617 in F3, 1,580 in A4, and 240 in F5 cells) and 1,399 that were downregulated (214 in F3, 865 in A4, and 320 in F5 cells). The expression profiling data of three hNSC lines with NURR1 overexpression have been submitted to the GEO (accession No. GSE58475).

In comparing the three cell lines, we found that genes showing regulation in a coordinated manner (i.e., in all three lines and in the same direction) were not highly represented. For example, we found no genes among any of the three cell lines that were within 30 orders of the greatest magnitude of change by NURR1. We inferred from this outcome that NURR1 may modulate many genes indirectly rather than directly (Supplementary Fig. S1B).

Direct targets of NURR1 in hNSC lines

For an understanding of the logic of NURR1-dependent neurogenesis, the newly identified downstream genes must be placed within an underlying regulatory hierarchy. The transcription factor Nurr1 is characterized by binding as a monomer to the NBRE sequence motif, as a heterodimer with retinoid X receptor (RXR) to DR5, or as a dimer to NurRE (Perlmann and Jansson, 1995). To filter genes that NURR1 targets directly, based on the presence of NURR1 consensus binding sites, we searched promoter sequences (retrieved from the University of California, Santa Cruz at http://genome.ucsc.edu/cgi-bin/hgTables) within 3 kb of the transcriptional start site (TSS) for NURR1 binding elements with no more than one base pair difference compared to the consensus sites (Fig. 1A). Of the putative NURR1 direct targets, using the Eisen lab software, we detected 187 sequences representing 180 genes, based on their NCBI UniGene designations, that were predicted to modulated by NURR1 overexpression (Supplementary Table S2). A hierarchical cluster diagram is shown in Fig. 1C, and Venn diagrams depicting the number of genes commonly and differentially up- or downregulated in the three hNSC lines are shown in Figs. 2A and 2B.

Based on the functions of the genes assigned by the Eisen lab software, NURR1 direct target genes were classified into categories (Fig. 1D). The classes of genes represented included nervous system development, cytoskeletal components, cell adhesion, DNA replication, regulation of progression through the cell cycle, transcription, and protein folding and secretion. In addition, we derived several important observations from this analysis. A significant proportion of the 180 differentially expressed genes were regulated in the same (coordinated) manner among the three hNSC lines, indicating transcriptional modulation that was independent of cell context (Fig. 1C, groups iv and v; Fig. 2C). Our microarray data also showed that several uncharacterized NURR1 targets, including chromosome 9 open reading frame 16 (C9orf16; highlighted in blue in Fig. 2D), are modulated in hNSC lines. Thus, NURR1 appears to modulate a number of characterized and uncharacterized genes, acting both as an activator and a repressor.

Genes regulated in a coordinated manner in all 3 lines

Next, to obtain a more restricted list of unambiguously regulated genes, we applied two stronger filters based on coordinated regulation of expression among the three cell lines (Fig. 2). First, we considered for further validation and quantitative analysis only those genes regulated in the same direction in all three NURR1-overexpressing hNSC lines (shown in Fig. 2C, Supplementary Tables S3 and S4). Of this group, 23 genes, including H2A.W histone (HIST3H2A), were upregulated, and 15 genes, including UDP-glucose ceramide glucosyltransferase (UGCG), were downregulated in all three lines. Almost all genes that were downregulated or upregulated in this coordinated manner showed a strong association with neuronal function, indicating a high stringency of our microarray data (Tables 1-3).

As a second step, we identified genes showing the greatest degree of change with NURR1 overexpression in each cell line (Fig. 2D). Upregulated genes that were within the 20th order of the greatest degree of change in all three lines were NMU, calponin 1 (CNN1), and C9orf16. The downregulated gene that met this criterion was G0S2.

Gene ontology (GO) and pathways of coordinately regulated genes

Because CNS development is still being investigated and novel functions and genes are always being added in this area, the GO analysis does not include the most recently published results. To identify in detail the function of each gene showing coordinated expression, we searched the possible functions of each candidate using UniProt (http://www.uniprot.org/), BioGPS, the Human Protein Atlas (http://www.proteinatlas.org/), Eukaryotic Linear Motif analysis, and PubMed. The results point to some striking features. In Tables 1-3, we present a detailed listing of these findings according to gene groups and pathways.

Overall, the analyses identified distinct categories of genes modulated by NURR1 in hNSCs, although they are not the well-known DA neuronal markers such as tyrosine hydroxylase (TH) and dopamine transporter (DAT). Among the gene sets showing coordinated upregulation by NURR1 in all three cell lines, we found a marked enrichment of synapse formation-associated dendritic morphogenesis categories as well as morphogenic movements of neuronal precursors during brain development, suggesting the importance of NURR1 in new synaptic formation and brain morphological development. Of the set of 38 genes included in this analysis, 13 of 23 upregulated genes (~56%) showed this enrichment, as did 1 of 15 downregulated genes (Table 1). Of this group, 10 genes, including sushi repeat containing protein X-linked 2 (SRPX2), localize in the synapse, primarily in postsynaptic elements (Table 1) (Sia et al., 2013; Soteros et al., 2018). In addition, fibronectin leucine-rich transmembrane protein 2 (FLRT2) was the most upregulated gene (3.5-fold) in the microarray data set for both HB1.F3 and HB1.F5 cells (Fig. 2D). Four genes, including FLRT2, regulate radial migration of pyramidal neurons and their tangential spread during embryonic brain development (Table 1) (Seiradake et al., 2014; Yamagishi et al., 2011). Of interest, C9orf16, the only uncharacterized gene, may function in synapse formation through trafficking of neurotransmitter-containing vesicles, based on its interaction with proteins in high-throughput proteomics analyses (UniProt site) (Huttlin et al., 2015; Marley and von Zastrow, 2010).

Moreover, nine genes were implicated in DNA replication and cell cycle categories. Eight of these were from a group of 15 downregulated genes (~53%) that were overrepresented among the coordinated NURR1-downregulated gene sets in the three cell lines. These results implicate NURR1 in cell cycle exit and in establishing post-mitotic DA neuroblasts during brain development (Table 2). In addition, nuclear and mitochondrially encoded genes, G0S2 and mitochondrially encoded NADH dehydrogenase 3 (MTND3), are modulated by NURR1 to increase oxidative phosphorylation (Lee et al., 2015; Yamada et al., 2013). mDA neurons are highly energy-demanding and have high requirements for ATP-generating oxidative phosphorylation. Our findings suggest that NURR1 is an important transcription factor for sustaining this distinguishing trait of DA neurons, consistent with previous results (Kadkhodaei et al., 2013). In this set of coordinated downregulated genes, we also detected leucine rich repeat containing 17 (LRRC17), UGCG, desmoplakin (DSP), and QKI, KH domain containing RNA binding (QKI), all of which occur in non-neuronal lineage differentiation categories, suggesting that NURR1 inhibits differentiation into non-neuronal lineages (Table 3).

Validation of NURR1 direct targets identified by microarray analysis

Of the 38 genes that were modulated in the same direction in the three hNSC lines (Fig. 2B, Supplementary Tables S3 and S4), 10 genes were chosen for validation by real time RT-PCR (Fig. 3). We selected eight genes involved in nervous system development: BAALC binder of MAP3K1 and KLF4 (BAALC) (Wang et al., 2005), calcium/calmodulin dependent protein kinase II inhibitor 1 (CAMK2N1) (Lucchesi et al., 2011; Saha et al., 2006), myosin light chain kinase (MYLK) (Li et al., 2016a), NMU (Kasper et al., 2016), C9orf16 (Huttlin et al., 2015; Marley and von Zastrow, 2010), prune homolog 2 with BCH domain (BMCC1) (Arama et al., 2012), FLRT2 (Seiradake et al., 2014; Yamagishi et al., 2011), and CNN1 (Blunk et al., 2014; Liu and Jin, 2016; Rami et al., 2006). The other two included genes, protein regulator of cytokinesis 1 (PRC1) (Polak et al., 2017) and G0S2 (Lee et al., 2015; Yamada et al., 2013), are involved in regulation of the cell cycle and in oxidative phosphorylation. Among the genes on this list, NMU, CNN1, C9orf16, and G0S2 were in the top 20 showing the greatest degree of change in all three lines. In all cases, we used real-time RT-PCR to confirm the direction of change in expression indicated by microarray (Figs. 3A-3C). In some cases, the array data suggested that NURR1 could modulate the same gene in opposite directions in different hNSC lines. Of these genes, we chose three for further analysis. Figure 3D shows real-time PCR data confirming this cell line–based directionality for decorin (DCN), pleiotrophin (PTN), and galanin and GMAP prepropeptide (GAL). In 30 of 39 (~77%) comparisons of gene expression using microarray and RT-PCR, results with the two methods showed agreement regarding the direction of change by cell line (Fig. 3).

NURR1 modulates target gene expression in a hNSC-specific manner

To discern which NURR1 target genes identified in this microarray analysis were involved only in hNSCs, we analyzed NURR1-associated gene expression in another cell type, human embryonic kidney HEK293 cells (Fig. 4A). In this cell line, expression of most target genes was not affected by increased levels of NURR1, except for NMU and BAALC, which showed 0.71- and –0.85-fold changes, respectively. This result contrasts with the direction of NURR1-associated gene expression change in hNSCs and suggests differential NURR1 functions in hNSCs versus HEK293 cells.

To functionally characterize the role of NURR1 in gene expression, we also inhibited NURR1 using siRNA (Fig. 4B). For this purpose, we focused on five genes (NMU, BAALC, CAMK2N1, FLRT2, and CNN1) that showed a consistent direction of change in expression levels on microarray and with RT-PCR (Fig. 3). Figure 4B demonstrates the successful inhibition of NURR1 using targeted siRNA compared to cells transfected with control siRNAs. The expression direction of the five genes was reversed in hNSCs expressing diminished levels of NURR1 compared to those overexpressing NURR1 (Fig. 4B). These data offer further evidence supporting the results of the real-time PCR analysis.

Transactivation of the human NMU and G0S2 promoters by NURR1

We next evaluated whether NURR1 may be involved in transcriptional regulation of human NMU and G0S2 genes. If NURR1 does regulate these genes during embryonic brain development, NURR1 binding elements may be evolutionarily conserved among species. We used the mVista (Frazer et al., 2004) program to identify evolutionarily conserved regions in the upstream sequences of the TSSs of human NMU and G0S2 that were upregulated or downregulated in the same direction in three cell lines (Fig. 5A). As noted above, NURR1 can bind as a monomer, a heterodimer, or dimer depending on the target. Two potential NURR1-binding elements (NBEs) were identified in each of the promoters: one NBRE and one NurRE in the NMU promoter, and two NBREs in the G0S2 promoter (Fig. 5A). Alignment of the human, rat, and mouse promoters revealed that NurRE in the human NMU promoter and NBRE2 in the human G0S2 promoter are perfectly conserved across all three species.

We performed EMSA and a ChIP assay to verify that NURR1 is recruited to these sites in in vitro and in vivo (Figs. 5B and 5C). EMSA showed that NURR1 could bind evolutionarily conserved NBEs in both promoters in vitro (Fig. 5B). In addition, competition experiments and supershift assay confirmed that NURR1 specifically bound NurRE in the NMU promoter and NBRE2 elements in the G0S2 promoter (Fig. 5B). For the ChIP assay, HB1.F3 cells overexpressing NURR1 were cross-linked to chromatin and immunoprecipitated with an antibody against NURR1, followed by PCR using primers corresponding to NurRE or NBRE2. As shown in Fig. 5C, NURR1 was recruited to the NBEs in both promoters.

To verify whether the predicted binding sites for NURR1 are functional, we cloned a 2102-bp fragment of the human genomic sequence upstream of the NMU translation start site into a luciferase reporter plasmid. Transient transfection showed that NURR1 overexpression increased transcription from this promoter by 7.1-, 13.0-, and 3.5-fold in F3, A4, and F5 cells, respectively (Fig. 5D). To investigate whether the evolutionarily conserved NurRE is important for NURR1-dependent activation of the human NMU promoter, we used site-directed mutagenesis to introduce point mutations in the promoter that disrupt these sequences. These mutations clearly reduced the ability of NURR1 to activate the promoter (approximately 66.6%, 56.8%, and 31.3% decreases in F3, A4, and F5 cells, respectively) (Fig. 5E). The NurRE element thus is essential for NURR1-dependent NMU expression in hNSCs. We also evaluated the transcriptional regulation of the G0S2 promoter by NURR1 (Fig. 5F). A 2082-bp fragment of the human G0S2 promoter was cloned into a luciferase reporter plasmid. As shown in Fig. 5F, NURR1 comprehensively repressed human G0S2 promoter activity in all three hNSC lines in a dose-dependent manner. Collectively, the effect of NURR1 on the promoter-driven reporter activity in the NMU and G0S2 genes paralleled changes in mRNA levels detected by microarray and RT-PCR.

Expression of Nurr1 direct target genes from the ventral midbrain of embryonic mice during in vivo mDA differentiation

To confirm the in vivo expression pattern of these gene sets during midbrain development, we performed quantitative PCR analysis of mRNA extracted from the mouse embryonic ventral midbrain (Fig. 6). Because Nurr1 expression becomes evident at E10.5 during mDA cell fate specification (Fig. 6A) (Perlmann and Wallen-Mackenzie, 2004), we prepared embryonic midbrain tissues from E9.5 to E13.5 for RT-PCR analysis (Figs. 6B and 6C). We detected an increase in expression of Nmu, C9orf16, dihydropyrimidinase-like 3 (Dpysl3), G0s2, and Baalc at successive intervals from E9.5 to E13.5 (Figs. 6B and 6C). In addition, Camk2n1 and Cnn1 expression became evident at E10.5 and was maintained thereafter (Figs. 6B and 6C). Figure 6D shows a summary of the timing of target gene expression. These results validate our microarray data and indicate that these genes are truly in vivo targets of Nurr1.

BioGPS analysis: CNS-specific expression of NURR1 direct target genes

To discern which Nurr1 target genes identified in our microarray analysis are biologically relevant, we searched each gene on the BioGPS website (http://biogps.org), using the GeneAtlas U133A gene expression (Wu et al., 2016). The most highly represented systems were CNS tissues, and other systems did not show a significant percentage of tissue-specific genes. As indicated in Table 4, 12 transcripts (31.6%) exhibited ≥ 1.5 fold increased expression in the CNS or fetal brain as compared to the mean expression values of all human tissues in the datasets. The gene symbols and CNS and fetal brain enrichment for these 12 genes are provided in Table 4, and additional details on the microarray probe ID and CNS versus average human expression levels are given in >Supplementary Table S5. Of these 12 genes, three (BAALC, CAMKIIN1, and QKI) were highly expressed, with an expression value >30 multiples of the median in the brain/CNS, in accordance with previously established stringency thresholds (Pennings et al., 2011). The other dominant genes were DPYSL3 and FLRT2, both of which were specifically expressed in fetal brains.

GWAS of NURR1 direct target genes

NURR1 is a Parkinson’s disease susceptibility gene because it is a key factor in the neurodevelopmental control of DA neurons in the substantia nigra pars compacta and ventral tegmental area. In addition, NURR1 mutations are associated with both schizophrenia and bipolar disorder (Buervenich et al., 2000; Vuillermot et al., 2011; Xing et al., 2006). Considering a possible role for NURR1 direct target genes in DA neurogenesis, we hypothesized that genetic alterations within the loci of these genes (including non-coding DNA sequences and regulatory elements) may be associated with susceptibility to schizophrenia, bipolar disorder, and depression, as well as Parkinson’s disease. For this reason, we searched for SNPs associated with susceptibilities to these conditions, focusing on the loci of NURR1 direct target genes from the previously identified variants of GWAS datasets derived using different chip platforms (Li et al., 2012; 2016b).

About half of the up- and downregulated targets showed associations with these diseases (Supplementary Table S6). Of these targets, we identified 11 SNPs with Parkinson’s disease susceptibilities at seven independent loci: rs1385331 (chr 3) for phosphatidylinositol-glycan biosynthesis class X (PIGX); rs10056132 (chr 5) and rs10041339 (chr 5) for DPYSL3; rs4874150 (chr 8) for tissue specific transplantation antigen P35B (TSTA3); rs2453998 (chr 8) for BAALC; rs3003602 (chr 9) and rs2502731 (chr 9) for C9orf16; rs1011711 (chr 9) for PRUNE; and rs5966709 (chr X), rs4828037 (chr X), and rs932437 (chr X) for SRPX2 (Table 5). All of the variations were found in non-coding DNA within ~50 kb of NURR1 target genes. Because a non-coding DNA sequence at the loci of direct target genes may encode information guiding the proper level, time, and place of gene expression during DA neurogenesis (Lee et al., 2020), we searched the genomic region for each variant associated with Parkinson’s disease, and annotated this region for sequence conservation and candidate enhancers based on H3K4me1 and H3K27ac histone modifications and DNase hypersensitive sites (Fig. 7). Of 11 SNPs with Parkinson’s disease susceptibilities, two variations in the C9orf16 locus (rs3003602 and rs2502731), two in DPYSL3 (rs10056132 and rs10041339), one in BAALC (rs2453998), and one in PIGX (rs2453998) directly overlapped with a H3K4me1 and a H3K27ac marker of regulatory element sequences (candidate enhancer). These variants were not evolutionarily conserved except in primates, and we did not detect any overlap between the other five polymorphisms and any regulatory feature. These GWAS data strongly support the inference that these polymorphisms are likely to disrupt an enhancer of C9orf16, DPYSL3, BAALC, and PIGX genes. In addition, there were many SNPs associated with schizophrenia, bipolar disorder, and depression, which are implicated in DA neuronal activity in the CNS (Supplementary Table S6).

Human brain development is very different from mouse brain development (Bernard et al., 2012; Bohland et al., 2010; La Manno et al., 2016; Lein et al., 2007; Myers et al., 2015; Ng et al., 2009). Previous reports have described gene expression profiling by the DA transcription factor Nurr1 in a rodent system (Jacobs et al., 2009a; Sousa et al., 2007) but did not confirm whether the modulated genes are direct targets of Nurr1. We addressed this question by identifying NURR1 direct target genes in human neural progenitors.

Our analyses have yielded several important findings. First, several of the 38 genes differentially and directly affected by NURR1 were regulated in the same (coordinated) manner among the three hNSC lines, indicating that NURR1 transcriptional regulation occurs independent of cell context. The results also suggest that NURR1 modulates many genes, acting both as an activator and a repressor. These findings contrast with earlier reports (Kadkhodaei et al., 2013; Panman et al., 2014) describing Nurr1 as functioning mainly as a transcriptional activator to regulate a battery of genes expressed in DA neurons. This discrepancy may trace to the fact that these authors analyzed Nurr1 ablation of already mature and postmitotic DA neurons.

Moreover, genes associated with synapse formation (6 genes), neuronal morphogenesis (4 genes), and cell migration during brain development (4 genes) were highly upregulated in the current study, whereas genes involved in cell cycle progression and DNA replication (9 genes) were downregulated. In the E10.5-E11.5 mouse embryo, a series of events is important for mesencephalic DA cell differentiation (Wallen et al., 1999). Specifically, DA progenitor cells migrate from the ventricular zone to the mantle zone, exit the cell cycle (become post-mitotic), and start to innervate the striatal target area. Simultaneously, Nurr1 expression is strongly upregulated. Our microarray data clearly fill in gaps from previous reports showing that in contrast to wild-type animals, in Nurr1–/– mutant mice, DA progenitor cells in the midbrain cannot migrate laterally, fail to innervate the striatal target area, and become apoptotic (Wallen et al., 1999). These authors did not directly demonstrate that the reason for these failures is a lapse in regulation of genes involved in migration, synaptogenesis of DA progenitor cells, and cell cycle exit. To date, a direct relationship of Nurr1 expression with other DA differentiation–related events has not been demonstrated. Our findings fill in these gaps, providing direct evidence that NURR1 selectively contributes to the expression of genes involved in the correct positioning of migrating neurons and in synapse formation with correct partners, as well as ensuring postmitotic status during mDA neuron development.

We also found that genes encoding proteins involved in regulating the DA neurotransmitter phenotype and expressed in the final stage of DA differentiation, including TH, vesicular monoamine transporter, dopa decarboxylase, and DAT, were not upregulated by NURR1 in our hNSC lines (Supplementary Table S2). This result is consistent with our earlier work (Kim et al., 2013) and with another previous report (Kadkhodaei et al., 2013). An explanation may lie in the timing of expression relative to the window of development we examined here. Our hNSC lines show gene expression profiling similar to what is seen from E9.5 to E10.5 (Supplementary Table S1), with in vivo Nurr1 expression becoming evident in the midbrain and DA precursors shifting to a postmitotic DA neuroblast state. These other well-known markers of DA neurons may not have shown expression changes because the timing of their action is outside of this window. Nevertheless, we did identify a large number of neurogenesis-related genes as direct targets of NURR1, all candidates that, to the best of our knowledge, have not yet been described previously as downstream of NURR1. The restricted modulation pattern of these genes reveals potential roles for NURR1 in brain development, signaling, metabolism, and disease. Taken together, our findings are consistent with published reports that in adult DA neurons, Nurr1 ablation is associated with a progressive pathology involving reduced striatal DA, impaired motor behaviors, and dystrophic axons and dendrites, as well as downregulation of oxidative phosphorylation–related genes in Parkinson’s disease (Kadkhodaei et al., 2013).

This research was supported by the Bio & Medical Technology Development Program of the National Research Foundation (NRF) funded by the Korean government (MSIT) (2015M3A9C6028956), and by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (2013064121).

S.M.K., M.W.K., and S.R.R. designed and performed the experiments, analyzed data, and wrote the manuscript. S.Y.C. carried out bioinformatic analyses from microarray data and statistical analyses. H.S.S. helped to conduct collection of mouse embryonic brain tissues and performed its gene expression analysis. Y.H.S. and D.G. performed immunohistochemistry in mouse embryonic brain. M.A.L. supervised the entire project, designed the experiments, and wrote the manuscript.

The authors have no potential conflicts of interest to disclose.

Genes associated with synapse formation, neuronal morphogenesis and cell migration during brain development

Genes associated with DNA replication, cell cycle progression, apoptosis and mitochondrion functions

Genes associated with differentiation into non-neuronal lineages

BioGPS gene expression analysis of Nurr1 direct target genes

Genetic Associations for Parkinson’s Disease at Nurr1 direct target genes