Neurotoxin-Induced Pathway Perturbation in Human Neuroblastoma SH-EP Cells

Jin Hwan Do

Abstract

The exact causes of cell death in Parkinson’s disease (PD) remain unknown despite extensive studies on PD.The identification of signaling and metabolic pathways involved in PD might provide insight into the molecular mechanisms underlying PD. The neurotoxin 1-methyl-4-phenylpyridinium (MPP+) induces cellular changes characteristic of PD, and MPP+-based models have been extensively used for PD studies. In this study, pathways that were significantly perturbed in MPP+-treated human neuroblastoma SH-EP cells were identified from genome-wide gene expression data for five time points (1.5, 3, 9, 12, and 24 h) after treatment. The mitogen-activated protein kinase (MAPK) signaling pathway and endoplasmic reticulum (ER) protein processing pathway showed significant perturbation at all time points. Perturbation of each of these pathways resulted in the common outcome of upregulation of DNA-damage-inducible transcript 3 (DDIT3). Genes involved in ER protein processing pathway included ubiquitin ligase complex genes and ER-associated degradation (ERAD)-related genes. Additionally, overexpression of DDIT3 might induce oxidative stress via glutathione depletion as a result of overexpression of CHAC1. This study suggests that upregulation of DDIT3 caused by perturbation of the MAPK signaling pathway and ER protein processing pathway might play a key role in MPP+-induced neuronal cell death. Moreover, the toxicity signal of MPP+ resulting from mitochondrial dysfunction through inhibition of complex I of the electron transport chain might feed back to the mitochondria via ER stress. This positive feedback could contribute to amplification of the death signal induced by MPP+.

Keywords: 1-methyl-4-phenylpyridinium, gene regulation, Parkinson’s disease, pathway perturbation, SH-EP cells

INTRODUCTION

Parkinson’s disease (PD) is a progressive neurological disorder that results primarily from the death of dopaminergic (DAergic) neurons in the substantia nigra. The factors that trigger cell death in PD are currently unknown although neuromelanin accumulation, mitochondrial dysfunction, oxidative stress, exposure to iron and other metals, mutations in the α-synuclein gene, trauma, and dysfunction of the ubiquitin-proteasome system have all been implicated in the pathogenesis of PD (Gálvez-Jiménez, 2007). Since the discovery that people who are intoxicated with 1-methyl-4-phenyl-1, 2, 3, 6-tetrahydropyridine (MPTP) develop a syndrome nearly identical to PD (Langston et al., 1983), MPTP has been used to generate PD models in non-human primates and mice. In the brain, MPTP is oxidized to 1-methyl-4-phenyl-2,3-dihydropyridinium (MPDP+) by monoamine oxidase B (MAOB) in glia and serotonergic neurons, which is then converted to 1-methyl-4-phenyl-pyridium (MPP+), an active metabolite of MPTP. MPP+ is taken up by DAergic neurons via dopamine and noradrenaline transporters, which results in inhibition of complex I of the mitochondrial electron transport chain and the formation of reactive oxygen species (ROS) (Lotharius and O’Malley, 2000; Nakamura et al., 2000), which in turn leads to cellular dysfunction and cell death (Nicotra and Parvez, 2002).

