Single-nuclei paired multiomic analysis of young, aged, and Parkinson’s disease human midbrain reveals age-associated glial changes and their contribution to Parkinson’s disease

Levi Adams, Min Kyung Song, Yoshiaki Tanaka, Yoon-seong Kim
doi: https://doi.org/10.1101/2022.01.18.22269350

Abstract

Age is the primary risk factor for Parkinson’s disease (PD), but how aging changes the expression and regulatory landscape of the brain remains unclear. Here, we present a single-nuclei multiomic study profiling shared gene expression and chromatin accessibility of young, aged and PD post-mortem midbrain samples. Combined multiomic analysis of midbrain nuclei along a pseudopathogenesis trajectory reveals all glial cell types are affected by age, but microglia and oligodendrocytes are further altered in PD. We present evidence for a novel age-associated oligodendrocyte subtype that appears during normal aging, characterized by elevated protein folding and chaperone-mediated autophagy pathways. Differential gene and protein-protein interaction analyses show that these functions are significantly compromised in PD. Peak-gene association from our paired data identifies differential cis-elements linked to regulation of differentially expressed genes between PD patients and neurologically healthy controls. Our study suggests a previously undescribed role for oligodendrocytes in aging and PD pathogenesis.

Parkinson’s disease (PD) is the second most common neurodegenerative disease and is estimated to affect more than ten million people globally1. PD prevalence increases from 41 cases per 100,000 individuals in the fourth decade of life to 425 per 100,000 in the sixth decade of life and to 1908 per 100,000 in the eighth decade of life2. During PD pathogenesis, dopaminergic neurons in the substantia nigra degenerate, resulting in neurological, cognitive and motor symptoms3. While many genetic components and environmental risk factors have been identified, aging is the primary risk factor for PD4.

Despite the connection between aging and PD, relatively little attention has been given to the changes in the aging midbrain and how these changes may predispose individuals to development of PD. Transcriptomic studies of mouse or fly brains have provided some insight into aging, but have not focused on the midbrain5, 6. Glaab and Schneider7 used microarray datasets to investigate shared pathways and network alterations in aging and PD, but there remains a lack of information on how alterations in gene expression during aging affect the different cell types in the midbrain and how they contribute to PD pathogenesis.

The field of single-cell transcriptomics and epigenomics has facilitated important advances in understanding how gene expression and chromatin accessibility contribute to neurodegeneration8–10. Recent single-cell expression studies of the midbrain have highlighted a potential role of oligodendrocytes in the midbrain in PD pathogenesis, supported by strong genome-wide association study (GWAS) and transcriptomic data11–14. This is a surprising finding because vulnerable dopaminergic neurons that are lost during PD pathogenesis are sparsely myelinated15. Despite these crucial insights, there remains a lack of cell type-specific high-resolution data on what distinguishes ‘healthy’ aging from neurodegeneration and an incomplete picture of how the epigenetic landscape changes during aging and PD.

To address this gap, we isolated nuclei from the substantia nigra of post-mortem midbrains of young and aged donors with no neurological disease as well as PD patients and compared the single-nuclei transcriptome and genome-wide chromatin accessibility simultaneously from each nucleus along aging and PD trajectories. This approach allows us to directly infer cis-acting elements that contribute to gene expression in the same cells and provides a more complete picture of how the aging process affects gene regulation and expression in distinct cell types in the midbrain. This multiomic analysis reveals a novel age-associated oligodendrocyte that may contribute to PD pathogenesis.

Results

Single-nuclei RNA-seq and ATAC-seq profiling of the human post-mortem midbrain

We obtained frozen, post-mortem midbrain samples of young (mean 24 years old) and aged (mean 75 years old) neurologically healthy donors and patients diagnosed with PD (mean 81 years old) (Supplementary Table 1). We isolated nuclei from the substantia nigra and performed paired single-nucleus RNA-sequencing (snRNA-seq) and single-nucleus Assay for Transposase Accessible Chromatin sequencing (snATAC-seq) on each nucleus (10X Genomics Single Cell Multiome ATAC + Gene Expression kit) (Fig. 1a, Supplementary Fig. 1a-c). After filtering to remove low-quality reads or potentially confounding signals caused by multiple nuclei being captured in a single bead (Methods), we retained 69,289 high-quality nuclei from 31 individuals (9 young donors, 8 aged donors, 14 PD patients), showing low rates of doublets based on gene expression and robust ATAC data quality (Supplementary Fig. 1d-g). Using these nuclei, we performed batch correction, variable gene and principal component analysis and UMAP dimensional reduction with Seurat v416 (Methods). We identified 23 separate clusters of nuclei in our snATAC-seq and snRNA-seq datasets (Extended Data Fig. 1a). Interrogation of gene expression patterns of known cell type markers10, 17, 18 allowed us to classify our nuclei into seven major cell types (Fig. 1b,d, Extended Data Fig. 1b): Neurons (N), Oligodendrocytes (ODC), Astrocytes (AS), Microglia (MG), Oligodendrocyte Precursor Cells (OPC), Endothelial Cells (EC), and peripheral immune cells/T-cells (T). As expected, clustering based on ATAC profiling generated the same major clustering patterns for these cells with distinctive chromatin accessibility profiles at key loci (Fig. 1c,e, Extended Data Fig. 1b). We found the vast majority of cells in the midbrain are ODC, followed by MG, OPC and AS (Fig. 1f, Supplementary Fig. 1h). Similar to recently published data14, we did note a statistical difference in MG and ODC populations between groups, with a significantly higher proportion of ODC (86.9%) and lower proportion of MG (3.8%) in the aged group compared to either the young (ODC, 76.5%; MG, 9.7%) or PD group (ODC, 74.0%; MG, 9.4%) (Supplementary Table 2).

Extended Data Fig. 1
Extended Data Fig. 1 Multiomic analysis of human midbrain.

a, UMAP visualization of single nuclei by RNA (left) and ATAC (right) profiles. Nuclei are colored by 23 joint clusters. b, Heatmap showing Spearman correlation of average RNA expression (left) and ATAC peaks profiles (right) by cell types for each individual. Top and second color bar represents groups of donors and cell types, respectively. c, Heatmap showing RNA expression and ATAC peak intensity of cell type-specific genes and their nearest peaks, respectively. ATAC peak intensity is calculated in the nearest peaks from TSS. Pearson correlation across cell types in each gene-peak pair is shown in the right panel (Cor=0.797±0.226). d, Enrichment of motifs in each annotated cell type. e-g, Dot plot displaying enrichment of gene expression modules for functionally distinct subpopulations for ODC (e), MG (f) and AS (g) clusters in the human midbrain.

Fig. 1
Fig. 1 Multiomic analysis of human midbrain.

a, Schematic of isolation and single-nuclei ATAC plus Gene Expression analysis for human post-mortem midbrain. b,c, UMAP visualization of single nuclei by RNA (b) and ATAC (c) profiles. Nuclei are colored by identified cell type (left) and donor type (right). d, Expression of cell type-specific genes in each annotated cell type cluster. e, Chromatin accessibility at cell type-specific genes in each annotated cell type cluster. f, Percentage of each cell type passing QC criteria identified in young, aged and PD. Student’s t-test shows ODC are significantly changed from young to aged (p-value, 0.002), and aged to PD (p-value, 0.023). MG are significantly changed from young to aged (p-value, 0.007) and aged to PD (p-value, 0.042).

