Genetic risk of Parkinson disease and progression:: An analysis of 13 longitudinal cohorts.

OBJECTIVE: To determine if any association between previously identified alleles that confer risk for Parkinson disease and variables measuring disease progression.
METHODS: We evaluated the association between 31 risk variants and variables measuring disease progression. A total of 23,423 visits by 4,307 patients of European ancestry from 13 longitudinal cohorts in Europe, North America, and Australia were analyzed.
RESULTS: We confirmed the importance of GBA on phenotypes. GBA variants were associated with the development of daytime sleepiness (p.N370S: hazard ratio [HR] 3.28 [1.69-6.34]) and possible REM sleep behavior (p.T408M: odds ratio 6.48 [2.04-20.60]). We also replicated previously reported associations of GBA variants with motor/cognitive declines. The other genotype-phenotype associations include an intergenic variant near LRRK2 and the faster development of motor symptom (Hoehn and Yahr scale 3.0 HR 1.33 [1.16-1.52] for the C allele of rs76904798) and an intronic variant in PMVK and the development of wearing-off effects (HR 1.66 [1.19-2.31] for the C allele of rs114138760). Age at onset was associated with TMEM175 variant p.M393T (-0.72 [-1.21 to -0.23] in years), the C allele of rs199347 (intronic region of GPNMB, 0.70 [0.27-1.14]), and G allele of rs1106180 (intronic region of CCDC62, 0.62 [0.21-1.03]).
CONCLUSIONS: This study provides evidence that alleles associated with Parkinson disease risk, in particular GBA variants, also contribute to the heterogeneity of multiple motor and nonmotor aspects. Accounting for genetic variability will be a useful factor in understanding disease course and in minimizing heterogeneity in clinical trials.

Parkinson disease is one of the most common neurodegenerative diseases, with an estimated lifetime risk as high as 1%–2%.[1] Parkinson disease is traditionally characterized by motor features such as bradykinesia, rigidity, and tremor. However, in addition to these motor symptoms, patients with Parkinson disease also develop nonmotor symptoms (NMSs), which include depression, cognitive decline, sleep abnormalities, reduced olfaction, and autonomic dysfunction.[2] Collectively, the combined spectrum of motor and NMSs more accurately reflects the multisystem nature of the disease. Patients with Parkinson disease may present with various combinations of symptoms and show differences in the rates of progression.[3] The application of modern molecular genetic approaches over the last decade has revealed a significant number of genetic risk loci for idiopathic Parkinson disease.[4-7] However, in comparison with case-control genome-wide association study (GWAS), analyzing how genetic factors influence clinical presentation and progression requires longitudinal cohorts with much more detailed observations. Such data are sparse, and individual cohorts are often small in size and quite varied, posing a challenge both in sample size and heterogeneity.

In an attempt to address these issues, we collected data from 13 distinct longitudinal Parkinson disease cohorts with detailed clinical data, including assessment of disease progression. We sought to determine whether Parkinson disease genetic risk factors, either in the form of known GWAS variants or an aggregate genetic risk score (GRS), are associated with changes in clinical progression and the disease features.

Methods

Study design and participants

A total of 13 Parkinson disease cohorts from North America, Europe, and Australia participated in the study. Nine were prospective observational cohorts and the rest were from randomized clinical trials. The observational cohorts were Drug Interaction with Genes in Parkinson's Disease (DIGPD), Harvard Biomarkers Study (HBS), Oslo Parkinson's Disease study (partly including retrospective data), The Norwegian ParkWest study (ParkWest), Parkinson's Disease Biomarker Program (PDBP), Parkinsonism: Incidence and Cognitive and Non-motor heterogeneity In CambridgeShire (PICNICS), Parkinson's Progression Markers Initiative (PPMI), Profiling Parkinson's disease study (ProPark), and the Morris K. Udall Centers for Parkinson's Research (Udall). The 4 cohorts from randomized clinical trials were Deprenyl and Tocopherol Antioxidative Therapy of Parkinsonism (DATATOP), NIH Exploratory Trials in Parkinson's Disease Large Simple Study 1, ParkFit study (ParkFit), and Parkinson Research Examination of CEP-1347 Trial with a subsequent prospective study (PreCEPT/PostCEPT). Information on these cohorts can be found in appendix e-1 (links.lww.com/NXG/A169). Subsets of participants from the cohorts who provided DNA and were nonrelated participants with PD, diagnosed at age 18 years or later, and of European ancestry were included in the study. Participants' information and genetic samples were obtained under appropriate written consent and with local institutional and ethical approvals.