In vitro models such as pheochromocytoma (PC12) cells and human neuroblastoma (SH-SY5Y or SH-EP) cells have been widely used in studies to understand the pathogenesis of PD. When these cells were treated with neurotoxins such as MPP+, 6-hydroxydopamine (6-OHDA), or rotenone, many aspects of the DAergic neuron death are similar to those observed in PD. Conn et al. (2003) explored the transcriptional response of SH-SY5Y cells to MPP+ using cDNA microarray and identified transcription factors and cell cycle proteins as major classes of differentially expressed proteins. DNA microarray analysis is capable of profiling the expression levels of a large number of genes simultaneously and is therefore a promising technology for the elucidation of gene interactions. However, in many cases DNA microarray experiments yield a list of differentially expressed genes (DEGs) and it is not easy to link these DEGs to phenotypes that are the results of thousands of complex interactions occurring in various metabolic and signaling pathways. This indicates that the complex interactions in living cells should be taken into account to fully understand DNA microarray data. Recently, many efforts have been made to identify which biological pathways are perturbed using gene expression profiles (Pepe and Grassi, 2014; Tarca et al., 2009). Information for each pathway can be obtained from pathway databases such as KEGG (Kanehisa et al., 2012) and Reactome (Croft et al., 2014), which describe metabolic pathways and gene signaling networks and offer the potential for more complex and useful analysis. Identifying which pathways are perturbed in PD is important for improving our understanding of the molecular mechanisms underlying PD. In this study, I identified pathways that were significantly perturbed after treatment with MPP+ using the list of DEGs and their expression values with the commercial software Pathway Guide 3.0 (http://www.advaitabio.com).

Perturbed pathways at 1.5, 3, 9, 12, and 24 h after MPP+ treatment were examined using DEGs identified from whole genome expression data measured at each time point and their expression values. Two KEGG pathways, the mitogen-activated protein kinase (MAPK) signaling pathway and endoplasmic reticulum (ER) protein processing pathway, showed significant perturbation persistently and both resulted in positive perturbation of DNA-Damage-Inducible Transcript 3 (DDIT3). This provides a plausible explanation for the strong upregulation of DDIT3 observed in both microarray and real-time quantitative PCR (qPCR) analyses. Information on pathways that are perturbed by MPP+ treatment might provide insight into MPP+-induced neuronal cell death and the pathogenesis of PD.

MATERIALS AND METHODS

Cell culture and MPP+ treatment

Human neuroblastoma SH-EP cells were kindly provided by Dr. Talia Hahn at Kaplan Medical Center (Rehovot, Israel). For treatment with MPP+ iodide (Sigma-Aldrich, USA), 1 × 106 cells were plated in 100 mm2 dishes (Corning, USA) in 10 ml Dulbecco’s modified Eagle’s Medium (DMEM; Sigma-Aldrich) with 10% fetal bovine serum (FBS), 100 units/ml penicillin, and 100 mg/ml streptomycin at 37°C with 5% CO2 and cultured for 3 days. Freshly prepared MPP+ toxin was added to the cultures to a concentration of 1.25 mM and incubation continued at 37°C for 0 (control), 1.5, 3, 9, 12, and 24 h. The experiment was performed two times.

Assessment of cell viability

Cell viability was measured using the quantitative colorimetric MTT assay, which reveals the mitochondrial activity of living cells, as described previously (Nanjo et al., 1996). Briefly, MTT dissolved in phosphate-buffered saline was added at the end of incubation to a final concentration of 0.5 mg/ml. After 4 h incubation at 37°C and 5% CO2, the supernatants were removed, and the formazan crystals that formed in the viable cells were measured at 550 nm using a microplate reader (Molecular Devices, USA).

RNA isolation and quality control

At each time point after MPP+ treatment, total RNA was extracted using TRIzol® (Invitrogen Life Technologies, USA) and purified using RNeasy columns (Qiagen, USA) according to the manufacturers’ protocol. Purified RNA samples with DNase treatment were quantified. All aliquots were stored at −80°C until use. For quality control, RNA purity and integrity were evaluated by denaturing gel electrophoresis, optical density comparison of the 260/280 ratio, and analysis using an Agilent 2100 Bioanalyzer (Agilent Technologies, USA).

Microarray expression profiling

Samples were collected at 0 h (control, before MPP+ treatment) and at 1.5, 3, 9, 12, and 24 h after MPP+ treatment. Total RNA was isolated from each sample, amplified, and purified using the Ambion Illumina RNA amplification kit (Ambion, USA) to generate biotinylated cRNA. Briefly, total RNA was reverse-transcribed to single stranded cDNA using a T7 oligo(dT) primer, converted into double-stranded cDNA, and purified. An in vitro transcription reaction was then carried out in the presence of biotinylated UTP and CTP to produce biotin-labeled cRNA from double stranded cDNA. After purification, the cRNA was quantified using a ND-1000 Spectrophotometer (NanoDrop, USA) and 750 ng of labeled cRNA was hybridized to each human HT-12 expression v.4 bead array for 16–18 h at 58°C, according to the manufacturer’s instructions (Illumina). Detection of the array signal was performed using fluorolink streptavidin-Cy3 (GE Healthcare Bio-Sciences, UK) according to the bead array manual. Arrays were scanned with an Illumina bead array reader confocal scanner according to the manufacturer’s instructions.

Normalization of microarray data

Probe signals exported from GenomeStudio software v2011.1 (Illumina) were checked for detection against negative controls with a GenomeStudio internal algorithm and missing values were introduced to replace signals under the detection limit. After the probe signals were log2 transformed, quantile normalization was performed with the function lumiN from the lumi bioconductor package (www.bioconductor.org). As these normalized values were presented in terms of probe IDs, among which multiple IDs might correspond to the same gene, the probe IDs targeting the same gene were averaged and conversion of the probe ID to the gene name was performed using the lumi bioconductor package. After removing the probe IDs with no gene name, 48,803 probe IDs were finally mapped to 32,421 genes with a unique gene name (i.e., gene symbol).

Selection of differentially expressed genes

Genes that were differentially expressed between MPP+ treated and untreated (control) samples were identified by significance analysis of microarrays (SAM) (Tusher et al., 2001), which assigns a scoredi to each gene i based on the change in gene expression relative to the standard deviation of repeated measurements. SAM performs permutations of the repeated measurements to estimate the percentage of genes identified by chance, termed the false discovery rate (FDR). As the number of MPP+-treated and untreated samples is two for each time point in this study, the total number of distinct permutations is 6 (C24 = 6) for each time point. TheΔ value, which represents the distance from the expected di computed from permutation, was used for threshold for significance. A SAM analysis was carried out using samr in the R package (www.r-project.org). Genes with > 1.5-fold change andΔ = 0.4 were considered significantly regulated.

Identification of pathways that are perturbed by MPP+ treatment

Significantly perturbed pathways in MPP+ treated human neuroblastoma SH-EP cells were identified using the commercial software Pathway Guide 3.0 (http://www.advaitabio.com/). This program explores significantly affected pathways based on two types of evidence: (i) over-representation of DEGs in a given pathway, and (ii) abnormal perturbation of that pathway, measured by propagating measured expression changes across the pathway topology. These two types of evidence are captured by two independent probability values, PNDE and PPERT. This method provide the possibility to account for system-level dependencies and interactions, as well as identify perturbations and modifications at the pathway by overcoming the limitation of the most currently available approaches such as over-representation analysis (Drăghici et al., 2003) and functional class scoring (Goeman et al., 2004), which consider only the set of genes on any given pathway and ignore their position in those pathways.

The first probability, PNDE= P(X ≥ NDEH0), captures the significance of the given pathway Pi as provided by an over-representation analysis (ORA) (e.g., Doniger et al., 2003; Pandey et al., 2004) of the number of DEGs (NDE) observed on the pathway. H0 stands for the null hypothesis: that the genes that appear as DEG on a given pathway are completely random. The PNDE value represents the probability of obtaining a number of DEGs on the given pathway at least as large as the observed one, NDE. The second probability, PPERT, is estimated based on the amount of perturbation measured in each pathway. To calculate PPERT for a given pathway, a gene perturbation factor (PF) is defined as described by Tarca et al. (2009):

$P F ( g i ) = Δ E ( g i ) + ∑ j = 1 n β i j P F ( g i ) N d s ( g j )$

where ΔE(gi) is the signed normalized expression change in the gene gi (log2 fold-change in this study) and the second term is the sum of the perturbation factors of the genegjdirectly upstream of the target gene gi, normalized by the number of downstream genes of each such gene Nds(gj). The absolute value of βij quantifies the strength of the interaction between genes gj and gi. In this study, the value of β was +1 for induction/activation, and −1 for repression/inhibition. A non-zero value of β was only assigned for genes that directly interacted with the gene gi according to the pathway description. Given an edge directed from gene A to gene B in the pathway, A is considered upstream of B or B is considered downstream of A. The net perturbation accumulation at the level of each gene, ACCg, is calculated as the difference in the perturbation factor PF of a gene and its log2 fold-change:

$A c c ( g i ) = P F ( g i ) − Δ E ( g i )$

The total net accumulated perturbation of a given pathway is calculated as the sum of all perturbation accumulations for all genes in the pathway: $t A = ∑ i A c c ( g i )$ . The probability, PPERT, is defined as the probability of observing a total accumulated perturbation of the pathway, TA, more extreme than tA just by chance:

$P P E R T = P ( T A ≥ t A | H 0 )$

Where H0 is the null hypothesis for the distribution of TA values. The null distribution of TA values is empirically determined using a bootstrap approach. In this study, the bootstrapping number of iterations is set to 10,000.

To identify significantly perturbed pathways in human neuroblastoma SH-EP cells after MPP+ treatment, the two types of evidence, PNDE and PPERT, are combined into one global probability, PG. The global probability PG is calculated as in Tarca et al. (2009):

$P G = C i − C i • ln ( C i )$

Where ci= PNDE(iPPERT(i). Because PG is a combined probability value, it can be used not only to rank the pathways but also to choose a desired type I error level. Since more than a hundred pathways (e.g., 146 pathways for KEGG database) are tested simultaneously, small PG values can occur by chance. Therefore, the significantly perturbed pathways were identified by controlling the FDR with the popular FDR algorithm (Benjamini and Yekutieli, 2001).

cDNA synthesis and real-time qPCR

cDNA was produced using the Superscript™II RT-PCR System (Invitrogen, Germany) according to the manufacturer’s recommendations for oligo(dT)20 primed cDNA-synthesis. cDNA synthesis was performed on 500 ng of total RNA at 42°C. Using this first-strand cDNA, the gene expression levels of the 21 test genes and one control gene were quantified using TaqMan technology on a ABI PRISM 7900HT Sequence Detection System (Applied Biosystems, USA) in 384-well microtiter plates using a final volume of 10 μl Gene-specific primers and probes were available as TaqMan Gene Expression Assays (Applied Biosystems, Table 1). Optimum reaction conditions were obtained with 5 μl of Universal Master Mix (Applied Biosystems, USA) containing dNTPs with UTP, MgCl2, reaction buffer and AmpliTaq Gold® DNA polymerase, 90 nM of primer(s) and 250 nM fluorescence-labeled TaqMan probe. Finally, 2 μl template cDNA was added to the reaction mixture. Glyceraldehyde-3-phosphate dehydrogenase (GAPDH) was used as the internal control. Amplifications were performed with a 10 min template denaturation step at 95°C, followed by 40 cycles at 95°C for 15 s and 60°C for 1 min. All samples were amplified in triplicate and data were analyzed with Sequence Detector software (Applied Biosystems). The relative mRNA levels were estimated by the 2ΔΔCT method (Livak and Schmittgen, 2001).

RESULTS

Effect of MPP+ concentration on cell viability

To investigate the neurotoxicity of MPP+ to SH-EP cells, cell viability was measured 48 h after MPP+ treatment. Greater than 50% cell death was observed at MPP+ concentrations of 1.25 mM or higher (Fig. 1A). This indicates that 1.25 mM MPP+ is sufficient to induce apoptosis in SH-EP cells. To capture the early cellular events triggered by MPP+ treatment, it is important to control the exposure time to MPP+. The time-dependent cell viability for 1.25 mM MPP+ is shown in Fig. 1B. Cell viability greater than approx imately 78% could be achieved for up to 24 h of treatment. Therefore, in this study I examined gene expression and pathway perturbation for up to 24 h after MPP+ treatment.

The viability of MPP+ treated human neuroblastoma SH-EP cells. (A) Cell viability after cultivation for 48 h with various MPP+ concentrations. (B) Time-dependent cell viability after treatment with 1.25 mM ...

Selection of DEGs from a time series microarray experiment

Gene expression for MPP+-treated SH-EP cells was measured at six time points between 0 and 24 h using the human HT-12 expression v.4 bead array that included 47,231 probes with well-established or provisional annotation. A total of 12 bead arrays were used for examination of gene expression at the six time points (0, 1.5, 3, 9, 12, and 24 h; two replicates for each time point). The bead-summary data were log2 transformed and normalized using quantile normalization (see “Materials and Methods”). After averaging of multiple probe IDs targeting the same gene, 48,803 probe IDs were finally mapped to 32,421 genes with a unique gene name, i.e., gene symbol. DEGs between MPP+-treated and untreated SH-EP cells (control) were identified for each time point by SAM analysis with R package samr (see “Materials and Methods”).

The proportion of DEGs that were downregulated increased with time of exposure to MPP+ (0, 20, 55, 57 and 53% respectively for 1.5, 3, 9, 12 and 24 h) (Fig. 2). The fact that no DEGs with down-regulation were detected at 1.5 h suggests that the effect of MPP+ treatment might be initialized by activation of MPP+ responsive genes. The numbers of up-regulated DEGs were 844, 803, 567, 797, and 1,679 respectively for 1.5, 3, 9, 12, and 24 h (Fig. 2). Of these, 55 genes were commonly up-regulated at all time points, as shown in Table 2. As the expression of these genes was robustly upregulated during the time course, they might play a core role in neuronal cell death induced by MPP+. The continued upregulation of genes such as brain-derived neurotrophic factor (BNDF) and fibroblast growth factor 2 (FGF2) is particularly interesting. Almeida et al. (2005) reported that BNDF protects neurons via transient activation of the Ras/MAPK pathway and the PI3-K/Akt pathway. In addition, it was reported that intrastriatal infusions of FGF2 induced recovery of striatal DAergic fibers and DAergic content in the mouse MPTP model, and increased protection in the neurotoxicity induced lesion of nigrostriatal DA system (Date et al., 1993). Continued upregulation of the BNDF and FGF2 genes might be a cell survival strategy against MPP+ toxicity. Extraordinary high expression was observed at two genes including MIR1974 and CLDN1. MIR1974 is considered non-canocical miRNA because it maps to the mitochondrial tRNA and rRNA genes (Bandiera et al., 2011), and its upregulation was also detected in MPP+ treated SH-SY5Y cells (Kimet et al., 2013). CLDN1 plays a major role in tight junction-specific obliteration of the intercellular space and it has been reported that its reexpression induces apoptosis in breat tumor spheroids (Hoevel et al., 2004). The percentage of DEGs for each time point compared with all genes on the array (i.e., 32,241 genes) was 2.6, 3.1, 3.88, 5.69, and 10.9% respectively for 1.5, 3, 9, 12, and 24 h.

Upregulated/downregulated DEGs at 1.5, 3, 9, 12, and 24 h after MPP+ treatment. The red and green bars represent the number of upregulated and downregulated genes, respectively.

Analysis of perturbed pathways in MPP+-treated SH-EP cells by Pathway Guide 3.0

Identifying the pathways that are significantly affected by MPP+ exposure is a crucial step in understanding the underlying molecular mechanisms of MPP+-induced cell death. To examine pathways that are significantly perturbed by MPP+ treatment, two independent probability values, PNDE and PPERT, were calculated for each pathway with incorporating parameters such as the normalized fold change of the DEGs, statistical significance of the set of pathway genes, and the topology of the signaling pathway using the commercial software Pathway Guide 3.0 (http://www.advaitabio.com; see “Materials and Methods” for details). These two types of evidence were combined into a global probability value, PG, which is used to rank the pathways and test the hypothesis that the pathway is significantly perturbed by MPP+ treatment. The perturbation of pathways at 1.5, 3, 9, 12, and 24 h after MPP+ treatment was examined for 146 well-characterized human gene signaling pathways in KEGG database. In Fig. 3, two types of p-value are illustrated in two-dimensional plots for each time point. In this figure, the horizontal axis is used to plot the −log PNDE (over-representation evidence) and the vertical axis is used to plot −log PPERT (perturbation evidence). Since several pathways are tested simultaneously, the significance should be considered after adjustment of the p value for multiple comparisons. Pathways above the oblique red line are significant at 5% after Bonferroni correction (Hochberg, 1988) whereas those above the blue line are significant at 5% after FDR correction (Benjamini and Yekutieli, 2001). The number of significantly perturbed pathways dramatically increased at 24 h. This suggests that it might be difficult to restore the normal state of cells after 24 h of MPP+ treatment.

Two-dimensional plots showing two types of evidence, PNDE and PPERT, for each time point (1.5, 3, 9, 12, and 24 h) after MPP+ treatment. The numbers beside the dots on ...

The top 10 pathways ranked by PG at each time point are shown in Table 3. When the threshold for significance after FDR correction was 5%, all pathways listed in Table 3 satisfied this criterion. In particular, the MAPK signaling pathway showed significant perturbation at all time points. Closer inspection of genes associated with this pathway at 9 h after MPP+ treatment showed that several genes encoding ligands, including BDNF, FGF2, and transforming growth factor beta 3 (TGFB3), and, voltage-dependent L-type calcium channel subunit beta-3 (CACNB3) were differentially expressed (Fig. 4A). As these genes are located at very upstream in the pathway and affect entry points controlling this pathway, their perturbations are widely propagated throughout the pathway. MAPK genes including MAP2K5, MAP3K6, MAP4K2, MAPK12, MAPK8IP1 and MAPK7 were downregulated whereas genes near the end of pathways such as activating transcription factor 4 (ATF4), avian myelocytomatosis viral oncogene homolog (MYC), and DDIT3 showed upregulation. Based on these observed expression values, the total perturbation of each gene and its propagation in the MAPK signaling pathway was estimated (see “Materials and Methods”) and is shown in Fig. 4B. The red arrows represent propagation of perturbation. Through propagation of the perturbation of upstream genes, downstream genes were either upregulated (MAPK8/9/10 and MAP3K4) or downregulated (MAP3K5). The MAPK signaling pathway was more severely perturbed at 24 h after MPP+ treatment (Fig. 5A). At this time point, many ligand-encoding genes such as nerve growth factor (NGF), BNDF, FGF2, interleukin 1, alpha (IL1A), and TGFB2/3 were perturbed, as well as genes encoding receptors including platelet-derived growth factor alpha (PDGFRA), Fas cell surface death receptor (FAS), and monocyte differentiation antigen CD14 (CD14), with upregulation of MAP2K4, MAP3K1, and MAP3K7. The propagation of IL1A perturbation via interleukin 1 receptor type I (IL1R1) induced upregulation of caspase 3 (CASP3), whereas perturbation of FAS and TGFB1 resulted in upregulation of death-domain associated protein (DAXX) (Fig. 5B). Perturbation of the MAPK signaling pathway for 9 and 24 h ultimately had the same result of upregulation of DDIT3, ATF4, and v-rel avian reticuloendothe-liosisviral oncogene homolog B (RELB). Upregulation of DDIT3 has been also reported in MPP+-treated human neuroblastoma SH-SY5Y cells (Conn et al., 2002).

Perturbation of the MAPK signaling pathway at 9 h after MPP+ treatment. (A) Perturbed (upregulated/downregulated) genes in the MAPK signaling pathway. (B) Propagation of perturbation in the MAPK signaling pathway ...
Perturbation of MAPK signaling pathway at 24 h after MPP+ treatment. (A) Perturbed (upregulated/downregulated) genes in the MAPK signaling pathway. (B) Propagation of perturbation in MAPK signaling pathway and total ...

The other pathway that was frequently perturbed was the ERprotein processing pathway (Table 3). Perturbation of genes involved in this pathway at 24 h is shown in Fig. 6A. Altered expression of ubiquitin ligase complex genes such as ring finger protein 5 (RNF5), VCP-interacting membrane protein (VIMP), DnaJ homolog, subfamily b, member 2 (DNAJB2), ubiquitin-conjugating enzyme E2G 2(UBE2G2), and F-box protein 6 (FBXO6) might be a response to the ER stress caused by misfolded proteins. In this pathway, the most strongly upregulated gene was DDIT3, with a log2 fold change of 4.4. The B-cell lymphoma 2 (BCL2) also showed upregulation (log2 fold change = 0.95), but is inhibited by DDIT3, therefore overall BCL2 is strongly negatively perturbed by the propagation of DDIT3 perturbation (Fig. 6B). Moreover, apoptosis induced through a cascade starting with Bcl-2 homologous antagonist/killer (BAK1) and ending with Caspase 12 (CASP12) might be inhibited because BCL2 inhibits BAK1. In addition, downregulation of ER-to-nucleus signaling 1 (ERN1) might inhibit apoptosis through a cascade starting from ERN1 and ending with MAPK8. DDIT3 was markedly positively perturbed by upstream genes including ATF4, nuclear factor erythroid 2-related factor 2 (NFE2L2), and X-box binding protein 1 (XBP1). Together, these data suggest that DDIT3 might play a key role in MPP+-induced neuronal cell death.

Perturbation of ER protein processing pathway at 24 h after MPP+ treatment. (A) Perturbed (upregulated/downregulated) genes in the ER protein processing pathway (B) Propagation of perturbation in the ER protein ...

Among the pathways listed in Table 3, MAPK signaling pathway and ER protein processing pathway showed larger values oftA than other pathways at 12 and 24 h after MPP+ treatment. In addition, these two pathways were significantly perturbed at all observed time points, implicating an important role in neuronal death induced by MPP+. The number of significantly perturbed pathways at a level of 5% after FDR correction was 26, 19, 13, 24, and 51 for 1.5, 3, 9, 12, and 24 h respectively. The dramatic increase in the number of perturbed pathways at 24 h indicates that the cytotoxic effect by MPP+ might be propagated broadly into several cellular regulation networks. Most of the significantly perturbed pathways showed a tendency to have perturbed genes at or near the entry point of the pathway, thus differential expression of the ligand or receptor might have a large impact on the downstream pathway.

Validation of microarray data by real-time qPCR

To increase the fidelity of the microarray data, all gene expression values were measured twice with two independent samples for each time point. However, the quality of gene expression data obtained from microarrays can vary greatly according to the platform and procedures used because microarray experiments measure the expression of tens of thousands of genes simultaneously. For validation of the microarray data, real-time qPCR was carried out for 21 genes randomly selected from significantly impacted pathways shown in Table 3 at two time points, 3 and 24 h after MPP+ treatment, with TaqMan Gene Expression Assays (Table 1). All real-time qPCR experiments were performed in triplicate for each sample. Although there were some differences in the absolute value of the fold change between microarray and real-time qPCR results, the pattern of regulation showed good agreement. That is, five genes that were upregulated genes in microarray data at 3 h-BNDF,SMAD4, CREB5, FGF2 and HIF1A-also showed upregulation in real-time qPCR data at 3 h with a fold-change threshold of 1.5 and 2, respectively, for microarray and real-time qPCR data (Table 4). In addition, genes that were downregulated in microarray data, such as DDIT4 and EGR1, also showed downregulation in real-time qPCR experiment. Similarly, the regulation pattern at 24 h showed good agreement between microarray and real-time qPCR data, but the fold change values for DDIT3 and IFNB1 in real-time qPCR were very high compared with those for microarray data. This might be due to the high signal detection sensitivity and accuracy of real-time qPCR. The dynamic range for most commercial microarray platforms is usually 3–4 orders of magnitude, whereas TaqMan-based real-time PCR can achieve a dynamic range of 6–8 orders of magnitude (Wang et al., 2006). Therefore, the fold-changes obtained from real-time qPCR were considered more reliable than those from microarray in terms of detection sensitivity and accuracy. Nonetheless, the reliability of gene expression data obtained from microarray is comparable to that of real-time qPCR data in terms of the pattern of regulation.

DISCUSSION

Pathways perturbed in MPP+-treated human neuroblastoma SH-EP cells were identified using genome-wide gene expression data at five time points (1.5, 3, 9, 12, and 24 h) after MPP+ treatment. Two types of evidence, PNDE and PPERT, were used to identify pathways that were significantly perturbed by MPP+ treatment. Two pathways, MAPK signaling pathway and ER protein processing pathways, were highly ranked at each time point (Table 3). Perturbation of each these two pathways caused the common effect of upregulation of DDIT3 (Figs. 4B and 5B), which was confirmed by real-time qPCR. Tabuchi et al. (2014) demonstrated that treatment of the cells with the siRNA for DDIT3 significantly elevated the cell viability in sodium floride-elicited cell death accompanying ER stress with the induction of DDIT. It was also reported that the knockdown of DDIT3 in neuroblastoma cells significanlty reduced glucose deprivation-triggered cell death (Kögel et al., 2006). These evidences support that MPP+ induces the apoptosis with ER stress as described in previous reports (Ghribi et al., 2003; Selvaraj et al., 2012). Upregulation of DDIT3 was previously reported in MPP+-treated SH-SY5Y cells, which were generated by sub-cloning of the neuroblastoma SK-N-H cell line (Conn et al., 2002). The SH-EP cells used in this study were also originated from the SK-N-H cell line, but are morphologically distinct from SH-SY5Y cells; SH-SY5Y cells are N-type cell line expressing neurotransmitters and various neuronal cell surface markers whereas SH-EP cells are an epithelial substrate-adherent Schwannian (S-type) cell line that expresses proteins characteristic of Schwann cells and lacks neuronal markers (Do et al., 2011). The upregulation of DDIT3 by MPP+ treatment in both SH-SY5Y and SH-EP cells suggests that DDIT3 might play a key role in MPP+-induced neuronal cell death. Conn et al. (2002) suggested that upregulation of DDIT3 in MPP+-treated SH-SY5Y cells might be a cellular mechanism that is distinct from mitochondrial impairment or oxidative stress. Our results suggest that DDIT3 could be strongly upregulated by perturbation of two pathways-MAPK signaling pathway and ER protein processing pathway.

Perturbation of the MAPK signaling pathway was mainly driven by differential expression of ligand-encoding genes such as NFG, BDNF, IL1A, and TGFB2, and receptor genes such as FAS and CD14. These genes were upregulated at 24 h after MPP+ treatment and their perturbations were propagated to downstream genes in the MAPK signaling pathway (Fig. 5B). This suggests that perturbation of the pathway entry point might be more influential than perturbation of downstream genes in the pathway. Perturbation of the MAPK signaling pathway ultimately resulted in upregulation of DDIT3 and ATF4. Perturbation of the ER protein processing pathway involved ubiquitin ligase complex genes as well as ER associated degradation (ERAD)-related genes (Fig. 6A). Absence of the ERAD proteins that remove misfolded proteins from the ER might result in the accumulation of misfolded proteins and activation of the unfolded protein response (UPR) by ER stress. Indeed, UPR pathway genes such as eukaryotic translation initiation factor 2 subunit 1 (EIF2S1), NFE2L2, and ATF4 were upregulated after MPP+ treatment although there was no change in the expression of eukaryotic translation initiation factor 2-alpha kinase 1 (EIF2AK1) (Fig. 6A). DDIT3, which plays a role in apoptosis during ER stress (Marciniak et al., 2004; McCullough et al., 2001; Oyadomari and Mori, 2004), might be strongly upregulated as a result of perturbation of upstream genes including ATF4, NFE2L2, and XBP1.

In summary, the very high expression of DDIT3 in MPP+-treated human neuroblastoma SH-EP cells might result from the perturbation of the MAPK signaling pathway and ER protein processing pathway. Overexpression of DDIT3 results in a decrease in BCL2 protein and overexpression of BCL2 blocks DDIT3-induced apoptosis (Gotoh et al., 2004; Matsumoto et al., 1996). In this study, BCL2 was upregulated (fold-change = 1.9) at 24 h after MPP+ treatment, but this expression value was 11 times lower than that of DDIT3. Thus, the blocking effect of BCL2 on DDIT3-induced apoptosis might be weak. DDIT3 might also lead to glutathione depletion by overexpression of cation transport regulator homolog 1 (CHAC1), which is under the regulation of DDIT3 during the UPR of ER stress (Mungrue et al., 2009). The CHAC1 family of proteins function as γ-glutamylcyclotransferases acting specifically to degrade glutathione (Kumar et al., 2012), and overexpression of CHAC1 (fold change = 4.4) were observed at 24 h after MPP+ treatment. As overexpression of DDIT3 leads to translocation of BAX protein from the cytosol to the mitochondria (Gotoh et al., 2004), the DDIT-mediated death signal is finally transmitted to the mitochondria (Oyadomari and Mori, 2004), which function as integrators and amplifiers of the death pathway. This suggests that the toxicity signal of MPP+ is initiated by mitochondrial dysfunction caused by inhibition of complex I of the electron transport chain (ETC) and might feed back to the mitochondria via ER stress. This positive feedback mechanism could contribute to amplification of the death signal.

Article information

Mol. Cells.Sep 30, 2014; 37(9): 672-684.
Published online 2014-9-18.
Department of Biomolecular and Chemical Engineering, DongYang University, Yeongju 750-711, Korea
*Correspondence: jinhwando@dyu.ac.kr
Received June 23, 2014; Accepted August 11, 2014.

References

• Almeida, RD, Manadas, BJ, Melo, CV, Gomes, JR, Mendes, CS, Grãos, MM, Carvalho, RF, Carvalho, AP, and Duarte, CB (2005). Neuroprotection by BDNF against glutamate-induced apoptotic cell death is mediated by ERK and PI3-kinase pathways. Cell Death Differ. 12, 1329-1343.
• Bandiera, S, Rüberg, S, Girard, M, Cagnard, N, Hanein, S, Chrétien, D, Munnich, A, Lyonnet, S, and Henrion-Caude, A (2011). . PLos One. 6, e20746.
• Benjamini, Y, and Yekutieli, D (2001). The control of the false discovery rate in multiple testing under dependency. Ann Stat. 29, 1165-1188.
• Conn, KJ, Gao, WW, Ullman, MD, McKeon-O’Malley, C, Eisenhaur, PB, Fine, RE, and Wells, JM (2002). Specific upregulation of GADD153/CHOP in 1-methyl-4-phenyl-pyridinium-treated SH-SY5Y cells. J Neurosci Res. 68, 755-760.
• Conn, KJ, Ullman, MD, Larned, MJ, Eisenhauer, PB, Fine, RE, and Wells, JM (2003). cDNA microarray analysis of changes in gene expression associated with MPP+ toxicity in SH-SY5Y cells. Neurochem Res. 28, 1873-1881.
• Croft, D, Mundo, AF, Haw, R, Milacic, M, Weiser, J, Wu, G, Caudy, M, Garapati, P, Gillespie, M, and Kamdar, MR (2014). The Reactome pathway knowledgebase. Nucleic Acids Res. 42, D472-D477.
• Date, I, Yoshimoto, Y, Imaoka, T, Miyoshi, Y, Gohda, Y, Furuta, T, Asari, S, and Ohmoto, T (1993). Enhanced recovery of the nigrostriatal dopaminergic system in MPTP-treated mice following intrastriatal injection of basic fibroblast growth factor in relation to aging. Brain Res. 621, 150-154.
• Do, JH, Kim, IS, Lee, JD, and Choi, D-K (2011). Comparison of genomic profiles in human neuroblastic SH-SY5Y and substrate-adherent SH-EP cells using array comparative genomic hybridization. BioChip J. 5, 165-174.
• Doniger, SW, Salomonis, N, Dahlquist, KD, Vranizan, K, Lawlor, SC, and Conklin, RR (2003). MAPPfinder: using Gene Ontology and GenMAPP to create a global gene expression profile from microarray data. Genome Biol. 4, R7.
• Drăghici, S, Khatri, P, Martins, RP, Ostermeier, GC, and Krawetz, SA (2003). Global functional profiling of gene expression. Genomics. 81, 98-104.
• Gálvez-Jiménez, N (2007). Parkinson’s disease. Neurobiology of Disease, Gilman, S, ed. (USA:Elsevier Academic Press), pp. 55.
• Ghribi, O, Herman, MM, Pramoonjago, P, and Savory, J (2003). MPP+ induces the endoplasmic reticulum stress response in rabbit brain involving activation of the ATF-6 and NF-kappaB signaling pathways. J Neuropathol Exp Neurol. 62, 1144-1153.
• Goeman, JJ, van de Geer, SA, de Kort, F, and van Houwelingen, HC (2004). A global test for groups of genes: testing association with a clinical outcome. Bioinformatics. 20, 93-99.
• Gotoh, T, Takeda, K, Oyadomari, S, and Mori, M (2004). hsp70- DnaJ chaperone pair prevents nitric oxide-, CHOP-induced apoptosis by inhibiting translocation of Bax to mitochondria. Cell Death Differ. 11, 390-402.
• Hochberg, Y (1988). A sharper Bonferroni procedure for multiple tests of significance. Biometrika. 75, 800-802.
• Hoevel, T, Macek, R, Swisshelm, K, and Kubbies, M (2004). Reexpression of the TJ protein CLDN1 induces apoptosis in breast tumor spheroids. Int J Cancer. 108, 374-383.
• Kanehisa, M, Goto, S, Sato, Y, Furumichi, M, and Tanabe, M (2012). KEGG for integration and interpretation of large-scale molecular data sets. Nucleic Acids Res. 40, D109-114.
• Kim, IS, Choi, D-K, and Do, JH (2013). Genome-wide temporal responses of human neuroblastoma SH-SY5Y cells to MPP+ neurotoxicity. BioChip J. 7, 247-257.
• Kögel, D, Svensson, B, Copanaki, E, Anguissola, S, Bonner, C, Thurow, N, Gudorf, D, Hetschko, H, Müller, T, and Peters, M (2006). Induction of transcription factor CEBP homology protein mediates hypoglycaemia-induced necrotic cell death in human neuroblastoma cells. J Neurichem. 99, 952-964.
• Kumar, A, Tikoo, S, Maity, S, Sengupta, S, Sengupta, S, Kaur, A, and Bachhawat, AK (2012). Mammalian proapoptotic factor ChaC1 and its homologues function as γ-glutamylcyclotransferases acting specifically on glutathione. EMBO Rep. 13, 1095-1101.
• Langston, JW, Ballard, P, and Irwin, I (1983). Chronic parkinsonism in human due to a product of meperidine-analog synthesis. Science. 219, 979-980.
• Livak, KJ, and Schmittgen, TD (2001). Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta DeltaC(T)) method. Methods. 25, 402-408.
• Lotharius, J, and O’Malley, KL (2000). The parkinsonism-inducing drug 1-methyl-4-phenylpyridinium triggers intracellular dopamine oxidation. A novel mechanism of toxicity. J Biol Chem. 275, 38581-38588.
• Marciniak, SJ, Yun, CY, Oyadomari, S, Novoa, I, Zhang, Y, Jungreis, R, Nagata, K, Harding, HP, and Ron, D (2004). CHOP induces death by promoting protein synthesis and oxidation in the stressed endoplasmic reticulum. Genes Dev. 18, 3066-3077.
• Matsumoto, M, Minami, M, Takeda, K, Sakao, Y, and Akira, S (1996). Ectopic expression of CHOP (GADD153) induces apoptosis in M1 myeloblastic leukemia cells. FEBS Lett. 395, 143-147.
• McCullough, KD, Martindale, JL, Klotz, LO, Aw, TY, and Holbrook, NJ (2001). Gadd153 sensitizes cells to endoplasmic reticulum stress by down-regulating Bcl2 and perturbing the cellular redox state. Mol Cell Biol. 21, 1249-1259.
• Mungrue, IN, Pagnon, J, Kohannim, O, Gargalovic, PS, and Lusis, AJ (2009). CHAC1/MGC4504 is a novel proapoptotic component of the unfolded protein response, downstream of the ATF4-ATF3-CHOP cascade. J Immunol. 182, 466-476.
• Nakamura, K, Bindokas, VP, Marks, JD, Wright, DA, Frim, DM, Miller, RJ, and Kang, UJ (2000). The selective toxicity of 1-methyl-4-phenylpyridinium to dopaminergic neurons: the role of mitochondria complex I and reactive oxygen species revisited. Mol Pharmacol. 58, 271-278.
• Nanjo, F, Goto, K, Seto, R, Suzuki, M, Sakai, M, and Hara, Y (1996). Scavenging effects of tea catechins and their derivatives on 1,1-diphenyl-2-picrylhydrazyl radical. Free Radic Biol Med. 21, 895-902.
• Nicotra, A, and Parvez, S (2002). Apoptotic molecules and MPTP-induced cell dealth. Neurotoxicol. Teratol. 24, 599-605.
• Oyadomari, S, and Mori, M (2004). Roles of CHOP/GADD153 in endoplasmic reticulum stress. Cell Death Differ. 11, 381-389.
• Pandey, R, Guru, RK, and Mount, DW (2004). Pathway miner: Extracting gene association networks from molecular pathways for predicting the biological significance of gene expression microarray data. Bioinformatics. 20, 2156-2158.
• Pepe, D, and Grassi, M (2014). Investigating perturbed pathway modules from gene expression data via structural equation models. BMC Bioinformatics. 15, 132.
• Selvaraj, S, Sun, Y, Watt, JA, Wang, S, Lei, S, Birnbaumer, L, and Singh, BB (2012). Neurotoxin-induced ER stress in mouse dopaminergic neurons involves downregulation of TRPC1 and inhibition of AKT/mTOR signaling. J Clin Invest. 122, 1354-1367.
• Tabuchi, Y, Yunoki, T, Hoshi, N, Suzuki, N, and Kondo, T (2014). Genes and gene networks involved in sodium fluoride-elicited cell death accompanying endoplasmic reticulum stress in oral epithelial cells. Int J Mol Sci. 15, 8959-8978.
• Tarca, AL, Draghici, S, Khatri, P, Hassan, SS, Mittal, P, Kim, JS, Kim, CJ, Kusanovic, JP, and Romero, R (2009). A novel signaling pathway impact analysis. Bioinformatics. 25, 75-82.
• Tusher, VG, Tibshirani, R, and Chu, G (2001). Significance analysis of microarrays applied to the ionizing radiation response. Proc Natl Acad Sci USA. 98, 5116-5121.
• Wang, Y, Barbacioru, C, Hyland, F, Xiao, W, Hunkapiller, KL, Blake, J, Chan, F, Gonzalez, C, Zhang, L, and Samaha, RR (2006). Large scale real-time PCR validation on gene expression measurements from two commercial long-oligonucleotide microarrays. BMC Genomics. 7, 59.

Figure 1

The viability of MPP+ treated human neuroblastoma SH-EP cells. (A) Cell viability after cultivation for 48 h with various MPP+ concentrations. (B) Time-dependent cell viability after treatment with 1.25 mM MPP+.

Figure 2

Upregulated/downregulated DEGs at 1.5, 3, 9, 12, and 24 h after MPP+ treatment. The red and green bars represent the number of upregulated and downregulated genes, respectively.

Figure 3

Two-dimensional plots showing two types of evidence, PNDE and PPERT, for each time point (1.5, 3, 9, 12, and 24 h) after MPP+ treatment. The numbers beside the dots on the plane represent the KEGG pathway ID. Pathways above the oblique red line are significant at 5% after Bonferroni correction whereas those above the blue line are significant at 5% after false discovery rate (FDR) correction.

Figure 4

Perturbation of the MAPK signaling pathway at 9 h after MPP+ treatment. (A) Perturbed (upregulated/downregulated) genes in the MAPK signaling pathway. (B) Propagation of perturbation in the MAPK signaling pathway and total perturbation of each gene. The red arrows represent propagation of perturbation via the signaling cascade.

Figure 5

Perturbation of MAPK signaling pathway at 24 h after MPP+ treatment. (A) Perturbed (upregulated/downregulated) genes in the MAPK signaling pathway. (B) Propagation of perturbation in MAPK signaling pathway and total perturbation of each gene. The red arrows represent propagation of perturbation via the signaling cascade.

Figure 6

Perturbation of ER protein processing pathway at 24 h after MPP+ treatment. (A) Perturbed (upregulated/downregulated) genes in the ER protein processing pathway (B) Propagation of perturbation in the ER protein processing pathway and total perturbation of each gene. The red arrows represent propagation of perturbation via the signaling cascade.

Table 1.

Target genes for real-timeqPCR and their corresponding TaqMan Gene Expression Assay IDs

Gene Assay ID* Description
ATF4 Hs00909569_g1 Activating transcription factor 4 (tax-responsive enhancer element B67) (ATF4), transcript variant 1, mRNA
BDNF Hs02718934_s1 Brain-derived neurotrophic factor (BDNF), transcript variant 3, mRNA
CALM1 Hs00300085_s1 Calmodulin 1 (phosphorylase kinase, delta) (CALM1), mRNA
CREB5 Hs00329596_s1 cAMP responsive element binding protein 5 (CREB5), transcript variant 1, mRNA
DDIT3 Hs00358796_g1 DNA-damage-inducible transcript 3 (DDIT3), mRNA
DDIT4 Hs01111686_g1 DNA-damage-inducible transcript 4 (DDIT4), mRNA
EGR1 Hs00152928_m1 Early growth response 1 (EGR1), mRNA
EPAS1 Hs01026142_m1 Endothelial PAS domain protein 1 (EPAS1), mRNA
FGF2 Hs00960934_m1 Fibroblast growth factor 2 (basic) (FGF2), mRNA
HIF1A Hs00936370_m1 Hypoxia-inducible factor 1, alpha subunit (basic helix-loop-helix transcription factor) (HIF1A), transcript variant 2, mRNA
IFNB1 Hs01077958_s1 Interferon, beta 1, fibroblast (IFNB1), mRNA
MAOA Hs02383327_s1 Monoamine oxidase A (MAOA), nuclear gene encoding mitochondrial protein, mRNA
NFATC4 Hs01113412_m1 Nuclear factor of activated T-cells, cytoplasmic, calcineurin-dependent 4 (NFATC4), mRNA
NFKBIE Hs00914563_g1 Nuclear factor of kappa light polypeptide gene enhancer in B-cells inhibitor, epsilon (NFKBIE), mRNA
NGF Hs00171458_m1 Nerve growth factor (beta polypeptide) (NGF), mRNA
NTN4 Hs01003502_m1 Netrin 4 (NTN4), mRNA
PIM1 Hs01065498_m1 pim-1 oncogene (PIM1), mRNA
PLXNB1 Hs00963524_m1 plexin B1 (PLXNB1), mRNA
PRKCD Hs00178914_m1 Protein kinase C, delta (PRKCD), transcript variant 1, mRNA
*TaqMan Gene Expression Assays (Applied Biosystems)

Table 2.

Commonly up-regulated genes at five time points (1.5, 3, 9, 12, and 24 h) after MPP+ treatment

Gene symbol Fold changes at 1.5, 3, 9, 12, 24 h Description
SNORA25 2.7, 2.7, 3.1, 3.1, 2.1 Small nucleolar RNA, H/ACA box 25 (SNORA25), small nucleolar RNA
BIRC2 2.2, 2.3, 2.3, 2.8, 2.0 Baculoviral IAP repeat-containing 2 (BIRC2), mRNA
ERRFI1 2.3, 2.3, 3.5, 3.5, 1.8 ERBB receptor feedback inhibitor 1 (ERRFI1), mRNA
FERMT2 2.4, 2.8, 2.2, 2.2, 1.5 Fermitin family homolog 2 (Drosophila) (FERMT2), mRNA
MAK16 2.1, 2.5, 3.9, 3.9, 3.5 MAK16 homolog (S. cerevisiae) (MAK16), mRNA
MTMR9 2.0, 1.7, 1.9, 1.9, 2.1 Myotubularin related protein 9 (MTMR9), mRNA
ANKRD50 1.9, 2.0, 2.5, 3.0, 3.2 Ankyrin repeat domain 50 (ANKRD50), mRNA
PTPN11 1.7, 1.8, 1.6, 1.6, 1.6 Protein tyrosine phosphatase, non-receptor type 11 (PTPN11), mRNA
GRPEL2 1.9, 2.4, 4.8, 5.9, 4.1 GrpE-like 2, mitochondrial (E. coli) (GRPEL2), nuclear gene encoding mitochondrial protein, mRNA
BDNF 1.9, 1.8, 2.7, 3.4, 3.2 Brain-derived neurotrophic factor (BDNF), transcript variant 3, mRNA
HIVEP1 1.7, 2.5, 2.2, 2.2, 2.0 Human immunodeficiency virus type I enhancer binding protein 1 (HIVEP1), mRNA
ZNF295 2.0, 1.9, 2.2, 2.3, 2.1 Zinc finger protein 295 (ZNF295), transcript variant 1, mRNA
OBFC2A 1.8, 1.9, 3.4, 3.8, 1.7 Oligonucleotide/oligosaccharide-binding fold containing 2A (OBFC2A), mRNA
DDX21 2.0, 2.2, 2.4, 2.4, 1.9 DEAD (Asp-Glu-Ala-Asp) box polypeptide 21 (DDX21), mRNA
RAB23 1.8, 1.9, 2.0, 2.3, 1.8 RAB23, member RAS oncogene family (RAB23), transcript variant 2, mRNA
IFRD1 1.6, 1.8, 3.2, 3.5, 3.4 Interferon-related developmental regulator 1 (IFRD1), transcript variant 1, mRNA
FGF2 2.1, 2.6, 4.9, 8.2, 7.5 Fibroblast growth factor 2 (basic) (FGF2), mRNA
ZNHIT6 1.6, 1.8, 2.1, 2.3, 2.0 Zinc finger, HIT type 6 (ZNHIT6), mRNA
BTAF1 1.8, 2.0, 1.9, 2.0, 1.7 BTAF1 RNA polymerase II, B-TFIID transcription factor-associated, 170kDa (Mot1 homolog, S. cerevisiae) (BTAF1), mRNA
FNDC3B 1.5, 1.9, 1.9, 2.1, 1.7 Fibronectin type III domain containing 3B (FNDC3B), transcript variant 2, mRNA
DLC1 3.7, 3.0, 3.5, 3.4, 1.9 Deleted in liver cancer 1 (DLC1), transcript variant 3, mRNA
LOC100131336 1.7, 1.8, 2.3, 2.5, 3.1 PREDICTED: misc_RNA (LOC100131336), miscRNA
MKI67IP 1.6, 1.8, 2.4, 2.5, 1.7 MKI67 (FHA domain) interacting nucleolarphosphoprotein (MKI67IP), mRNA
NAV3 1.6, 2.7, 3.0, 4.1, 2.1 Neuron navigator 3 (NAV3), mRNA
NEDD4 1.6, 1.7, 2.0, 2.2, 1.6 Neural precursor cell expressed, developmentally down-regulated 4 (NEDD4), transcript variant 1, mRNA
OTUD4 1.7, 1.9, 1.7, 1.8, 2.1 OTU domain containing 4 (OTUD4), transcript variant 1, mRNA
SPTY2D1 1.6, 1.7, 1.7, 1.6, 1.8 SPT2, Suppressor of Ty, domain containing 1 (S. cerevisiae) (SPTY2D1), mRNA
ZNF286C 1.6, 1.7, 2.6, 2.5, 3.0 Zinc finger 286C pseudogene (ZNF286C), non-coding RNA
SELI 1.6, 1.7, 2.0, 2.0, 2.4 Selenoprotein I (SELI), mRNA
CLDN12 1.7, 1.8, 2.3, 2.4, 2.1 Claudin 12 (CLDN12), mRNA
LOC203547 1.6, 1.5, 1.6, 1.6, 1.5 Hypothetical protein LOC203547 (LOC203547), mRNA
KPNA4 1.6, 1.8, 2.1, 2.2, 2.7 Karyopherin alpha 4 (importin alpha 3) (KPNA4), mRNA
C5orf5 1.7, 1.8, 2.8, 2.9, 3.0 Chromosome 5 open reading frame 5 (C5orf5), mRNA
ZNF23 1.6, 1.8, 3.0, 3.2, 2.2 Zinc finger protein 23 (KOX 16) (ZNF23), mRNA
SMG1 1.6, 1.6, 1.9, 1.9, 2.9 PI-3-kinase-related kinase SMG-1 (SMG1), mRNA
POGZ 1.6, 1.6, 1.7, 1.9, 1.5 Pogo transposable element with ZNF domain (POGZ), transcript variant 3, mRNA
CLK1 1.9, 2.0, 2.5, 2.4, 4.3 CDC-like kinase 1 (CLK1), mRNA
PDCD1LG2 1.6, 1.7, 2.6, 2.6, 2.3 Programed cell death 1 ligand 2 (PDCD1LG2), mRNA
PRPF38B 1.6, 1.7, 2.3, 2.8, 3.9 PRP38 pre-mRNA processing factor 38 (yeast) domain containing B (PRPF38B), mRNA
GREM1 2.1, 3.2, 3.1, 3.3, 3.1 Gremlin 1, cysteine knot superfamily, homolog (Xenopuslaevis) (GREM1), mRNA
RSBN1 1.6, 1.8, 2.1, 2.2, 2.6 Round spermatid basic protein 1 (RSBN1), mRNA
FAM175B 1.6, 1.8, 2.0, 2.0, 2.2 Family with sequence similarity 175, member B (FAM175B), mRNA
C14orf138 1.5, 1.6, 1.9, 1.8, 1.8 Chromosome 14 open reading frame 138 (C14orf138), transcript variant 2, mRNA
KLHDC5 1.5, 1.7, 1.7, 1.5, 1.6 Kelch domain containing 5 (KLHDC5), mRNA
EPRS 1.6, 1.5, 1.7, 1.9, 1.9 Glutamyl-prolyl-tRNAsynthetase (EPRS), mRNA
C1orf71 1.7, 2.4, 3.4, 4.6, 3.1 Chromosome 1 open reading frame 71 (C1orf71), mRNA
HNRPDL 1.5, 1.6, 1.9, 2.1, 2.4 Heterogeneous nuclear ribonucleoprotein D-like (HNRPDL), transcript variant 3, transcribed RNA
C1orf124 1.5, 1.6, 1.9, 2.0, 1.6 Chromosome 1 open reading frame 124 (C1orf124), transcript variant 1, mRNA
RCAN1 2.7, 1.8, 3.4, 4.4, 4.9 Regulator of calcineurin 1 (RCAN1), transcript variant 3, mRNA
DGKD 1.7, 1.7, 1.6, 1.8, 1.7 Diacylglycerol kinase, delta 130kDa (DGKD), transcript variant 2, mRNA
MIR1974 11.0, 28.7, 81.1, 61.4, 45.7 microRNA 1974 (MIR1974), microRNA
CLDN1 1.7, 1.9, 6.6, 14.0, 33.4 Claudin 1 (CLDN1), mRNA
ZNF26 1.6, 1.8, 2.6, 2.8, 2.2 Zinc finger protein 26 (ZNF26), mRNA
SCHIP1 1.9, 2.6, 2.5, 2.7, 1.7 Schwannomin interacting protein 1 (SCHIP1), mRNA

Table 3.

Top 10 ranked KEGG pathways by PG for each time point (1.5, 3, 9, 12, and 24 h) after MPP+ treatment

Time (h) KEGG pathway KEGG ID tA PNDE PPERT PG PG,FDR*
1.5 TNF signaling pathway 04668 −0.2272 0.0000 0.8210 0.0000 0.0024
NOD-like receptor signaling pathway 04621 −1.1113 0.0000 0.2443 0.0001 0.0042
p53 signaling pathway 04115 0.0198 0.0000 0.9838 0.0002 0.0069
Protein processing in endoplasmic reticulum 04141 0.6061 0.0000 0.5293 0.0003 0.0090
MAPK signaling pathway 04010 0.8414 0.0001 0.4008 0.0004 0.0093
NF-kappa B signaling pathway 04064 −1.2155 0.0002 0.2270 0.0004 0.0093
Cell cycle 04110 −0.8842 0.0001 0.3758 0.0005 0.0099
Pathways in cancer 05200 0.4030 0.0001 0.6860 0.0006 0.0099
Proteoglycans in cancer 05205 0.6451 0.0001 0.5214 0.0007 0.0099
Small cell lung cancer 05222 −0.3325 0.0001 0.7427 0.0008 0.0099
3 TNF signaling pathway 04668 −1.1952 0.0000 0.2133 0.0000 0.0003
Pathways in cancer 05200 1.0934 0.0000 0.2735 0.0000 0.0003
MAPK signaling pathway 04010 1.2421 0.0000 0.1937 0.0000 0.0003
Proteoglycans in cancer 05205 3.0568 0.0003 0.0059 0.0000 0.0008
NOD-like receptor signaling pathway 04621 −1.1165 0.0000 0.2224 0.0000 0.0008
Hepatitis B 05161 2.2398 0.0002 0.0239 0.0001 0.0016
Small cell lung cancer 05222 −0.9612 0.0001 0.3417 0.0003 0.0055
PI3K-Akt signaling pathway 04151 1.1016 0.0001 0.2859 0.0003 0.0055
Colorectal cancer 05210 4.8055 0.0446 0.0013 0.0006 0.0092
RNA degradation 03018 0.0734 0.0001 0.8701 0.0007 0.0098
9 MAPK signaling pathway 04010 1.1397 0.0000 0.2430 0.0000 0.0011
Neurotrophin signaling pathway 04722 −0.1552 0.0000 0.8715 0.0001 0.0049
Phosphatidylinositol signaling system 04070 2.7042 0.0014 0.0206 0.0003 0.0133
Colorectal cancer 05210 2.8492 0.0032 0.0103 0.0004 0.0133
Protein processing in endoplasmic reticulum 04141 5.5431 0.1303 0.0007 0.0009 0.0269
Cell cycle 04110 −1.0684 0.0005 0.2623 0.0012 0.0284
Dorso-ventral axis formation 04320 −0.9378 0.0003 0.5636 0.0016 0.0328
Axon guidance 04360 1.0965 0.0008 0.2446 0.0018 0.0329
TNF signaling pathway 04668 −1.2830 0.0015 0.1618 0.0023 0.0355
PPAR signaling pathway 03320 1.0002 0.0020 0.1502 0.0028 0.0355
12 MAPK signaling pathway 04010 4.2762 0.0000 0.0004 0.0000 0.0000
Cell cycle 04110 −1.0880 0.0000 0.2481 0.0000 0.0008
Protein processing in endoplasmic reticulum 04141 4.9870 0.0039 0.0012 0.0001 0.0030
PI3K-Akt signaling pathway 04151 2.2328 0.0006 0.0328 0.0002 0.0081
Hippo signaling pathway 04390 −1.1356 0.0001 0.2330 0.0003 0.0088
Neurotrophin signaling pathway 04722 −0.6262 0.0001 0.5220 0.0004 0.0093
Colorectal cancer 05210 3.3836 0.0095 0.0047 0.0005 0.0093
HTLV-I infection 05166 −1.0903 0.0002 0.2638 0.0005 0.0093
Chronic myeloid leukemia 05220 −1.8518 0.0010 0.0533 0.0006 0.0094
Alcoholism 05034 0.5809 0.0001 0.5562 0.0008 0.0116
24 Protein processing in endoplasmic reticulum 04141 5.4592 0.0000 0.0002 0.0000 0.0000
MAPK signaling pathway 04010 3.2812 0.0000 0.0039 0.0000 0.0001
Pathogenic Escherichia coli infection 05130 1.6970 0.0000 0.0712 0.0000 0.0001
Measles 05162 −0.8260 0.0000 0.3883 0.0000 0.0002
HTLV-I infection 05166 −0.8087 0.0000 0.4108 0.0000 0.0002
Pathways in cancer 05200 0.2457 0.0000 0.8038 0.0000 0.0003
Cell cycle 04110 0.0753 0.0000 0.9428 0.0000 0.0003
Neurotrophin signaling pathway 04722 −0.8145 0.0000 0.3944 0.0000 0.0004
RNA transport 03013 −0.8182 0.0000 0.4355 0.0000 0.0004
Vibrio cholerae infection 05110 −1.7378 0.0000 0.0667 0.0000 0.0006
*p-value after FDR correction

Table 4.

Average fold changes for 21 genes obtained from microarray and real-time qPCR at 3 and 24 h after MPP+ treatment

Gene Microarray
qPCR
3 h 24 h 3 h 24 h
BDNF 1.836668 3.399638 3.084422 5.443998
MAOA 0.80517 0.617155 1.135504 1.024083
NFKBIE 1.014369 2.670331 1.406393 5.589296
NTN4 1.131652 1.600827 1.643001 9.767572
PLXNB1 0.733329 0.434851 1.23799 0.63581
CALM1 0.876405 0.466762 1.411602 0.671131
CREB5 1.62003 2.31009 4.179509 11.5301
EPAS1 0.847435 0.712773 1.597704 1.756455
FGF2 2.624728 2.464603 3.001465 14.27415
NGF 1.122245 2.026697 1.445598 6.5251
PRKCD 0.685482 0.496559 0.840508 0.578611
ATF4 1.281656 2.86627 1.857463 5.157481
DDIT3 1.112806 16.66675 1.690754 89.51143
DDIT4 0.369263 2.743092 0.460413 7.120606
EGR1 0.199824 0.306295 0.106826 0.252088
HIF1A 1.647623 2.612308 3.129647 5.819923
IFNB1 0.952799 8.165385 0.844694 49.18001
NFATC4 0.733197 0.596169 0.846941 0.829895
PIM1 0.605094 0.755351 0.773425 1.344746