High Throughput Sequencing Identifies MicroRNAs Mediating α-Synuclein Toxicity by Targeting Neuroactive-Ligand Receptor Interaction Pathway in Early Stage of Drosophila Parkinson's Disease Model

Parkinson’s disease (PD) is a prevalent neurodegenerative disorder with pathological features including death of dopaminergic neurons in the substantia nigra and intraneuronal accumulations of Lewy bodies. As the main component of Lewy bodies, α-synuclein is implicated in PD pathogenesis by aggregation into insoluble filaments. However, the detailed mechanisms underlying α-synuclein induced neurotoxicity in PD are still elusive. MicroRNAs are ~20nt small RNA molecules that fine-tune gene expression at posttranscriptional level. A plethora of miRNAs have been found to be dysregulated in the brain and blood cells of PD patients. Nevertheless, the detailed mechanisms and their in vivo functions in PD still need further investigation. By using Drosophila PD model expressing α-synuclein A30P, we examined brain miRNA expression with high-throughput small RNA sequencing technology. We found that five miRNAs (dme-miR-133-3p, dme-miR-137-3p, dme-miR-13b-3p, dme-miR-932-5p, dme-miR-1008-5p) were upregulated in PD flies. Among them, miR-13b, miR-133, miR-137 are brain enriched and highly conserved from Drosophila to humans. KEGG pathway analysis using DIANA miR-Path demonstrated that neuroactive-ligand receptor interaction pathway was most likely affected by these miRNAs. Interestingly, miR-137 was predicted to regulate most of the identified targets in this pathway, including dopamine receptor (DopR, D2R), γ-aminobutyric acid (GABA) receptor (GABA-B-R1, GABA-B-R3) and N-methyl-D-aspartate (NMDA) receptor (Nmdar2). The validation experiments showed that the expression of miR-137 and its targets was negatively correlated in PD flies. Further experiments using luciferase reporter assay confirmed that miR-137 could act on specific sites in 3’ UTR region of D2R, Nmdar2 and GABA-B-R3, which downregulated significantly in PD flies. Collectively, our findings indicate that α-synuclein could induce the dysregulation of miRNAs, which target neuroactive ligand-receptor interaction pathway in vivo. We believe it will help us further understand the contribution of miRNAs to α-synuclein neurotoxicity and provide new insights into the pathogenesis driving PD.


Introduction
Parkinson's disease (PD) is the second most prevalent neurodegenerative disorder affecting the elderly population [1]. Its predominant pathological features are death of dopaminergic (DA) neurons in the substantia nigra pars compacta and intraneuronal accumulations of Lewy bodies [2]. As the main component of Lewy bodies, α-synuclein contributes to PD by aggregation into insoluble filaments. Multiplication of α-synuclein or mutations such as A53T, A30P and E46K were found in familial forms PD patients [3][4][5].However, the detailed mechanisms underlying α-synuclein induced neurotoxicity in PD still need further investigation.
PD animal models have been established by ectopic expression of human α-synuclein in yeast, Caenorhabditis elegans, Drosophila melanogaster, rat, mouse, and non-human primates [6][7][8][9][10][11]. Drosophila models have been widely used to study neurodegenerative diseases including Alzheimer's disease (AD), Huntington's disease (HD) and PD [8,12,13]. In addition to the advantages of short lifespan and convenience for genetic manipulation, Drosophila conceives complicated central and peripheral nervous systems which are analogous to those of human. Panneuronal expression of human wild type and mutant α-synuclein (A53T and A30P) demonstrate adult onset PD pathological features including DA neuronal loss, decreased dopamine level, impaired locomotive ability and shortened lifespan [8,14,15]. Drosophila models provide efficient tools for screening genes participate in PD and potential drugs against PD.
MicroRNAs are~20nt small RNA molecules that fine-tune gene expression at posttranscriptional level [16]. They usually bind to 3'UTR of target mRNA and lead to translational inhibition or target degradation. It is estimated that more than half of human genes are regulated by miRNAs and the regulatory mechanisms are highly conserved among invertebrates and vertebrates. Since the discovery in 1990s, miRNAs have been found to exert essential roles in development, homeostasis and diseases. A plethora of miRNAs have been found to be dysregulated in the brain and blood of PD patients [17][18][19][20]. However, the underlying mechanisms and their functions in PD are still elusive.
In the present study, we examined the expression of miRNAs in a PD Drosophila model expressing α-synuclein by high throughput small RNA sequencing technology. We found that five miRNAs (dme-miR-133-3p, dme-miR-137-3p, dme-miR-13b-3p, dme-miR-932-5p, dme-miR-1008-5p) were upregulated in PD flies. Among them, miR-13b, miR-133, miR-137 are brain enriched and highly conserved from Drosophila to Homo sapiens. Validation experiment using qRT-PCR confirmed that these miRNAs were elevated in PD flies. KEGG pathway analysis indicated that neuroactive-ligand receptor interaction pathway was most likely affected by these miRNAs. Further studies showed miR-137 targeted multiple molecules in this pathway as predicted, including dopamine receptor (DopR, D2R), GABA receptor (GABA-B-R1, GABA-B-R3) and NMDA receptor (Nmdar2). The mRNA levels of these molecules were significantly decreased in PD flies. Our findings indicated that α-synuclein could induce the dysregulation of miRNAs, which target neuroactive ligand receptor interaction pathway in vivo.

