Association between Residential Greenness and Incidence of Parkinson’s Disease: A Population-Based Cohort Study in South Korea

It is widely known that exposure to residential greenness is beneficial for health. However, few studies have analyzed the association between greenery and Parkinson’s disease (PD). We selected 313,355 participants who matched the inclusion criteria from the National Health Insurance Service-National Sample Cohort, followed up from 2007 to 2015. Residential greenness, represented by the normalized difference vegetation index (NDVI), was obtained from satellite measurements. We estimated hazard ratios of PD associated with a 0.1-unit increase in long-term greenness exposure at the district level for the previous 1 year of each year until a censoring/event occurred, using time-varying Cox proportional hazard models, adjusted for individual- and area-level characteristics. During the 2,745,389 person-years of follow-up, 2621(0.8%) participants developed PD. Exposure to higher levels of residential greenness was found to be associated with a decreased risk of PD incidence (21% per 0.1-unit increase, 95% confidence interval (CI): 0.74–0.84). In subgroup analyses, stronger protective effects were observed in participants aged over 50 years, females, overweight/obese participants, non-urban residents, non-smokers, alcoholics, and those with comorbidities. Long-term exposure to greenness was beneficial to incident PD, and our findings could aid in the development of public-health strategies.


Introduction
Parkinson's disease (PD) is the second most prevalent neurodegenerative disease after Alzheimer's disease (AD), posing a great socio-economic burden to our aging society [1]. Globally, 6.1 million individuals had PD in 2016, which was much higher than 2.5 million PD cases in 1990, and disability-adjusted life-years and deaths due to PD were 3.2 million and 211,000 in 2016, respectively [2]. Given the rapid increase in the aging population worldwide, the global burden of PD is expected to increase drastically in the coming decades [2,3]. Moreover, many motor and non-motor processes are affected by PD, including tremor, bradykinesia, depression, and cognitive dysfunction [4,5], all of which contribute to patients' poor quality of life [6] and even suicide [7]. Despite its severity, currently there is no cure for PD, and establishing basic preventative treatment measures can reduce the disease's growing burden.
Considering that 90% of PD cases have no distinct genetic cause [8], many previous studies have identified behavioral or environmental protective factors for PD, including smoking, caffeine consumption, use of ibuprofen, and physical activity [9], though only physical activity could be a reasonable intervention. Recently, residential greenness has been widely examined as an emerging environmental factor that reduces various health comorbidities, including neurological diseases [10][11][12][13]. Good residential greenness has been associated with reduced prevalence of depression [13], slower cognitive decline [14], increased cortical thickness [15], reduced mortality from neurological diseases [11], and reduced risk of AD and non-AD dementia [12]. This recent evidence supports the hypothesis that residential greenness could be a neuroprotective factor for PD. To the best of our knowledge, there have been only two epidemiological studies assessing the association between residential greenness and the risk of PD: One study reported a significant protective effect of residential greenness for incident PD using 47,516 cohort participants in Ningbo, China [16], and the other study reported a marginal protective effect based on 634,432 participants in Metro Vancouver, Canada [17].
To contribute to the limited evidence on the protective role of residential greenness in PD, we conducted a national-level cohort study to investigate the association between residential greenness and the incidence of PD using a nationwide population-based cohort of South Korean citizens.

Data Source and Study Population
Data for this nationwide population-based cohort study were obtained from the National Health Insurance Service-National Sample Cohort (NHIS-NSC, ver. 2.0). The NHIS is a universal insurance provider that was launched in 2000 as a single-payer health insurance system for all citizens of South Korea [18]. The NHIS established the National Health Insurance Database (NHID), which contains medical records of all citizens, including sociodemographic information (sex, age, insurance premium, and residential addresses), healthcare utilization (outpatient visits or admission), medication prescriptions, and national health examination results. From the 2006 NHID, the NHIS created the NHIS-NSC of one million individuals, who were randomly and proportionally sampled based on sex, age, residential area, and insurance premium from the total 48 million population of South Korea. The NHIS-NSC includes individuals' information for 2002-2015 merged retrospectively and prospectively.
We selected 329,491 participants who had undergone a national health examination at least once between 2002 and 2006 from one million individuals in the NHIS-NSC because data on important PD-related factors-including smoking status, physical activity, alcohol consumption, and body mass index-were obtained from the national health examination [9]. In addition, we excluded participants who had died before 2007 (n = 0) and were diagnosed with PD before 2007 (n = 680), those with incomplete address information (n = 1579), and those with incomplete information for any covariates (n = 13,877). Finally, 313,355 participants were included in this cohort study ( Figure S1).