Genotyping SNPs and calculation of GRS

Oslo samples were genotyped on the Illumina Infinium OmniExpress array, DIGPD samples were genotyped by Illumina Multi-Ethnic Genotyping Array, and all other samples were genotyped on the NeuroX array.[8] The quality control process of variant calling included GenTrain score <0.7, minor allele frequency (MAF) >0.05 (for sample quality control but not in our analysis of rare risk factors), and Hardy-Weinberg equilibrium test statistic >10−6. Sample-specific quality control included a sample call rate of >0.95, confirmation of sex through genotyping, homozygosity quantified by F within ± 3 SD from the population mean, European ancestry confirmed by principal-components analysis with 1000 Genomes data as the reference, and genetic relatedness of any 2 individuals <0.125. Detailed information regarding NeuroX and the quality control process has been described previously.[9] In the present study, we investigated 31 single nucleotide polymorphisms (SNPs) previously shown to be significantly associated with Parkinson disease.[10-12] In addition, we also calculated a GRS for each participant based on these variants. The scores were transformed into Z-scores within each cohort and treated as an exposure, with effect estimates based on 1 SD change from the population mean. The list of 31 SNPs and the GRS calculation method are provided in table e-1 (links.lww.com/NXG/A170).

Furthermore, principal components (PCs) were created for each data set from genotypes using PLINK. For the PC calculation, variants were filtered for MAF (>0.05), genotype missingness (<0.05), and Hardy-Weinberg equilibrium (p ≥ 10−5). The remaining variants were pruned (using a 50-kb window, with a 5 SNP shift per window and r threshold of 0.5), and PCs were calculated using the pruned variants.

Measurements

The following clinical measurements and binomial outcomes were recorded longitudinally (table e-2 links.lww.com/NXG/A171): total and subscores of the Unified Parkinson's Disease Rating Scale (UPDRS) or the Movement Disorder Society revised UPDRS version (MDS-UPDRS); modified Hoehn and Yahr scales (HY); modified Schwab and England Activities of Daily Living Scale; and scores for the Mini-Mental State Examination (MMSE), The SCales for Outcomes in PArkinson's disease (SCOPA)-Cognition, and Montreal Cognitive Assessment (MoCA). Each was treated as a continuous outcome. For the UPDRS and MDS-UPDRS scores specifically, we took Z-scores of the total and subscores (except for part 4 at baseline) to compare the original and revised UPDRS versions. The conversion was applied to the scores for all subsequent visits. For UPDRS part 4, most participants had very low scores or 0 at baseline, so we normalized across all observations within each cohort. We also analyzed binomial outcomes. If we had access to the raw data, we used common cutoff values, which had been tested and reported specificity of 85% or more in patients' population. The binomial outcomes include existence of family history (1st-degree relative. 1st- and 2nd-degree relatives in HBS, PreCEPT, ProPark, and Udall), hyposmia (University of Pennsylvania Smell Identification Test < 21,[13] or answering “yes” to question 2 in the NMS questionnaire), cognitive impairment (SCOPA-Cognition < 23, MMSE < 27, or MoCA < 24,[14,15] or diagnosed with The Diagnostic and Statistical Manual of Mental Disorders -IV criteria for dementia), wearing off (UPDRS/MDS-UPDRS part 4 off time >0 or physician's diagnosis), dyskinesia (UPDRS/MDS-UPDRS part 4 dyskinesia time >0 or physician's diagnosis), depression (Beck Depression Inventory > 14 [PICNICS used 9 instead of 14], Hamilton Depression Rating Scale > 9, Geriatric Depression Scale [GRS] > 5,[16] or physician's diagnosis), constipation (MDS-UPDRS part 1 item 11 > 0, or answering “yes” to question 5 in the NMS questionnaire), excessive daytime sleepiness (Epworth sleepiness scale > 9,[17] insomnia [MDS-UPDRS part 1, item 7 > 0], rapid eye movement sleep behavior disorder [RBD]) (answered “yes” to question 1 on the Mayo Sleep Questionnaire [MSQ],[18] or RBD screening questionnaire [RBDSQ >5],[19] restless legs syndrome [RLS]) (answered “yes” to MSQ question 3,[20] or RLS diagnosis positive by RBDSQ), and the progression to HY ≥ 3 (HY3, representing moderate to severe disease). The individual definitions of these binomial outcomes are summarized in table e-2 (links.lww.com/NXG/A171). Age, sex, years of education, age at motor symptom onset, and whether the patient was treated with levodopa or dopamine agonists at each visit were also recorded for adjustments.

