Residential Greenness and Cardiovascular Disease Incidence, Readmission, and Mortality.

BACKGROUND: Living in greener areas of cities was linked to increased physical activity levels, improved mental well-being, and lowered harmful environmental exposures, all of which may affect human health. However, whether living in greener areas may be associated with lower risk of cardiovascular disease incidence, progression, and premature mortality is unclear.
OBJECTIVES: We conducted a cohort study to examine the associations between residential green spaces and the incidence of acute myocardial infarction (AMI) and heart failure (HF), post-AMI and HF hospital readmissions, and mortality.
METHODS: We simultaneously followed four large population-based cohorts in Ontario, Canada, including the entire adult population, adults free of AMI and HF, and survivors of AMI or HF from 2000 to 2014. We estimated residential exposure to green spaces using satellite-derived observations and ascertained health outcomes using validated disease registries. We estimated the associations using spatial random-effects Cox proportional hazards models. We conducted various sensitivity analyses, including further adjusting for property values and performing exploratory mediation analysis.
RESULTS: Each interquartile range increase in residential greenness was associated with a 7% [95% confidence interval (CI): 4%, 9%] decrease in incident AMI and a 6% (95% CI: 4%, 7%) decrease in incident HF. Residential greenness was linked to a ∼10% decrease in cardiovascular mortality in both adults free of AMI and HF and the entire adult population. These associations remained consistent in sensitivity analyses and were accentuated among younger adults. Additionally, we estimated that the decreases in AMI and HF incidence associated with residential greenness explained ∼53% of the protective association between residential greenness and cardiovascular mortality. Conversely, residential greenness was not associated with any delay in readmission or mortality among AMI and HF patients.
CONCLUSIONS: Living in urban areas with more green spaces was associated with improved cardiovascular health in people free of AMI and HF but not among individuals who have already developed these conditions. https://doi.org/10.1289/EHP6161.

Introduction

Globally, cardiovascular disease (CVD) is a major public health concern, affecting people and resulting in 18 million deaths each year (WHO 2013a). Although CVD mortality has declined in many high-income countries during the past 20 y owing to advances in public health and medical care, the decreases have plateaued in recent years (Roth et al. 2017). Meanwhile, only small improvements in cardiovascular mortality have occurred in low-income regions (Roth et al. 2017). Thus, addressing the continuing threat of CVD to public health has been a global priority (WHO 2013b).

Over the last few decades, a growing body of literature has suggested that living in greener areas of cities may increase opportunities for improving physical activity levels (Fong et al. 2018; Markevych et al. 2017), a key target for individuals to prevent the onset and complications of CVD (Myers 2003). More recent studies further suggested that living in greener, natural environments may contribute to increased social interactions, improved mental well-being, and reduced exposure to harmful environmental stressors such as noise, heat, and air pollution (Fong et al. 2018; Markevych et al. 2017), all of which are likely to be conducive to human health. Indeed, a number of panel studies have found positive physiological responses, such as improved heart rate (Lanki et al. 2017; Ottosson and Grahn 2005), reduced stress (Aspinall et al. 2015; Gidlow et al. 2016; Honold et al. 2016), and enhanced parasympathetic nerve activity and restoration (Lee et al. 2011; Park et al. 2007) in participants after visiting parks and other green spaces.

To date, however, whether increasing urban greenness may be associated with reduced risk of cardiovascular morbidity and mortality is not well understood. A small number of population-based studies have been conducted to examine the relationship between residential exposure to urban greenness and cardiovascular mortality (Crouse et al. 2017; James et al. 2016; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017); some reported a protective relationship (Crouse et al. 2017; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017), but others did not find an association (James et al. 2016). Evidence is sparser about the associations between greenness and cardiovascular incidence and progression. The few existing studies were limited by a lack of comprehensive CVD outcome data and the spatially resolved estimates of exposure to greenness, small sample size, and cross-sectional study design (Donovan et al. 2015; Pereira et al. 2012; Seo et al. 2019; Tamosiunas et al. 2014). With continuing urbanization globally and the ever-growing burden of CVD, understanding of the influence of urban greenness on cardiovascular morbidity and mortality, especially how it may unfold along different stages of cardiovascular health, can have a broad public health implication. We thus conducted a population-based cohort study to examine the relationships between residential greenness and the incidence of two leading CVDs [acute myocardial infarction (AMI) and heart failure (HF)], their progression, and ultimately death in Ontario, the most populous province in Canada.

Methods

Study Design and Population

We constructed four population-based cohorts in Ontario, including: a) the entire adult population of Ontario (referred to full cohort); b) the entire adult population of Ontario who had no history of physician-diagnosed AMI and HF before cohort inception (referred to as incidence cohort); and c) survivors of AMI (AMI cohort) or HF (HF cohort) (Figure S1). The full cohort was used to facilitate comparison of this study with previous population-based studies that investigated urban greenness and mortality (Crouse et al. 2017; James et al. 2015; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017). We restricted the four cohorts to urban residents who were 35–100 years of age at baseline, lived in common postal-code areas shared by cohorts, and were followed up simultaneously over time. The four cohorts consisted of participants at different stages of cardiovascular health, but they were exposed to similar green environments. This design allowed us to better disentangle the potentially different impacts of greenness on various cardiovascular outcomes.

To establish the full cohort comprising all Ontario adults, we used data from the Ontario Population Health and Environment Cohort (ONPHEC), a large population-based comprising all Ontario residents who, on 1 April 1996 and onward, were registered with provincial health insurance and lived in Ontario for (Chen et al. 2016a). ONPHEC was constructed through record linkage of population-based health administrative databases (e.g., hospitalizations, emergency department visits, physician office visits) (Chen et al. 2016a). From the full cohort, we further restricted to those who were free of physician-diagnosed AMI and HF at baseline, yielding the incidence cohort. The follow-up of both cohorts extended from 1 January 2000 through 31 December 2014.