Primary Outcome
The primary outcome of this study was the incidence of PD. PD was identified based on the International Classification of Diseases, 10th Revision (ICD-10) code G20, and medication prescriptions. The incidence of PD was defined as the first healthcare utilization (admission or outpatient visit) due to PD as a primary or secondary diagnosis for 2007-2015, simultaneously accompanied by the prescription of anti-PD medications such as carbidopa and levodopa.

Residential Greenness
Exposure to residential greenness was represented by the satellite-derived NDVI, which has been widely used to estimate the effects of vegetation on mortality [19] and morbidity [20]. NDVI obtained by the Moderate Resolution Imaging Spectroradiometer (MODIS) from aboard Terra satellites has been produced at 16-day intervals at a 250 m spatial resolution since 2000 [21]. NDVI was calculated as the land surface reflectance using near-infrared radiation (NIR) and visible red, because green vegetation absorbs red wavelengths but reflects NIR compared to non-green vegetation [22]. Therefore, given the mathematical formula (NIR − red)/(NIR + red), negative values close to −1 indicate soil and area without vegetation, while positive values close to 1 indicate a high density of green vegetation. We determined the daily NDVI for 16 days between 2000 and 2017 by 282 spatial units named Si-Gun-Gu, which is the practical administrative district unit in Korea, and assigned annual average individual exposure of NDVI from the previous 1-5 years at baseline to censoring date by considering the annually updated participants' residence. Exposure to greenness within the buffer area could not be measured because of a lack of data, although a previous study considered exposure to greenness in buffer areas having radii of 300 m and 1000 m around a residential address [16].

Statistical Analyses
Time-varying Cox proportional hazard regression was conducted to assess the hazard ratios (HRs) and 95% confidence intervals (CIs) of PD incidence associated with annually updated exposure to NDVI from the previous 1 year at baseline. The time-varying model reflects the variability in exposure that changes over time. Survival time was measured in years from baseline (1 January 2007) until PD incidence date, mortality date, or end of follow-up (31 December 2015). Figure 1 shows the distribution of residential greenness on July 2006 (A) and 2014 (B). Thus, the time-varying model was suitable to reflect the distribution of green areas that change over time. We applied the NDVI as a continuous and categorical variable. Based on the quartile of previous 1-year NDVI exposure at baseline, we estimated the HRs of each quartile compared to the lowest quartile as reference and tested the linear trend using the ordinal rank for each quartile. For continuous NDVI, the HRs by a 0.1 increase in NDVI were estimated. (admission or outpatient visit) due to PD as a primary or secondary diagnosis for 2007-2015, simultaneously accompanied by the prescription of anti-PD medications such as carbidopa and levodopa.

Residential Greenness
Exposure to residential greenness was represented by the satellite-derived NDVI, which has been widely used to estimate the effects of vegetation on mortality [19] and morbidity [20]. NDVI obtained by the Moderate Resolution Imaging Spectroradiometer (MODIS) from aboard Terra satellites has been produced at 16-day intervals at a 250 m spatial resolution since 2000 [21]. NDVI was calculated as the land surface reflectance using near-infrared radiation (NIR) and visible red, because green vegetation absorbs red wavelengths but reflects NIR compared to non-green vegetation [22]. Therefore, given the mathematical formula (NIR − red)/(NIR + red), negative values close to −1 indicate soil and area without vegetation, while positive values close to 1 indicate a high density of green vegetation. We determined the daily NDVI for 16 days between 2000 and 2017 by 282 spatial units named Si-Gun-Gu, which is the practical administrative district unit in Korea, and assigned annual average individual exposure of NDVI from the previous 1-5 years at baseline to censoring date by considering the annually updated participants' residence. Exposure to greenness within the buffer area could not be measured because of a lack of data, although a previous study considered exposure to greenness in buffer areas having radii of 300 m and 1000 m around a residential address [16].