Statistical analysis

Cohort-level analysis

We analyzed the association between exposures and outcomes using appropriate additive models. Covariates of interest were not available for all cohorts; therefore, the model specifications were slightly different between cohorts (detailed in table e-3, links.lww.com/NXG/A172). Briefly, the associations between an SNP/GRS and age at onset were analyzed by linear regression modeling adjusting for population stratification (PC1 and PC2). The association between family history of Parkinson disease and SNP/GRS was analyzed with a logistic regression model adjusting for PC1/2. For continuous variables, linear regression modeling adjusting for sex, education, PC1/2, age at onset, years from diagnosis, family history, and treatment status was applied. For those who had multiple observations, random intercept was added to adjust for repeated measurements of the same individual. For binomial outcomes, the logistic regression at baseline observation was applied using the same covariates as the continuous models. Those that were negative at baseline were further analyzed by a Cox regression with the same covariates but with treatment status as a time-varying covariate. Observations with missing variables were excluded from the analyses.

Meta-analysis

We applied inverse weighting (precision method) for each combination of outcome-predictor association and combined the estimates from the 13 different cohorts in a fixed effect model. Multiple test correction for SNPs was controlled with an overall false discovery rate (FDR) of 0.05 per outcome being considered significant. Similarly, multiple testing of outcomes for GRS was corrected with an FDR of 0.05, but across all traits. In addition, as a test of homogeneity, I2 indices and forest plots were used for quantitative assessment. As a sensitivity analysis, we conducted up to 13 iterations of the meta-analyses for the 12 cohorts excluding each cohort per iteration. This analysis provides information regarding heterogeneity of the cohorts and how one specific cohort exclusion affects the results. The range of estimates and maximum p values for the iterations were included. Finally, we conducted the 13-cohort meta-analysis in a random effects model with restricted maximum likelihood estimation using the same multiple testing correction.

All the above analyses were conducted with PLINK version 1.9, and R version 3.4.4 (64-bit). Statistical tests were all 2 sided.

Data availability

Qualified investigators can request raw data through the organizations' homepages (PDBP: pdbp.ninds.nih.gov/, PPMI: ppmi-info.org/) or collaboration.

Results