Fly stocks and maintenance
The elav-C155 and UAS-α-synuclein flies were obtained from Bloomington Stock Center (Indiana University, USA). Flies were raised in standard yeast agar food at 25°C with a 12/12 hours light/dark cycle. After backcrossing with w 118 flies for 6 generations, elav-C155 virgin flies were crossed with w 118 or UAS-α-synuclein A30P males. The F1 generation offspring expressed α-synuclein in panneuronal manner and were used for further experiments.

Lifespan analyses
Two days after the eclosion, mated males and females were discriminated and transferred to different vials. Each vial contained 10 flies and at least 100 in total for each group. The vials were changed 3 times a week and deaths were recorded. Data was presented as survival curves and analysis was performed using log-rank tests to compare between groups.

Climbing assay
In order to characterize behavior defects in PD flies, climbing assay was performed as described previously [8]. Briefly, twenty male flies were transferred into an empty plastic vial and gently tapped to the bottom. The numbers of flies that could climb to the top (above 8cm) or remained at bottom in 18 seconds were recorded. The climbing assay was performed at least 3 times for every vial at each time point.

High throughput sequencing for miRNAs
Total RNA of each sample (three biological repeats for PD and control fly heads) was used to prepare the miRNA sequencing library through following steps: 1) 3'-adapter ligation with T4 RNA ligase; 2) 5'-adapter ligation with T4 RNA ligase; 3) cDNA synthesis with RT primer; 4) PCR amplification; 5) extraction and purification of~125-145 bp PCR amplified fragments (correspond to~15-35 nt small RNAs) from the PAGE gel. After the completed libraries were quantified with Agilent 2100 Bioanalyzer, the DNA fragments in the libraries were denatured with 0.1M NaOH to generate single-stranded DNA molecules, captured on Illumina flow cells, amplified in situ and finally sequenced for 36 cycles on Illumina HiSeq2000 according to the manufacturer's instruction. Raw sequences were generated as clean reads from Illumina HiSeq by real-time base calling and quality filtering. Subsequently, the 3' adapter sequence was trimmed from the clean reads and the reads with lengths shorter than 15 nt were discarded. As the 5'-adaptor was also used as the sequencing primer site, the 5'-adaptor sequence is not present in the sequencing reads. The trimmed reads (length 15 nt) were aligned to the fly pre-miRNA in miRBase 21, using novoalign software. The miRNA expression levels were measured and normalized as transcripts per million of total aligned miRNA reads (TPM). When comparing profile differences two groups of samples (PD and Control), the "fold change" (i.e. the ratio of the group averages) and p-value were caculated. miRNAs having fold changes 1.2, P-value 0.05 or fold change 2.0 were selected as the differentially expressed miRNAs.
qRT-PCR for miRNA CAGTGCGTGTCGTGGAGTCGGCAATTGCACTGGATACGACCTACGT-3' (dme-miR-137-3p); 5'-GTCGTATCCAGTGCGTGTCGTGGAGTCGGCAATTGCACTGGATACGA CACTCGT-3'(dme-miR-13b-3p). qPCR reaction was performed with 2×SYBR Green PCR Master Mix (Arraystar) and ViiA 7 Real-time PCR System (Applied Biosystems) with the program: 95°C for 10 min to denature DNA templates, followed by 40 cycles of 95°C for 15 s, 60°C for 60s.PCR primers for miRNAs and U6 were listed in S1 Table. miRNA Targets Prediction and Pathway Enrichment Analysis The target genes of differentially expressed miRNAs were predicted by miRanda-mirSVR algorithm and then subjected to GO analysis using Database for Annotation, Visualization and Integrated Discovery (DAVID) (count cutoff 10, EASE 0.01). The significantly affected GO terms (p<0.05) in biological process, cellular component and molecular function were identified.
DIANA-miRPath is an efficient tool for analyzing the combinatorial effect of microRNAs on target pathways. We uploaded the dysregulated miRNAs and predicted potential target pathways using DIANA-microT-CDS algorithm. The significantly influenced pathways (p<0.05) were identified.