Genes with cell type-specific expression patterns also showed differential ATAC peaks nearest to the transcription start site (TSS), which implies unique promoter accessibility for each cell type (Extended Data Fig. 1c), and distinctive enrichment of cell type-specific transcription factor binding motifs (Extended Data Fig. 1d). Within each cell type, we analyzed differentially expressed genes (DEG) between young and aged groups, as well as aged and PD groups (Supplementary Table 3). Several genes that have been previously associated with neurological disorders, such as NEAT1, FKBP5 and SLC38A2, were differentially expressed in multiple cell types (Supplementary Fig. 2a,b, Methods). We identified cell type-specific changes in gene expression across groups (Supplementary Fig. 2c,d,f) that showed potential alterations in neurological function and metabolic pathways (Supplementary Fig. 2e,g). The majority of differentially expressed genes were cell-type specific (Supplementary Fig. 3a). Integration of chromatin accessibility data showed that the accessibility of the nearest peak to the TSS of differentially expressed genes generally did not change between young and aged (Supplementary Fig. 3b,c) or aged and PD (Supplementary Fig. 3d,e) within the same cell type. Only approximately 10% of putative promoters for genes differentially expressed during aging or PD pathogenesis had correlation values over 0.85, as opposed to cell-type specific genes which strongly correlate with promoter peaks (Extended Data Fig. 1c), suggesting that the expression of DEG between groups might be controlled by other distal regulatory elements, not by a promoter.

We noticed an inverse correlation between number of nuclei and number of DEG for each group, with the larger groups having fewer DEG. This is a common effect in transcriptomics known as Simpson’s Paradox (thoroughly discussed by Cole Trapnell19, 20) where analyzing large, heterogeneous datasets can lead to erroneous conclusions from the presence of unlike groups. When we examine the gene expression signatures of previously defined functional subtypes21–26 (Methods) with our ODC, MG and AS clusters, we saw functionally distinct cluster patterns. In ODC, we identified clusters corresponding to newly formed ODC, myelin forming ODC and mature ODC (Extended Data Fig. 1e, Supplementary Fig. 4a,d). As previously reported, these clusters showed minimal overlap23, 26. We also noticed clusters with high expression of genes relating to neuronal and synaptic support (such as NRXN3 and NFASC), and clusters with increased expression of ODC-neuron adhesion markers (such as STMN1 and HAPLN2). We found three clusters of microglia including homeostatic MG, aging MG and a small population of Stage1 MG (TREM2-independent) disease-associated MG (Extended Data Fig. 1f, Supplementary Fig. 4b,e)21, 24. Three clusters were found in astrocytes that included signatures for disease-associated astrocytes (previously identified in Alzheimer’s disease), A1 reactive astrocytes, and homeostatic GFAP-low astrocytes (Extended Data Fig. 1g, Supplementary Fig. 4c,f)22, 25.

Multiomic analysis of peak-gene associations

The paired single-nucleus gene expression and ATAC-seq data from the same nuclei provides a unique advantage because it allows us to infer the relationship between chromatin accessibility and gene expression in cis. We utilized the analytical framework developed by Ma et al. to analyze paired snRNA-seq and snATAC-seq data8 to generate peak-gene associations (Fig. 2a). Correlated associations were calculated in the chromatin regions within ±500 kb of the TSS of annotated genes where there is co-variation between chromatin accessibility and gene expression. In each cell type of each cohort, an average of 193,732 significantly peak-gene associations (5.3 peaks per gene) were identified. Notably, we found most peak-gene associations (66.6-92.2%) were exclusive to young, aged or PD cohorts (Supplementary Fig. 5a), and a minority of associations (0.274-10.6%) were shared between groups, which was consistent in all cell types (Supplementary Fig. 5b).

Fig. 2
Fig. 2 Analysis of gene-peak connections in human midbrain.

a, Schematic of gene-peak analysis between samples. b, Comparison of the number of connected peaks in each gene between young and aged (left), and aged and PD (right) midbrain. Some genes show altered peak-connections from young in both aged and PD (red dots). Some genes show altered peak-connections only in PD (blue dots). c,d, Gene-peak connection plots for (c) NEAT1 and (d) RASGRF1 between young, aged and PD midbrain. Top panel: normalized ATAC read count distribution. Bottom panel: arcs show ATAC peaks significantly correlated with gene expression. Arc height represents statistical significance.

More than half of shared peak-gene associations (50.9-56.6%) were near the TSS (within 5 kb) and unique peaks were linked in more distal regions, reinforcing the idea that expression changes between groups are likely driven by more distal regions such as enhancers (Supplementary Fig. 5c). In agreement with previous data8, we also observed that not all peaks were connected to a gene (91,177 out of 210,609 peaks were connected to at least one gene), and there was an average of 6,812 connecting to more than 10 peaks (Supplementary Fig. 5d-e). Comparing peak-gene associations in ODC between groups identified a subset of genes whose correlated associations were significantly changed during aging and PD progression (Fig. 2b). Increased peak-gene associations in both aged and PD groups compared to young (Fig. 2b, red dots) may indicate aging-specific changes, and notably include the previously discussed aging marker, NEAT1. Conversely, some associations were unchanged between young and aged, but increased only in the PD cohort, including genes previously linked with longevity such as RASGRF127 and MSRA28 (Fig. 2b, blue dots). Peak-gene associations for NEAT1 were largely shared in older (aged and PD) groups compared to young samples (Fig. 2c, Extended Data Fig. 2a), while for RASGRF1, similar patterns were observed in normal control donors (young and aged) with much higher associations compared to PD (Fig. 2d, Extended Data Fig. 2b). Analysis of cis-regulatory motifs of the 43 correlated peaks for NEAT1 specific to the older groups (aged and PD) showed enrichment of binding motifs for transcription factors related to aging such as EGR1/229 (Extended Data Fig. 2a, c). For RASGRF1, we found 61 significant shared peaks in the control groups (young and aged), and enriched motifs for NRF2 and ASCL1, which are related to brain health and neurogenesis30, 31. These motifs were not utilized in PD samples (Extended Data Fig. 2b, d).

Extended Data Fig. 2
Extended Data Fig. 2 Changes in gene-peak connections in young, aged, and PD for NEAT1 and RASGRF1. a,b,

Venn diagram of the number of associated peaks with NEAT1 (a) and RASGRF1 (b) for young, aged and PD midbrain. c,d, Heatmap showing enrichment of TF binding motifs in associated peaks with NEAT1 (c) and RASGRF1 (d).

Defining pseudopathogenesis trajectory in oligodendrocytes

Next, we analyzed each major cell type (ODC, MG, OPC and AS) individually to explore potential changes in gene expression and epigenetic dynamics during the aging and disease process. To highlight heterogeneity and subtle changes within each cell type, we re-clustered each type of cell individually (Methods). Surprisingly, we noticed that for both ATAC and RNA data, nuclei from normal and control groups clustered more closely than PD nuclei (Fig. 3a) although they did not cluster by individual donors (Supplementary Fig. 6a). To identify potential transitions in gene expression and chromatin accessibility between clusters, we performed pseudotime analysis using Monocle332. Pseudotime analysis generated a striking trajectory in both snRNA and snATAC clusters moving from young to aged and then to PD patients in a way we termed pseudopathogenesis (Fig. 3a). Taking advantage of the paired multiomic nature of our data, we combined the snRNA and snATAC pseudopathogenesis scores and corrected for error to establish a combined pseudopathogenesis (cPP) trajectory (Fig. 3b). This cPP score takes into account the contribution of both gene expression and epigenomic status of each individual nucleus. Comparison of the cPP scores between groups exhibited significant increases from young to aged, and from aged to PD (Fig. 3c). Next, we investigated differentially expressed genes and peaks across the cPP trajectory and identified 299 increased and 474 decreased genes, and 912 increased and 1590 decreased peaks (Fig. 3d, Extended Data Fig. 3a, Supplementary Fig. 6b, Supplementary Table 4). Notably, aging or PD relevant genes that we identified in peak-gene association analysis, such as NEAT1 and RASGRF1 (Fig. 2c,d), were also correlated with cPP, supporting the ability of this analysis method to capture physiologically relevant changes (Extended Data Fig. 3b). Gene ontology analysis showed that pathways related to Response to Unfolded Protein, Chaperone-Mediated Autophagy (CMA) and Negative Regulation of Cell Death were increased while Myelination, Receptor Clustering and Regulation of Membrane Potential were decreased with advancing pseudopathogenesis (Fig. 3d). Of note, we observed an increase in genes related to protein translation and ribosomal function (RPS family genes) that have been previously reported to be increased during aging6 (Supplementary Table 4). Combined analysis of gene expression modules for Myelination and Neuronal and Synaptic Support exhibited decreases along the cPP trajectory while gene modules for Unfolded Protein Response and CMA increased with cPP. We noticed that almost 80% of ODC had low cPP scores of less than 7 (99.8% of young, 90% of aged and 50% of PD), and gene modules for normal functions like myelination and synaptic support showed little change at lower cPP scores. To examine this in greater detail, we divided the nuclei into nine groups of equal size based on cPP score (Supplementary Fig. 6c) and examined combined module scores. Notably, a sharp increase in protein folding and CMA pathways were only observed in the top two cPP groups, corresponding with a decline in normal ODC functions such as myelination (Extended Data Fig. 3c).