A total of 23,423 visits by 4,307 patients with a median follow-up period of 2.97 years (quartile range of [1.63–4.94] years) were eligible for the analysis. The baseline characteristics of the cohorts are shown in table 1. The mean ages at onset varied from 54 to 69 years; the average disease durations at cohort entry ranged from less than 1 to 10 years, and the mean observation periods were between 1.2 and 6.8 years. All DATATOP, ParkWest, PPMI, and PreCEPT participants were dopaminergic therapy naive at baseline; patients in the other cohorts were not. In the primary analysis of 13 cohorts, 17 associations were identified as significant after FDR correction (table 2, and more information in table e-4, links.lww.com/NXG/A173). Overwhelmingly, 10 were associated with GBA variants. In particular, GBA p.E365K (rs2230288) was associated with 2.37- (1.53–3.66) (95% CI) fold higher odds of having cognitive impairment at baseline (p = 1.09 × 10−4) and 2.78- (1.88–4.11) fold higher hazard ratio (HR) of developing cognitive impairment during follow-up among those who were negative for cognitive impairment at baseline (p = 2.97 × 10−7). This SNP was also associated with a higher mean on the HY at 0.10 (0.04–0.16) (p = 1.53 × 10−3), but the test of homogeneity was rejected (p = 0.017, I2 = 48.9%). In addition, it was associated with the development of an RBD among those who did not have the disorder at baseline. Other GBA mutations, p.N370S (rs767763715) and p.T408M (rs75548401), were both associated with a higher HR of reaching HY3 (4.59 [2.60–8.10] for p.N370S [p = 1.58 × 10−7] and 1.93 [1.34–2.78] for p.T408M [p = 4.40 × 10−4]). GBA p.N370N was also associated with a higher risk of developing wearing-off, dyskinesia, and daytime sleepiness. p.T408M was associated with a 6.48 (2.04–20.60) times higher odds ratio (OR) of having an RBD symptom at baseline (p = 1.53 × 10−3).

Table 1

Summary characteristics of 13 cohorts

Table 2

Meta-analysis for 13 cohorts and the results of sensitivity analysis

Two LRRK2 variants in our 31 SNPs of interest were significantly associated with outcomes. LRRK2 p.G2019S (rs34637584) was associated with higher odds of having a family history of Parkinson disease (OR 3.54 [1.72–7.29], p = 6.06 × 10−4), and the T allele of rs76904798 (intergenic at the 5′ end of LRRK2) was associated with a higher HR of reaching HY3 (HR 1.33 [1.16–1.52] for the T allele, p = 5.27 × 10−5).

Age at onset was inversely associated with the Z value of the GRS (−0.60 [−0.89 to −0.31] years per +1 SD, p = 5.33 × 10−5). Moreover, it was associated with rs34311866 (TMEM175 p.M393T), the C allele of rs199347 (intronic region of GPNMB), and the G allele of rs1106180 (intronic region of CCDC62).

The majority (14/17) of associations showed good accord across cohorts (I2 < 50%), and the forest plots (figures 1–3) also illustrate this qualitatively. Furthermore, up to 13 iterations of the leave-one-out analysis assessed 15 associations of which outcomes were measured in more than 2 cohorts and showed a small range of betas. The maximum p value of 13 iterations was less than 0.05 for all associations except for rs114138769 (intron of PMVK) and rs76763715 (GBA p.N370S) for wearing-off. A meta-analysis with a random effect model also detected 9 associations after the same FDR correction, although the model is more conservative than a fixed model.

Forest plots for GBA (p.N370S and p.T408M) variants and symptoms of Parkinson disease

DATATOP = Deprenyl and Tocopherol Antioxidative Therapy of Parkinsonism; DIGPD = Drug Interaction with Genes in Parkinson's Disease; HBS = Harvard Biomarkers Study; NET-PD_LS1 = NIH Exploratory Trials in Parkinson's Disease Large Simple Study 1; Oslo = Oslo PD study; ParkFit = ParkFit study; ParkWest = the Norwegian ParkWest study; PDBP = Parkinson's Disease Biomarker Program; PICNICS = Parkinsonism: Incidence and Cognitive and Non-motor heterogeneity In CambridgeShire; PPMI = Parkinson's Progression Markers Initiative; PreCEPT/PostCEPT = Parkinson Research Examination of CEP-1347 Trial with a subsequent prospective study; ProPark = Profiling Parkinson's Disease study; Udall = Morris K. Udall Centers for Parkinson's Research. * Indicates Beta in a Cox model; ** indicates Beta in a logistic model at baseline.