Validation of target mRNA Expression
In order to validate the expression of predicted targets for dysregulated miRNAs, qRT-PCR was performed according to previously reported methods as mentioned before. PCR primers for mRNAs were listed in S2 Table.

Luciferase Reporter Assay
The 3'UTR fragments flanking miR-137 targeting sites of Nmdar2, D2R and GABA-B-R3 were cloned from Drosophila cDNA library and inserted into pGL3-promoter vectors respectively. Each of these vectors was co-transfected with Renilla plasmid pRL-TK and dme-miR-137-3p mimics (Genepharma, Shanghai) into HEK 293 cells in 12-well plates using Lipofactamine 2000 (Invitrogen). 24 hours post-transfection, luciferase activity was measured with Dual-Glo (Promega) according to the manufacturer's instructions. PCR primers for amplification were listed in S3 Table. Mutant constructs were made by site-directed mutagenesis to replace seed sequence with BglII cleavage site.

Statistics
Log-rank tests were performed to compare lifespan between groups. For other experiments, the significance of the difference was analyzed with Student's t test using GraphPad Prism software, and p <0.05 were considered statistically significant.

Characterization of PD Drosophila models
We establised PD fly models according to literatures [8,14,15]. Briefly, elav-Gal4 flies were crossedwith UAS-α synuclein (A30P) to ectopically express human α synuclein in nervous systems. As reported previously, we found that PD flies exhibited shorter lifespan ( Fig 1A) and impaired locomotive ability ( Fig 1B) compared with control flies. These results indicated that Drosophila models successfully developed adult-onset PD like phenotype in age dependent manner. Climbing ability of PD Drosophila was comparable with control flies at day 10 post eclosion ( Fig 1B), which was consistent with the results from Feany et al [8]. At this time point, loss of dopaminergic cells in PD Drosophila was also detected. Therefore, we chose day 10 flies post eclosion as early PD stage to investigate miRNA expression profiles.

Deep sequencing data analysis and verification
The miRNA samples from heads of control and PD flies were sequenced using Illumina HiSeq2000 platform. The total numbers of the reads at the sequencing data processing stages are listed for each sample (3 biological repeats for control and PD flies) in Table 1. The majority of small RNAs were 20-22nt which were the typical length for miRNAs (Figs 2 and S1). The reads can be divided into several groups (miRNAs, tRNAs, rRNAs, sRNAs, snRNAs, other ncRNAs). As shown with pie charts in Figs 3 and S2, the majority was miRNAs (86.1%-90.8%).

Functional annotations for targets of differentially expressed miRNAs
As four of the dysregulated miRNAs in PD flies including dme-miR-133-3p, dme-miR-137-3p, dme-miR-13b-3p and dme-miR-932-5p were brain enriched, we predicted targets of them and then submit to DAVID for Gene Ontology analysis (Fig 6 and S7 Table). GO enrichment analysis revealed that the target genes were functionally enriched in neuron related biological process (neurodifferentiation, neuron development, neuron projection development, neuron projection morphogenesis). In addition, cell component analysis showed that these targets were enriched in the membrane proteins.
DIANA miRPath is a powerful tool to analyze the combinational effects of miRNA on signaling pathways [21]. Using this method, we found dysregulated miRNAs significantly affect four pathways, of which neuroacitve-ligand receptor interaction was most significant (Fig 7  and Table 3). Four of the total dysregulated miRNA could target 8 genes of this pathway. As   Fig 1. α-synuclein-induced locomotion defects and shortened lifespan. (A) Expression of A30P α-synuclein specifically in the nervous system shortened lifespan. Survival curves were compared using the log-rank test (P<0.01 between elav>w 1118 and elav>α-synuclein A30P flies). (B) There is no difference for climbing ability for flies expressing A30P α-synuclein and control genotype at day 10. In contrast, PD flies showed significant age dependent locomotive impairments at days 30 Table 4, miR-137-3p potentially targeted Nmdar2 (receptor for N-Acetylaspartyl glutamate and Glutamate, L-asparate, L-cysteic acid, L-homocysteic acid), (Relaxin receptor) and AR2 (Galanin receptor) were predicted to be targeted by miR133-3p and miR-13b-3p respectively. In addition, miR-932-5p was proposed to act on AlstR (Galanin receptor) and GABA-B-R1 (GABA receptor). These results indicated that dysregulation of miRNAs potentially lead to interruption of neuroactive-ligand receptor signaling pathway and contributed to α-synuclein induced PD flies.

