Front. Neurosci.Frontiers in NeuroscienceFront. Neurosci.1662-453XFrontiers Media S.A.10.3389/fnins.2021.733501NeuroscienceOriginal ResearchTranscriptomic Signatures Associated With Regional Cortical Thickness Changes in Parkinson’s DiseaseKeoArlin12DzyubachykOleh3van der GrondJeroen3van HiltenJacobus J.4ReindersMarcel J. T.125*†MahfouzAhmed125*†1Leiden Computational Biology Center, Leiden University Medical Center, Leiden, Netherlands2Delft Bioinformatics Lab, Delft University of Technology, Delft, Netherlands3Department of Radiology, Leiden University Medical Center, Leiden, Netherlands4Department of Neurology, Leiden University Medical Center, Leiden, Netherlands5Department of Human Genetics, Leiden University Medical Center, Leiden, Netherlands
Edited by: Jiajia Zhu, First Affiliated Hospital of Anhui Medical University, China
Reviewed by: Jolanta Dorszewska, Poznan University of Medical Sciences, Poland; Huaigui Liu, Tianjin Medical University General Hospital, China; Roger Jacob Mullins, National Institutes of Health (NIH), United States
*Correspondence: Marcel J. T. Reinders, m.j.t.reinders@tudelft.nlAhmed Mahfouz, a.mahfouz@lumc.nl
†These authors have contributed equally to this work and share last authorship
This article was submitted to Brain Imaging Methods, a section of the journal Frontiers in Neuroscience
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
Cortical atrophy is a common manifestation in Parkinson’s disease (PD), particularly in advanced stages of the disease. To elucidate the molecular underpinnings of cortical thickness changes in PD, we performed an integrated analysis of brain-wide healthy transcriptomic data from the Allen Human Brain Atlas and patterns of cortical thickness based on T1-weighted anatomical MRI data of 149 PD patients and 369 controls. For this purpose, we used partial least squares regression to identify gene expression patterns correlated with cortical thickness changes. In addition, we identified gene expression patterns underlying the relationship between cortical thickness and clinical domains of PD. Our results show that genes whose expression in the healthy brain is associated with cortical thickness changes in PD are enriched in biological pathways related to sumoylation, regulation of mitotic cell cycle, mitochondrial translation, DNA damage responses, and ER-Golgi traffic. The associated pathways were highly related to each other and all belong to cellular maintenance mechanisms. The expression of genes within most pathways was negatively correlated with cortical thickness changes, showing higher expression in regions associated with decreased cortical thickness (atrophy). On the other hand, sumoylation pathways were positively correlated with cortical thickness changes, showing higher expression in regions with increased cortical thickness (hypertrophy). Our findings suggest that alterations in the balanced interplay of these mechanisms play a role in changes of cortical thickness in PD and possibly influence motor and cognitive functions.
cortical thicknessneurodegenerative diseasesneuroimaging dataimaging-geneticsgene expression analysis1272117126Stichting voor de Technische Wetenschappen10.13039/501100003958Introduction
Parkinson’s disease (PD) is a neurodegenerative disorder characterized by a progressive loss of dopaminergic and non-dopaminergic neurons in the brain and peripheral and autonomic nervous system (Hirsch et al., 2012). Cortical atrophy occurs during the later disease stages and has been associated with cognitive decline, including executive, attentional, memory, and visuospatial deficits (Aarsland et al., 2017; Wilson et al., 2019). Although MRI studies of patient brains have tried to link regional cortical atrophy to clinical features of the disease (Rosenberg-Katz et al., 2016; Wang et al., 2016; Chen et al., 2017; Li et al., 2018; Zheng et al., 2019), little is known about the pathobiology that underlies the selective cortical vulnerability in PD.
Analyzing the transcriptome in vulnerable cortical regions may aid in better understanding the underlying molecular mechanisms of atrophy in PD. Although gene expression data of human post-mortem PD brains is available, most findings relate to studies that focused only on one or few coarse brain regions (Oerton and Bender, 2017). To perform whole brain analysis of both gene expression and imaging data, studies turn to the Allen Human Brain Atlas (AHBA), a high resolution gene expression atlas covering the entire brain of six adult donors without any history of neurological disorders (Hawrylycz et al., 2015; Arnatkevičiūtė et al., 2019). The AHBA has been combined with functional MRI data of PD patients and revealed that the regional expression of MAPT, but not SNCA, correlated with the loss of regional connectivity (Rittman et al., 2016). Using a similar approach, correlations were identified between a cortical atrophy pattern and the regional expression of 17 genes implicated in PD (Freeze et al., 2018). Although both studies used spatial transcriptomics to explore gene expression across the whole brain, they only analyzed the expression of a limited set of genes that are of interest to PD, e.g., those that are known as genetic risk factors.
To investigate the relationship between high dimensional genome-wide expression patterns and imaging data, multivariate analysis methods are required. Partial least squares (PLS) regression has been used to perform simultaneous analysis of brain-wide gene expression from the AHBA and neuroimaging data of adolescents, healthy adults, and Huntington’s disease patients (Vértes et al., 2016; Whitaker et al., 2016a; McColgan et al., 2018). The PLS approach allows the linking of multiple predictor variables (genes) and multiple response variables (imaging features) and deals with multicollinearity by projecting variables to a smaller set of components that are maximally correlated between both datasets. Thus, PLS is an attractive model to identify gene expression patterns associated with imaging features.
Here, we exploited PLS regression to find transcriptomic signatures that are related to changes in cortical thickness (CT) in PD. MRI data was obtained from patients and age-matched controls to find CT changes across all cortical regions. Gene expression samples from healthy donors in the AHBA were anatomically mapped to the cortical regions to find brain-wide gene expression patterns predictive of the CT changes observed in PD patients. In addition, we assessed the relationships between CT and clinical scores in PD patients and used a second PLS model to find expression patterns associated with these relationships across all cortical regions. With these models we address three research questions: (1) Which cortical regions show CT changes in PD, (2) Which genes and biological pathways show expression patterns associated with these regional changes, and (3) Which molecular mechanisms underlie the relationships between CT and clinical scores in PD. To answer these questions, we explored the whole transcriptome in cortical regions of the healthy brain to find expression signatures predictive of imaging features in PD.
Materials and MethodsMRI Data Acquisition
MRI images of 149 PD patients (mean age = 64.8 years; 65.7% male) were obtained from a cross-sectional cohort study and is part of the “PROfiling PARKinson’s disease” (PROPARK) study (de Schipper et al., 2017; Supplementary Table 1). PD patients were recruited from the outpatient clinic for Movement Disorders of the Department of Neurology of the Leiden University Medical Center and nearby university and regional hospitals. All participants fulfilled the United Kingdom Parkinson’s Disease Society Brain Bank criteria for idiopathic PD (Gibb and Lees, 1988); written consent was obtained from all participants. The Medical Ethics Committee of the LUMC approved the study. Three-dimensional T1-weighted anatomical images were acquired on a 3 Tesla MRI scanner (Philips Achieva, Best, Netherlands) using a standard 32-channel whole-head coil. Acquisition parameters were: repetition time = 9.8 ms, echo time = 4.6 ms, flip angle = 8°, field of view 220 × 174 × 156 mm, 130 slices with a slice thickness of 1.2 mm with no gap between slices, resulting in a voxel size of 1.15 mm × 1.15 mm × 1.20 mm.
Three−dimensional T1−weighted images from 369 controls (mean age = 65.7 years; 48.1% male) were acquired in a different cohort (Altmann-schneider et al., 2012), where all imaging was performed on a whole body 3 Tesla MRI scanner (Philips Medical Systems, Best, Netherlands), using the following imaging parameters: TR = 9.7 ms, TE = 4.6 ms, FA = 8°, FOV = 224 × 177 × 168 mm. The anatomical images covered the entire brain with no gap between slices resulting in a nominal voxel size of 1.17 × 1.17 × 1.4 mm. Acquisition time was approximately 5 min.
Cortical Thickness Changes in Segmented Cortical Regions
Cortical thickness in cortical regions of PD patients and controls was determined using cortical parcellation implemented in FreeSurfer version 5.3.0 (Fischl and Dale, 2000). The FreeSurfer algorithm automatically parcellates the cortex and assigns a neuroanatomical label to each location on a cortical surface model based on probabilistic information. The parcellation scheme of the Desikan–Killiany atlas was used to divide the cortex into 34 regions per hemisphere (Desikan et al., 2006).
To assess CT changes between patients (149) and controls (369), a two-tailed t-test assuming unequal variances was applied in SPSS Statistics version 23. P-values were corrected for multiple testing across 68 cortical regions using the Benjamini-Hochberg (BH) method. A two-tailed t-test was also used to assess CT differences between the left and right hemisphere for each one of the 34 cortical regions, with P-values being BH-corrected across the 34 cortical regions.
Clinical Scores
All patients underwent standardized assessments, and an evaluation of demographic and clinical characteristics (de Schipper et al., 2017). MDS-UPDRS is a clinical rating scale consisting of four parts: (I) Non-motor Experiences of Daily Living; (II) Motor Experiences of Daily Living; (III) Motor Examination; and (IV) Motor Complications (Goetz et al., 2008). UPDRSTOTSCR is the total score of all four parts. The SENS-PD scale is a composite score of non-dopaminergic symptoms (van der Heeden et al., 2016), LED is the levodopa equivalent dose (Tomlinson et al., 2010), and MMSE is the mini-mental state examination (Folstein et al., 1975).
Relationship Between Cortical Thickness and Clinical Scores
We used CT data and clinical scores from 149 PD patients to determine the relationships between CT and clinical domains. We selected nine clinical features with numeric (non-nominal) values for which scores were available for 82–123 patients: AGEONSET, SENSPDSC, MDS_UPDRS_3, MMSE, LED, MDS_UPDRS_1, MDS_UPDRS_2, MDS_UPDRS_4, and UPDRSTOTSCR (Supplementary Figure 1).
The correlation between CT and the scores of each clinical feature individually was determined across patients by applying linear regression. To obtain maximum correlation, separate linear regression models were used for each combination of a region and clinical feature:
CTi=α+β1Kj+β2Age+ε
where CTi is the CT of one region i across patients, K_j is the score of one clinical feature j across patients. Age is taken into account to correct for the age of patients. α is the background term, β1 is the regression coefficient for K_j, β2 is the regression coefficients for Age, and ε is the residual. The regression coefficient β1 was used to determine the relationship between CT and clinical domain scores, and assessed for statistical significance where P-values were BH-corrected for 34 regions and nine clinical features (t-test, H0: β1 = 0, P < 0.05).
Mapping Transcriptomic Data to Cortical Regions
We downloaded normalized gene expression data from the Allen Human Brain Atlas (AHBA1), a microarray data set of 3,702 anatomical brain regions from six non-neurological individuals [5 males and 1 female, mean age 42, range 24–57 years (Hawrylycz et al., 2015)]. Preprocessing steps are described in Supplementary Methods. To analyze the transcriptome in the cortical regions, we used the mapping of AHBA samples to cortical regions in neuroimaging data proposed in Arnatkevičiūtė et al. (2019), where they applied Freesurfer on T1 MRIs of the six donors in the AHBA to segment the cortical regions according to the Desikan-Killiany atlas. AHBA samples were mapped to 34 cortical regions from the left hemisphere, since for only two out of six brains samples were collected from both hemispheres and for four brains they only sampled from the left hemisphere. By only analyzing the left hemisphere, we assumed that there are small to no differences in gene expression between the left and right hemisphere (Hawrylycz et al., 2015). Samples were assigned to a segmented cortical region when their MNI coordinates corresponds to a voxel within a parcel, including samples that are up to 2 mm away from any voxel in the parcel. In total 1,284 samples from the AHBA were assigned to the 34 cortical regions.
Partial Least Squares (PLS) Model-1 and Model-2
We used PLS regression (R-package pls 2.7) to find gene expression patterns across the 34 cortical regions that are predictive of gray matter atrophy and possibly their relationship to scores of nine clinical domains (Supplementary Methods). PLS regression and principal component analysis regression are both methods where the original measurements are projected to latent variables to study the data in reduced dimensions (Figure 1A). PLS, however, projects variables from each dataset to latent variables such that they are maximally correlated between two datasets X and Y (Figure 1B). In this study, the predictor X is a gene expression matrix of 34 regions (n) in the left hemisphere and all 20,017 genes (m) and is used to predict imaging variables (p) in the same set of 34 cortical regions. For each cortical region and each gene, expression levels were averaged across samples that fall within that cortical region and then averaged across the six donors from the AHBA, such that the input matrix of predictor variables contains one expression value for every gene per cortical region. We implemented two PLS models (Figure 1C): one single-response PLS model, model-1, to predict CT changes, measured as the t-statistics of ΔCT between PD patients and controls, and one multi-response PLS model, model-2, to predict the correlation between CT and clinical scores in PD patients, measured as the t-statistics of the coefficients β1 in Eq. 1.
Principal of partial least squares regression (PLS). (A) Principal component analysis (PCA) and PLS project measurements to a new latent space. Unlike PCA, PLS tries to find a latent space that is maximally correlated with another measurement y from dataset Y on the same samples. (B) The first latent component t1 of dataset X is maximally correlated with the first latent factor u1 of dataset Y. T and U scores determine the outer relations of individual datasets in the model. The coefficient β determines the inner relation between both datasets X and Y in the model (more details in Supplementary Methods). (C) In PLS model-1, we used regional gene expression as input to predict the regional t-statistics of ΔCT. Given the PLS model, R in Eqs. 5 and 6 in Supplementary Methods is used as gene weights. In PLS model-2, we used the same input to predict the t-statistics of correlation coefficients β1 of clinical features from Eq. 1.
Pathway Enrichment
Pathway enrichment analysis was done using gene set enrichment analysis (GSEA) and 2,225 pathways from the Reactome database in ReactomePA R-package version 1.28. Genes were ranked based on their weights to each PLS component; R in Eqs. 5 and 6 in Supplementary Methods. Pathways were significant when the FDR-adjusted P < 0.05.
Data and Code Availability
Transcriptomic data from the AHBA is available at http://human.brain-map.org/. All scripts were run in R version 4 and can be found online at https://github.com/arlinkeo/pd_pls.
ResultsCortical Thickness Changes Between Parkinson’s Disease Patients and Controls
We analyzed CT changes between PD patients and healthy controls (ΔCT) as a measure for gray matter loss (Figure 1C). Each of the 68 cortical regions from both hemispheres was assessed, for which ΔCT was statistically significant in 10 cortical regions (t-test, BH-corrected P < 0.05; Figure 2A and Supplementary Table 2). The lateral occipital cortex showed decreased CT in patients compared to controls in both the left hemisphere and right hemisphere. The left caudal anterior cingulate, right isthmus cingulate, and right pericalcarine also showed decreased CT in patients. Cortical regions with increased CT in patients included the pars opercularis from both the left hemisphere and right hemisphere, the right rostral middle frontal cortex, right temporal pole, and right superior temporal cortex. In general, we observed more decreased CT (atrophy) in caudal regions of the cortex compared to rostral regions that showed increased CT (hypertrophy).
t-Statistics of cortical thickness changes (ΔCT) across cortical regions. (A) CT was assessed between PD patients and healthy controls. Higher t-statistics (red) indicate a larger CT in controls compared to the CT in patients and thus corresponds to cortical atrophy. (B) CT in the left hemisphere compared to the right hemisphere in PD patients. Higher t-statistics (red) indicate a larger CT in the right hemisphere compared to the left hemisphere and thus corresponds to cortical atrophy in the left hemisphere. P-values are BH-corrected and significant regions (P < 0.05) are labeled.
Cortical Thickness Changes Between Hemispheres in Parkinson’s Disease
Clinical symptoms appear asymmetrical at disease onset with the left hemisphere being more susceptible to degeneration than the right (Claassen et al., 2016). To assess whether this asymmetry is reflected also in the observed atrophy patterns, we compared the CT between the left and right hemisphere for each of the 34 cortical regions in PD patients. We found six cortical regions that showed significant hemispheric differences (BH-corrected P < 0.05; Figure 2B and Supplementary Table 3). For five out of six significant regions, CT was indeed smaller in the left hemisphere compared to the right: banks of superior temporal sulcus, entorhinal cortex, temporal pole, medial orbitofrontal cortex, and lateral occipital cortex. For the lateral orbitofrontal cortex, the CT was larger in the left hemisphere compared to the right.
Gene Expression Patterns Predictive of Cortical Thickness Changes in Parkinson’s Disease Patients
To identify the molecular mechanisms underlying CT changes in PD, we integrated the imaging features with brain-wide gene expression profiles from the AHBA (Figure 1C). Using PLS model-1 (see section “Materials and Methods”), the expression of all 20,017 genes in 34 brain regions from the left hemisphere was used as predictor variables and we used the t-statistics of ΔCT between PD patients and controls in the 34 regions (Supplementary Table 2) as a single response variable. The number of AHBA samples varied between 0 and 92 for each one of the six brain donors and 34 cortical regions (Supplementary Table 4).
The PLS components that explain maximum covariance between the input space and the response variable are derived from successively deflated predictor and response matrices. Hence, the first component of the predictor matrix, component-1, has maximum covariance with the first component of the response matrix, and the second component of the predictor matrix, component-2, has maximum covariance with the second component of the response matrix, etc. Since PLS model-1 has a single response variable, component-1 of the response matrix is equal to a scaled version of the single response variable. As such, we only examined PLS component-1 of the predictor matrix (additional checking with leave-one-out cross-validation showed that the optimal number of components is indeed one, Supplementary Figure 2).
The scores of PLS component-1 of the predictor variables (genes) showed a caudal-to-rostral expression pattern (Figure 3A) that was correlated with CT changes in PD brains (Figure 3B), i.e., gene expression of PLS component-1 was high in caudal regions associated with atrophy and low in regions associated with hypertrophy. The Pearson correlation between the PLS component-1 scores of the predictor variables (gene expression) and the response variable (t-statistics of ΔCT) was 0.58, and explained 20.5% of the variance in gene expression and 34.2% of the variance in CT changes. Cortical atrophy was highest in the lateral occipital cortex and related to high PLS component-1 scores. The pericalcarine showed the highest PLS component-1 score. These results showed that the expression profiles of a weighted combination of genes can be predictive of CT changes in PD.
Transcriptional signatures to predict t-statistics of ΔCT between PD patients and controls in PLS model-1. (A) PLS component-1 scores of predictor variables (gene expression) visualized in cortical regions (lateral and medial view of the left hemisphere). (B) Regression fit of the latent predictor variable, PLS component-1 scores, with the single response variable, CT changes in PD measured as the t-statistics of ΔCT between PD patients (149) and controls (369) across the 34 cortical regions. (C) Mean expression of genes in the top 30 significant pathways (rows) across cortical regions (columns). A complete heatmap with all significant pathways is given in Supplementary Figure 4. The correlation between transcriptomic signatures and CT changes in PD across cortical regions is predicted by the gene weights for PLS component-1 shown in boxplots for each pathway where the median weight is either negative or positive. Negatively correlated pathways show high gene expression in regions with low t-statistics of ΔCT and gene expression decreases in regions with higher t-statistics of ΔCT. In our analysis, negative t-statistics correspond to increased CT (cortical hypertrophy) and positive t-statistics of ΔCT correspond to decreased CT (cortical atrophy). Positively correlated pathways show low expression in regions with low t-statistics of ΔCT and expression increases in regions with higher t-statistics of ΔCT. * = APC:Cdc20 mediated degradation of cell cycle proteins prior to satisfaction of the cell cycle checkpoint. ** = APC/C:Cdh1 mediated degradation of Cdc20 and other APC/C:Cdh1 targeted proteins in late mitosis/early G1.
Functionality of Genes Predictive of Cortical Thickness Changes
A PLS component of the predictor variables is a linear combination of weighted gene expression. We used the gene weights of PLS component-1 to perform GSEA analysis and revealed significant enrichment of 90 pathways, which were among others involved in DNA damage checkpoints, stabilization of p53, regulation of apoptosis, mitochondrial translation, and SUMOylation of chromatin organization proteins (Supplementary Table 5). High overlap of genes between the enriched pathways suggested that these functional processes are highly related to each other (Supplementary Figure 3).
Significant pathways are either positively or negatively correlated with CT changes based on the median weight of genes within pathways. Out of the 90 pathways that were significantly enriched, three pathways were positively correlated with the t-statistic of ΔCT. These included SUMOylation of chromatin organization proteins, signaling by cytosolic FGFR1 fusion mutants, and class C/3 (Metabotropic glutamate/pheromone receptors). Higher mean expression of genes within these three pathways is related to cortical atrophy (higher t-statistics of ΔCT); as apparent in the lateral occipital cortex (Figure 3C and Supplementary Figure 4). The positive correlation also indicates that a lower expression of these pathways is related to cortical hypertrophy (lower t-statistics of ΔCT). We found 87 negatively correlated pathways (median gene weight < 0). These pathways seem to play a role in the mitochondrial regulation of mitosis as we found pathways for mitochondrial translation, the regulation of mitotic cell cycle, p53-(in)dependent DNA damage checkpoints, and the degradation of mitotic proteins, such as cyclins A, and D. In general, the mean expression of genes in the negatively correlated pathways was high in cortical regions that showed hypertrophy, such as the pars opercularis or the entorhinal cortex.
Relationships Between Clinical Scores and Cortical Thickness
Next, we set to understand the relationship between CT in 34 cortical regions and clinical scores of PD patients. Linear regression was used to predict clinical scores from CT across patients and obtain regression coefficients, β1, for each cortical region and clinical domain (Eq. 1). We assessed the t-statistics of the regression coefficients instead of the coefficients β1 themselves (H0: β1 = 0) (Figure 4). Negative t-statistics showed that most combinations of cortical regions and clinical features are negatively correlated. For all clinical features, higher scores also indicate more severe symptoms, except for MMSE scores where lower scores indicate more severe symptoms, and thus showed positive relationships with CT. In most regions, age at onset (AGEONSET) also showed positive relationships with CT, indicating that age at onset has an effect on the loss of CT. While these general interpretations apply to most cortical regions, some regions showed different relationships with CT. For example, CT in the rostral anterior cingulate is negatively related to age at onset, and positively related to MDS-UPDRS 4 scores.
Relationship between clinical scores and CT across PD patients. Linear regression was used to predict clinical scores from CT across at most 123 PD patients (Supplementary Figure 1). Separate models were used for each clinical feature (row) and cortical region (column) to obtain regression coefficients, see Eq. 1. The heatmap shows the two-sided t-statistics of the regression coefficient when tested for H0: β1 = 0. Regions (columns) are clustered based on complete linkage of the Euclidean distance of the t-statistics of β1.
Genes Predictive of Relationships Between Clinical Scores and Cortical Thickness
With PLS model-2, we examined gene expression patterns that are predictive of the relationship between CT and clinical scores measured as t-statistics of the correlation coefficients β1in Eq. 1 (Figure 1C). We selected the first two PLS components for further analysis, which explained 36% of the variance of the predictor variables and 37% of the variance of the response variables (Supplementary Figure 5). PLS component-1 scores of the predictor variables showed a ventral-to-dorsal gene expression pattern (Figure 5A) that is correlated with the PLS component-1 scores of the response variables (Pearson r = 0.76, Figure 5B). The dorsal regions include the postcentral gyrus which is part of the primary somatosensory cortex. PLS component-2 scores of the predictor variables showed a caudal-to-rostral gene expression pattern (Pearson r = 0.56, Figure 6A) that is correlated with the PLS component-2 scores of the response variables (Figure 6B). Moreover, we assessed PLS component-3 (Pearson r = 0.76 between the predictors and response variables), which additionally explained 9% variance of the predictor variables and 11% variance of the response variables. However, further analysis revealed there were no enriched pathways for component-3 limiting the functional interpretation of this component.
Transcriptional signatures of PLS component-1 in PLS model-2 predictive of the relationship between cortical thickness (CT) and clinical scores. (A) PLS scores for PLS component-1 of the predictor variables (gene expression) and (B) its correlation with PLS component-1 of the response variables (t-statistic of β1 in Eq. 1). Axes show the percentage of explained variance for each component; r indicates the Pearson correlation. (C) Mean expression across cortical regions (columns) of genes in the top 30 significant pathways (rows). A complete heatmap with all significant pathways is given in Supplementary Figure 8. * = Respiratory electron transport, ATP synthesis by chemiosmotic coupling, and heat production by uncoupling proteins.
Transcriptional signatures of PLS component-2 in PLS model-2 predictive of the relationship between cortical thickness (CT) and clinical scores. (A) PLS scores for PLS component-2 of the predictor variables (gene expression) and (B) its correlation with PLS component-1 of the response variables (t-statistic of β1 in Eq. 1). Axes show the percentage of explained variance for each component; r indicates the Pearson correlation. (C) Mean expression across cortical regions (columns) of genes in the top 30 significant pathways (rows). A complete heatmap with all significant pathways is given in Supplementary Figure 9.
Partial least squares component-1 and component-2 of the predictor variables showed 144 and 230 significantly enriched pathways, respectively, with 54 overlapping pathways between the two components (Supplementary Tables 6, 7). Both components showed a cluster of related pathways involved in anterograde and retrograde transport between Golgi and endoplasmic reticulum (ER), and asparagine N-linked glycosylation (Supplementary Figures 6, 7). Other pathways that overlapped between the two components included macroautophagy, mitochondrial translation, mitochondrial biogenesis, mitochondrial protein import, DNA damage/telomere stress induced senescence, oxidative stress induced senescence, and protein localization.
Furthermore, PLS component-1 showed enrichment of pathways involved in tRNA and rRNA processing in the nucleus and mitochondrion, voltage-gated potassium channels, uptake and actions of bacterial toxins, and interleukin signaling. PLS component-2 showed strong enrichment of neutrophil degranulation, DNA replication, p53-(in)dependent DNA damage response, and chaperonin-mediated protein folding, and tubulin folding. Notably, the gene expression pattern of PLS component-2 was also associated with several sumoylation pathways and pathways involved in mitotic cell cycles and the degradation of mitotic proteins (Supplementary Table 7).
The enriched pathways for PLS component-1 and component-2 either showed negative or positive median gene weights that inform about the sign of the correlation between genes within a pathway and the PLS component score of the response variables (Figures 5C, 6C and Supplementary Figures 8, 9). For example, the expression of genes within pathways relating to mitochondrial processes increases for higher PLS component-1 scores of the response variables. We further assessed PLS component-1 and component-2 scores of the predictor variables and their correlation with each individual response variable, which are the clinical features and their relationship with CT in PD patients (Figure 7). The rostral-to-dorsal expression pattern of PLS component-1 is highly predictive of the relationship between CT and MMSE score in patients (Pearson’s r = 0.71). Thus, pathways associated with PLS component-1 may play an important role in cognitive circuits, which seems to be apparent based on their expression in the postcentral gyrus, but also the entorhinal cortex. PLS component-2 scores showed low correlations with the clinical features and their relation with CT across cortical regions, and suggests weak associations between the expression patterns of PLS component-2 and the response variables.
Correlations between PLS model-2 component-1 and component-2 scores of the predictor variables and individual response variables. Each plot shows the correlation between the predictor variables of gene expression (x-axis) and the response variables which are the relationships between CT and scores of a clinical feature across cortical regions (y-axis). On top of each plot, the Pearson correlation and the Y-loadings (Q in Eqs. 4 and 8) are shown; both values tell something about the sign (−/+) and magnitude (high/low) of the correlation. Each point or sample is one of the 34 cortical regions. Regions are labeled for those with minimum or maximum value along one of the axes.
Discussion
To examine the selective vulnerability of brain regions to PD, we explored PLS regression to find correlations between gene expression signatures across the healthy brain and cortical thinning patterns in PD brains. PLS regression is a suited method to identify relationships between gene expression and neuroimaging data, especially when the number of predictor variables (genes) are highly interdependent or multi-collinear, which is the case for gene expression data. This was shown before by earlier studies (Vértes et al., 2016; Whitaker et al., 2016b; McColgan et al., 2018) and we show here that PLS can also be used to identify molecular signatures that are relevant to PD, namely, we found genes that participate in different pathways involved in cellular maintenance mechanisms. By mapping genes to latent components, genes with similar expression patterns are grouped and can be explored to better understand their collective behavior and their relationship with structural changes in the PD brain. In our study, we found a caudal-to-rostral gene expression pattern that was correlated with CT changes in PD (PLS model-1); cortical atrophy was found in caudal regions while rostral regions showed cortical hypertrophy. This transcriptional signature was highly enriched for genes in biological pathways associated with mitochondrial translation and mitotic cell cycle regulation. We also found a ventral-to-dorsal and caudal-to-rostral gene expression pattern that was correlated with the relationship between CT and clinical domains of PD (PLS model-2). Both transcriptional signatures were associated with similar pathways, including macroautophagy and Golgi-ER trafficking, and may be involved in the effect of CT on clinical scores, namely MMSE scores for cognitive assessment.
The CT analyses between disease conditions and hemispheres in patients revealed cortical regions that are susceptible to atrophy. Cortical atrophy in PD commonly occurs asymmetrical, with a preference for the left hemisphere, particularly in the early disease stages (Brück et al., 2004; Mak et al., 2014; Pereira et al., 2014; Claassen et al., 2016). Here, we showed that five out of six regions with significant CT changes between hemispheres, indeed revealed more atrophy in the left hemisphere. Two cortical regions that showed significant changes between patients and controls, also showed changes between the left and right hemisphere. Our findings are in line with those of a previous study showing that cortical atrophy in PD most prominently affects the lateral occipital cortex, particularly in the left hemisphere (Freeze et al., 2018). The temporal pole showed hypertrophy in patients compared to controls, which was only significant in the right hemisphere. However, our analysis between hemispheres of PD brains suggests that the left temporal pole is more susceptible to CT loss than the right hemisphere. The remaining regions that were susceptible to CT changes showed atrophy in either the left or right hemisphere; however, differences between hemispheres in patients could not be confirmed. All 10 regions that were different between patients and controls, except the pericalcarine, were earlier identified as part of two structural covariance networks that were related to gray matter atrophy in the same PD dataset as in this study (de Schipper et al., 2017). Overall, we observed atrophy in caudal regions, which earlier has been associated with late stage PD (Claassen et al., 2016).
With our findings of the PLS models we interpret gene expression patterns of the healthy brain in relation to imaging features observed in PD. The six adult donors of the AHBA had no known neuropsychiatric or neuropathological history (Hawrylycz et al., 2012), however, it is unknown whether these individuals could have developed neurodegenerative diseases later in life. The observed spatial gene expression patterns reflect the physiological conditions in the adult healthy brain and are informative of important molecular mechanisms that are vulnerable in PD. The biological pathways found for PLS model-1 were closely related as they shared many similar genes. These interrelated pathways suggest a strong functional relationship between molecular processes involving mitotic cell cycle, mitochondrial translation, transport between ER and Golgi, DNA damage checkpoints, and sumoylation. We found that differential regulation of these molecular processes across the brain was associated with CT changes observed in PD. Similar pathways were found in PLS model-2 with multiple response variables corresponding to the relationships between CT and nine clinical domain scores in PD.
There is evidence that impaired cell cycle control plays a role in the pathogenesis of neurodegenerative diseases. In healthy conditions, differentiated neuronal cells become quiescent cells that cannot re-enter the cell cycle, however, in neurodegenerative diseases they are reactivated which is associated with increased cell death (Bonda et al., 2018). Cell cycle checkpoints are controlled by cyclins that guide the cell from one phase to the next phase and its expression can induce cell cycle re-initiation (Walton et al., 2019). Here, we found that regional expression of pathways associated with the degradation of cell cycle proteins in healthy conditions were negatively correlated with CT changes in PD, i.e., higher expression was associated with cortical hypertrophy in rostral regions such as the pars opercularis and temporal pole. Reversely, we observed low expression of protein degradation pathways in caudal regions that were associated with atrophy, and therefore suggests that regions with low expression are more vulnerable to improper degradation of cell cycle proteins leading to cell cycle initiation. This indicates that regions with low expression of such essential pathways are predisposed to neurodegeneration.
We found that the expression of several pathways associated with DNA replication and p53-(in)dependent DNA damage responses and checkpoints were correlated with CT changes. DNA replication during the S-phase may control the survival of post-mitotic cells by DNA repair mechanisms or apoptosis followed by DNA damage, which seems to be the case in neurodegenerative diseases (Tokarz et al., 2016). Furthermore, DNA damage response signaling can be modulated by tumor suppressor p53 and may also contribute to apoptosis in aging and age-related neurodegenerative disorders (Mohammadzadeh et al., 2019). These pathways showed similar expression patterns as those associated with the mitotic cell cycle, and therefore a lower expression of these DNA damage response pathways in caudal regions is related to cortical atrophy in PD.
Similar caudal-to-rostral expression patterns were found for pathways associated with mitochondrial translation. Increased risk for PD has been associated with mutations in SNCA, PARK2 (parkin), PINK1, DJ-1, and LRRK2 which have been linked to mitochondrial function and oxidative stress (Yan et al., 2013). PINK1 and parkin mediates clearance of damaged mitochondria by mitophagy and may therefore influence mitotic cell cycle progression (Sarraf et al., 2019). PINK1 also regulates both retrograde and anterograde axonal transport of mitochondria via axonal microtubules (Liu et al., 2012). The interaction between PINK1 and parkin is likely involved in mitochondrial quality control mechanisms, where anterograde transport of damaged mitochondria is reduced and retrograde transport is enhanced for elimination by mitophagy in the neuronal cell body (Lionaki et al., 2015).
A cluster of pathways involved in ER-Golgi traffic were found enriched for PLS model-2 component-1 and component-2, and involved both ER-to-Golgi anterograde and Golgi-to-ER retrograde transport. Component-1 showed a ventral-to dorsal gene expression pattern that was associated with higher correlations between CT and clinical scores, namely, the mental state of PD patients and the performance of motor functions. The pathways involved in ER-Golgi traffic were notably high expressed in the postcentral gyrus which contains the somatosensory cortex that is known for its role in processing sensory information and the regulation of emotion (Kropf et al., 2019). Our results suggest that genes in ER-Golgi traffic pathways are important for cognitive functions controlled by the postcentral gyrus. Genes involved in ER-Golgi vesicle trafficking have the ability to modify α-synuclein toxicity in yeast (Cooper et al., 2006). Moreover, fragmentation to the Golgi apparatus has been associated with the accumulation of aberrant proteins in neurodegenerative diseases, including α-synuclein (Fan et al., 2008). A study in yeast models has showed that α-synuclein expression modulates ER stress signaling response and inhibits viral infections and viral replication (Beatman et al., 2016). We found several pathways associated to HIV and influenza infections that were correlated to the relationship between CT and clinical scores. Another pathway that shared overlapping genes with those involved in ER-Golgi traffic was asparagine N-linked glycosylation, which is a biochemical linkage important for the structure and function of proteins. The N-glycosylated proteins are synthesized essentially in the ER and Golgi through sequential reactions and aberrant glycolysation of proteins may lead to inflammation and mitochondrial dysfunction in PD and consequently to a cellular overload of dysfunctional proteins (Videira and Castro-Caldas, 2018).
We found that the expression of genes involved in sumoylation of chromatin organization proteins was correlated with CT changes, i.e., higher expression within caudal brain regions, such as the pericalcarine and the lateral occipital cortex, was associated with greater atrophy in PD. Therefore, higher activity of sumoylation events may play a role in the regional vulnerability to neurodegeneration observed in PD. On the other hand, lower expression of these pathways, such as in the pars opercularis, was associated with hypertrophy in rostral regions, suggesting that lower expression of sumoylation pathways has a protective effect. Additionally, the higher expression of sumoylation pathways was associated with higher correlations between CT and clinical scores as projected by PLS component-2 in model-2. Sumoylation involves small ubiquitin-like modifier (SUMO) proteins that increase in response to cellular stress, such as DNA damage and oxidative stress, and can promote α-synuclein aggregation and Lewy body formation (Bologna and Ferrari, 2013; Eckermann, 2013; Rott et al., 2017). Several proteins associated with inherited forms of PD are targets modified by SUMO regulating mitochondrial processes, these include α-synuclein, DJ-1, and parkin (Guerra de Souza et al., 2016). Sumoylation has been associated with several diseases, including cancers, cardiac diseases, and neurodegenerative diseases (Yang et al., 2017). In cancer, sumoylation mediates cell cycle progression and plays an essential role during mitosis (Eifler and Vertegaal, 2015). SUMO seems to promote cell death mediated by the p53 tumor suppressor protein, which may be responsible for the cell death of dopaminergic neurons in PD (Eckermann, 2013). Our findings are in support of these hypotheses, and further suggest that sumoylation is important in specific cortical regions that are atrophic in PD, such as the lateral occipital cortex.
Spatial gene expression data from PD brains are limited in the number of brain donors and brain regions, which is mainly due to the limited availability of well-defined post-mortem PD patients. Therefore, we used healthy gene expression from the AHBA to perform unbiased whole brain and whole transcriptome analysis. Gene expression for all the six healthy adult donors in AHBA was only available for the left hemisphere. Therefore, this study was restricted to the analysis of the left hemisphere when combining gene expression with MRI data. Furthermore, it is generally assumed that gene expression changes with age, however, due to the limited number of brain donors in the AHBA, age-related differences in gene expression were not taken into account. In addition, MRI data from the patient and control groups were collected from different studies in separate cohorts and were age-matched, but the difference in the percentage of men in was not taken into account. In addition, the different scanner parameter settings were used in both studies. However, both datasets were processed with FreeSurfer which is a widely used tool to reliably measure thickness of gray matter in the cerebral cortex and was shown to be robust to variations in scanner platforms, sequence parameters, scan sessions, scanner manufacturer, and field strength (Fischl, 2012). Brain volumetric measurements by FreeSurfer have also been shown to be reproducible between different scanners in multiple sclerosis (Guo et al., 2019) and in vivo assessments of cortical thickness from MRI are similar to histological examinations of cortical thickness (Scholtens et al., 2015). Furthermore, an Alzheimer’s disease study showed that FreeSurfer competed with manual measurements and encourages the use of FreeSurfer in clinical practice (Clerx et al., 2015). Finally, to determine whether genes and pathways truly have predictive power of imaging features, both PLS models need to be validated with an independent imaging cohort of PD patients.
Imaging cohorts of PD patients are generally quite heterogeneous because PD is a complex disorder with a wide spectrum of symptoms that vary substantially across patients. To better understand the different forms of PD, previous neuroimaging-genetics studies have grouped PD patients based on the presence of a genetic mutation associated with PD, e.g., LRRK2 and GBA (van der Vegt et al., 2009; Weingarten et al., 2015), however, PD diagnosis cannot be confirmed based on genetic mutations. It is nowadays based on clinical observations, but true diagnosis can only be confirmed by pathological examination when patients are diseased. Therefore, it should be noted that patients with different forms of PD cannot be clearly distinguished based on clinical manifestations, genetic overlap or neuroimaging findings.
The 34 brain regions defined by the Desikan-Killiany atlas consist of different volumes and also differ in the number of gene expression samples that fall within a brain region. Since PLS requires the same number of samples for the predictor and response datasets, the transcriptomic and neuroimaging data was processed such that both datasets had an equal number of samples, which are the 34 brain regions. For the transcriptomic data this meant that the expression for one brain region was based on the average expression of all samples that fall within the brain region. For the neuroimaging data, the CT reported by FreeSurfer is the average CT for a brain region given its volume. Because volume and sample size can affect these estimates, the average gene expression and average CT, the sample size can also affect the correlations predicted by the PLS model. Finally, our PLS models also do not account for the number of subjects used in this study. Future studies may improve in applying machine learning models that are better fitted to the data to find statistical associations that are more relevant to the disease being studied.
Conclusion
We set out to find biological explanations for the selective regional vulnerability in PD. For this purpose, we applied PLS to assess the healthy transcriptome across the whole brain and find correlations with cortical thickness changes in PD, which can be observed as atrophy and hypertrophy patterns in neuroimaging data. Previous PD studies analyzed gene expression in only few brain regions due to the limited availability of PD donors, however, we made use of the AHBA to study the healthy transcriptome across the whole brain at a high resolution. We found genes that point toward pathways involved in cellular maintenance mechanisms that are well known in PD and other neurodegenerative diseases, but here we show that these pathways are differently regulated across brain regions. More specifically, sumoylation pathways showed opposite expression patterns across the brain compared to pathways associated with the regulation of mitotic cell cycle, p53-(in)dependent DNA damage response, mitochondrial translation, and ER-Golgi trafficking (Figure 8). In addition, multiple genes and biological pathways identified in this study have been associated to PD before, however, their relationship with cortical thickness and clinical features was previously not known. Also, similar pathways were identified that were associated with the severity of clinical symptoms in PD, which could be a consequence of cortical atrophy or hypertrophy. All identified pathways were highly interconnected as shown by the number of shared genes and suggest a balanced interplay between sumoylation events and the other molecular mechanisms that seem to be important in controlling CT in different cortical regions. We believe that these particular pathways are interesting for further research to better understand the shared molecular mechanisms between the multiple pathways that are involved in PD progression. With our multivariate PLS approach we were able to combine multiple data modalities to provide meaningful new insights into the selective vulnerability of brain regions to PD.
Schematic overview of the balance between biological pathways and their influence on CT across cortical brain regions. The big arrow indicates the caudal-to-rostral (red-to-blue) or rostral-to-caudal (blue-to-red) change in CT across cortical brain regions of PD patients with red indicating decreased CT (atrophy) in caudal regions and blue indicating increased CT (hypertrophy) in rostral regions. Genes within pathways associated with sumoylation showed that the expression of these genes within the pathways increases from rostral to caudal regions. Other biological pathways that were correlated with CT changes in PD included regulation of mitotic cell cycle, mitochondrial translation, DNA damage responses, and ER-Golgi traffic, and the involved genes showed decreasing expression patterns from rostral to caudal regions (or increasing from caudal to rostral regions). All enriched pathways shared many common genes and were generally associated with cellular maintenance mechanisms. Literature studies suggest that these biological pathways may be involved in the pathobiology of PD through their interaction with genetic risk variants.
Data Availability Statement
Publicly available datasets were analyzed in this study. This data can be found here: http://human.brain-map.org.
Ethics Statement
The studies involving human participants were reviewed and approved by the Medical Ethics Committee of Leiden University Medical Center. The patients/participants provided their written informed consent to participate in this study.
Author Contributions
AK, MR, and AM designed the study. JG provided the neuroimaging data. OD processed the images. AK performed the data analysis. AK, JH, MR, and AM interpreted the data and wrote the manuscript with input from all authors. AM and MR supervised the overall project. All authors read and approved final manuscript.
Conflict of Interest
The authors declare that this study received funding from AbbVie, Hoffmann-La-Roche, and Lundbeck. The funders were not involved in the study design, collection, analysis, interpretation of data, the writing of this article or the decision to submit it for publication.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Funding
This research received funding from the Netherlands Technology Foundation (STW), as part of the STW Project 12721 (Genes in Space). OD received funding from The Dutch Research Council (NWO) project 17126 (3DOmics). JH received grants from Alkemade-Keuls Foundation; Stichting Parkinson Fonds (Optimist Study); The Netherlands Organisation for Health Research and Development (#40-46000-98-101); The Netherlands Organisation for Scientific Research (#628.004.001); Hersenstichting; AbbVie; Hoffmann-La-Roche; Lundbeck; and Centre for Human Drug Research outside the submitted work.
We would like to thank Laura J. de Schipper for discussions on atrophy patterns in PD.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnins.2021.733501/full#supplementary-material
ReferencesAarslandD.CreeseB.PolitisM.ChaudhuriK. R.FfytcheD. H.WeintraubD. (2017). Cognitive decline in Parkinson disease.Nat. Rev. Neurol.13217–231. 10.1038/nrneurol.2017.2728257128Altmann-schneiderI.de CraenA. J. M.SlagboomP. E.WestendorpR. G. J.van BuchemA.MaierA. B. (2012). Brain tissue volumes in familial longevity: the Leiden longevity study.Aging Cell11933–939. 10.1111/j.1474-9726.2012.00868.x22882385ArnatkevičiūtėA.FulcherB. D.FornitoA. (2019). A practical guide to linking brain-wide gene expression and neuroimaging data.Neuroimage189353–367. 10.1016/j.neuroimage.2019.01.01130648605BeatmanE. L.MasseyA.ShivesK. D.BurrackK. S.ChamanianM.MorrisonT. E. (2016). Alpha-synuclein expression restricts RNA viral infections in the brain.J. Virol.902767–2782. 10.1128/JVI.02949-15.EditorBolognaS.FerrariS. (2013). It takes two to tango: ubiquitin and SUMO in the DNA damage response.Front. Genet.4:106. 10.3389/fgene.2013.0010623781231BondaD. J.CasadesusG.ZhuX.SmithM. A. (2018). Review: cell cycle aberrations and neurodegeneration.Neuropathol. Appl. Neurobiol.36157–163. 10.1111/j.1365-2990.2010.01064.x.ReviewBrückA.KurkiT.KaasinenV.VahlbergT.RinneJ. O. (2004). Hippocampal and prefrontal atrophy in patients with early non-demented Parkinson’s disease is related to cognitive impairment.J. Neurol. Neurosurg. Psychiatry751467–1469. 10.1136/jnnp.2003.03123715377698ChenB.WangS.SunW.ShangX.LiuH.LiuG. (2017). Functional and structural changes in gray matter of Parkinson’s disease patients with mild cognitive impairment.Eur. J. Radiol.9316–23. 10.1016/j.ejrad.2017.05.01828668411ClaassenD. O.McDonellK. E.DonahueM.RawalS.WylieS. A.NeimatJ. S. (2016). Cortical asymmetry in Parkinson’s disease: early susceptibility of the left hemisphere.Brain Behav.6:e00573. 10.1002/brb3.57328031997ClerxL.GronenschildE. H. B. M.EchavarriC.VerheyF.AaltenP.JacobsH. I. L. (2015). Can FreeSurfer compete with manual volumetric measurements in Alzheimer’s disease?Curr. Alzheimer Res.12358–367. 10.2174/156720501266615032417481325817249CooperA. A.GitlerA. D.CashikarA.HaynesC. M.HillK. J.BhullarB. (2006). a-synuclein blocks ER-Golgi traffic and Rab1 rescues neuron loss in Parkinson’s models.Science313324–328. 10.1126/science.112946216794039de SchipperL. J.van der GrondJ.MarinusJ.HenselmansJ. M. L.van HiltenJ. J. (2017). Loss of integrity and atrophy in cingulate structural covariance networks in Parkinson’s disease.Neuroimage Clin.15587–593. 10.1016/j.nicl.2017.05.01228652971DesikanR. S.SégonneF.FischlB.QuinnB. T.DickersonB. C.BlackerD. (2006). An automated labeling system for subdividing the human cerebral cortex on MRI scans into gyral based regions of interest.Neuroimage31968–980. 10.1016/j.neuroimage.2006.01.02116530430EckermannK. (2013). SUMO and Parkinson’s disease.Neuromol. Med.15737–759. 10.1007/s12017-013-8259-523979994EiflerK.VertegaalA. C. O. (2015). SUMOylation-mediated regulation of cell cycle progression and cancer SUMO: a ubiquitin-like modifier that regulates nuclear processes.Trends Biochem. Sci.40779–793. 10.1016/j.tibs.2015.09.00626601932FanJ.HuZ.ZengL.LuW.TangX.ZhangJ. (2008). Golgi apparatus and neurodegenerative diseases.Int. J. Dev. Neurosci.26523–534. 10.1016/j.ijdevneu.2008.05.00618599251FischlB. (2012). FreeSurfer.Neuroimage62774–781. 10.1016/j.neuroimage.2012.01.02122248573FischlB.DaleA. M. (2000). Measuring the thickness of the human cerebral cortex from magnetic resonance images.Proc. Natl. Acad. Sci. U.S.A.9711050–11055. 10.1073/pnas.20003379710984517FolsteinM. F.FolsteinS. E.McHughP. R. (1975). “Mini-mental state” a practical method for grading the cognitive stte of patients for the clinician.J. Psychiatr. Res.12189–198. 10.1016/0022-3956(75)90026-6FreezeB. S.AcostaD.PandyaS.ZhaoY.RajA. (2018). Regional expression of genes mediating trans-synaptic alpha-synuclein transfer predicts regional atrophy in Parkinson disease.Neuroimage Clin.18456–466. 10.1016/j.nicl.2018.01.00929868450GibbW. R.LeesA. J. (1988). The relevance of the Lewy body to the pathogenesis of idiopathic Parkinson’s disease.J. Neurol. Neurosurg. Psychiatry51745–752.GoetzC. G.TilleyB. C.ShaftmanS. R.StebbinsG. T.FahnS.Martinez-MartinP. (2008). Movement disorder society-sponsored revision of the unified Parkinson’s disease rating scale (MDS-UPDRS): scale presentation and clinimetric testing results.Mov. Disord.232129–2170. 10.1002/mds.2234019025984Guerra de SouzaA. C.PredigerR. D.CimarostiH. (2016). SUMO-regulated mitochondrial function in Parkinson’s disease.J. Neurochem.137673–686. 10.1111/jnc.1359926932327GuoC.FerreiraD.FinkK.WestmanE.GranbergT. (2019). Repeatability and reproducibility of FreeSurfer, FSL-SIENAX and SPM brain volumetric measurements and the effect of lesion filling in multiple sclerosis.Eur. Radiol.291355–1364. 10.1007/s00330-018-5710-x30242503HawrylyczM.MillerJ. A.MenonV.FengD.DolbeareT.Guillozet-BongaartsA. L. (2015). Canonical genetic signatures of the adult human brain.Nat. Neurosci.181832–1844. 10.1038/nn.417126571460HawrylyczM. J.LeinE. S.Guillozet-bongaartsA. L.ShenE. H.NgL.MillerJ. A. (2012). An anatomically comprehensive atlas of the adult human brain transcriptome.Nature489391–399. 10.1038/nature1140522996553HirschE. C.VyasS.HunotS. (2012). Neuroinflammation in Parkinson’s disease.Park. Relat. Disord.18S210–S212. 10.1016/s1353-8020(11)70065-7KropfE.SyanS. K.MinuzziL.FreyB. N. (2019). From anatomy to function: the role of the somatosensory cortex in emotional regulation.Braz. J. Psychiatry41261–269. 10.1590/1516-4446-2018-018330540029LiX.XingY.Martin-BastidaA.PicciniP.AuerD. P. (2018). Patterns of grey matter loss associated with motor subscores in early Parkinson’s disease.Neuroimage Clin.17498–504. 10.1016/j.nicl.2017.11.00929201638LionakiE.MarkakiM.PalikarasK.TavernarakisN. (2015). Mitochondria, autophagy and age-associated neurodegenerative diseases: new insights into a complex interplay.Biochim. Biophys. Acta18471412–1423. 10.1016/j.bbabio.2015.04.01025917894LiuS.SawadaT.LeeS.YuW.SilverioG.AlapattP. (2012). Parkinson’s disease-associated kinase PINK1 regulates miro protein level and axonal transport of mitochondria.PLoS Genet.8:e1002537. 10.1371/journal.pgen.100253722396657MakE.ZhouJ.TanL. C. S.AuW. L.SitohY. Y.KandiahN. (2014). Cognitive deficits in mild Parkinson’s disease are associated with distinct areas of grey matter atrophy.J. Neurol. Neurosurg. Psychiatry85576–580. 10.1136/jnnp-2013-30580524133286McColganP.GregoryS.SeunarineK. K.RaziA.PapoutsiM.JohnsonE. (2018). Brain regions showing white matter loss in Huntington’s disease are enriched for synaptic and metabolic genes.Biol. Psychiatry83456–465. 10.1016/j.biopsych.2017.10.01929174593MohammadzadehA.Mirza-aghazadeh-attariM.HallajS.MajidiniaM. (2019). Crosstalk between P53 and DNA damage response in ageing.DNA Repair808–15. 10.1016/j.dnarep.2019.05.00431176959OertonE.BenderA. (2017). Concordance analysis of microarray studies identifies representative gene expression changes in Parkinson’s disease: a comparison of 33 human and animal studies.BMC Neurol.17:58. 10.1186/s12883-017-0838-x28335819PereiraJ. B.SvenningssonP.WeintraubD.BrønnickK.LebedevA.WestmanE. (2014). Initial cognitive decline is associated with cortical thinning in early Parkinson disease.Neurology822017–2025. 10.1212/WNL.000000000000048324808018RittmanT.RubinovM.VértesP. E.PatelA. X.GinestetC. E.GhoshB. C. P. (2016). Regional expression of the MAPT gene is associated with loss of hubs in brain networks and cognitive impairment in Parkinson disease and progressive supranuclear palsy.Neurobiol. Aging48153–160. 10.1016/j.neurobiolaging.2016.09.00127697694Rosenberg-KatzK.HermanT.JacobY.KliperE.GiladiN.HausdorffJ. M. (2016). Subcortical volumes differ in Parkinson’s disease motor subtypes: new insights into the pathophysiology of disparate symptoms.Front. Hum. Neurosci.10:356. 10.3389/fnhum.2016.0035627462214RottR.SzargelR.ShaniV.HamzaH.SavyonM.AbdF. (2017). SUMOylation and ubiquitination reciprocally regulate α-synuclein degradation and pathological aggregation.Proc. Natl. Acad. Sci. U.S.A.11413176–13181. 10.1073/pnas.170435111429180403SarrafS. A.SiderisD. P.GiagtzoglouN.NiL.KankelM. W.SenA. (2019). PINK1/Parkin influences cell cycle by sequestering TBK1 at damaged mitochondria, inhibiting mitosis.Cell Rep.29225–235. 10.1016/j.celrep.2019.08.08531577952ScholtensL. H.de ReusM. A.van den HeuvelM. P. (2015). Linking contemporary high resolution magnetic resonance imaging to the von economo legacy: a study on the comparison of MRI cortical thickness and histological measurements of cortical structure.Hum. Brain Mapp.363038–3046. 10.1002/hbm.2282625988402TokarzP.KaarnirantaK.BlasiakJ. (2016). Role of the cell cycle re-initiation in DNA damage response of post-mitotic cells and its implication in the pathogenesis of neurodegenerative diseases.Aging Dis.19131–140. 10.1089/rej.2015.171726214710TomlinsonC. L.StoweR.PatelS.RickC.GrayR.ClarkeC. E. (2010). Systematic review of levodopa dose equivalency reporting in Parkinson’s disease.Mov. Disord.252649–2653. 10.1002/mds.2342921069833van der HeedenJ. F.MarinusJ.Martinez-MartinP.van HiltenJ. J. (2016). Evaluation of severity of predominantly non-dopaminergic symptoms in Parkinson’s disease: the SENS-PD scale.Park. Relat. Disord.2539–44. 10.1016/j.parkreldis.2016.02.01626922004van der VegtJ. P. M.Van NuenenB. F. L.BloemB. R.KleinC.SiebnerH. R. (2009). Imaging the impact of genes on Parkinson’s disease.Neuroscience164191–204. 10.1016/j.neuroscience.2009.01.05519409223VértesP. E.RittmanT.WhitakerK. J.Romero-GarciaR.VášaF.KitzbichlerM. G. (2016). Gene transcription profiles associated with inter-modular hubs and connection distance in human functional magnetic resonance imaging networks.Philos. Trans. R. Soc. B Biol. Sci.371:20150362. 10.1098/rstb.2015.036227574314VideiraP. A. Q.Castro-CaldasM. (2018). Linking glycation and glycosylation with inflammation and mitochondrial dysfunction in Parkinson’s disease.Front. Neurosci.12:381. 10.3389/fnins.2018.0038129930494WaltonC. C.ZhangW.Patiño-parradoI.Barrio-alonsoE.GarridoJ.FradeJ. M. (2019). Primary neurons can enter M-phase.Sci. Rep.9:4594. 10.1038/s41598-019-40462-430872646WangM.JiangS.YuanY.ZhangL.DingJ.WangJ. (2016). Alterations of functional and structural connectivity of freezing of gait in Parkinson’s disease.J. Neurol.2631583–1592. 10.1007/s00415-016-8174-427230857WeingartenC. P.SundmanM. H.HickeyP.ChenN. K. (2015). Neuroimaging of Parkinson’s disease: expanding views.Neurosci. Biobehav. Rev.5916–52. 10.1016/j.neubiorev.2015.09.00726409344WhitakerK. J.VértesP. E.Romero-garciaR.MoutoussisM.PrabhuG. (2016a). Adolescence is associated with genomically patterned consolidation of the hubs of the human brain connectome.Proc. Natl. Acad. Sci. U.S.A.1139105–9110. 10.1073/pnas.160174511327457931WhitakerK. J.VértesP. E.Romero-GarciaR.VášaF.MoutoussisM.PrabhuG. (2016b). Adolescence is associated with genomically patterned consolidation of the hubs of the human brain connectome.Proc. Natl. Acad. Sci. U.S.A.1139105–9110.WilsonH.NiccoliniF.PellicanoC.PolitisM. (2019). Cortical thinning across Parkinson’s disease stages and clinical correlates.J. Neurol. Sci.39831–38. 10.1016/j.jns.2019.01.02030682518YanM. H.WangX.ZhuX. (2013). Mitochondrial defects and oxidative stress in Alzheimer disease and Parkinson disease.Free Radic. Biol. Med.6290–101. 10.1016/j.freeradbiomed.2012.11.01423200807YangY.HeY.WangX.HeG.ZhangP.ZhuH. (2017). Protein SUMOylation modification and its associations with disease.Open Biol.7:170167.ZhengD.ChenC.SongW. C.YiZ. Q.ZhaoP. W.ZhongJ. G. (2019). Regional gray matter reductions associated with mild cognitive impairment in Parkinson’s disease: a meta-analysis of voxel-based morphometry studies.Behav. Brain Res.371:111973. 10.1016/j.bbr.2019.11197331128163