Extended Data Fig. 3
Extended Data Fig. 3 cPP analysis of ODC.

a, Example genes that have correlated increase (HSP90AA1) and decrease (CTNNA3) of expression with cPP trajectory. b, Plot of expression of NEAT1 and RASGRF1 across cPP trajectory. NEAT1 expression is correlated with increasing cPP score; RASGRF1 expression is inversely correlated with cPP score (Cor > 0.1 or < −0.1). c, Dot plot of gene expression modules for functional cellular subtypes in the midbrain for ODC over cPP. Groups contain equal number of nuclei, which are ordered by increasing cPP. Violin plots show expression for example genes from each of the modules.

Fig. 3
Fig. 3 Establishment of pseudopathogenesis trajectory in ODC.

a, UMAP plot of ODC nuclei colored by young, aged and PD donor (left), and heatmap of RNA- and ATAC-based pseudopathogenesis trajectory (right). b, Schematic of calculation of combined pseudopathogenesis (cPP) score from transcriptome and ATAC profiles. Cor = −0.419 (p<2.2e-16). c, Significant difference in cPP scores of individual ODC nuclei from young, aged and PD midbrain. (One-way ANOVA, p-value for Y/A < 1e-199, p-value for A/P < 1e-199). d, Heatmap showing genes correlated with cPP trajectory. X-axis represents individual cells sorted by cPPs. Y-axis of heatmap represents positively (upper)- and negatively (bottom)-correlated genes. Representative genes and significant GO terms are shown in right panel (Cor > 0.1 or < −0.1). e, Change in gene expression modules across cPP trajectory in ODC. Top panel shows individual ODC nuclei along cPP scores and donor group. X-axis shows cPP score. Y-axis is combined expression score for all genes in the expression module.

Pseudopathogenesis analysis shows distinct changes in microglia

We applied the same analytical framework to microglia and generated a cPP score for each nucleus (Fig. 4a,b, Supplementary Fig. 7a). We found significant increases in cPP between young, aged and PD samples for MG (Fig. 4c). While 894 genes were increased as cPP score increased, 254 genes were decreased (Fig 4d, Supplementary Table 4). Gene ontology analysis of cPP-relevant genes showed a loss of cell adhesion and chemotaxis during cPP increase and elevated immune activation and cytokine-mediated signaling pathway (Fig. 4d). Combined gene module analysis showed a decrease in homeostatic gene signatures and increases in aging and stage1 (TREM2 independent) disease-associated MG during cPP increase (Fig. 4e) as well as in three equal-sized groups based on cPP (Fig. 4f). We also found distinctive ATAC peak changes with cPP progression (Supplementary Fig. 7b). We analyzed astrocyte clusters by applying the same strategy (Extended Data Fig. 4a,b, Supplementary Fig. 8a). We found that there is a statistically significant increase in cPP score from young to aged samples, but no difference between aged and PD nuclei (Extended Data Fig. 4c). This suggests that, in the midbrain, despite changes in gene expression for AS during the aging process, they may not play a major role in PD pathogenesis. Gene module analysis indicated increases in A1-reactive AS signatures, and disease-associated AS signatures near the middle of the trajectory, with a concomitant reduction of GFAP-low module signatures across cPP, as well in three equal-sized groups based on cPP (Extended Data Fig. 4d). Gene ontology analysis of trajectory-relevant gene expression (Supplementary Table 4) showed generalized increases in apoptosis resistance and CMA pathways, but a reduction in neuronal support (Supplementary Fig. 8b). OPC exhibited similar results to AS (Supplementary Fig. 9a,b,10a, Supplementary Table 4) with a significant increase in cPP scores from young to aged nuclei, but no difference in PD samples (Supplementary Fig. 9c). Gene ontology analysis of trajectory-relevant genes showed a similar result as astrocytes with an increase in general stress-resistance markers during cPP progression, but a reduction in differentiation and synaptic signaling pathways (Supplementary Fig. 10b).

Extended Data Fig. 4
Extended Data Fig. 4 Establishment of pseudopathogenesis trajectory in AS.

a, UMAP plot of AS nuclei colored by young, aged and PD donor. b, UMAP plot of AS nuclei colored by cPP. c, cPP scores of individual AS nuclei from young, aged and PD midbrain are significantly changed over aging but not disease state. (One-way ANOVA, p-value for Y/A = 0.002, p-value for A/P = 0.24). d, Gene expression modules across AS cPP trajectory (left). Top panel shows individual nuclei AS along cPP scores and donor group. X-axis shows cPP score. Y-axis is combined expression level for all genes in the expression module. Dot plot of AS splits into three equal groups based on cPP displaying enrichment of gene expression modules (middle). Violin plots show expression for example genes from each of the modules (right).

Fig. 4
Fig. 4 Establishment of cPP trajectory in MG.

a, UMAP plot of MG nuclei colored by young, aged and PD donors. b, UMAP plot of MG nuclei colored by cPP. c, Significant difference in cPP scores of individual MG nuclei from young, aged and PD midbrain. (One-way ANOVA, p-value for Y/A < 7e-130, p-value for A/P <9.8e-15). d, Heatmap showing genes correlated with cPP trajectory. X-axis represents individual cells sorted by cPPs. Y-axis represents positively (upper)- and negatively (bottom)-correlated genes. Representative genes and significant gene ontology terms are shown in right panel (Cor > 0.1 or < −0.1). e, Gene expression modules across MG cPP trajectory. Top panel shows individual MG nuclei along cPP scores and donor group. X-axis shows cPP score. Y-axis is combined expression score for all genes in the expression module. f, Dot plot of MG split into three equal groups based on cPP displaying enrichment of gene expression modules. Violin plots show expression of example genes from each of the modules.

In-depth analysis of high cPP oligodendrocytes shows differences between PD and aged nuclei

Pseudopathogenesis trajectory analysis showed that the ODC belonging to lower cPP groups (below 7) encapsulated most normal functional subtypes previously described, while those functional pathways were largely lost in higher cPP groups (Extended Data Fig. 3c). The nuclei in the two highest cPP groups originated from every aged and PD donor with only 0.4% of nuclei coming from a single young donor.

