Parkinson’s disease (PD) is the second most common neurodegenerative disease, presenting with motor system disfunction including slow movements, tremor, gait and balance disturbances, as well as non-motor system symptoms including loss of smell, cognitive decline, and dementia [1]. PD affects approximately 1% of the population above the age of 60 [2] and displays heritability estimates of approximately 30% [3, 4]. Individuals carrying known PD-associated genetic variants can display a wide range of disease presentations including variable age of onset, symptoms, and progression rates. In addition to genetic risk factors, specific environmental exposures such as organochlorine pesticides, including rotenone and paraquat [5], cocaine and amphetamines [6], and 1-methyl-4-phenyl-1,2,3,6-tetrahydropyridine [7], have been identified as major environmental contributors to parkinsonism.

The presence of Lewy bodies and the progressive loss of dopaminergic neurons (DaN) in the substantia nigra pars compacta (SNpc) are the main pathological features of PD [8, 9]. Additionally, multiple biological pathways encompassing mitochondrial dysfunction, neuroinflammation, and protein turnover have been implicated in the disease process and are believed to intersect with the accumulation of α-synuclein—the primary component of Lewy bodies—and the elevated degeneration susceptibility of DaNs in the SNpc. While the discovery of the loss of DaNs in PD has led to attempted treatments such as dopamine substitution and deep brain stimulation, the disease remains incurable, as available therapeutics can only target specific symptoms or minimally slow disease progression [10]. Indeed, elucidating the complex interactions between perturbed pathways and the neuroanatomical changes seen in the parkinsonian brain will advance the development of therapeutics in PD.

One promising approach for characterizing pathological interactions in PD is transcriptome analysis, which allows for global differential gene expression profiling. By comparing samples from individuals with disease to samples from healthy controls, we can learn how specific disease events influence the dynamic nature of gene expression or perhaps how these events lead to pathogenesis [11]. After identifying statistically significant differentially expressed genes or transcripts, we can apply gene set enrichment analysis (GSEA) to learn whether the transcriptional products are implicated in common biological pathways. Understanding how these pathways become disrupted, whether via hypo- or hyper-activity, can guide us in implementing the appropriate therapeutic strategy for halting disease progression [12].

In the early 1990s, the ability to interrogate RNA sequences on a large scale was facilitated by the sequential introductions of the expressed sequence tag (EST) and the serial analysis of gene expression (SAGE) methods [13]. However, microarray (hybridization-based) transcriptome profiling emerged in 1995 and superseded EST and SAGE given its rapidity and affordability, allowing researchers to employ large-scale studies to generate quantitative gene expression changes in disease states [14]. Despite advancements, microarray-based expression analyses are dependent on cross-hybridization, and their potential is restricted by a requirement for a priori knowledge of genetic sequences and a limited capacity to detect and quantify low- and high-abundance transcripts [15]. Moreover, sub-optimally designed probes and incorrect probe annotations reduce reproducibility and cross-platform consistency [16].

Next-generation RNA sequencing (RNA-seq) is a more efficient, less-expensive technology capable of transcriptome-wide analyses [11]. RNA-seq approaches begin with a cell suspension comprised of dissociated cells from the tissue of interest. Entire RNA molecules are used to generate complementary DNA (cDNA) sequences, which are fragmented and sequenced by massively parallel deep sequencing. The sequences are then aligned to a reference genome to produce a tissue-wide gene expression profile (Fig. 1). Referred to as total or bulk RNA-seq, this method involves RNA being extracted from all cells within a tissue (typically blood, skin, or muscle) and mixed together, resulting in an average expressional profile [17]. However, in the absence of single-cell isolation, an inherent limitation of bulk RNA-seq stems from sequencing tissues comprised of many cell types. As such, biological interpretations are impeded by the lack of cell-type specific gene expression.