To establish the AMI and HF cohorts, we used data from the Enhanced Feedback For Effective Cardiac Treatment (EFFECT) study, a large randomized trial in Ontario (Tu et al. 2009). This EFFECT study comprised all patients admitted to 86 Ontario hospitals with AMI or HF during two periods (April 1999 to March 2001 and April 2004 to March 2005) (see Supplemental Material for more detail). An important feature of the EFFECT study is that trained nurses extracted detailed clinical information (e.g., laboratory tests, treatment, and medical history) from medical charts of all study participants. Because long-term outcomes after AMI or HF were strongly influenced by clinical severity and quality of care, using data from the EFFECT study allowed for better controlling for these potential confounders, which were otherwise unavailable from using routine health administrative databases. Like previous studies (Berglind et al. 2009; Chen et al. 2016b; Tonne and Wilkinson 2013), we restricted AMI and HF patients to those who had survived for at least 28 d after index hospital discharge (baseline). Thus, for the AMI and HF cohorts, the follow-up extended from the 29th day after the index hospital discharge until 31 December 2014.

All four cohorts were linked using unique encoded personal identifiers and analyzed at ICES. Participants living in residential postal-code areas not shared by all four cohorts were excluded. The use of data in this analysis was authorized under section 45 of Ontario’s Personal Health Information Protection Act, which does not require review by a research ethics board.

Outcomes

The study outcomes of interest were the incidence of AMI and HF [AMI: International Classification of Diseases, Ninth Revision (ICD-9) code 410 and ICD-10 code I21-22; HF: ICD-9 code 428 and ICD-10 code I50], first readmission after the index hospital discharge from any cardiovascular cause (ICD-9: 401-459 and ICD-10: I10-I99), and death from any cardiovascular cause (ICD-9: 401-459 and ICD-10: I10-I99). To assess potential bias because of unmeasured confounding and other errors that may lead to spurious inference, we also ascertained nonaccidental mortality (ICD-9: and ICD-10: A00-R99) as a positive control outcome in this study (Bowe et al. 2018) because nonaccidental mortality has been regularly investigated by previous studies of urban greenness (Crouse et al. 2017; James et al. 2015; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017).

To ascertain incident AMI and HF, we used record linkage to disease registries of these two conditions in Ontario. The two disease databases have been previously validated using chart reviews and shown to have high sensitivity (85%–89%) and specificity (93%–97%) (Austin et al. 2002; Yeung et al. 2012). We also used record linkage to the Ontario Registrar General’s Death database to ascertain deaths and the Canadian Institute for Health Information’s hospital discharge abstracts for hospital readmissions. Prevalent cases of AMI and HF diagnosed before the baseline in 2000 were excluded from the incidence cohort, but for other three cohorts, they were retained as comorbidities.

Green Spaces

Exposure to urban green spaces was estimated using the satellite-derived Normalized Difference Vegetation Index (NDVI). Calculated as a function of vegetative absorption of infrared light and reflection of near-infrared light, NDVI is an objective measure of vegetative greenness that has been widely used in epidemiological studies of urban green spaces (Crouse et al. 2017; James et al. 2016; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017). NDVI values range from to 1, with negative values corresponding to bodies of water and positive values corresponding to green vegetation.

We obtained NDVI measures using all cloud-free images from the Landsat 5 Thematic Mapper for 1996–2012 (available until 2012). We removed values for water because they would artificially lower NDVI levels (Crouse et al. 2017; Hystad et al. 2014; Markevych et al. 2017). Following the methods shown in previous studies (Crouse et al. 2017; Hystad et al. 2014; James et al. 2016), we calculated NDVI annual measures for a buffer area around each residential postal code across Ontario over the periods of 1996–2000, 2001–2005, 2006–2010, and 2011–2012. This overcame the limitation of relatively low temporal resolution of the Landsat images, especially under cloud-free conditions in Ontario. As a sensitivity analysis, we also estimated annual averages of NDVI for , , and buffer areas (Crouse et al. 2017). We noted that annual averages of NDVI and maximum green season NDVI in Ontario were highly correlated (Pearson’s correlation coefficient , depending on buffer size).

To account for residential mobility and long-term changes in urban greenness, we assigned the closest annual average measures of NDVI according to participants’ yearly postal-code residences during the follow-up. For example, a participant’s exposure for 2000 would be estimated using annual average measures of NDVI for the period 1996–2000. Participants’ annual postal codes were obtained from the Ontario Registered Persons database, a registry of all Ontario residents.

Covariables

For all cohorts, we considered a priori age at baseline, sex, and seven time-varying variables for neighborhood-level sociodemographic status (SES): a) percentage of population age with less than high school education; b) percentage of recent immigrants; c) unemployment rate; d) income quintile, a measure of relative household income accounting for household size; e) deprivation based on the Ontario Marginalization Index that quantifies the degree of marginalization in health and social well-being (Matheson et al. 2012); f) population density (individuals per square kilometer); and g) average property values, using Canadian census dissemination area data from the closest census year (Table S1). A dissemination area (with 400–700 people) is the smallest census geographic area for which census data are disseminated. In addition, we derived area-level active commuting and transit use using census tract data, because active commuting is known to improve fitness and cardiovascular health (Celis-Morales et al. 2017). Additionally, we created variables for residential proximity to the nearest acute-care hospital (in meters) and access to primary care based on density of family physicians (per 10,000 individuals, at the census subdivision level) (see Supplemental Material for descriptions about census geographic areas). Furthermore, to control for regional differences in mortality and hospitalization, we created a dichotomous variable classifying Ontario into the Greater Toronto Area, a densely populated urban mega-region, and all other areas. Toronto tends to differ from other areas in Ontario with respect to socioeconomic and demographic characteristics, health care, and mortality rate (Chen et al. 2016b).

