Alpha-synuclein overexpression induces epigenomic dysregulation of glutamate signaling and locomotor pathways

Abstract

Parkinson's disease (PD) is a neurological disorder with complex interindividual etiology that is becoming increasingly prevalent worldwide. Elevated alpha-synuclein levels can increase risk of PD and may influence epigenetic regulation of PD pathways. Here, we report genome-wide DNA methylation and hydroxymethylation alterations associated with overexpression of two PD-linked alpha-synuclein variants (wild-type and A30P) in LUHMES cells differentiated to dopaminergic neurons. Alpha-synuclein altered DNA methylation at thousands of CpGs and DNA hydroxymethylation at hundreds of CpGs in both genotypes, primarily in locomotor behavior and glutamate signaling pathway genes. In some cases, epigenetic changes were associated with transcription. SMITE network analysis incorporating H3K4me1 ChIP-seq to score DNA methylation and hydroxymethylation changes across promoters, enhancers, and gene bodies confirmed epigenetic and transcriptional deregulation of glutamate signaling modules in both genotypes. Our results identify distinct and shared impacts of alpha-synuclein variants on the epigenome, and associate alpha-synuclein with the epigenetic etiology of PD.

Figure 1

aSyn overexpression was associated with decreased DNAm at the SNCA 5′ UTR. (A) Gene-level SNCA mRNA expression in nRPKM. p-values: DESeq2 Wald test with Benjamini–Hochberg adjustment. WT aSyn vs. control and A30P aSyn vs. control padj < 2.2 × 10⁻³⁰⁸. (B) Transcript-level SNCA mRNA expression in TPM. See Supplementary Material, Table S1 for statistics. (C) Top: hg19 coordinates are shown for the SNCA gene, with Ensembl transcripts below. Boxes represent exons, lines represent introns. Bottom: DNAm beta values are shown for each sample, colored by genotype. Only EPIC probes specific to endogenous SNCA are shown. Control cells: n = 7 biological replicates; WT aSyn cells: n = 8; A30P aSyn cells: n = 8. (D) hg19 coordinates, DNAm beta values, and transcription factor binding sites for the SNCA 5′ UTR/intron 1/TSS upstream region. p-values: limma empirical Bayes moderated t test with Benjamini–Hochberg adjustment. Bottom: ENCODE transcription factor ChIP-seq clusters.

Figure 2

Overexpression of WT and A30P mutant aSyn altered the DNA methylome of dopaminergic neurons. (A) Volcano plots comparing DNAm patterns between control (n = 7 biological replicates) and WT aSyn (n = 8) LUHMES cells, control and A30P aSyn (n = 8) LUHMES cells, and WT aSyn and A30P aSyn LUHMES cells. Colored points passed thresholds of delta beta ≥0.05 and padj ≤ 0.05 (limma empirical Bayes moderated t test with Benjamini–Hochberg adjustment). (B) Number of differentially methylated probes unique to each comparison and probes shared between comparisons (padj < 0.001, 10 000 permutations). (C) Relative enrichment/depletion of differentially methylated sites across genomic contexts, permuted against array background (10 000 iterations). (D) Relative enrichment/depletion of differentially methylated sites by relation to CpG island, permuted against array background (10 000 iterations).

Figure 3

WT and A30P aSyn overexpression were correlated with increased DNAhm levels in dopaminergic neurons. (A) Volcano plots comparing DNAhm patterns between control (n = 7 biological replicates) and WT aSyn (n = 8) LUHMES cells, control and A30P aSyn (n = 8) LUHMES cells, and WT aSyn and A30P aSyn LUHMES cells. Colored points passed thresholds of delta beta ≥0.05 and padj ≤ 0.05 (limma empirical Bayes moderated t test with Benjamini–Hochberg adjustment). (B) Number of differentially hydroxymethylated probes unique to each comparison and probes shared between comparisons. All padj = 1 (10 000 permutations). (C) Relative enrichment/depletion of differentially hydroxymethylated sites across genomic contexts, permuted against array background (10 000 iterations). (D) Relative enrichment/depletion of differentially hydroxymethylated sites by relation to CpG island, permuted against array background (10 000 iterations).