Bulk and single-cell/single-nucleus RNA sequencing workflows. Samples are acquired from disease-relevant post-mortem tissue. Single-cell/nuclei isolation in the single-cell and single-nucleus RNA sequencing (scnRNA-seq) workflow is achieved by first dissociating the tissue, followed by mRNA barcoding via microfluidics, bead vortexing, SPLiT-seq, laser capture microdissection (LCM), or fluorescence activated cell sorting (FACS) to tag the transcripts from each cell with cell-specific barcodes. In bulk RNA-seq, RNA must be extracted from the cells or tissue to be analyzed prior to library preparation. Complementary DNA (cDNA) library synthesis allows RNA-seq experiments to be carried out on technically mature commercial instruments designed for DNA-based sequencing. Sequencing on next-generation sequencing platforms provides “reads,” defined as strings of sequence data. Quality control of sequencing data ensures the absence of intrinsic biases associated with library preparation, including nucleotide composition bias and GC bias. The number of reads mapped to each gene is used to generate gene expression profiles for each sample or cell for bulk and scnRNA-seq, respectively. Computational cell-type clustering in the scnRNA-seq workflow groups cells with similar expression profiles to identify distinct cellular species. The resulting expression data from bulk RNA-seq and scnRNA-seq can be used for tissue-level or cell-level biological analyses, respectively

Single-cell RNA-seq (scRNA-seq) was introduced in 2009 by Tang et al. to capture an entire transcriptome library of one cell [18]. Further, in 2011, Islam et al. applied scRNA-seq to highly multiplexed samples comprising 96 cells to generate large-scale single-cell expression profiles [19]. In 2015, Drop-seq introduced the first scalable method to obtain single-cell transcriptomes on thousands of cells by separating individual cells into aqueous droplets, tagging each cell’s RNA with a unique molecular barcode, and sequencing all cells simultaneously [20]. In addition to microfluidic technologies (such as Drop-seq), bead vortexing, in situ barcoding (SPLiT-seq), laser capture microdissection (LCM), and fluorescent activated cell sorting (FACS) can also be used to isolate single cells prior to sequencing. By sequencing individual cells, scRNA-seq mitigates the limitations of bulk RNA-seq by providing transcriptional data at a single-cell resolution, capturing cell-to-cell variability in gene expression, and revealing the cell-type heterogeneity of a tissue. Accordingly, scRNA-seq is particularly important for studying neurodegenerative disease as it facilitates the exploration of specific cell types that are most vulnerable to degeneration and can inform mechanisms by which cell populations interact to promote disease pathology.

One limitation of scRNA-seq is that it may provide a biased representation of captured cell types, as certain cellular species are more vulnerable to the dissociation process [21]. For example, neuronal cells from the human neocortex dissociate poorly compared to non-neuronal cell types and are therefore underrepresented in cell suspensions [22]. One way to overcome this limitation is to harness single-nucleus RNA-seq (snRNA-seq), which is similar to scRNA-seq as it is used to sequence transcriptomes of individual cells; however, only nuclear RNA is captured. Given that nuclei are more resistant to mechanical stress, snRNA-seq ensures that vulnerable cellular species are adequately represented in the final dataset, and, unlike single cells, they can be isolated from frozen post-mortem tissues [23,24,25]. Furthermore, it has been shown that nuclei yield expression profiles resembling entire cells, thus making snRNA-seq an adequate substitution for scRNA-seq in most applications [24].

scRNA-seq and snRNA-seq (scnRNA-seq) workflows resemble bulk RNA-seq; however, single cells or nuclei must be isolated by microfluidics, bead vortexing, in situ barcoding, FACS, or LCM prior to cDNA synthesis and sequencing. After cells or nuclei have been sequenced, the results are embedded onto gene expression space and computationally clustered to identify transcriptionally unique populations that can be annotated based on the expression of established marker genes, enabling cell-type specific analyses downstream (Fig. 1). Importantly, the current application of scnRNA-seq assays to cohort-scale analyses is restricted by relatively high monetary cost and stringent single-cell isolation processes.