Due to the widespread presence of these nuclei in both aging cohorts and their reduction of canonical ODC functional gene signatures, we reasoned that these high cPP nuclei may represent a previously undescribed age-associated oligodendrocyte (AAODC). Interestingly, comparison of DEG across each cPP group revealed little change between aged and PD ODC in low cPP groups. However, the high-cPP AAODC showed a dramatic increase in DEG between non-disease controls and PD patients (Fig. 5a), as well as changes in associated gene ontology pathways (Extended Data Fig. 5a). This interesting finding indicates that the majority of ODC in the midbrain may maintain their normal functions during aging and even in disease, reinforcing the importance of subpopulation-focused approaches. Comparison of DEG between aged and PD AAODC found 57 increased genes in aged nuclei and 43 genes in PD (Supplementary Table 5). Among those DEG, 27 genes in aged and 22 in PD overlapped with the 299 genes identified as the high cPP genes in ODC (Fig. 3d, Fig. 5b). We compared both gene ontology and protein-protein interactions (https://string-db.org/) for the cPP-relevant DEG in the AAODC (Extended Data Fig. 5b,c) and noticed that, in particular, genes upregulated in aged control donors compared to PD were largely related to protein folding and chaperone mediated autophagy, such as HSP90AA1, HSPH1, DNAJB1 and ST13 (Fig. 5c). In PD, we noted increases in ferritin pathway genes, such as FTH1, and genes related to aggregated protein response such as CRYAB and CLU (Fig. 5c). Peak-gene association analyses for DEG between aged and PD in AAODC revealed significant changes in gene-peak connections between aged and PD in AAODC (Fig. 5d,e). We examined the transcription factor motifs in these correlated peaks and found that the SOX family motif was enriched in genes that were upregulated in aged nuclei, while the ATF family motif was enriched in PD.

Extended Data Fig. 5
Extended Data Fig. 5 Cellular pathways and GO of DEG between aged and PD AAODC.

a, Enrichment of GO terms increased in PD (top panel) and aged (bottom panel) across cPP group. b, STRING protein-protein interaction network for AAODC genes increased in aged group. Significant GO terms are shown in the lower panel. c, STRING Protein-protein interaction network for AAODC genes increased in PD donors. Significant GO terms are shown in the lower panel.

Fig. 5
Fig. 5 Analysis of AAODC.

a, Differentially expressed genes between aged and PD in ODC increase at high cPP. b, Venn diagram of cPP relevant genes compared with differentially expressed genes between AAODC for aged and PD. c, Dot plot showing differentially expressed genes and percentage expression between AAODC for aged and PD. d,e, Comparison of significant gene-peak connections for AAODC from aged and PD donors. Top panel: Normalized ATAC read count distribution. Bottom panel: arcs show ATAC peaks significantly correlated with gene expression. Arc height represents significance. Significant gene-peak connections for HSP90AA1 and PTGES1 (d), and FTH1 and PDE1A (e). Locations of transcription factor binding motifs for SOX family (d) and ATF family (e) are marked.

Taken together, these results suggest that in the midbrain, a subset of ODC with distinct expression signatures marked by high levels of genes related to protein folding and aggregation appear during the aging process, and this ODC subset loses this crucial signature in PD.

Discussion

In this study, we assessed paired snRNA-seq and snATAC-seq patterns in the SN region of the midbrain of young, aged controls and PD patients to elucidate how aging processes may predispose an individual to PD. By adding data from young control donors to the current single-cell studies of the brain, our study may significantly extend our understanding of the aging process and neurodegeneration at single-cell resolution. Leveraging paired analysis of RNA expression and chromatin accessibility from the same cells provides a unique opportunity to identify cellular pathways and their contribution to disease processes.

We found there are dramatic and widespread changes in links between expressed genes and discreet chromatin regions during aging and PD. Trajectory analysis revealed all cells change over the course of aging, while there appear to be significant further alterations in oligodendrocyte and microglia populations in PD. In particular, we identified an oligodendrocyte subtype with a novel gene expression signature that only appears in aged individuals (AAODC). Comparison of AAODC from aged individuals without neurological disease and PD patients revealed that AAODC may play important roles in resistance to cellular stresses incurred during aging, such as protein folding and chaperone-mediated autophagy, and that this crucial function is lost in PD.

Each brain region has a unique cell-type composition. Our data agrees with previous data11 indicating that ODC is the major cell type in the midbrain, representing approximately 75% of all cells, which has led to increasing recognition of ODC as a potential major player in PD11–14. After ODC, MG, OPC and AS are the next most common cell types; however, we failed to faithfully recover neuronal populations that were reported in previous studies in the midbrain14, potentially due to their increased sensitivity to nuclear permeabilization for ATAC-seq preparation.

We adapted the peak-gene association analysis developed by Ma et al.8 to link distal peaks to genes in cis based on co-variation in chromatin accessibility and gene expression, which is a crucial step forward in unraveling the complex cellular heterogeneity of the brain. The recent work by Morabito et al.33 elegantly combined separate ATAC and RNA expression data to form correlative relationships between two similar nuclei. Peak-gene association studies further advance this technique by analyzing paired data from the same nuclei. The critical importance of this approach is highlighted by our unexpected findings showing surprisingly little variation in ATAC peaks during aging or PD regardless of changes in gene expression. Despite this, peak-gene analysis revealed widespread alterations in gene expression-peak relationships and allowed us to decipher cis-motifs that are associated gene expression in the midbrain during aging and PD pathogenesis. For example, we identified the aging-associated lncRNA, NEAT1, is linked to EGR1/2 motif-containing peaks, a downstream NEAT1 effector34, even though there were no noticeable changes in chromatin accessibility in the regions containing this motif during aging.

A more in-depth study of these cis-regulatory mechanisms underlying functional alterations in cellular subpopulations is necessary to understand fully the enormously complex changes that occur during aging and neurodegeneration.

Adapting pseudotime to infer pathogenic changes is a critical step toward expanding our understanding of how aging may alter the midbrain and how those alterations predispose an individual to PD. Our pseudopathogenesis analysis revealed that the majority of oligodendrocytes in the midbrain maintain their normal functions during aging and PD, but a small subset with distinct signatures appears during aging. AAODC exhibit highly specific changes between control and PD samples, and these changes are not reflected in the broader ODC population. Identification of differential activation of pathways involving protein folding/aggregation and chaperone-mediated autophagy between normal subjects and PD individuals in AAODC adds to the growing body of evidence that oligodendrocytes may play key roles in PD pathogenesis outside of their traditional myelination role. Multiomic analysis also indicates that microglia in the midbrain undergo significant population shifts during aging into PD and show an increase in markers previously associated with aging and neurodegenerative disease. Notably, our analysis did not identify a strong correlation of midbrain astrocytes or OPC with PD development, although this does not preclude them from having important roles in other brain regions in PD.

There is increasing evidence suggesting oligodendrocytes may support neurons outside of the typical myelination functions35–37. Our findings open the possibility that AAODC in particular may play a role in protecting neurons from normal age-induced stressors. Whether the loss of such functions in AAODC leads to dopaminergic degeneration, or these changes are secondary to the disease process, remains to be elucidated. It has been well described that oligodendrocytes are able to take up protein aggregates from the extracellular milieu that are released from neurons38. It is plausible that AAODC could support neuronal health by refolding or degrading misfolded or aggregated proteins during the normal aging process. This is particularly intriguing because misfolded alpha-synuclein plays discrete roles in the pathogenesis of PD and multiple system atrophy (MSA). MSA is characterized by accumulation of alpha-synuclein selectively in ODC39, 40, and although there is no strong evidence of alpha-synuclein accumulation in ODC in PD, there remain significant gaps in our understanding of the complex mechanisms of alpha-synuclein aggregation in different disease conditions. In PD, AAODC with enhanced function in handling misfolded proteins released from neurons may play an important role in supporting neuronal health, although the mechanistic details remain unclear. Further investigation of this subset of ODC in MSA and PD might provide new insights into differences in alpha-synuclein pathology between these two synucleinopathies.

Our findings that ODC and MG are altered in aging and PD add to the body of evidence of the recently highlighted roles of these cells in PD development11–14. Specific PD subpopulations identified in this study may yield mechanistic insight into how these cells are implicated in dopaminergic degeneration and may lead to novel therapeutic targets. In future studies, a more nuanced gradient of aging samples may be required to understand the continuum of cell states and how the populations shift during aging. It will also be interesting to see what findings of similar studies in different brain regions will uncover, as this is essential to building a brain-wide atlas of human aging and neurodegeneration.

Data Availability

Processed data for the samples presented in this study is available in the Gene Expression Omnibus (GEO) database under accession number GSE193688. Raw data will be made available through the dbGap portal, and specific access information will be added once the deposit request has been finalized.

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE193688

Methods

Nuclei isolation and sequencing

Human tissue samples were obtained through the NIH Neurobiobank from the Human Brain and Spinal Fluid Resource Center (VA Greater Los Angeles Healthcare System). To capture the SN region of the midbrain precisely, we dissected approximately 3 x 3 x 5 mm sections using postmortem photos provided (Supplementary Fig.1a), yielding tissue blocks ranging from 50-75 mg. We processed the tissue in a dounce homogenizer and then centrifuged it through 1.5 M sucrose to isolate the nuclei. With this method, we are able to reliably isolate > 1×106 intact nuclei from ∼50 mg postmortem midbrain tissue. Centrifugation through an iodixanol gradient as previous described9 produced high quality nuclei with low mitochondrial contamination, but the actual yield of nuclei was very low (2.5-5×104 nuclei from ∼50 mg postmortem midbrain tissue), so we did not continue with this method as the SN is a small region and tissue amount is limited. We stained the nuclei with 7AAD and sorted to purify nuclei from mitochondria. After sorting, we permeabilized according to current 10XGenomics protocols (Protocol CG000375 Rev A: 10 mMTris-HCl (pH 7.4), 10 mM NaCl, 3 mM MgCl2, 0.1% Tween-20, 0.01% Digitonin, 1% BSA, 1 mM DTT, 1 U/uL RNase Inhibitor) for 2 min on ice). Nuclei were washed once before counting and proceeding with 10XGenomics protocols for transposition, nuclei isolation and barcoding, and library preparation exactly as written (10XGenomics Protocol Chromium Next GEM Single Cell Multiome ATAC + Gene Expression CG000338 Rev A).