Forest plots for GBA (p.E365K) variants and symptoms of Parkinson disease

DATATOP = Deprenyl and Tocopherol Antioxidative Therapy of Parkinsonism; DIGPD = Drug Interaction with Genes in Parkinson's Disease; HBS = Harvard Biomarkers Study; NET-PD_LS1 = NIH Exploratory Trials in Parkinson's Disease Large Simple Study 1; Oslo = Oslo PD study; ParkFit = ParkFit study; ParkWest = the Norwegian ParkWest study; PDBP = Parkinson's Disease Biomarker Program; PICNICS = Parkinsonism: Incidence and Cognitive and Non-motor heterogeneity In CambridgeShire; PPMI = Parkinson's Progression Markers Initiative; PreCEPT/PostCEPT = Parkinson Research Examination of CEP-1347 Trial with a subsequent prospective study; ProPark = Profiling Parkinson's Disease study; Udall = Morris K. Udall Centers for Parkinson's Research. * Indicates Beta in a Cox model; ** indicates Beta in a logistic model at baseline; *** indicates Beta in a linear mixed model.

Forest plots for non-GBA risk variants/genetic risk score and symptoms or features of Parkinson disease

DATATOP = Deprenyl and Tocopherol Antioxidative Therapy of Parkinsonism; DIGPD = Drug Interaction with Genes in Parkinson's Disease; HBS = Harvard Biomarkers Study; NET-PD_LS1 = NIH Exploratory Trials in Parkinson's Disease Large Simple Study 1; Oslo = oslo PD study; ParkFit = ParkFit study; ParkWest = the Norwegian ParkWest study; PDBP = Parkinson's Disease Biomarker Program; PICNICS = Parkinsonism: Incidence and Cognitive and Non-motor heterogeneity In CambridgeShire; PPMI = Parkinson's Progression Markers Initiative; PreCEPT/PostCEPT = Parkinson Research Examination of CEP-1347 Trial with a subsequent prospective study; ProPark = Profiling Parkinson's Disease study; Udall = Morris K. Udall Centers for Parkinson's Research. * Indicates Beta in a Cox model; ** indicates Beta in a logistic model at baseline; *** indicates Beta in a linear mixed model.

Discussion

We conducted a meta-analysis with 13 longitudinal patient cohorts and identified multiple associations between genotypes and clinical phenotypic characteristics, including progression rates. Among these, GBA coding variants showed clear associations with the rate of cognitive decline (binomial outcome or UPDRS part 1 score) and motor symptom progression (HY, HY3), consistent with previous studies.[12,21-25]

In addition, we found associations between GBA variants and RBD and daytime sleepiness. A previous cross-sectional study with 120 Ashkenazi-Jewish patients reported a higher frequency of RBDSQ-detected RBD symptoms in GBA variant carriers.[26] Our finding suggests that GBA is associated not only with baseline clinical presentation but also with disease progression.