Statistical Analyses
Time-varying Cox proportional hazard regression was conducted to assess the hazard ratios (HRs) and 95% confidence intervals (CIs) of PD incidence associated with annually updated exposure to NDVI from the previous 1 year at baseline. The time-varying model reflects the variability in exposure that changes over time. Survival time was measured in years from baseline (1 January 2007) until PD incidence date, mortality date, or end of follow-up (31 December 2015). Figure 1 shows the distribution of residential greenness on July 2006 (A) and 2014 (B). Thus, the time-varying model was suitable to reflect the distribution of green areas that change over time. We applied the NDVI as a continuous and categorical variable. Based on the quartile of previous 1-year NDVI exposure at baseline, we estimated the HRs of each quartile compared to the lowest quartile as reference and tested the linear trend using the ordinal rank for each quartile. For continuous NDVI, the HRs by a 0.1 increase in NDVI were estimated. We explored the effects of continuous and categorized NDVI on PD incidence in several models. Model 1 is a crude model adjusted for age and sex. Model 2 was additionally adjusted for individual characteristics such as residential province, household income (three levels: low, medium, and high), insurance type (four levels: self-insured, dependent on self-insured, employer-insured, and dependent on employer-insured citizens), alcohol consumption (three levels: rarely, sometimes, and frequently), smoking status (three levels: never, former, and current), exercise (three levels: rarely, sometimes, and frequently), body mass index (BMI), fasting glucose, total cholesterol, frequency of healthcare utilization, history of traumatic brain injury and mental disorders, Charlson comorbidity index (CCI), and medication history (ibuprofen, antipsychotics, and hormone). Model 3 is a fully adjusted model further adjusted by regional characteristics such as percentage of highly educated (high school graduation or over) people in the 2000 census, percentage of elderly (≥65 years) people in the 2000 census, gross regional domestic product in 2005, type of area (three levels: urban, suburban, rural), population density in 2008, and mean concentration of particles <10 µm (PM 10 ). Personal exposure to the annual average of PM 10 was estimated in the district (Si-Gun-Gu) unit by inverse distance weighting (IDW), which considers the square of the inverse distance between monitoring stations as weight [23]. IDW is one of the most widely used spatial interpolation methods that predicts the unmeasured values using the nearest points to the location, and is suitable because monitoring stations are non-equally distributed nationwide.
To investigate more affected subgroups in the association, we conducted stratified analyses with exposure to greenness as continuous variable by key confounding factors [1] including sex, age (<50, ≥50 years), BMI (<25, ≥25 kg/m 2 ), residential province (seven metropolitan cities and nine provinces), smoking status (never, ever), alcohol consumption (rarely, sometimes/frequently), and comorbid states (CCI 0, 1-17). In addition, several sensitivity analyses were conducted to test the robustness of the long-term effects of NDVI in a continuous term on PD in the main model. First, we assessed the effects of different exposure windows of NDVI from the previous 2 to 5 years. Second, we applied conventional Cox proportional hazard regression models with baseline NDVI exposure from the previous 1 year to 5 years. Third, we restricted the study population to non-movers during the study period. Finally, we adjusted for nitrogen dioxide (NO 2 ), sulfur dioxide (SO 2 ), and carbon monoxide (CO), instead of PM 10 , to examine possible confounding effects of various air pollutants.
All analyses were conducted using SAS Enterprise Guide version 7.2 (SAS Institute Inc., Cary, NC, USA) and R Studio version 1.0.136 (packages survival and survminer).

Results
Among 313,355 participants over 2,745,389 person-years of follow-up, we observed 2621 incident PD patients between 2007 and 2015 ( Table 1). The mean (standard deviation; SD) age of the total population was 48.9 (14.2) years, and 35.6% of participants were in the 35-49 age group. More than half (54.1%) were male and less than half (48.5%) were employer-insured citizens, most of them rarely consumed alcohol (71.5%), never smoked (66.0%), and rarely exercised (81.5%). Compared to the study participants in the lowest quartile of greenness, those in the highest quartile of greenness were older, more likely to consume alcohol, exercise frequently, and had more incident PD events and comorbid diseases, and lived in provinces with distinct regional characteristics (higher proportion of the elderly, lower proportion of highly educated people, low gross regional domestic product, and population density). Table 2 presents the summary statistics of the NDVI exposure and air pollution concentrations of the study participants for the first study year, 2007. The mean (SD) of NDVI exposure for the previous 1 year was 0.41 (0.11), and those for the previous 2-5 years were similar, ranging from 0.40 (0.11) to 0.42 (0.12). The 1-year average concentrations of PM 10 , SO 2 , NO 2 , and CO were 59.25 µg/m 3 , 5.76 ppb, 27.06 ppb, and 0.82 ppm, respectively.   Table 3 shows the estimated HRs (95% CIs) of PD incidence associated with exposure to continuous and categorized NDVI in the previous 1 year. In Model 3, where individualand area-level covariates were fully adjusted, a 0.1 unit increase in NDVI was significantly associated with a 21% decrease in the risk of incident PD (95% CI, 0.74-0.84). In the categorized NDVI analysis, exposure to the highest quartile of NDVI compared to the lowest quartile was associated with a 36% decreased risk of incident PD (95% CI, 0.52-0.77), showing a linear trend (p for trend = 0.01). In Model 1 and Model 2, where individual and/or area-level covariates were not adjusted, the association between NDVI and PD was negative but not significant. Table 3. Association between exposure to normalized difference vegetation index and incidence of Parkinson's disease in 1-year time-varying Cox hazards models. In stratified analyses, we identified more affected subgroups per a 0.1-unit increase in NDVI (Table 4)   Our main findings were that the significant inverse associations between residential greenness and incident PD were robust to several sensitivity analyses. The analyses with longer exposure windows of NDVI yielded consistent results with the main findings (Table S1), while the conventional time-constant Cox analysis (satisfied proportional hazards assumption; p > 0.05) resulted in decreased effect estimates (HR 0.92-0.93) ( Table  S2). The analysis restricted to 173,661 non-movers showed similar results to the main results (Table S3). Finally, the analyses adjusting for SO 2 , NO 2 , and CO instead of PM 10 also showed consistent results with the primary findings (Table S4).