In addition, we obtained the annual mean concentrations of fine particulate matter () and nitrogen dioxide () at participants’ residential postal codes in each year during the follow-up using the satellite-based estimates of (van Donkelaar et al. 2015) and a national land-use regression (LUR) model for (Hystad et al. 2011). Briefly, estimates of ground-level concentrations of were derived from satellite retrievals of aerosol optical depth, in conjunction with outputs from a global atmospheric chemistry transport model (GEOS-Chem CTM) (van Donkelaar et al. 2015). The estimates were further calibrated using the information on urban land cover, elevation, and aerosol composition using a geographically weighted regression, thus producing annual mean concentration of ground-level at a spatial resolution of for each year between 1998 and 2012. These annual estimates of have been shown to closely agree with ground measurements at fixed-site monitoring stations across North America (, ) (van Donkelaar et al. 2015). Similarly, to derive exposure to , we made use of a national LUR model developed using measurements of at the fixed-site stations of Environment Canada’s National Air Pollution Surveillance Network (Hystad et al. 2011). This model was constructed by regressing measured annual mean concentrations of in Canada against an array of predictors (e.g., satellite estimates of , area of industrial land use, road length) to capture background and regional variations of (Hystad et al. 2011). The estimates were then augmented by incorporating local-scale variations of caused by vehicle emissions by applying spatially varying multipliers that represented distance–decay gradient in . The LUR model explained 73% of the variability in annual 2006 measurements of , with a root mean square error of 2.9 ppb (Hystad et al. 2011). The resulting LUR estimates were available for each year between 1996 (5 y before baseline) and 2014, after applying temporal adjustment, as done previously (Chen et al. 2017). Similarly, we applied the temporal calibration to produce annual exposure surfaces between 1995 and 1997 and between 2013 and 2014, respectively.

Furthermore, we ascertained the presence of comorbidities at baseline, including hypertension, diabetes, stroke, asthma, chronic obstructive pulmonary disease (COPD), and dementia, using validated chronic disease databases in Ontario (see Table S2 for ICD-9 and ICD-10 codes). These databases were created using hospital discharge abstracts from the Canadian Institute for Health Information, physician service claims from the Ontario Health Insurance Plan database, and emergency department records from National Ambulatory Care Reporting System (NACRS), and they were validated using chart reviews, which showed sensitivity of 72%–89% and specificity of 76%–100%, depending on the disease (Austin et al. 2002; Gershon et al. 2009a, 2009b; Hux et al. 2002; Jaakkimainen et al. 2016; Tu et al. 2008, 2013; Yeung et al. 2012). In addition, we ascertained prior history of cancer using the Ontario Cancer Registry, which captures diagnostic and demographic information on all Ontario residents who were diagnosed with cancer.

For the AMI and HF cohorts, we further retrieved information from medical charts about in-hospital care, clinical severity, cardiovascular medications (e.g., statins) at discharge, family history of coronary artery disease (for AMI patients only), and history of atrial fibrillation (for HF only). To characterize the clinical severity of AMI patients, we assessed the type of AMI [ST-elevated MI (STEMI) vs. non-STEMI] and derived the Global Registry of Acute Cardiac Events (GRACE) risk score based on age, heart rate, systolic blood pressure, and several other prognostic variables (Bradshaw et al. 2006). Similarly, for HF patients, we assessed the type of HF (ischemic etiology or not) (Chun et al. 2012) and derived the EFFECT-HF mortality risk score based on age, respiratory rate, serum sodium level, and other prognostic variables (Lee et al. 2003). To characterize in-hospital care for both patient cohorts, we obtained data on the length of hospital stay (days), attending physician specialty (general/cardiologist/internist), and hospital type (teaching/community/small). Additionally, we extracted information from medical charts on smoking status (current smoker or not), employment (employed or not), marital status, hyperlipidemia, pulmonary edema, and any history of percutaneous coronary intervention. Last, we obtained information about coronary revascularization during follow-up through data linkage to hospital discharge abstracts and physician service claims.

Statistical Analysis

We used two-level spatial random-effects Cox proportional hazards models to estimate the association between residential green space and each outcome separately. Briefly, a multilevel random-effects Cox model can be expressed as: , where represents a random effect associated with the th cluster (Austin 2017). It entails that all individuals within the same cluster share the same random effect, and it has been extended to allow for a nested set of clusters (Therneau 2020). Like previous studies (Chen et al. 2016b; Crouse et al. 2012), we defined a first cluster level by census divisions (equivalent to counties) and a second level defined by census tracts within census divisions. We assumed that any two census divisions were independent, as were census tracts within each census division. This approach accounted for spatial clustering in the health status of participants from the same communities.

We developed separate models for each cohort and outcome. In all models, time variable was survival time from the baseline. Study participants who did not have the outcome of interest were censored at the time of death, ineligible for provincial health insurance (i.e., migration out of Ontario) (Table S3), or end of study (31 December 2014). We also allowed the baseline hazard function to differ by region and by the period of EFFECT enrollment (for the AMI and HF cohorts only). To assess chronic greenness exposure, we estimated 5-y running average exposure for the full and incidence cohorts. Conversely, for the AMI and HF cohorts, we modeled time-weighted exposure from baseline (i.e., hospital discharge) until outcome to assess postdischarge exposure, with weights for each patient defined by the time spent at each residence. We included participants with nonmissing information on exposure and covariates, except for marital status and employment status, for which we created a separate category of missing values to avoid losing substantial statistical power.

We grouped similar covariates together and added them incrementally to obtain fully adjusted Cox models. Our final models included age, sex, and time-varying area-level SES variables (i.e., education, income, unemployment, percent of immigrants, and population density). For the AMI and HF cohorts, we further adjusted the final models for clinical severity, in-hospital care, medications at discharge, smoking, and individual-level SES variables (marital and employment status). For all models, we tested the proportional hazards assumption by adding the cross-product of each variable with the natural logarithm of the time variable, but we did not find any violation of this assumption (). Hazards ratios () were scaled to an interquartile increase (IQR) in NDVI for all models.

Sensitivity Analyses

To evaluate the robustness of our effect estimates, we performed sensitivity analyses by further adjusting for exposure to and , comorbidities, and access to primary care; adjusting for proximity to acute-care hospitals using a natural spline with three degrees of freedom; and additionally adjusting for neighborhood deprivation. We estimated exposure to each air pollutant using the same exposure time window as that of greenness for each outcome. In addition, we adjusted for neighborhood-level active commuting and average property values, and restricting to participants who were followed for due to a concern that short follow-up may limit our ability to detect the effects. As well, to assess whether HRs were sensitive to the buffer size of NDVI, we repeated analysis using alternative NDVI measures with , , and buffer areas, respectively. Furthermore, we carried out a series of sensitivity analyses without restricting to common postal code areas shared by the cohorts.