An association between GBA and daytime sleepiness has been rarely documented. One study reported an association between sleep problems (as assessed by the Parkinson's Disease Sleep Scale) and GBA.[27] However, this scale is a combined measure of daytime sleepiness and other aspects of sleep problems.

Finally, a GBA variant (p.N370S) was also associated with treatment-related complications of wearing-off and dyskinesia. Two studies have reported the association of GBA variants with these complications, with 1 positive and 1 negative result.[28,29] The negative result may be due to insufficient power with only 19 patients with GBA mutations.

Overall, our study provides a distinct clinical profile of patients with GBA variants compared with those without. We note that with 63 carriers for p.N370S, 166 for p.T408M, and 217 for p.E365K, we have a reasonable power, but the number is yet not enough. And this may affect the results in seemingly different magnitudes of associations and the association for different traits per variants (e.g., motor complications with p.N370S and cognitive impairment with p.E365K). Another possible explanation is that although the effects are associated with the same gene, the biological activity or molecular mechanism could be different. Such an example has already been reported for LRRK2 p.G2019S and p.G2385R.[30]

Aside from GBA variants, the associations between close intergenic (5′_end) variant of LRRK2, rs76904798, and the faster development of motor symptom, and the intronic region variant of PMVK, rs114138760, and the development of wearing-off, were significant. This variant is 4.3 kb upstream from the 5′ end of LRRK2 and reported to be associated with LRRK2 gene expression changes in recent blood cis-expression quantitative trait loci (eQTL) study from the eQTLGen Consortium.[31] In contrast, we did not find an association between rs34637584, LRRK2 coding mutation (p.G2019S) and motor progression. The p.G2019 variant is a rare variant (MAF 0.5% in our study), and our sample size was not adequate barring an extremely large effect size. The intronic region variant of PMVK, rs114138760, and the development of wearing-off was another finding. The biological effect of PMVK on PD has not been reported, but the variant is also located at close proximity of the GBA-SYT11 locus, so it is possible that its association was through a similar mechanism as GBA. Including the results of cross-sectional analysis, the associations of age at onset with rs34311866 (TMEM175, p.M393T), rs199347 (intron of GPNMB), and rs11060180 (intron of CCDC62) were found. TMEM175 has been reported to impair lysosomal and mitochondrial function and increase α-synuclein aggregation,[32] although no functional data for this missense variant were studied. Of interest, the variant has recently been reported in another study as being associated with the age at onset.[33] rs199347 is an eQTL increasing the brain expression of GPNMB,[34] suggesting a causal link. Regarding rs1160180, no functional data are available in this locus.

We also evaluated the association between genetic risk variants and clinical outcomes by 2-step meta-analysis. This analysis is exploratory, and we acknowledge that this is biased toward the null due to power issues when partitioning studies randomly. However, we believe that it is helpful to assess the rigorousness of the associations we found in the primary analysis and to explore potential missed associations.

A strength of the current study was its design, incorporating multiple distinct independent Parkinson disease cohorts with longitudinal follow-ups. Although the cohorts contained patients at different disease stages, and some of the definition of outcomes were not identical, we analyzed each cohort separately and combined the results. Thus, the significant findings are consistent and applicable to the wider Parkinson disease populations. The forest plots showed that most of the estimates agree with each other despite the relative differences in the cohort characteristics. Another strength is the size of the study. The total number of genotyped and phenotyped patients with Parkinson disease (N = 4,307) is one of the largest to date for an investigation of disease progression.

The limitations of our study were as follows. First, we only included patients of European ancestry. It is uncertain whether the associations in the current study are also applicable to people from different ethnic backgrounds and further research is needed. Second, the current analysis could not distinguish causality, only basic associations. Different approaches, such as molecular-level assessment and Mendelian randomization, are crucial. Third, interaction effects between genes and other factors are another important research target not addressed in this report because of power constraints. For example, gene-by-smoking interactions for Parkinson disease were indicated recently[35] and highlight the importance of correctly modeling gene-environment interactions. Finally, compared with the typical GWAS analysis (which includes tens of thousands of cases), the number of participants was small, and the outcomes of interest were not as simple or easily defined as with case-control distinctions in GWAS. Acknowledging the limitations, the list of associations provided here is valuable as a foundation for further studies and as an example that illustrates the potential of efforts to define the genetic basis of variability in presentation and course. Accounting for this variability, even in part, has the potential to positively affect etiology-based clinical trials by reducing variability between placebo and treatment groups and by providing better predictions of expected individual progression.