Discussion
In this nationwide population-based cohort study, we investigated the association between long-term exposure to residential greenness and incident PD events after adjusting for demographic, lifestyle, physical, laboratory, and regional characteristics, as well as exposure to ambient air pollution. We found that exposure to higher levels of residential greenness significantly decreased the risk of PD. Stronger inverse associations were observed in participants aged over 50 years, females, overweight/obese participants, non-urban residents, never-smokers, alcohol consumers, and those with comorbidities. Our primary findings were consistent with those of various sensitivity analyses To the best of our knowledge, our study is the first national-level cohort study on residential greenness and incident PD based on 313,355 nationwide-representative cohort participants. To date, only two epidemiological studies have investigated the association between long-term greenness exposure and PD [16,17]. Using a prospective cohort of 47,516 participants aged over 18 years in Ningbo, China, between 2015 and 2018, a previous study found a significant protective effect of residential greenness on PD incidence, which corresponds to our findings [16]. The effect estimates from the previous study were also similar to our findings in both the linear NDVI (HR per IQR [0.25] increment: 0.80, 95% CI: 0.65-0.98) and the categorized NDVI analyses (HR for Q4 [>0.51] vs. Q1 [≤0.26] of greenness within radii of 300 m around residence: 0.57, 95% CI: 0.37-0.85), although our effect estimates were slightly higher. A possible reason for the slight difference may be the different spatial resolution used to define residential greenness; our study used districtlevel exposure based on a 250 m buffer, whereas the previous study used a 300 m buffer around residences. Additionally, differences in several factors-including vegetation status , study area, and study populationmight also contribute to the slightly higher effect estimates in our study. Another cohort study of 634,432 participants aged 45-85 years in Metro Vancouver, Canada between 1999 and 2003 reported a marginal protective effect of greenness in 100 m radius buffer on PD incidence (HR per IQR [0.11] increment: 0.96, 95% CI: 0.92-1.01), which was consistent with our findings [17]. Due to restricted data availability, the 1999-2002 greenness average was assigned as the baseline (1999) exposure and a time-constant Cox model was utilized, which could have resulted in increased exposure measurement error. Moreover, the annual greenness exposure in the previous study ranged from −0.11 to 0.78, while that in our study ranged from 0.16 to 0.64.
Recent studies have studied the neurological effects of residential greenness and reported significant neuroprotective effects, including reduced depression [13], slower cognitive decline [14], increased cortical thickness [15], reduced mortality from neurological diseases [11], and reduced risk of AD and non-AD dementia [12]. Although this study provides evidence about the neurological effect of residential greenness based on PD incidence, there is still a dearth of information, requiring further studies to elucidate the effects of greenness on neurological disorders.
The biological pathways of residential greenness in neurological disorders, including PD, remains unclear; however, greenness exposure would be associated with a lower chance of developing PD. Greenness may encourage higher levels of moderate or vigorous physical activities [24], such as walking or jogging, which promotes neural growth factors, contributing to the survival of dopaminergic neurons [25]. In addition, green infrastructure could mitigate the impact of atmospheric particles, which causes neurotoxicity and neurodegeneration through oxidative stress by dry deposition onto the leaves of vegetation [26,27]. Given the limited evidence on the biological mechanisms, future studies focusing on the biological mechanisms of nervous system diseases associated with greenness exposure should be conducted.
Our findings of the stronger inverse associations in older participants (≥50 years), overweight/obese participants, alcoholics, and those with comorbidities can be explained by the susceptibility of these subgroups, implying that exposure to greenness could be more advantageous to sensitive subgroups. In the sex-specific analysis, we observed a stronger protective effect in females, which is consistent with the findings of a previous PD study in Ningbo, China [16]. Studies on other health outcomes besides PD found inconsistent results regarding the sex-specific effects of greenness; some studies reported a stronger effect in females [28][29][30], while several studies reported a stronger association among males [31,32]. The difference might be explained by different behavioral patterns in which males are more likely to be more active than females in green spaces [33], but this could not be verified in this study. Further studies on different behavioral patterns, as well as biological differences among males and females, could help to elucidate the effects of greenness on health, including PD.
The stronger protective effect of greenness in never smokers than ever smokers in our study is consistent with a recent study from Ningbo, China. [16]. The study found that the highest quartile of greenness exposure compared to the lowest quartile was significantly associated with a lower risk of incident PD in never smokers (HR, 0.53; 95% CI: 0.34-0.83) but not in ever smokers [16]. Similarly, a longitudinal cohort study on the association between residential greenness and mortality reported a stronger association in never smokers (HR: 1.09, 95% CI 1.04-1.14) than in smokers (HR: 1.03, 0.91-1.17) [34]. In contrast, smokers benefitted more from the effects of greenness than never smokers with decreased blood pressure [35] and the prevalence of chronic kidney disease [36]. Therefore, both smokers and never-smokers can benefit from exposure to residential greenness.
We found a noticeable protective effect of greenery on PD in provinces than in metropolitan cities. Aside from Parkinson's disease, for which no previous research has been done, various studies have looked into the relationship between urbanity and cancer mortality [37], general health [38], hypertension [39], and cardiovascular mortality [40]. Research on cancer mortality [37] and general health [38] indicated a stronger effect of residential greenness in urban areas than in rural regions, whereas studies on hypertension [39] and cardiovascular mortality [40] found a stronger effect of residential greenness in rural areas. The possible urban-rural differences can be explained by several factors, including the amount of residential greenness [41], accessibility to green space [42], types of vegetation, and structure of green space [43], which need to be validated in further research.
Our study has several limitations. First, our exposure assessment based on districtlevel residential addresses may not reflect accurate individual exposure and lead to exposure misclassification, despite the fact that we obtained high spatial resolution greenness measurements from satellites. However, district-level residential greenness assessments may be more accurate in capturing participants' travel from home to work or other activities. Moreover, the NDVI captured by the satellite could not indicate vegetation type [44] or green structure [45] which could be possible confounders in the association between greenness and incident PD; however, NDVI could be an indicator of objectively measured vegetation density. Second, other potential confounders-such as occupation, education level, genetic factors, and caffeine intake-could not be adjusted due to limited data availability, although we extensively adjusted various sociodemographic and laboratory data as well as area-level characteristics. Lastly, the outcome ascertainment using diagnostic codes might have caused outcome misclassification; however, the possibility of misclassification is low because, in South Korea, the government pays 90% of medical expenses for PD, and PD patients must obtain confirmatory diagnosis of PD using imaging examinations to be eligible for these benefits.
Despite several limitations, our research has a number of advantages. First, we conducted a national-level cohort analysis based on the government's largest and most trustworthy national claims database, which contains actual healthcare utilization of Korean citizens. Second, important individual-level covariates were addressed in this study, such as medication history, lifestyle (smoking/drinking/exercise), laboratory findings, and concomitant disorders. Third, we employed spatial interpolation of air pollution to overcome the limitations of monitoring stations concentrated in urban areas to adjust for the primary confounding factor, long-term air pollution exposure. Finally, the results of various stratified studies by important features could aid in the identification of more benefited subpopulations and the development of effective health policies.

Conclusions
With no clear treatment for PD, it is important to identify the factors which protect against developing PD. Our study provides scientific evidence that living in higher levels of residential greenness could be beneficial for lowering the risk of PD incidence. With the rapid urbanization and rapid aging of population worldwide, our findings should be useful for policy makers and city planners to establish appropriate public health policies and interventions.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/ijerph19063491/s1, Figure S1: Flow chart of study population; Table S1. Hazard ratios of long-term exposure to normalized difference vegetation index and on Parkinson's disease in time-dependent Cox proportional hazards model and Cox proportional hazards model using various exposure windows; Table S2. Hazard ratios of long-term exposure to normalized difference vegetation index and on Parkinson's disease in conventional Cox proportional hazards model using various exposure windows; Table S3. The association between 0.1 increase of normalized difference vegetation index and incidence of Parkinson's disease among non-movers; Table S4. Hazard ratios of long-term exposure to normalized difference vegetation index and on Parkinson's disease according to adjustment of different air pollutants