For the AMI and HF cohorts, we additionally adjusted for coronary revascularization and excluded because of their limited exposure patients who were discharged from a hospital and subsequently died within the same winter (November to March). For the full and incidence cohorts, we further assessed the potential influence of unmeasured individual-level sociodemographic and behavioral factors, especially education and smoking, on our effect estimates. To do this, we used a method to mathematically adjust the HRs for education and smoking while simultaneously adjusting for all variables available in the model (i.e., age, sex, and SES) (Shin et al. 2014). Briefly, this method requires spatial associations between the unmeasured (i.e., smoking and education) and observed variables (e.g., age, sex, SES) from an auxiliary data set. Following previous Canadian studies using this method (Chen et al. 2017, Crouse et al. 2015), we obtained the relationships using data from Ontario respondents to the 2000/2001, 2003, and 2005 cycles of the Canadian Community Health Survey. These population-based surveys are conducted routinely by Statistics Canada and cover nearly all household residents age 12 y or older in Ontario and other provinces (Statistics Canada 2019). This information along with estimated associations between the unmeasured variables and incident AMI and HF and mortality from the literature were used to estimate their impact on the HRs. Based on the strength of evidence from recent systematic reviews of these outcomes (Table S4), we adjusted for smoking and education for the HRs for green spaces with incident AMI and HF and mortality from any cardiovascular and nonaccidental causes, respectively.

Further Analyses

In addition, we conducted a post hoc mediation analysis to examine the proportion of the association between residential green space and cardiovascular mortality that was mediated through the development of AMI and HF using the incidence cohort, following VanderWeele (2011). The mediator was defined as the occurrence of either AMI or HF during the follow-up. The outcome and mediator models were fitted using the random-effect Cox models adjusting for all available covariates, as in the primary analysis. We assumed no unmeasured exposure–outcome, mediator–outcome, and exposure–mediator confounding, as well as no exposure–mediator interaction (VanderWeele 2011). We computed proportion mediated and their 95% CIs using bootstrapping procedures, as done previously (VanderWeele 2016).

Furthermore, we examined whether age, sex, income quintile, disease severity (for the AMI and HF cohorts only), and smoking (for the AMI cohort only) may modify the HRs, by assessing whether the interaction term that was the cross-product of each variable with NDVI value was statistically significant ().

We used R statistical package (version 3.3.1) with the Coxme library for conducting statistical analysis and ArcGIS (version 10.4) (ESRI) for mapping.

Results

Study Participants

Among the entire population of long-term adult residents of Ontario, we made the following restrictions: a) those living in postal-code areas common to all four study cohorts with and nonmissing data on area-level SES variables and and , yielding the full cohort ( adults); b) further restricting to those free of both AMI and HF (incidence cohort, ); and c) those discharged alive with AMI or HF in the EFFECT study (AMI cohort: 10,790 and HF cohort: 10,676) (Figure S1 and Table S5). Most study participants lived in a narrow band above the Canada–U.S. border, at approximately the 49th parallel, following the pattern of the majority of Canadians (see a map of study area in Figure S2). In addition, of all four cohorts, the distribution of the size of their postcodes ranged from 0.0019 to (25th percentile: ; 75th percentile: ; and median: ).

Of the incidence and full cohorts, the mean age was at baseline, 47% were male, 30% had hypertension, and 10% had COPD (Table 1, Table S6, and as a comparison, baseline description of residents excluded from the full cohort is presented in Table S7). In comparison, the AMI and HF patients were considerably older and had more chronic conditions at baseline (Table 2). Of the AMI cohort, 4,419 patients were readmitted to hospitals for any cardiovascular cause and nearly half of the cohort died during the follow-up (mean follow-up: 9.4 y). Of the HF cohort, 5,482 were readmitted, and most of the cohort died during the study (mean follow-up: 5.4 y).

Table 1

Baseline characteristics of two population-based cohorts of urban residents in Ontario, Canada in 2000 (percent or count, otherwise specified).

Baseline characteristics	Full cohorta (n=1,362,723)		Incidence cohorta (n=1,290,288)
	Subject count	%	Subject count	%
Age at entry [y (mean±SD)]	56.3±14.8	—	55.3±14.3	—
Sex
Male	637,754	47	602,564	47
Female	724,969	53	687,724	53
Preexisting comorbidityb
Diabetes	124,008	9	103,223	8
Hypertension	402,004	30	350,959	27
Acute myocardial infarction	27,255	2	—	—
Heart failure	54,509	4	—	—
Stroke	38,156	3	28,386	2
Asthma	29,980	2	27,096	2
Chronic obstructive pulmonary disease	138,998	10	114,836	9
Dementia	25,892	2	19,354	2
Cancer	83,127	6	73,547	6
Area-level risk factorsc
Income quintile
Lowest	246,244	18	229,155	18
Lower middle	273,907	20	257,671	20
Middle	283,447	21	268,896	21
Upper middle	276,360	20	263,864	20
Upper	281,947	21	270,057	21
% of recent immigrants	—	4	—	4
% population age >15 y with under high school education	—	35	—	34
% population age ≥15 y without employment	—	9	—	9

Note: —, no data; AMI, acute myocardial infarction; HF, heart failure; SD, standard deviation.

Full cohort comprised all urban residents ages 35–100 y in Ontario who lived in the same postal-code areas as in two patient cohorts with AMI or HF (see Table 2). Incidence cohort comprised all individuals in the full cohort who were free of any AMI and HF.

During the past 10 y before cohort inception.

From the 2001 Canadian Census, at the census dissemination area level.

Table 2

Baseline characteristics of two cohorts of patients diagnosed with acute myocardial infarction or heart failure, respectively, enrolled in the EFFECT study, and lived in an urban area across Ontario, Canada (percent or ).

Baseline characteristics	AMI cohort	HF cohort
	(n=10,790)	(n=10,676)