Libraries were sequenced at GeneWiz using Illumina Novaseq S4 flowcells. Our actual average sequencing depth was 1.8E8 reads/sample for gene expression libraries and 1.75E8 reads/sample for ATAC libraries. After sequencing, three samples (one each of Young, Aged, and PD) showed strong ATAC results but poor RNA read quality indicative of RNA degradation in the original sample. The nuclei from these three donors were excluded from downstream analysis.

Authors’ note on nuclear isolation

This isolation method was sufficient to obtain high quality nuclei with minimal blebbing and worked well in our initial pilot study. However, after in-depth data analysis, we noted that barcodes for nuclei with neuronal expression patterns consistently showed unusually high RNA read count, suggesting that these nuclei were clumping together and forming multiplets, and thus, these neurons were excluded by our QC filtering (see below). In nuclear isolations this is commonly caused by a partial loss of membrane integrity that allows genomic DNA to leak out and cause nuclei to stick together. We feel it is likely the digitonin permeabilization step may have been too harsh for these samples, and that the larger neuronal nuclei were more sensitive to digitonin. Previous publications have also noted the effect of differential responses in cell types to different isolation protocols6. In the months since we performed the initial isolation, multiple alternative protocols for isolation of frozen post-mortem tissue have been posted and shared online in unofficial forums that omit digitonin for this reason. In future studies, we will also adjust our isolation protocols accordingly.

Quality control of single-nuclei RNA+ATAC multiome