In this review, we leverage the transcriptional landscape of the Parkinsonian brain, which has been thoroughly investigated using both bulk and single-cell approaches, as a case study for evaluating the utility of bulk and scnRNA-seq in neurodegeneration. By comparing and integrating the findings derived from short-read next-generation sequencing methods, we evaluate whether transcriptomic analyses at different resolutions provide compatible results and demonstrate the complementarity of both assays. Finally, we evaluate the current feasibility of bulk and scnRNA-seq methods through the lens of the PD transcriptomic literature to illustrate the necessity of both technologies for achieving a holistic insight into the mechanism by which gene expression promotes neurodegenerative disease.

Herein, we compare and integrate findings from bulk and scnRNA-seq studies that have explored the post-mortem frontal cortex, prefrontal cortex, midbrain, substantia nigra (SN), ventral tegmental area (VTA), substantia innominata, and hypothalamus in the context of PD (Fig. 2). We focus on the loss of DaNs in the SNpc and highlight key pathways that overlap across RNA-seq studies including synaptic dysfunction, inflammation, mitochondrial dysfunction, and the unfolded protein response (UPR). A comprehensive analysis of these RNA-seq studies is shown in Table 1. To the best of our knowledge, we have included all RNA-seq analyses of human post-mortem brain tissue focussed on PD.

RNA sequencing (RNA-seq) studies on Parkinson’s disease (PD). RNA-seq studies performed on human post-mortem brain tissue to study PD, which were reviewed herein. Studies are categorized by RNA-seq technology, and brain tissue. Abbreviations: scnRNA-seq, single-cell and single-nucleus RNA sequencing

The progressive loss of DaNs in the SNpc is one of the main pathological features of PD [9]; 50% of DaNs were found to be lost through histological analysis of post-mortem tissue [26]. Interestingly, some DaNs in the SNpc appear resistant to degeneration and survive into the late stages of the disease, while DaNs in the VTA (adjacent to the substantia nigra) of PD patients show a much lower magnitude of degeneration compared to those in the SNpc [27,28,29,30,31,32]. Yet, the molecular basis leading to the vulnerability of a subset of DaNs in the SNpc has remained elusive.

The resolution of bulk RNA-seq is limited to providing a whole-tissue-level understanding of the transcriptome, permitting the comparison of SNpc to VTA but not subtypes within one region. One exciting application for scnRNA-seq assays in PD is the ability to characterize the transcriptional landscape of the distinct cellular subgroups of DaNs that are particularly vulnerable or resistant to degeneration within the SNpc. Transcriptional data at a single-cell resolution can then inform disease-related changes in cell-type proportions and the presence of transcriptionally distinct, disease-specific cell types. For example, upon selective enrichment of SNpc samples for neurons, snRNA-seq corroborated that DaNs from the SNpc showed the largest decline in PD samples compared to non-PD controls by identifying one highly vulnerable, transcriptionally distinct DaN subtype that was selectively depleted in PD samples and two transcriptionally distinct DaN subtypes that were selectively enriched in PD samples [33], implying that these subtypes are resistant to PD-related cell death. Consistent with patterns of DaN loss in PD [34], single-molecule fluorescence in situ hybridization (FISH) revealed that the highly vulnerable DaN subtype was confined to the ventral tier of the SNpc and that the two resistant DaN subtypes were confined to the dorsal tier [33]. Importantly, it should be noted that distinct cellular populations identified by scnRNA-seq warrant further validation using different modalities, such as FISH, histology, or spatial multiomics to support the findings.