Demographic characteristics
Age at entry [y (mean±SD)]	67.3±13.3	75.9±11.2
Male	65	50
Marital status
Married	63	46
Single, separated, divorced, or widowed	33	49
Unknown	4	5
Employment status
Employed	27	6
Not employed	66	90
Unknown	7	4
Clinical severity
Acute pulmonary edema	5	20
Patients with acute myocardial infarction
GRACE risk score	142±36	—a
ST elevation myocardial infarction (STEMI)	43	—
Patients with heart failure
EFFECT-heart failure mortality risk score	—	87±25
Ischemic etiology	—	60
In-hospital care
Length of stay (d)	7.6±6.9	8.1±7.9
Specialty of attending physician
Cardiology	49	27
Internal medicine	28	32
General practice	23	41
Characteristics of hospitals
Teaching	23	24
Community	76	75
Small	1	1
Cardiovascular medication at dischargeb
Statins	54	16
Aspirin	86	40
ACE inhibitor	62	64
Beta-blockers	74	28
Cardiac risk factors and history
Current smoker	32	12
Family history of coronary artery disease	33	—b
Preexisting comorbiditiesc
Diabetes	25	36
Hypertension	52	59
Hyperlipidemia	39	30
Atrial fibrillation	—	33
Acute myocardial infarction	23	38
Stroke	8	17
Previous percutaneous coronary intervention	4	5
Asthma	4	7
Chronic obstructive pulmonary disease	23	42
Dementia	4	8
Cancer	12	18
Area-level risk factorsd
Income quintile
Lowest	20	24
Lower middle	22	22
Middle	21	20
Middle upper	19	17
Highest	18	16
% of recent immigrants	4	4
% population age >15 y with less than high school education	32	29
% population age >15 y without employment	8	7

Note: EFFECT cohorts comprised all patients who were discharged alive with AMI or HF from one of 86 hospitals across Ontario in two periods, 1999–2000 and 2004–2005, respectively. In the present study, we restricted EFFECT cohorts to urban residents ages 35–100 y. AMI, acute myocardial infarction; EFFECT, Enhanced Feedback For Effective Cardiac Treatment study; GRACE, Global Registry of Acute Coronary Syndromes; HF, heart failure; SD, standard deviation.

Not applicable unless specified otherwise. For HF cohort, information for the family history of coronary artery disease was unavailable.

For HF cohort, only data from the first phase of EFFECT study (1999–2000) were collected.

During the past 10 y before cohort inception.

From the 2001 Canadian Census, at the census dissemination area level.

Associations of Urban Greenness with Cardiovascular Incidence and Mortality

From 2000 to 2014, we identified 58,553 and 134,655 incident cases of AMI and HF, respectively (Figure 1). In addition, a total of 277,174 nonaccidental deaths occurred in the incidence cohort and 330,560 in the full cohort, of which were attributed to CVD.

Figure 1.

Associations between green spaces within of study participants’ residences and cardiovascular mortality, nonaccidental mortality, and cardiovascular-related incidence and hospital readmission in four population-based cohorts of urban residents in Ontario, Canada. Note that hazard ratios were scaled to an interquartile increase in NDVI (full and incidence cohorts: ; AMI cohort: ; HF cohort: ). For the full and incidence cohorts, the models adjusted for age, sex, region (lived or not in the Greater Toronto Area), area-level unemployment, percent less than high school education, percent recent immigrants, and household income (quintiles), and population density. For the AMI and HF cohorts, the models further adjusted for clinical severity, in-hospital care, medications at discharge, smoking, and individual-level SES variables. Note: AMI, acute myocardial infarction; CVD, cardiovascular disease; HF, heart failure.

In the fully adjusted models, we observed a of 0.93 [95% confidence interval (CI): 0.91, 0.96] and 0.94 (95% CI: 0.93, 0.96) for the associations between residential greenness and incident AMI and HF, respectively (Figure 1). Increased amounts of greenness were also associated with a of 0.90 (95% CI: 0.88, 0.92) for cardiovascular mortality in the incidence cohort. A similar protective association was observed for the full cohort. In addition, higher levels of greenness were linked to lower nonaccidental mortality in both cohorts (incidence cohort: , 95% CI: 0.88, 0.91; full cohort: , 95% CI: 0.88, 0.90).

Conversely, there was no association between residential greenness and hospital readmissions in AMI and HF patients ( to 1.02, depending on the cohort) (Figure 1). We did not find any strong evidence associating urban greenness with mortality in these two patient cohorts (e.g., cardiovascular death: ; 95% CI: 0.93, 1.05 for AMI cohort and ; 95% CI: 0.94, 1.04 for HF cohort).

The associations between greenness exposure and cardiovascular incidence, readmissions, and mortality were insensitive to further adjustments for air pollution, comorbidities (e.g., diabetes, hyperlipidemia, COPD, and cancer), deprivation, and proximity to hospitals (Table 3). Using the NDVI with alternative radiuses of , , and yielded similar results as those using a radius (Table 3 and Table S8). As well, the magnitude of the associations was not materially altered after further restricting to participants who were followed up for and adjusting for neighborhood-level estimates of property values and active commuting (Table 3). Furthermore, the estimated associations without restricting to common postal code areas shared by cohorts were broadly consistent with the main analysis, although the association with AMI incidence was somewhat attenuated (Table S9).

Table 3

Sensitivity analysis for associations between green spaces and cardiovascular incidence, readmission, and mortality in four population-based cohorts of urban residents in Ontario, Canada, 2000–2014.

Outcome	+NO_2 and PM_2.5c	+Comorbiditiesc	+Proximity to hospitalsc	+Deprivationc	+Active commutingc	+Average property valuec	Restricted to those survivingc >5 y	Within a 500-m buffer
	HR (95% CI)	HR (95% CI)	HR (95% CI)	HR (95% CI)	HR (95% CI)	HR (95% CI)	HR (95% CI)	HR (95% CI)