Figure 4

Differentially methylated genes in both genotypes were enriched for locomotor behavior and glutamate signaling functions. (A) GO biological process terms enriched at multifunctionality-corrected, multiple test-corrected hypergeometric p-value (CorrectedMFPvalue) < 0.05 in genes shared between control (n = 7 biological replicates) vs. WT aSyn (n = 8) and control vs. A30P aSyn (n = 8) DNAm analyses. (B) Heat map showing beta values for all probes annotated to the ‘locomotory behavior’ pathway and differentially methylated in at least one comparison. Row labels (left to right): probes that passed significance thresholds (absolute delta beta ≥0.05 and padj ≤ 0.05, limma empirical Bayes moderated t test with Benjamini–Hochberg adjustment) in control vs. WT aSyn comparison, probes that passed significance thresholds in control vs. A30P aSyn comparison. (C) Representative examples of 353 total locomotory behavior-related CpG sites differentially methylated only in WT aSyn cells (DAB1 padj = 1.91e⁻⁴, CRH padj = 9.34e⁻⁵). (D) Representative examples of 22 total locomotory behavior-related CpG sites differentially methylated only in A30P aSyn cells (CYTH2 padj = 1e⁻⁵, FXYD1 padj = 0.028). (E) Representative examples of 67 total locomotory behavior-related CpG sites differentially methylated in both genotypes (PCDH15 WT aSyn padj = 7.8e⁻⁶, PCDH15 A30P aSyn padj = 8.61e⁻⁵, PAPP2 WT aSyn padj = 1.68e⁻⁹, PAPP2 A30P aSyn padj = 3.55e⁻⁷). p-values: limma empirical Bayes moderated t test with Benjamini–Hochberg adjustment.

Figure 5

Modules from SMITE analysis were driven by gene expression, with contributions from DNAm and/or DNAhm. (A) Top: Names and ranks of 18 SMITE modules in control vs. WT aSyn LUHMES cells. Bottom: Distribution of significant (p ≤ 0.05) and non-significant (p > 0.05) SMITE combined p-values (Fisher’s method) for each score category, across 136 genes with at least one significant category in control vs. WT aSyn analysis. (B) Top: Names and ranks of 24 SMITE modules in control vs. A30P aSyn LUHMES cells. Bottom: Distribution of significant (p ≤ 0.05) and non-significant (p > 0.05) SMITE combined p-values (Fisher’s method) for each score category, across 225 genes with at least one significant category in control vs. A30P aSyn analysis. DNA(h)m data: n = 7 control LUHMES cells and n = 8 WT aSyn and A30P aSyn LUHMES cells, expression data: n = 3/group (9), H3K4me1 ChIP-seq data: n = 3/group. All samples are distinct biological replicates.

Figure 6

GRIK2 was differentially methylated, hydroxymethylated, and expressed in LUHMES cells overexpressing either WT or A30P aSyn. (A) GRIK2 protein–protein interaction network showing differential DNAm, DNAhm, and expression in control vs. WT aSyn cells for the underlying genes. (B) GRIK2:DLG3 protein–protein interaction network showing differential DNAm, DNAhm, and expression in control vs. A30P cells for the underlying genes. (C) Relative expression of GRIK2 mRNA in each genotype, normalized relative to control cells (9). WT aSyn padj = 1.46e⁻⁶, A30P aSyn padj = 2.71e⁻¹⁵ (DESeq2 Wald test with Benjamini–Hochberg adjustment). Error bars: SE. (D) Top: UCSC hg19 coordinates are shown, with GRIK2 transcripts below. Middle: LUHMES cell H3K4me1 consensus peaks are shown in blue. Bottom: DNAm beta values are shown for each sample, colored by genotype. (E) DNAm at the GRIK2 TSS200 region. (F) DNAhm levels for EPIC array probes across the GRIK2 gene. p-values: limma empirical Bayes moderated t test with Benjamini–Hochberg adjustment. DNA(h)m data: n = 7 control LUHMES cells and n = 8 WT aSyn and A30P aSyn LUHMES cells, expression data: n = 3/group (9), H3K4me1 ChIP-seq data: n = 3/group. All samples are distinct biological replicates.