MPTP-induced PD motor deficits
We evaluated the exercise capacity of MPTP-induced mice by pole test and rotarod test. Fig. 1A, B showed the time required for the movement of the mice in the pole test and rotarod test. In the pole test, MPTP-induced mice were used longer time than the control group (Fig. 1a). Compared with the control group, the MPTP-induced mice had poor coordination and the rotarod test time was significantly shortened (Fig. 1b). This indicates that MPTP-induced PD mice were successfully constructed.
Summary of RNA-Seq data sets
To understand the transcriptome information of different brain regions in PD mice, 24 libraries were constructed in the CC, HP, ST, and CB regions. Subsequently, the libraries were sequenced on the Illumina HiSeq 2500 sequencing platform. The raw reads of each sample ranged from 27.4 to 63.2 million (Additional file 1: Table S1). After quality control and filtering, clean reads were obtained. We compared the clean reads to the reference genome with Hisat2. The reads for each gene were calculated using HTSeq-count. Principal component analysis (PCA) was performed based on the DESeq2 method to analyze whether the same brain region samples were together. The results showed that the four brain regions of the control group were divided into different clusters (Fig. 2a). Compared with the CC and HP regions, the transcriptional characteristics of the CB region were significantly different. In addition, the transcriptome characteristics of the ST were significantly different from those of the CC and HP regions. However, the transcriptome characteristics of the CC and HP regions were not significantly different. The four brain regions of the model group were divided into different clusters (Fig. 2b). The difference between the CB region and the other three brain regions was the most obvious.
Differentially expression gene in different brain regions
The transcriptome characteristics of the model group and control group mice were analyzed. DEG analysis explained the difference between the transcription in the model group and control group. According to the values of｜fold change｜ > 2 and P < 0.05, DEGs were obtained using the DESeq2 R package.
Cluster analysis was performed using the correlation distance and hierarchical algorithm to show gene expression patterns in different brain regions of the control group and model group (Fig. 3a and b). In the same pathological state, the gene expression in the CB was significantly different from that in the CC, HP, and ST regions.
To analyze the effects of gene expression in different brain regions on PD, the DEGs in the model and control groups in the same brain region were compared. As shown in Fig. 3c, 110 DEGs were obtained in the CC region, of which 64 genes were upregulated and 46 genes were downregulated. A total of 179 DEGs were obtained in the HP region, of which 99 genes were upregulated and 80 genes were downregulated. A total of 521 DEGs were obtained in the ST region, of which 485 genes were upregulated and 36 genes were downregulated. There were 96 DEGs in the CB region, of which 67 were upregulated and 29 were downregulated. Fig. 3d shows the overlap of DEGs in the four brain regions. Of the four brain regions, only two genes (Gxylt1 and C920006O11 Rik) were co-expressed.
Differentially expressed genes in HP and ST regions of Parkinson mice
To further understand the transcript information in the HP and ST regions of mice with PD, RNA sequencing analysis was undertaken. Based on the results, unbiased hierarchical clustering analysis was performed for all genes in the HP region. The data of all genes were visualized using a heat map (Fig. 4a). For data quality control, we performed PCA (Fig. 4b). The results showed that the gene expression in the HP region of the model group was different from that in the HP region of the control group. A total of 179 DEGs were found in the HP region in the model group. Additional file 1: Table S2 lists the top 10 upregulated genes and top 10 downregulated genes with the most significant differences in the HP region. In addition, hierarchical clustering analysis of the DEGs showed that there were differences in gene expression in the ST region of the PD model mice (Fig. 4c). Additional file 1: Table S3 lists the top 10 upregulated genes and top 10 downregulated genes with the most significant differences in the ST region.
Compared with the control group, there were 110 (64 upregulated, 46 downregulated) and 96 (67 upregulated, 29 downregulated) DEGs in the CC and CB regions of the model group mice, respectively. Additional file 1: Table S4 and Table S5 listed the top 10 up-regulated genes and the top 10 down-regulated genes in CC and CB regions, respectively. A heat map of the top 100 genes differentially expressed between the control group and model group in the CB region is provided in Additional file 2: Figure S1.
Gene ontology enrichment analysis
Gene ontology enrichment analysis of the DEGs in HP and ST regions
The biological functions of DEGs were identified by GO analysis. The GO analysis database can classify DEGs, including “biological processes,” “cellular components,” and “molecular functions.” The GO term corrected by P < 0.05 was defined as significant enrichment of DEGs. In Fig. 5, the significant enrichment of GO terms in the three main categories (cellular component, biological process, and molecular function) is shown. In the HP region, the main five subcategories of biological process were “DNA methylation or demethylation,” “posttranscriptional regulation of gene expression,” “rhythmic process,” “protein localization to organelle,” and “regulation of DNA metabolic process,” and the two subcategories for molecular function were “guanyl-nucleotide exchange factor activity” and “protein serine/threonine kinase activity” (Fig. 5a). In the ST region, the main five subcategories of cell component were “post synapse,” “dendritic tree,” “presynapse,” “perikaryon”, and “endoplasmic reticulum subcompartment.” The five subcategories for biological processes were “positive regulation of nervous system development,” “synapse organization,” “regulation of transmembrane transport,” “regulation of exocytosis,” and “regulation of intracellular transport”; and five subcategories for molecular function, which are “protein domain specific binding,” “protein kinase binding,” “protein heterodimerization activity,” “GABA receptor binding,” and “guanyl-nucleotide binding” (Fig. 5d). Fig. 5b, 5c and 5e, 5f show the top 20 clusters of significantly enriched terms in the HP and ST regions, respectively.
Gene ontology enrichment analysis of the DEGs in CC and CB regions
GO enrichment analysis was performed on DEGs in the CC and CB regions to determine the effects of PD on these two brain regions. The results showed that 10 GO terms were enriched in the CC region (Fig. 6a), in which the cell components were significantly enriched in two GO terms (“axon” and “filopodium”), the biological process was significantly enriched in six GO terms (“negative regulation of protein kinase activity,” “regulation of neuron death,” “negative regulation of cell cycle,” “ribosome biogenesis,” “chromatin organization,” and “organophosphate biosynthetic process”), and the molecular function was significantly enriched in two GO terms (“myosin binding” and “lipase activity”). A total of 14 GO terms were enriched in the CB region (Fig. 6b), in which the cellular components were significantly enriched in two GO terms, the biological process was significantly enriched in nine GO terms, and the molecular function was significantly enriched in three GO terms. In the CC region, the top five significantly enriched GO terms were “axon,” “associative learning,” “myosin binding,” “negative regulation of protein kinase activity,” and “regulation of neuron death” (Fig. 6c). In the CB region, the top five significantly enriched GO terms were “negative regulation of transmembrane transport,” “translation factor activity, RNA binding,” “histone lysine methylation,” “response to glucose,” and “dioxygenase activity” (Fig. 6d).
KEGG enrichment analysis
KEGG pathway analysis of DEGs in the HP and ST regions of the model group helped identify the biological pathways related to DEGs. The P < 0.05 pathway was defined as the DEG enrichment pathway. A total of 179 and 521 DEGs in the HP and ST regions were enriched in 122 and 213 KEGG pathways, respectively. Table 1 shows the pathways and related genes that were significantly enriched in the HP region. In the ST region, the top 10 significantly enriched signaling pathways and related genes are shown in Table 2. There were three signaling pathways (“dopaminergic synapse,” “glutamatergic synapse,” and “adrenergic signaling in cardiomyocytes”), which was consistent with the study by Fu et al. (2019). However, only a few KEGG pathways were enriched in the CC and CB regions. Only the “focal adhesion” and “regulation of actin cytoskeleton” pathways were significantly enriched in the CC region (Table 3), and only the “adrenergic signaling in cardiomyocytes” and “RNA transport” pathways were significantly enriched in the CB region (Table 3).
There were 427 proteins in the PPI network and the nodes represented proteins (Fig. 7a). The most interacting proteins in the PPI network were LRRK2, DA receptor D2 (DRD2), protein kinase A catalytic subunits (PRKACA), PPP2R5C, insulin-like growth factor 1 (IGF-1), PIK3R1, guanine nucleotide binding protein alpha inhibiting 1 (GNAI1), and guanine nucleotide binding protein alpha inhibiting 3 (GNAI3). A total of 15 proteins that may be associated with PD were extracted from the PPI network and the protein network interactions were reconstructed (Fig. 7b).