Both single-nuclei RNA and ATAC-seq reads were spontaneously mapped to GRCh38 human reference transcriptome and genome, respectively, by Cell Ranger ARC (v1.0.1) with “--mempercore=8 --localcores=12” (https://support.10xgenomics.com/single-cell-multiome-atac-gex/software/downloads/latest). Quality control (QC) of single-nuclei transcriptome and chromatin profiles was conducted by Seurat (v4.0.0) and Signac (v1.1.1) in R package (v4.0.2), respectively41, 42. Briefly, we removed cells with 1) less than 100 or more than 7,000 detected genes, 2) less than 500 or more than 20,000 reads and 3) more than 5% mitochondria-derived reads (damaged/dead cells or doublet) in the RNA datasets. We also evaluated the doublet frequency by counting cells expressing both STMN2 and AQP4 (Fig S1d), which are usually exclusively expressed in cortical neuron and astrocyte, respectively43, 44. If the ratio of STMN2+AQP4+ cells to total number of STMN2+AQP4-, STMN2-AQP4+ and STMN2+AQP4+ cells was more than 3%, the libraries were excluded from subsequent analyses. For QC of ATAC datasets, we calculated four measurements: 1) ratio of mononucleosomal fragments to nucleosome-free fragments, 2) enrichment score in transcription start site (TSS), 3) total number of ATAC-seq reads per cell and 4) ratio of reads mapped to “blacklist regions” that are genomic regions inducing aberrant mapping and artificial signals45. We filtered out cells with 1) more than 4 mononucleosome:nucleosome-free ratio, 2) less than 2 TSS enrichment score, 3) less than 1000 and more than 60,000 ATAC reads, 4) more than 2% of reads in the blacklist regions and 5) less than 1,000 or more than 25,000 ATAC peaks. Our initial isolation yielded 82735 nuclei. After QC we retained 69,289 high-quality nuclei that met the above criteria of RNA and ATAC for subsequent analyses (Fig. S1d-g).

Integration of single-cell gene expression profiles

Gene expression profiles were harmonized using Seurat (v4.0.0)3 as we have done previously17, 18, 43, 44, 46. In each snRNA-seq library, raw UMI count was normalized to total UMI count. Highly variable genes were then identified by variance stabilizing transformation with 0.3 loess span and automatic setting of clip.max value. The top 2,000 variable genes were used to identify cell pairs anchoring different snRNA-seq libraries using 20 dimensions of canonical correlation analysis (CCA). After scaling gene expression values across all integrated cells, we performed dimensional reduction using principal component analysis (PCA). For the visualization, we further projected single cells into two-dimensional Uniform Manifold Approximation and Projection (UMAP) space from 1st and 20th PCs. Graph-based clustering was then implemented with shared nearest neighbor method from 1st and 20th PCs and 0.8 resolution value. Differentially-expressed genes (DEGs) in each cluster were identified with more than 1.25-fold change and p<0.05 by two-sided unpaired T test. Gene Ontology (GO) analysis was performed to DEGs by GOstats Bioconductor package (v2.56.0)47. False discovery rate was adjusted by p.adjust function in R with “method=”BH”” option.

Cluster annotation

Cell types were assigned in each cluster “island” with uniquely-expressed genes18, 46. First, we assigned neuronal clusters with expression of STMN2 and TBR1. Oligodendrocyte and astrocyte clusters were classified by myelination markers (MBP and MOG) and astrocyte-specific proteins (GFAP and AQP4) expression. OLIG1+OLIG2+ clusters without any neuronal and astrocyte markers were defined as oligodendrocyte precursor cells (OPCs). Endothelial cells were annotated by substantial expression of FLT1, VWF and PDGFRB. We annotated clusters with GPR34, TREM2 and C1QC without any OPC and astrocyte markers. Since T-cell infiltration into substantia nigra is a hallmark of PD pathology, we also defined T-cell clusters with their specific marker expression (CD3E and CD8A).

Integration of single-cell chromatin profiles

Integrative analysis of single-cell chromatin profiles was performed by Seurat (v4.0.0), Signac (v1.1.1), GenomicRanges (v1.42.0) and Harmony (v1.0) R packages41, 42, 48, 49. First, ATAC peaks from all ATAC libraries were merged by reduce function in GenomicRanges. Huge (>10,000bp) and tiny (<20bp) combined peaks were removed from subsequent analyses. ATAC reads were recounted in the combined peaks by using the FeatureMatrix function in Signac. After merging the peak x cell matrices from all ATAC libraries, the ATAC read counts were normalized by term frequency inverse document frequency. Latent sematic indexing (LSI) was then computed from the merged count matrix by the singular value decomposition method. To minimize the batch and technical difference across libraries, we ran the Harmony algorithm using LSI embeddings with “project.dim=F” option. For data visualization, all cells were embedded into two-dimensional UMAP space from 2nd and 30th Harmony-corrected LSI. The 1st index of LSI may represent technical variations, such as sequencing depth and was not used for UMAP projection.

Analysis of differential gene expression and chromatin accessible regions

Gene expression and chromatin accessibility profiles were compared between PD and healthy aged groups, between healthy aged and young groups or across cell types. The differentially-expressed genes and differential open chromatin regions were identified in each cluster by 1.25-fold change and p<0.05 unpaired T test. Cellular events and functions related to the differentially-expressed genes and differential open chromatin regions were analyzed by the enrichment of Gene Ontology using GOstats (v2.56.0) as described above47.

Genomic distribution of all ATAC peaks and differential open chromatin regions was identified by annnotatePeak.pl script in HOMER suite (v4.11.1) with default parameter50. Analysis of previous association of differentially expressed genes was done manually6, 51–74.

Analysis of relationship between gene expression and ATAC peaks

To investigate the relationship between regulatory sites and differential gene expression, Spearman correlation was calculated using the normalized gene expression values and ATAC read counts in each gene-peak pair8. We chose all ATAC peaks within ± 500kbp of TSSs of genes for the correlation analysis. In each peak, we generated 100 permutated background peaks with the same accessibility. Spearman correlation was also calculated to the background peaks to estimate the background correlation distribution. We assumed that the background correlation follows Gaussian distribution. Thus, we determined p-value of the observed Spearman correlation using pnorm function in R with mean and standard deviation of the background correlation distribution. We defined significant gene-peak association with p<0.05.

Subgrouping and pseudopathogenesis analysis in each cell type

In each cell type, individual cells were re-scaled in both RNA and ATA profiles by SCTransform (v0.3.2). Dimensionality reduction was then implemented by PCA and UMAP embedding by Seurat as described above. Subsequently, we performed cell trajectory analysis using Monocle3 (v0.2.3.0) to estimate pathological stages of individual cells32. Briefly, Seurat object was converted into Monocle3 “cell_data_set” format by SeuratWrapper (v0.3.0). Clustering was performed by Leiden community detection. Principal graph of cell trajectory was constructed using learn_graph function with default parameter and used for pseudotime calculation by choosing cells in young group-specific clusters as root cells. In each cell i, two pathological measurements, RNA- and ATAC-derived pseudopathogenesis (PPRNA,i and PPATAC,i), were combined by dividing their variances of all cells (σ2RNA and σ2ATAC) as follows: Embedded Image Here, we called the combined values with corrected errors PseudoPathogenesis (cPPi).

To identify differential genes and peaks along pseudopathogenesis, we calculated Spearman correlation with the scaled gene expression and peak intensity values. We chose genes and peaks with > 0.1 or < − 0.1 for subsequent analyses. Peaks associated with at least one pseudopathogenesis-correlated genes were used for subsequent motif enrichment analysis.

Functional subtype gene modules

To identify potential subsets of nuclei with functional differences, we used literature-derived gene sets that have been previously reported21, 22, 24–26. For each set of genes, we used the AddModuleScore function of Seurat to establish a combined expression score for each nucleus for all the genes in the list, and added it as metadata. The gene sets for each oligodendrocyte module are as follows: newly formed ODC (TCF7L2, CASR, CEMIP2, ITPR2); myelin forming ODC (MAL, MOG, PLP1, OPALIN, SERINC5, CTPS1); mature ODC (KLK6, APOD, SLC5A11, PDE1A); Synaptic Support ODC (NFASC, NRXN3, CNTNAP2, ANK3); ODC-Neuron Adhesion Markers (HAPLN2, STMN1, MAP1B, SEMA5A, EPHB2, S100B, PRKCA). After cPP analysis, two additional modules were identified based on cPP-relevant gene expression and GO analysis: Protein Folding (HSP90AA1, HSPA1A, HSPA1B, CRYAB, FKBP5) and Chaperone Mediated Autophagy (LAMP2, ST13, DNAJB1, STIP1, HSPA8, BAG3). The gene sets for microglia modules are as follows: homeostatic MG (HEXB, CST3, CX3CR1, CTSD, CSF1R, CTSS, SPARC, TMSB4X, P2RY12, C1QA, C1QB); stage 1 TREM2 independent DAM MG (TYROBP, CTSB, APOE, B2M, FTH1); stage 2 TREM2 dependent DAM MG (TREM2, AXL, CST7, CTSL, LPL, CD9, CSF1, ITGAX, CLEC7A, LILRB4, TIMP2); and Aging MG (IL15, CLEC2B, DOCK5). The gene sets for Astrocyte modules are as follows: Disease-Associated AS (GFAP, CSTB, VIM, OSMR, GSN, GGTA1P); A1 Reactive AS (SERPING1, GBP2, FKBP5, PSMB8, SRGN, AMIGO2); A2 Reactive AS (CLCF1, TGM1, PTX3, S100A10, SPHK1, CD109, EMP1, SLC10A6, TM4SF1, CD14), and GFAP-low AS (LUZP2, SLC7A10, MFGE8).

Motif enrichment analysis

Transcription factor binding motifs in ATAC peaks were identified using the HOMER (v4.11) suite. Briefly, we create BED format files for specific sets of ATAC peaks and ran findMotifsGenome.pl script with “-size given” option to evaluate statistical significance of the motif enrichment. -log10(FDR) value was used as a score of motif enrichment in the set of ATAC peaks.

Funding

Study was supported by NIH 1R01-NS100919 and 1R01-NS101461 (to Y.S.K) and by startup funds from Centre de recherche de l’Hôpital Maisonneuve-Rosemont, Université de Montréal (to Y.T), and Fonds de recherche du Québec - Santé (FRQS) (to Y.T), and transition grant from Cole Foundation (to Y.T).

Author contributions

L.A. and Y.K designed the study. Nuclear isolation and library preparation were performed by L.A. and M.S. Computational analysis was performed by Y.T. and L.A. L.A., Y.K., and Y.T. wrote the manuscript; all authors contributed to review and revision of the manuscript.

Human Samples

Human midbrain samples were obtained through NIH NeuroBioBank requests (The Human Brain and Spinal Fluid Resource Center, University of Maryland Brain and Tissue Bank, Harvard Brain Tissue Resource Center, University of Miami Brain Endowment Brank). This project utilized de-identified post-mortem brain samples, so is not considered to meet federal definitions for IRB jurisdiction, and falls outside the purview of the Rutgers IRB committee.

Competing interests

There are no competing interests among the authors.

Data and materials availability

Processed data for the samples presented in this study is available in the Gene Expression Omnibus (GEO) database under accession number GSE193688. Raw data will be made available through the dbGap portal, and specific access information will be added once the deposit request has been finalized.

Acknowledgement

Authors gratefully acknowledge NIH Neurobiobank for providing all postmortem brain samples. Authors also thank Mr. Andrew Knott for language editing.

References

  1. 1.↵
    Parkinson’s Foundation. Understanding Parkinson’s/Statistics. <;https://www.parkinson.org/Understanding-Parkinsons/Statistics>. (2020).
  2. 2.↵
    Pringsheim, T., Jette, N., Frolkis, A. & Steeves, T.D. The prevalence of Parkinson’s disease: a systematic review and meta-analysis. Mov Disord 29, 1583–90 (2014).
    OpenUrlCrossRefPubMed
  3. 3.↵
    Dauer, W. & Przedborski, S. Parkinson’s disease: mechanisms and models. Neuron 39, 889–909 (2003).
    OpenUrlCrossRefPubMedWeb of Science
  4. 4.↵
    Hindle, J.V. Ageing, neurodegeneration and Parkinson’s disease. Age Ageing 39, 156–61 (2010).
    OpenUrlCrossRefPubMedWeb of Science
  5. 5.↵
    Davie, K. et al. A Single-Cell Transcriptome Atlas of the Aging Drosophila Brain. Cell 174, 982–998 e20 (2018).
    OpenUrlCrossRefPubMed
  6. 6.↵
    Ximerakis, M. et al. Single-cell transcriptomic profiling of the aging mouse brain. Nat Neurosci 22, 1696–1708 (2019).
    OpenUrlCrossRef
  7. 7.↵
    Glaab, E. & Schneider, R. Comparative pathway and network analysis of brain transcriptome changes during adult aging and in Parkinson’s disease. Neurobiol Dis 74, 1–13 (2015).
    OpenUrl
  8. 8.↵
    Ma, S. et al. Chromatin Potential Identified by Shared Single-Cell Profiling of RNA and Chromatin. Cell 183, 1103–1116 e20 (2020).
    OpenUrlCrossRef
  9. 9.↵
    Mathys, H. et al. Single-cell transcriptomic analysis of Alzheimer’s disease. Nature 570, 332–337 (2019).
    OpenUrl
  10. 10.↵
    Nott, A. et al. Brain cell type-specific enhancer-promoter interactome maps and disease-risk association. Science 366, 1134–1139 (2019).
    OpenUrlAbstract/FREE Full Text
  11. 11.↵
    Agarwal, D. et al. A single-cell atlas of the human substantia nigra reveals cell-specific pathways associated with neurological disorders. Nat Commun 11, 4183 (2020).
    OpenUrlCrossRefPubMed
  12. 12.
    Bryois, J. et al. Genetic identification of cell types underlying brain complex traits yields insights into the etiology of Parkinson’s disease. Nat Genet 52, 482–493 (2020).
    OpenUrl
  13. 13.
    Reynolds, R.H. et al. Moving beyond neurons: the role of cell type-specific gene regulation in Parkinson’s disease heritability. NPJ Parkinsons Dis 5, 6 (2019).
  14. 14.↵
    Smajic, S. et al. Single-cell sequencing of human midbrain reveals glial activation and a Parkinson-specific neuronal state. Brain (2021).
  15. 15.↵
    Braak, H. & Del Tredici, K. Poor and protracted myelination as a contributory factor to neurodegenerative disorders. Neurobiol Aging 25, 19–23 (2004).
    OpenUrlCrossRefPubMedWeb of Science
  16. 16.↵
    Hao, Y. et al. Integrated analysis of multimodal single-cell data. Cell 184, 3573–3587 e29 (2021).
    OpenUrlCrossRefPubMed
  17. 17.↵
    Cakir, B. et al. Engineering of human brain organoids with a functional vascular-like system. Nat Methods 16, 1169–1175 (2019).
    OpenUrl
  18. 18.↵
    Tanaka, Y., Cakir, B., Xiang, Y., Sullivan, G.J. & Park, I.H. Synthetic Analyses of Single-Cell Transcriptomes from Multiple Brain Organoids and Fetal Brain. Cell Rep 30, 1682–1689 e3 (2020).
    OpenUrlCrossRef
  19. 19.↵
    Simpson, E.H. The Interpretation of Interaction in Contingency Tables. Journal of the Royal Statistical Society Series B-Statistical Methodology 13, 238–241 (1951).
    OpenUrl
  20. 20.↵
    Trapnell, C. Defining cell types and states with single-cell genomics. Genome Res 25, 1491–8 (2015).
    OpenUrlAbstract/FREE Full Text
  21. 21.↵
    Galatro, T.F. et al. Transcriptomic analysis of purified human cortical microglia reveals age-associated changes. Nat Neurosci 20, 1162–1171 (2017).
    OpenUrlCrossRefPubMed
  22. 22.↵
    Habib, N. et al. Disease-associated astrocytes in Alzheimer’s disease and aging. Nat Neurosci 23, 701–706 (2020).
    OpenUrl
  23. 23.↵
    Jakel, S. et al. Altered human oligodendrocyte heterogeneity in multiple sclerosis. Nature 566, 543–547 (2019).
    OpenUrlCrossRefPubMed
  24. 24.↵
    Keren-Shaul, H. et al. A Unique Microglia Type Associated with Restricting Development of Alzheimer’s Disease. Cell 169, 1276–1290 e17 (2017).
    OpenUrlCrossRefPubMed
  25. 25.↵
    Liddelow, S.A. et al. Neurotoxic reactive astrocytes are induced by activated microglia. Nature 541, 481–487 (2017).
    OpenUrlCrossRefPubMed
  26. 26.↵
    Marques, S. et al. Oligodendrocyte heterogeneity in the mouse juvenile and adult central nervous system. Science 352, 1326–1329 (2016).
    OpenUrlAbstract/FREE Full Text
  27. 27.↵
    Borras, C. et al. RasGrf1 deficiency delays aging in mice. Aging (Albany NY) 3, 262–76 (2011).
    OpenUrl
  28. 28.↵
    Salmon, A.B. et al. Effects of transgenic methionine sulfoxide reductase A (MsrA) expression on lifespan and age-dependent changes in metabolic function in mice. Redox Biol 10, 251–256 (2016).
    OpenUrl
  29. 29.↵
    Duclot, F. & Kabbaj, M. The Role of Early Growth Response 1 (EGR1) in Brain Plasticity and Neuropsychiatric Disorders. Front Behav Neurosci 11, 35 (2017).
  30. 30.↵
    Kim, E.J., Ables, J.L., Dickel, L.K., Eisch, A.J. & Johnson, J.E. Ascl1 (Mash1) defines cells with long-term neurogenic potential in subgranular and subventricular zones in adult mouse brain. PLoS One 6, e18472 (2011).
  31. 31.↵
    Sandberg, M., Patil, J., D’Angelo, B., Weber, S.G. & Mallard, C. NRF2-regulation in brain health and disease: implication of cerebral inflammation. Neuropharmacology 79, 298–306 (2014).
    OpenUrlCrossRefPubMed
  32. 32.↵
    Cao, J. et al. The single-cell transcriptional landscape of mammalian organogenesis. Nature 566, 496–502 (2019).
    OpenUrlCrossRefPubMed
  33. 33.↵
    Morabito, S. et al. Single-nucleus chromatin accessibility and transcriptomic characterization of Alzheimer’s disease. Nat Genet 53, 1143–1155 (2021).
    OpenUrl
  34. 34.↵
    Pistoni, M. et al. Long Noncoding RNA NEAT1 Acts as a Molecular Switch for BRD4 Transcriptional Activity and Mediates Repression of BRD4/WDR5 Target Genes. Mol Cancer Res 19, 799–811 (2021).
    OpenUrlAbstract/FREE Full Text
  35. 35.↵
    Bernstein, H.G., Keilhoff, G., Dobrowolny, H., Guest, P.C. & Steiner, J. Perineuronal oligodendrocytes in health and disease: the journey so far. Rev Neurosci 31, 89–99 (2019).
    OpenUrl
  36. 36.
    Philips, T. & Rothstein, J.D. Oligodendroglia: metabolic supporters of neurons. J Clin Invest 127, 3271–3280 (2017).
    OpenUrlCrossRefPubMed
  37. 37.↵
    Simons, M. & Nave, K.A. Oligodendrocytes: Myelination and Axonal Support. Cold Spring Harb Perspect Biol 8, a020479 (2015).
  38. 38.↵
    Pukass, K. & Richter-Landsberg, C. Oxidative stress promotes uptake, accumulation, and oligomerization of extracellular alpha-synuclein in oligodendrocytes. J Mol Neurosci 52, 339–52 (2014).
    OpenUrl
  39. 39.↵
    Peng, C. et al. Cellular milieu imparts distinct pathological alpha-synuclein strains in alpha-synucleinopathies. Nature 557, 558–563 (2018).
    OpenUrlCrossRefPubMed
  40. 40.↵
    Tu, P.H. et al. Glial cytoplasmic inclusions in white matter oligodendrocytes of multiple system atrophy brains contain insoluble alpha-synuclein. Ann Neurol 44, 415–22 (1998).
    OpenUrlCrossRefPubMedWeb of Science

Methods References

  1. 41.↵
    Stuart, T. et al. Comprehensive Integration of Single-Cell Data. Cell 177, 1888–1902.e21 (2019).
    OpenUrlCrossRefPubMed
  2. 42.↵
    Stuart, T., Srivastava, A., Lareau, C. & Satija, R. Multimodal single-cell chromatin analysis with Signac. BioRxiv (2020).
  3. 43.↵
    Xiang, Y. et al. hESC-Derived Thalamic Organoids Form Reciprocal Projections When Fused with Cortical Organoids. Cell Stem Cell 24, 487–497.e7 (2019).
    OpenUrlCrossRefPubMed
  4. 44.↵
    Xiang, Y. et al. Fusion of Regionally Specified hPSC-Derived Organoids Models Human Brain Development and Interneuron Migration. Cell Stem Cell 21, 383–398.e7 (2017).
    OpenUrlCrossRefPubMed
  5. 45.↵
    Amemiya, H.M., Kundaje, A. & Boyle, A.P. The ENCODE Blacklist: Identification of Problematic Regions of the Genome. Sci Rep 9, 9354 (2019).
    OpenUrlCrossRefPubMed
  6. 46.↵
    Xiang, Y. et al. Dysregulation of BRD4 Function Underlies the Functional Abnormalities of MeCP2 Mutant Neurons. Mol Cell 79, 84–98.e9 (2020).
    OpenUrl
  7. 47.↵
    Falcon, S. & Gentleman, R. Using GOstats to test gene lists for GO term association. Bioinformatics 23, 257–8 (2007).
    OpenUrlCrossRefPubMedWeb of Science
  8. 48.↵
    Korsunsky, I. et al. Fast, sensitive and accurate integration of single-cell data with Harmony. Nat Methods 16, 1289–1296 (2019).
    OpenUrlPubMed
  9. 49.↵
    Lawrence, M. et al. Software for computing and annotating genomic ranges. PLoS Comput Biol 9, e1003118 (2013).
  10. 50.↵
    Heinz, S. et al. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities. Mol Cell 38, 576–89 (2010).
    OpenUrlCrossRefPubMedWeb of Science
  11. 51.↵
    Al-Dalahmah, O. et al. Single-nucleus RNA-seq identifies Huntington disease astrocyte states. Acta Neuropathol Commun 8, 19 (2020).
  12. 52.
    Baldo, B. et al. A screen for enhancers of clearance identifies huntingtin as a heat shock protein 90 (Hsp90) client protein. J Biol Chem 287, 1406–14 (2012).
    OpenUrlAbstract/FREE Full Text
  13. 53.
    Barry, G. et al. The long non-coding RNA NEAT1 is responsive to neuronal activity and is associated with hyperexcitability states. Sci Rep 7, 40127 (2017).
  14. 54.
    Barry, G., Guennewig, B., Fung, S., Kaczorowski, D. & Weickert, C.S. Long Non-Coding RNA Expression during Aging in the Human Subependymal Zone. Front Neurol 6, 45 (2015).
  15. 55.
    Binder, E.B. et al. Polymorphisms in FKBP5 are associated with increased recurrence of depressive episodes and rapid response to antidepressant treatment. Nat Genet 36, 1319–25 (2004).
    OpenUrlCrossRefPubMedWeb of Science
  16. 56.
    Blair, L.J. et al. Accelerated neurodegeneration through chaperone-mediated oligomerization of tau. J Clin Invest 123, 4158–69 (2013).
    OpenUrlCrossRefPubMedWeb of Science
  17. 57.
    Chang, L.C., et al. A conserved BDNF, glutamate- and GABA-enriched gene module related to human depression identified by coexpression meta-analysis and DNA variant genome-wide association studies. PLoS One 9, e90980 (2014).
  18. 58.
    Evers, M., Cunningham-Rundles, C. & Hollander, E. Heat shock protein 90 antibodies in autism. Mol Psychiatry 7 Suppl 2, S26–8 (2002).
    OpenUrlCrossRefPubMed
  19. 59.
    Karis, K. et al. Altered Expression Profile of IgLON Family of Neural Cell Adhesion Molecules in the Dorsolateral Prefrontal Cortex of Schizophrenic Patients. Front Mol Neurosci 11, 8 (2018).
  20. 60.
    Katsel, P. et al. The expression of long noncoding RNA NEAT1 is reduced in schizophrenia and modulates oligodendrocytes transcription. NPJ Schizophr 5, 3 (2019).
  21. 61.
    Kim, J.J. et al. Identification of antibodies to heat shock proteins 90 kDa and 70 kDa in patients with schizophrenia. Schizophr Res 52, 127–35 (2001).
    OpenUrlCrossRefPubMedWeb of Science
  22. 62.
    Labadorf, A. et al. RNA Sequence Analysis of Human Huntington Disease Brain Reveals an Extensive Increase in Inflammatory and Developmental Gene Expression. PLoS One 10, e0143563 (2015).
  23. 63.
    Lee, C.H. et al. Heat shock protein 90 and its cochaperone, p23, are markedly increased in the aged gerbil hippocampus. Exp Gerontol 46, 768–72 (2011).
    OpenUrlPubMed
  24. 64.
    Liu, Y. & Lu, Z. Long non-coding RNA NEAT1 mediates the toxic of Parkinson’s disease induced by MPTP/MPP+ via regulation of gene expression. Clin Exp Pharmacol Physiol 45, 841–848 (2018).
    OpenUrl
  25. 65.
    Maccarrone, G. et al. Psychiatric patient stratification using biosignatures based on cerebrospinal fluid protein expression clusters. J Psychiatr Res 47, 1572–80 (2013).
    OpenUrl
  26. 66.
    Mariani, E. et al. Meta-Analysis of Parkinson’s Disease Transcriptome Data Using TRAM Software: Whole Substantia Nigra Tissue and Single Dopamine Neuron Differential Gene Expression. PLoS One 11, e0161567 (2016).
  27. 67.
    Patel, H. et al. Transcriptomic analysis of probable asymptomatic and symptomatic alzheimer brains. Brain Behav Immun 80, 644–656 (2019).
    OpenUrl
  28. 68.
    Patel, N., Crider, A., Pandya, C.D., Ahmed, A.O. & Pillai, A. Altered mRNA Levels of Glucocorticoid Receptor, Mineralocorticoid Receptor, and Co-Chaperones (FKBP5 and PTGES3) in the Middle Frontal Gyrus of Autism Spectrum Disorder Subjects. Mol Neurobiol 53, 2090–9 (2016).
    OpenUrl
  29. 69.
    Sha, L. et al. Pharmacologic inhibition of Hsp90 to prevent GLT-1 degradation as an effective therapy for epilepsy. J Exp Med 214, 547–563 (2017).
    OpenUrlAbstract/FREE Full Text
  30. 70.
    Sinclair, D., Fillman, S.G., Webster, M.J. & Weickert, C.S. Dysregulation of glucocorticoid receptor co-factors FKBP5, BAG1 and PTGES3 in prefrontal cortex in psychotic illness. Sci Rep 3, 3539 (2013).
    OpenUrlPubMed
  31. 71.
    Sunwoo, J.S. et al. Altered Expression of the Long Noncoding RNA NEAT1 in Huntington’s Disease. Mol Neurobiol 54, 1577–1586 (2017).
    OpenUrlCrossRefPubMed
  32. 72.
    Uryu, K. et al. Convergence of heat shock protein 90 with ubiquitin in filamentous alpha-synuclein inclusions of alpha-synucleinopathies. Am J Pathol 168, 947–61 (2006).
    OpenUrlCrossRefPubMedWeb of Science
  33. 73.
    Woo, J.A. et al. Enhanced tau pathology via RanBP9 and Hsp90/Hsc70 chaperone complexes. Hum Mol Genet 26, 3973–3988 (2017).
    OpenUrlCrossRef
  34. 74.↵
    Zhao, M.Y. et al. The long-non-coding RNA NEAT1 is a novel target for Alzheimer’s disease progression via miR-124/BACE1 axis. Neurol Res 41, 489–497 (2019).
    OpenUrl
Back to top
Posted January 24, 2022.
Download PDF

Supplementary Material

Data/Code
Email
Single-nuclei paired multiomic analysis of young, aged, and Parkinson’s disease human midbrain reveals age-associated glial changes and their contribution to Parkinson’s disease
Levi Adams, Min Kyung Song, Yoshiaki Tanaka, Yoon-seong Kim
medRxiv 2022.01.18.22269350; doi: https://doi.org/10.1101/2022.01.18.22269350
Twitter logo Facebook logo LinkedIn logo Mendeley logo
Citation Tools

Subject Area