Full cohorta
Mortality
CVD	0.89 (0.87, 0.91)	0.91 (0.89, 0.92)	0.90 (0.89, 0.92)	0.91 (0.89, 0.92)	0.89 (0.87, 0.91)	0.91 (0.89, 0.92)	0.91 (0.89, 0.93)	0.89 (0.87, 0.90)
NAC	0.90 (0.88, 0.91)	0.91 (0.90, 0.92)	0.90 (0.89, 0.91)	0.90 (0.89, 0.91)	0.90 (0.88, 0.91)	0.90 (0.89, 0.91)	0.91 (0.90, 0.92)	0.89 (0.88, 0.90)
Incidence cohorta
Mortality
CVD	0.89 (0.87, 0.91)	0.91 (0.90, 0.93)	0.90 (0.88, 0.92)	0.91 (0.89, 0.92)	0.89 (0.87, 0.91)	0.91 (0.89, 0.92)	0.91 (0.89, 0.93)	0.89 (0.88, 0.91)
NAC	0.90 (0.88, 0.91)	0.91 (0.90, 0.92)	0.90 (0.89, 0.91)	0.91 (0.89, 0.92)	0.90 (0.88, 0.91)	0.90 (0.89, 0.91)	0.91 (0.90, 0.92)	0.89 (0.88, 0.90)
Incidence
AMI	0.95 (0.92, 0.97)	0.94 (0.92, 0.96)	0.96 (0.94, 0.98)	0.93 (0.91, 0.95)	0.93 (0.91, 0.95)	0.93 (0.91, 0.96)	0.93 (0.91, 0.95)	0.91 (0.89, 0.93)
HF	0.95 (0.94, 0.97)	0.95 (0.93, 0.96)	0.95 (0.93, 0.96)	0.92 (0.91, 0.94)	0.94 (0.92, 0.95)	0.94 (0.93, 0.96)	0.94 (0.92, 0.95)	0.93 (0.92, 0.95)
AMI cohortb
Readmission (mortality)	1.03 (0.98, 1.07)	1.02 (0.98, 1.06)	1.02 (0.98, 1.06)	0.99 (0.94, 1.05)	1.03 (0.98, 1.07)	1.02 (0.98, 1.06)	1.03 (0.98, 1.08)	1.01 (0.96, 1.05)
CVD	0.99 (0.92, 1.05)	0.99 (0.93, 1.05)	0.99 (0.93, 1.05)	0.98 (0.91, 1.05)	0.99 (0.93, 1.05)	0.99 (0.94, 1.05)	0.96 (0.87, 1.04)	0.96 (0.90, 1.02)
NAC	0.99 (0.94, 1.03)	1.00 (0.96, 1.04)	1.00 (0.95, 1.04)	0.98 (0.93, 1.03)	0.99 (0.95, 1.04)	1.00 (0.95, 1.04)	0.98 (0.92, 1.04)	0.99 (0.94, 1.03)
HF cohortb
Readmission (mortality)	1.00 (0.95, 1.06)	0.99 (0.94, 1.04)	0.99 (0.95, 1.04)	0.99 (0.94, 1.04)	1.00 (0.95, 1.04)	0.99 (0.95, 1.03)	0.99 (0.92, 1.06)	0.99 (0.95, 1.04)
CVD	0.99 (0.94, 1.05)	0.98 (0.93, 1.03)	0.98 (0.93, 1.04)	1.01 (0.96, 1.06)	0.99 (0.94, 1.03)	0.99 (0.94, 1.03)	0.99 (0.90, 1.08)	0.98 (0.93, 1.03)
NAC	0.99 (0.95, 1.03)	0.99 (0.95, 1.03)	0.98 (0.94, 1.02)	0.99 (0.95, 1.03)	0.98 (0.95, 1.02)	0.99 (0.96, 1.03)	0.96 (0.90, 1.03)	0.99 (0.95, 1.03)

Note: Each covariate was added individually to the main analysis presented in Figure 1. AMI, acute myocardial infarction; CI, confidence interval; CVD, cardiovascular disease; HF, heart failure; HR, hazard ratio; IQR, interquartile range; NAC, nonaccidental.

Two-level nested, spatial random-effects Cox proportional hazards model (level one: census division, level two: census tract). Hazard ratios were scaled to an interquartile increase in NDVI (full and incidence cohorts: ). The fully-adjusted model included age, sex, region (lived or not in the Greater Toronto Area), area-level unemployment, percent less than high school education, percent recent immigrants, and household income (quintiles), and population density. The number of events are: (CVD death) and (nonaccidental death) for full cohort; (CVD death), (nonaccidental death), (AMI incidence), and (HF incidence) for incidence cohort.

For the AMI and HF cohorts, the models further adjusted for clinical severity, in-hospital care, medications at discharge, smoking, and individual-level SES variables (AMI cohort: ; HF cohort: ). The number of events are: (CVD death), (nonaccidental death), (CVD readmission) for AMI cohort; (CVD death), (nonaccidental death), (CVD readmission) for HF cohort.

Within of study participants’ postal-code residences.

For the full and incidence cohorts, further adjustment for smoking, education, and access to primary care did not appreciably alter the associations with incident AMI and HF and mortality (Table 4). For the AMI and HF cohorts, the null associations remained similarly unchanged with additional control for access to primary care and coronary revascularization and exclusion of patients who were discharged and subsequently died within the same winter (Table S10).

Table 4

Sensitivity analysis for associations between urban green spaces within of study participants’ postal code residence and cardiovascular incidence, readmission, and mortality in four population-based cohorts of urban residents in Ontario, Canada, 2000–2014.

Outcome	Number of events	Further adjusted for smoking and educationc	Further adjusted for access to primary cared
		Hazard ratio (95% CI)	Hazard ratio (95% CI)
Full cohorta
Mortality
CVD	114,208	0.91 (0.89, 0.94)	0.89 (0.87, 0.91)
Nonaccidental	330,560	0.91 (0.89, 0.93)	0.89 (0.88, 0.90)
Incidence cohorta
Mortality
CVD	88,263	0.92 (0.89, 0.95)	0.90 (0.88, 0.92)
Nonaccidental	277,174	0.91 (0.89, 0.93)	0.89 (0.88, 0.90)
Incidence
AMI	58,553	0.95 (0.92, 0.98)	0.94 (0.92, 0.96)
HF	134,655	0.96 (0.93, 0.98)	0.95 (0.93, 0.96)
AMI cohortb
Readmission: CVD	4,419	—	1.03 (0.98, 1.07)
Mortality
CVD	2,788	—	0.99 (0.93, 1.04)
Nonaccidental	5,463	—	0.99 (0.95, 1.03)
HF cohortb
Readmission: CVD	5,482	—	1.00 (0.95, 1.05)
Mortality
CVD	4,981	—	0.98 (0.93, 1.03)
Nonaccidental	9,151	—	0.98 (0.94, 1.02)