Although the molecular basis responsible for the selective vulnerability of a subset of DaNs in the PD SNpc remains elusive, snRNA-seq of post-mortem midbrain tissues has provided notable insights. The vulnerable, ventrally located DaN subtype—selectively depleted in PD samples—from the SNpc were most strongly enriched for the expression of PD risk genes identified by genome-wide association studies (GWAS), including SNCA, MAPT, GAK, WNT3, and IGSF9B [33, 35], suggesting that the presence of these genetic risk factors, and possibly their altered activity in these vulnerable neurons, comprises an inherent risk that influences their survival in PD. These same DaNs were enriched for pathways related to regulation of neuron death and WNT signaling [33]. Similarly, in a separate study, a transcriptionally distinct neuronal subtype was predominantly found in the midbrain of PD samples compared to non-PD controls, which showed elevated expression of CADPS2 and TIAM1 [36]. CADPS2 has been previously linked to catecholamine uptake, while TIAM1 is involved in WNT/DVl/RAC1 signalling [37]. The WNT signalling pathway is critical for development and maintenance of DaNs [38]; thus, the convergence of both snRNA-seq analyses for aberrant WNT signalling in vulnerable DaNs emphasizes this pathway as an important contributor to the selective degeneration of neuronal cells in the PD midbrain. Indeed, there is evidence showing that the loss of function of parkin, a PD-associated protein, can lead to increased WNT signaling and consequent post-mitotic DaN death [39]. Despite the apparent importance of the WNT signaling pathway in the PD disease process, bulk RNA-seq analyses of post-mortem PD brains have yet to find evidence for this pathological feature, likely owing to the relative sparsity of DaNs that results in expressional signals being diminished upon averaging among neighboring cells comprising the tissue.

Altogether, these findings demonstrate the value of scnRNA-seq for quantifying disease-related changes in cell-type abundances, characterizing cellular heterogeneity, and nominating aberrant pathways in disease pathology that are restricted to sparse cellular species. In contrast to scnRNA-seq, it is unsurprising that bulk RNA-seq is insufficient for detecting transcriptional changes associated with rare cellular species; however, this does not discount the validity of the assay since the tissue-level resolution is critical for detecting the most dominant pathways in neurodegenerative disease, as discussed below.

Beyond the loss of DaNs in the SNpc, intracytoplasmic inclusions known as Lewy bodies are a hallmark of PD pathology [40]. Αlpha-synuclein is the major protein component of Lewy bodies and forms aggregates in a prion-like manner in the disease state [41]. In PD, α-synuclein aggregates are associated with compromised neurotransmission, correlated with greater cognitive decline, and are believed to be a primary cause of DaN degeneration [40, 42,43,44]. In a normal physiological state, α-synuclein regulates neurotransmitter release, synaptic function, and plasticity, making the synapse a critical area in PD [44]. Furthermore, synaptic dysfunction related to α-synuclein aggregates may lead to PD [40]. Synaptic dysfunction in PD comprises a valuable case study for demonstrating the compatibility of bulk and scnRNA-seq, whereby transcriptional analyses at a bulk-tissue resolution are sufficient to detect the dominant pathway and define important avenues for further investigation, while scnRNA-seq can inform the cell types likely responsible for the signals observed at a bulk resolution.

Two bulk RNA-seq analyses of post-mortem midbrain tissue [45, 46] have reported changes in pathways associated with synaptic dysfunction in PD. In turn, scnRNA-seq has made multiple important contributions that corroborate the findings from bulk tissue. First, scRNA-seq analysis of LCM DaNs from the SNpc and VTA from non-PD controls revealed upregulation of pathways promoting synapse integrity in the VTA [47], suggesting that the maintenance of synaptic function is integral to the elevated resistance of DaNs in the VTA compared to DaNs in the SNpc. Second, an analysis of oligodendrocytes from the PD midbrain identified differentially expressed genes (DEG) in synaptic transmission pathways compared to non-PD controls, which coincided with a reduced fraction of myelinating oligodendrocytes in PD samples [36]. Third, an oligodendrocyte subtype, whose transcripts were enriched for processes linked to synapse organization, was found to be depleted in PD samples compared to non-PD controls [48]. Indeed, oligodendrocytes are critical in the maintenance of synaptic function in a healthy state [49], yet the snRNA-seq data suggest that a subset of oligodendrocytes lose their preservative functions in PD pathology. Interestingly, snRNA-seq revealed that oligodendrocyte and oligodendrocyte precursor cells from PD samples expressed genes that were significantly enriched within PD GWAS risk loci [50]; thus, it is possible that these genetic risk factors pose inherent risks that influence the actions of oligodendrocytes in PD pathology. It must be noted, however, that these GWAS enrichment findings derive from a pre-print manuscript which has not yet been peer-reviewed and is subject to changes following potential updates to the analyses. Nonetheless, upon amalgamating the RNA-seq evidence for synaptic dysfunction in the PD brain, it is tempting to speculate that oligodendrocyte dynamics may contribute to the elevated susceptibility of DaNs located in the SNpc. Yet, it remains to be elucidated whether a mere depletion of oligodendrocytes, a transition towards a pathogenic state, or a combination of both possibly contribute to the loss of DaNs in the PD brain.