The mRNA levels of predicted targets were downregulated in PD flies
We examined the transcriptional levels of miR-137 targets in neuroactive ligand-receptor interaction pathway. Within five predicted targets, Nmdar2, GABA-B-R3, GABA-B-R1 and D2R were confirmed to be downregulated in PD flies (Fig 8). Particularly, the NMDA receptor Nmdar2 and GABA receptor GABA-B-R3 were most significant. Interestingly, hsa-miR-137-3p was also predicted to target KEGG pathways including Glutamatergic synapse (hsa04724) (p = 0.001749507) and GABAergic synapse (hsa04727) (p = 0.007160067) by DIANA miRPath analysis. GABA-B receptor (GABRA1, GABRA6, GABBR2) and NMDA receptor (GRIN2A) were identified as hsa-miR-137-3p targets (Table 5). Our results were consistent with previous reports that PD was associated with neuroactive ligand-receptor interaction pathway [22] and miR-137 could regulate synaptogenesis and neuronal transmission [23]. The regulatory effects of miR-137 on GRIN2A expression have been confirmed in human neuronal like SH-SY5Y cells [23]. Luciferase reporter assay showed that miR-137 could target GRIN2A directly in Rats [24]. The regulatory mechanisms seemed to be highly conserved from Drosophila to humans. In order to further confirm dme-miR-137 could directly regulate targets in neuroactiveligand receptor interaction pathway, we constructed luciferase reporter plasmids carrying Nmdar2, D2R and GABA-B-R3 3'UTR fragments containing miR-137 binding sites (Fig 9A).
Our study using high throughput sequencing of miRNAs identified miR-13b, miR-133, miR-137, miR-932 and miR-1008 consistently upregulated in early stage PD flies. Among the dysregulated miRNAs, miR-13b, miR-133 and miR-137 were highly conserved from Drosophila to H. sapiens and their expression was validated by qRT-PCR. MiR-13b's human homologue is miR-499 [57] that expressed in brain region and its polymorphism is associated with ischemic stroke [58]. Previously, we found miR-13b was also upregulated in adult onset AD flies [59]. These results indicate that miR-13b/miR-499 play important roles in pathogenesis of brain insults. MiR-133a and miR-133b are human orthologs of dme-miR-133 and enriched in human brain. Exosomes containing miR-133b from mesenchymal stem cells (MSCs) regulate neurite outgrowth of neural cells [60]. Morphine regulates dopaminergic neuron differentiation via miR-133b [61]. In addition to its physical functions, miR-133b is essential for functional recovery after spinal cord injury in adult zebrafish [62]. By targeting Pitx3, miR-133b was found to regulate the maturation and function of midbrain dopaminergic neurons, contributing to PD pathogenesis [55]. MiR-137 is also a highly conserved miRNA and exerts  important roles in neuronal development and diseases. By regulating expression of nuclear receptor tailless (TLX) and lysine-specific demethylase 1 (LSD1) in neural stem cells, miR-137 controls the dynamics between neural stem cell proliferation and differentiation during neural development [63]. MiR-137 could also regulate neuronal maturation by targeting ubiquitin ligase mind bomb-1 [64]. Recently, it was reported that miR-137 and its seed-similar fly homologue miR-1000 regulated vesicular glutamate transporter (VGlut) expression and fine-tune excitatory synaptic transmission [65]. In addition, miR-137 also plays important roles in brain disorders. MiR-137 is associated with intellectual disability [66]. miR-137 is also proved to be associated with schizophrenia susceptibility, which usually accompanied with PD [67][68][69]. The mechanistic studies reveal that miR-137 regulates gene sets involved in synaptogenesis and neuronal transmission as well as glucocorticoid receptor-dependent signalling network, contributing to etiology of schizophrenia [23,70]. In another neurodegenerative disorder Alzheimer's disease, miR-137 is found to be associated with serine palmitoyltransferase (SPT) and amyloid β (Aβ) levels [71].
The reason for α-synuclein induced miRNA dysregulation in vitro could be explained in various mechanisms. Firstly, α-synuclein overexpression and aggregation in neuronal cells may influence signaling pathways and transcription factors that mediate miRNA expression. α-synuclein expression could influence signaling pathways including IRS-1/insulin/Akt, mTOR/S6K, MAPK, p53, GSK3β, PKC, synaptic transmission, ubiquitin protein pathway (UPS) and the autophagy pathway [72][73][74][75][76][77][78][79]. These pathways could further stimulate transcription factors and lead to miRNA dysregulation. ɑ synuclein could increase the activity of transcription factors including NRF2, NFAT, MEF2C-PGC1α, CREB, NF-κB, p53, Nurr1, and FOXP1 [74,[80][81][82][83][84][85][86]. We analyzed the promoter region of dme-miR-137 (5kb upstream of predme-miR-137) using AliBaba2.1 based on TRANSFAC 4.0 and found three CREB binding sites as well as six NF-κB binding sites. In addition, CREB and NF-κB were also predicted to bind to hsa-miR-137 promoter, indicating the regulatory mechanisms were highly conserved. Taken together, ɑ synuclein may induce miR-137 expression by transcription factor CREB and NF-κB. Second, α-synuclein overexpression and aggregation in neuronal cells may stimulate cells to release different factors that induce miRNA expression. These factors include brainderived neurotrophic factor, glial cell line-derived neurotrophic factor, reactive oxygen species, nitric oxide and other factors [87][88][89][90]. These factors may act on other cells and actviate relevant signaling pathways as well as downstream transcription factors and induce miRNA expression. Third, the extracellular α-synuclein aggregates may act on neurons to regulate miRNA expression. It was reported that exogenous α-synuclein fibrils induce could activate singlaing pathways including PI3/Akt, calpain-dependent CDK5, LKB1/AMPK/Raptor, leading to synaptic dysfunction and neuron death [91][92][93][94]. Extracellular alpha-synuclein may also induce miRNA expression in vitro. Detailed experiments are required to clarify this problem In order to elucidate which signaling pathways potentially affected by these dysregulated miRNAs in PD flies, DIANA-miRPath analysis was performed and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway neuroactive-ligand receptor interaction was identified. Consistent with our findings, Huang et al. reported that when applying to a genome-wide association study (GWAS) dataset for Parkinson disease, extended Bayesian lasso (EBLasso) identified three significant pathways including the neuroactive-ligand receptor interaction, the primary bile acid biosynthesis pathway, and the mitogen-activated protein kinase (MAPK) signaling pathway [95]. Our validation experiments showed that downregulations of NMDA results showed that dme-miR-137-3p could inhibit the luciferase activities for all these vectors does dependently as compared with miR-negative control. (C) The inhibitory effects were abolished when all the miR-137 targeting sites within the amplified sequences in D2R 3'UTR were mutated. (* p<0.05, ** p<0.01).
doi:10.1371/journal.pone.0137432.g009 miRNAs Mediate α-Synuclein Toxicity in PD Drosophila receptor (Nmdar2) and GABA receptor (GABA-B-R3) were most significant. NMDA receptor GRIN2A was also predicted to be targeted by miR-137 in Homo sapiens, which have been validated in human SY-SH5Y cells [23]. Luciferase reporter assay showed that miR-137 could target GRIN2A directly in Rats [24], suggesting the regulatory mechanisms seemed to be highly conserved from Drosophila to humans. Interestingly, Genome-wide gene-environment study identifies glutamate receptor gene GRIN2A as a Parkinson's disease modifier gene via interaction with coffee [95]. Activation of GABAB receptors within the substantia nigra pars reticulata (SNr), but not the globus pallidus (GP), reverses reserpine-induced akinesia in rats. The success of intracerebroventricular injection of baclofen suggests a potential for GABAB receptor agonists in the treatment Parkinson's disease [96]. Hillman R et al reported that GABA rescue the loss of climbing activity in this PD fly models [40]. More specifically, GABA(B) agonists baclofen and the allosteric agonists CG 7930 and GS 39783 could also ameliorate locomotive defects, which diminished when flies are cofed with the GABA(B) receptor antagonist 2-hydroxysaclofen. In contrast, GABA(A) receptor agonist muscimol has no effect. This result indicated the important roles for neuroactive-ligand pathways in PD. Next step, we will use genetic manipulations and pharmacological methods to clarify the role miRNA-targets axis we identified within this pathway in PD.

Conclusions
Our findings indicated that α-synuclein could induce the dysregulation of highly conserved and brain enriched miRNAs, which target neuroactive ligand-receptor interaction pathway in vivo. We believe it will contribute to understanding miRNA functions in mediating α-synuclein toxicity and provide new insights into the pathogenesis driving PD.