Note: —, no data; AMI, acute myocardial infarction; CI, confidence interval; CVD: cardiovascular disease; HF, heart failure.

Two-level nested, spatial random-effects Cox proportional hazards model (level one: census division, level two: census tract). Hazard ratios were scaled to an interquartile increase in NDVI (full and incidence cohorts: ; AMI cohort: ; HF cohort: ). The fully adjusted model included age, sex, region (lived or not in the Greater Toronto Area), area-level unemployment, percent less than high school education, percent recent immigrants, and household income (quintiles), and population density.

For the AMI and HF cohorts, the models further adjusted for clinical severity, in-hospital care, medications at discharge, smoking, and individual-level SES variables (AMI cohort: ; HF cohort: ).

Indirect adjustment of smoking was conducted using data from the 2000/2001, 2003, and 2005 cycles of Canadian Community Health Survey. This was not applicable to AMI and HF cohort.

Access to primary care was derived using the density of family physicians. This sensitivity analysis was based on the main analysis presented in Figure 1.

Using mediation analysis, we estimated that the decreases in AMI and HF incidence associated with urban greenness explained 52.7% (95% CI: 37.6%–67.3%) of the protective association between urban greenness and cardiovascular mortality. In the subgroup analyses, we observed stronger, more protective associations between NDVI and incident AMI and HF and mortality in younger adults (e.g., AMI incidence: ; 95% CI: 0.86, 0.93 in younger adults vs. ; 95% CI: 0.94, 1.03 in older adults) (Figure 2 and Tables S11). The protective association of NDVI with cardiovascular mortality was also more pronounced in those with higher income (e.g., ; 95% CI: 0.80, 0.94 vs. ; 95% CI: 0.95, 1.04 between highest and lowest income). For AMI and HF patients, there was no strong evidence for effect modification of greenness with factors examined (Figures 2 and Tables S12–S14).

Figure 2.

Exploratory stratified analysis of associations between green spaces within of study participants’ residences and incidence of AMI and HF, cardiovascular mortality, and cardiovascular readmission among four population-based cohorts of urban residents in Ontario, Canada, according to age, sex, and income quintile. Note that hazard ratios were scaled to an interquartile increase in NDVI (full and incidence cohorts: ; AMI cohort: ; HF cohort: ). The number of events are: (CVD death), and (nonaccidental death) for full cohort; (CVD death), (nonaccidental death), (CVD death), and (nonaccidental death) for incidence cohort; (cardiovascular death), (nonaccidental death), and (cardiovascular readmission) for AMI cohort; and (cardiovascular death), (nonaccidental death), (cardiovascular readmission) for HF cohort. Note: AMI, acute myocardial infarction; CVD, cardiovascular disease; HF, heart failure.

Discussion

In this large cohort study in Ontario, living in urban areas with more green spaces was associated with lower risks of developing AMI and HF and dying from any cardiovascular cause. Each IQR increase in NDVI was associated with a 6% reduction in HF incidence, a 7% reduction in AMI incidence, and an 11% reduction in cardiovascular mortality among individuals free of AMI and HF. The protective associations were robust with various sensitivity analyses and were further accentuated among younger adults. In addition, of the protective association between greenness and cardiovascular mortality was found to be mediated through the decreases in AMI and HF incidence associated with greenness. On the other hand, living in these same areas was not associated with the risk of hospital readmission or mortality in AMI and HF patients.

We conducted, for the first time, a study to simultaneously examine the influence of urban greenness on a wide array of cardiovascular outcomes, from the development of two major CVDs, to hospital readmission, and ultimately to premature death. Our findings of an inverse relationship between green spaces and cardiovascular mortality are consistent with several population-based cohort studies (Crouse et al. 2017; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017). For example, in a cohort of urban dwellers who participated in the 2001 Canadian Census (Crouse et al. 2017), each IQR increment in NDVI (0.15) was associated with an 8% and 9% in the reduction of mortality from cardiovascular causes and from nonaccidental causes respectively. A similar protective association between NDVI and mortality was reported in 4.2 million adults in Switzerland (e.g., cardiovascular mortality: and nonaccidental mortality: per ) (Vienneau et al. 2017). On the other hand, in the Nurses’ Health Study, NDVI was not found to be associated with cardiovascular mortality (James et al. 2016). This inconsistency may be attributable to differences in population characteristics (the general population’s characteristics vs. female nurses’ characteristics), misclassification in exposure (resulting from inherent differences in the nature of greenness represented by NDVI), or chance. It is noteworthy that the effect size of our results closely agrees with those in the previous studies reporting an association with greenness (Crouse et al. 2017; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017) and remained consistent in all sensitivity analyses examined.

In contrast, the associations between green space and the incidence of AMI and HF have not been previously reported, although the plausibility of a relationship with cardiovascular incidence was suggested (Donovan et al. 2015; Pereira et al. 2012; Seo et al. 2019; Tamosiunas et al. 2014). Among the few studies conducted, a cross-sectional study in Perth, Australia, associated residential green space with a 7% decrease in the odds of self-reported heart disease per IQR increase in exposure (Pereira et al. 2012). More recently, a seven-city study in South Korea reported that regions with greater coverage of greenness tended to have lower incidence rates of cardiovascular diseases (Seo et al. 2019). Using the spatially resolved NDVI, we further found that higher levels of greenness in the vicinity of participants’ residences were related to a 7% decrease in incident AMI and a 6% decrease in HF. The decreases in AMI and HF incidence associated with residential greenness explained an estimated 53% of the protective association between greenness and cardiovascular mortality (under the assumptions of no unmeasured exposure–outcome, mediator–outcome, and exposure–mediator confounding, as well as no exposure–mediator interaction), implying that the possible cardiovascular benefits of green spaces may operate, to a large part, through its influence on the development of major CVDs. Living in greener areas of cities may contribute to reducing the risks of AMI and HF through promoting physical activity (Fong et al. 2018; Markevych et al. 2017). Exposure to urban greenness has also been shown to increase social engagement, promote psychological restoration, reduce stress, and buffer harmful environmental exposures (e.g., air pollution, noise, and heat), all of which are likely beneficial to cardiovascular system (Fong et al. 2018; Markevych et al. 2017). Given the importance of developing primary prevention strategies for CVD, future research should seek to examine the associations of urban greenness with the risk of developing other major CVDs.