In the prefrontal cortex, bulk RNA-seq analyses identified transcriptomic signatures related to synaptic processes in PD samples; however, these findings were attenuated upon correction for cellularity [51]. Nonetheless, snRNA-seq revealed that neuronal cells in the PD prefrontal cortex may undergo compensatory synaptogenesis as the most significantly upregulated pathways were associated with synaptic assembly [52]. While, these observations were unique to a pre-print manuscript that has not undergone peer review, the findings at a bulk resolution provide support to the conclusions. Furthermore, expressional changes related to postsynaptic pathways in excitatory neurons showed the strongest correlation with Lewy pathology score—measured by phospho-S129-α-synuclein immunohistochemistry [52]. These findings reflect the idea that alterations to the neuronal architecture in the prefrontal cortex and the accumulation of α-synuclein aggregates are associated with disrupted synaptic processes in PD pathogenesis [40]. Indeed, the convergence on synaptic-related pathways by both technologies supports its pathological role in the PD prefrontal cortex but raises a limitation with applying cell-type corrections to bulk RNA-seq data as it may discard relevant biological signals. Accordingly, the importance of validating the findings from bulk RNA-seq with the improved resolution of scnRNA-seq emerges, especially when attempting to conduct analytical procedures that address the cellularity of the transcriptional profiles.

Although neuroinflammation in the pathophysiology of PD has garnered support through epidemiological and genetic studies, many questions concerning its role in the PD disease process remain unsolved [53]. Among them is what cell types contribute to the innate immune response [53], which is primed for investigation by scnRNA-seq assays given the single-cell resolution. Further, a major limitation of bulk RNA-seq assays for investigating neuroinflammation in neurodegenerative disease is establishing whether altered expression derives from changes in cell-type proportion, regulatory events that affect gene expression, or a combination of both. While cell-type deconvolution tools may help elucidate this ambiguity, issues surrounding the accuracy of current deconvolution algorithms have been raised [54]; thus, these findings often warrant further validation by different modalities, including scnRNA-seq, which can be used to quantify both changes in cellular composition and gene expression between phenotypes.

Bulk RNA-seq analyses of the midbrain [45], prefrontal cortex [51], and frontal cortex [55] have all identified transcriptional evidence for inflammatory pathways in PD samples. The strongest signal was derived from the midbrain, as bulk RNA-seq disclosed that approximately half of the DEGs between PD samples and non-PD controls were related to the immune-inflammatory response [45]. In this analysis, cell-type deconvolution revealed equal distributions of cellular species according to disease status, suggesting that the immune-inflammatory signal resulted from regulatory events affecting gene expression, as opposed to changes in the number of cells expressing the genes. In contrast, an activated inflammatory response in the PD prefrontal cortex was inferred by bulk RNA-seq through significant increases in microglia and oligodendrocyte marker gene profiles compared to non-PD controls, without any differential gene expression to suggest an elevated inflammatory response in PD [51]. In the frontal cortex, Cappelletti et al. leveraged bulk RNA-seq data from 84 samples, allowing them to maintain power upon stratifying the PD samples into three groups according to their degree of Lewy body pathology—based on the Braak Lewy body stage—which was used as a proxy for disease progression [