Although residential greenness was associated with reduced cardiovascular incidence and mortality in people free of AMI and HF, we did not find any association between greenness and any delay in disease progression or mortality in patients with these conditions. This finding is in contrast to a randomized experimental study in Kaunas, Lithuania, which showed patients () with coronary artery disease walking in the park experienced favorable changes in heart rate and diastolic blood pressure relative to those walking in an urban street (Grazuleviciene et al. 2015). Improved survival was also noted in a small cohort of stroke survivors living in greener areas (Wilker et al. 2014). These previous studies, however, differ substantially from ours, especially in exposure assignment and cohort characteristics. The lack of associations in the AMI and HF cohorts may be due to these factors or to our efforts to reduce potential confounding by disease severity and quality of care. In addition, unlike other subpopulations, AMI and HF patients were often limited in their mobility after coronary events, constraining their interactions with outdoor green environments; thus they might not be amenable to the (active) impacts of greenness. Indeed, 66% of the AMI cohort and 90% of the HF cohort were not employed, which might be indicative of their reduced mobility or severe disease status, which would limit their ability to use green spaces in comparison with individuals in the full and incidence cohorts. Understanding why AMI and HF patients did not benefit from urban greenness and identifying opportunities to improve the benefits of greenness for these patients merit future attention.

The findings of this study should be interpreted in light of some limitations. First, we were unable to identify undiagnosed incident cases of AMI and HF. However, the HRs remained unchanged when we adjusted for access to health care and deprivation. Given universal health care in Ontario, incomplete diagnosis may lead to an underestimation of the true effect because this measurement error was likely independent of greenness exposure. Second, we lacked information on individual-level socioeconomic and behavioral variables (e.g., smoking, education, and employment status) for the full and incidence cohorts. Similarly, although we directly adjusted for smoking (current smoker or not) with the AMI and HF cohorts, we were unable to fully account for the other important aspects of smoking (e.g., intensity) because of data availability. To further account for potential confounding by these factors, we adjusted for various indicators for neighborhood SES and comorbidities. Because neighborhood SES is strongly associated with individual socioeconomic and behavioral variables (Gan et al. 2012; Janssen et al. 2006), and comorbidities and cardiovascular health share some common behavioral factors, adjusting for these variables would reduce the influence of these unmeasured variables on the HRs. Additionally, we adjusted indirectly for smoking and education for the full and incidence cohorts and found similar results. Furthermore, our observed greenness–mortality association in the full cohort closely agrees with that of previous studies using similar general populations (Crouse et al. 2017; Vienneau et al. 2017; Villeneuve et al. 2012; Wang et al. 2017). Despite these efforts, we acknowledge that our inability to fully account for the potential impact of individual-level covariates such as smoking on the HRs is a limitation of this study. Third, although the NDVI is an objective measure of greenness, it does not describe how people use green spaces around their homes, nor does it separate parks from green-belt areas, which can be assumed to have different impacts on cardiovascular health. Due to the inherent imprecision of NDVI as urban greenness variable, we were unable to further elucidate which types of natural environments (e.g., forests, public parks, or private yards) may be more beneficial to cardiovascular health. Additionally, we used annual NDVI levels to represent long-term exposure to average NDVI levels. Some previous studies used summer or maximum NDVI levels (Vienneau et al. 2017). We found a high correlation between maximum NDVI and annual NDVI levels (), but our mean and IQR NDVI levels may be lower because we used NDVI measures throughout the year. However, this approach likely better reflects actual long-term NDVI exposure. Furthermore, we excluded negative values of NDVI in our analysis and thus did not take into account of the full range of NDVI. Future research to evaluate the negative values of NDVI (often related to water bodies) in association with cardiovascular effects is needed. Last, we were able to ascertain exposure to greenness at participants’ residential postal codes (rather than street addresses), which inevitably led to spatial misalignment that could further result in exposure misclassification and for only a limited time period (as opposed to over an individual’s lifetime among people free of AMI and HF) because of data availability. These limitations further lead to an underestimation of associations between greenness and cardiovascular risks in this study.

Our study has several strengths that merit highlighting, including its large size and population-based representation of individuals who were either at risk or lived with AMI and HF. By simultaneously examining the impact of urban greenness on subpopulations at different stages of cardiovascular health but with similar exposures to green environments, this study provided a refined understanding of green space and cardiovascular health. In addition, we obtained extensive individual-level information including clinical data and in-hospital care for the AMI and HF patients, which allowed for better control for clinical severity and quality of care that are known to strongly influence post-AMI and HF survival. Aspects of our analytic approach further reduced concerns about confounding, such as the use of spatial random-effect Cox models, which accounted for the possibility that spatial dependencies in cardiovascular risk may not be completely explained by included variables. Furthermore, our study benefited from having detailed information on medical and residential history, allowing us to adjust for important cardiovascular risk factors (e.g., diabetes and hypertension) and account for the influence of residential mobility on exposure.

Conclusions

In this large population-based study, we found that living in urban areas with more green space was associated with lower risks of developing AMI and HF and dying from cardiovascular causes in people free of these outcomes. Among people who have already developed AMI or HF, however, increased exposure to green spaces was not associated with any delay in disease progression or mortality. With growing urban populations worldwide, this study adds weight to previous studies that suggested that increasing urban greenness could have a broad public health implication, especially for improving cardiovascular health, and that future research is necessary to understand why certain subpopulations (e.g., AMI and HF patients) did not benefit from urban greenness and identify opportunities to improve the benefits of greenness for these subpopulations.

Click here for additional data file.

Click here for additional data file.