Exposure to Residential Greenness as a Predictor of Cause-Specific Mortality and Stroke Incidence in the Rome Longitudinal Study.

BACKGROUND: Living in areas with higher levels of surrounding greenness and access to urban green areas have been associated with beneficial health outcomes. Some studies suggested a beneficial influence on mortality, but the evidence is still controversial.
OBJECTIVES: We used longitudinal data from a large cohort to estimate associations of two measures of residential greenness exposure with cause-specific mortality and stroke incidence.
METHODS: We studied a population-based cohort of 1,263,721 residents in Rome aged [Formula: see text], followed from 2001 to 2013. As greenness exposure, we utilized the leaf area index (LAI), which expresses the tree canopy as the leaf area per unit ground surface area, and the normalized difference vegetation index (NDVI) within 300- and [Formula: see text] buffers around home addresses. We estimated the association between the two measures of residential greenness and the outcomes using Cox models, after controlling for relevant individual covariates and contextual characteristics, and explored potential mediation by air pollution [fine particulate matter with aerodynamic diameter [Formula: see text] [Formula: see text] and [Formula: see text]] and road traffic noise.
RESULTS: We observed 198,704 deaths from nonaccidental causes, 81,269 from cardiovascular diseases [CVDs; 29,654 from ischemic heart disease (IHD)], 18,090 from cerebrovascular diseases, and 29,033 incident cases of stroke. Residential greenness, expressed as interquartile range (IQR) increase in LAI within [Formula: see text], was inversely associated with stroke incidence {hazard ratio (HR) 0.977 [95% confidence interval (CI): 0.961, 0.994]} and mortality for nonaccidental [HR 0.988 (95% CI: 0.981, 0.994)], cardiovascular [HR 0.984 (95% CI: 0.974, 0.994)] and cerebrovascular diseases [HR 0.964 (95% CI: 0.943, 0.985)]. Similar results were obtained using NDVI with 300- or [Formula: see text] buffers.
CONCLUSIONS: Living in greener areas was associated with better health outcomes in our study, which could be partly due to reduced exposure to environmental hazards. Further research is required to understand the underlying mechanisms. https://doi.org/10.1289/EHP2854.

Introduction

Urban settings are often characterized by artificial environments and buildings, high levels of stress, traffic, noise, air pollution, and heat, and the opportunities to have contact with nature are rare. According to the biophilia hypothesis, human beings have developed an innate affinity for natural environments and other forms of life through the evolution processes (Wilson 1984). Hence, it has been proposed that exposure to greenness and access to green spaces for urbanites may have a positive influence on their health through a variety of mechanisms.

Some studies have suggested that associations between exposure to green spaces and beneficial health outcomes may be mediated by increasing physical activity and a reduction in obesity and sedentary lifestyle (Lachowycz and Jones 2011; Richardson et al. 2013). Other studies highlighted potential beneficial effects on mental health. For example, green areas may provide space for social interaction and thus may benefit mental health by lowering psychological stress and depression (Alcock et al. 2014; Gascon et al. 2015a; Taylor et al. 2015). Still others have hypothesized beneficial effects due to a reduction of the exposure to urban environmental hazards, like air pollution and noise (Klingberg et al. 2017; Manes et al. 2016; Nowak et al. 2006).

An increasing number of studies have explored associations between greenness and nonaccidental and cause-specific mortality, although availability of individual-level data and exposure assessment procedures of environmental factors vary among studies. A 2015 systematic review and meta-analysis reported a significant reduction in cardiovascular disease (CVD) and all-cause mortality when using high vs. low categories of exposure to greenness (Gascon et al. 2015b). More recently, a nationwide prospective cohort study in the United States found evidence of a reduction in nonaccidental mortality among women participating in the U.S. Nurses’ Health Study (James et al. 2016). Two large cohort studies in Canada and Switzerland also found significant associations with reduced all-cause mortality and mortality due to cardiovascular, cerebrovascular, and respiratory diseases, as well as a community-based cohort study of elderly individuals living in Hong Kong (Crouse et al. 2017; Vienneau et al. 2017; Wang et al. 2017).

In this study, we aimed to investigate the association between residential greenness and all-cause and cause-specific mortality in a prospective cohort study of over a million subjects from the Rome Longitudinal Study (Cesaroni et al. 2013). Evidence of an association between residential greenness and cerebrovascular mortality led us to perform an additional investigation of the association between residential greenness and stroke incidence. Finally, we explored the mediating effect of fine particulate matter with aerodynamic diameter , , and noise, and studied effect modification by individual- and area-level characteristics [i.e., age group, sex, educational level, occupational status, and socioeconomic position (SEP)] to find possible differences in susceptibility.

Materials and Methods

Study Design and Population

Rome is the capital and main metropolitan city in Italy, and it covers an area of about . The population of the city, according to the 2001 Italian census (ISTAT 2001), was about 2.5 million inhabitants, most of whom were concentrated within the Main Ring Road (Figure 1). The Rome Longitudinal Study is the administrative 2001 census fixed cohort of Rome, is part of the National Statistical Program for the years 2017–2019 (www.sistan.it), and was approved by the Italian Data Protection Authority. It is described in detail elsewhere (Cesaroni et al. 2010, 2013). Briefly, information on sex, age, birthplace, marital status, education, and occupation collected from census are linked with the Regional Mortality Register and the Hospital Discharge Register through a personal identifying code, under strict controls to protect individual privacy. Causes of death from the Mortality Register are coded according to the Ninth Revision of International Classification of Diseases (ICD-9) (WHO 1977). The diagnoses of discharges from public and private hospitals (up to six) are coded according to the Ninth Revision of International Classification of Diseases Clinical Modification (ICD-9-CM).

Figure 1.

LAI values on 17 July 2015 for Rome municipality. Resolution: . Note: LAI, leaf area index, of leaf surface over of soil area (dimensionless quantity). The image was elaborated from a Landsat 8-Operational Land Imager/Thermal Infrared Sensor remote sensing image downloaded from the U.S. Geological Survey database using the ArcGIS software. See “Materials and Methods” for further details.

In the present study, we selected 1,265,058 subjects aged 30 or more who were resident in Rome from at least 5 y and who were not living in institutions at enrollment (prisons, hospitals, or nursing homes). We excluded subjects without a personal identifying code (1,332) and those with missing information on any of the covariates (only five subjects), obtaining a study population of 1,263,721 subjects.

Follow-up started on the census reference day (21 October 2001) and ended on the date of death, migration from Rome, or 31 December 2013, whichever came first.

For each subject, the baseline address of residence was geocoded with the ArcGIS 10 software (ESRI), using the World Geodetic System (1984) with the Universal Transverse Mercator zone 33 North Projection (WGS84 UTM 33N). To locate the address on the map, we used the Italian road network (Tele Atlas, Italy). The process resulted in 64% of addresses geocoded with an exact match and 36% with an approximation with scores .

Outcomes

The selected mortality outcomes included nonaccidental causes (ICD-9 codes ), CVDs (ICD-9: 390–459), ischemic heart disease - IHD (ICD-9: 410–414), cerebrovascular diseases (ICD-9: 430–438), and respiratory diseases (ICD-9: 460–519).

In addition, we analyzed the association with incident cases of stroke. We identified cases of stroke following the protocol used in the European Study of Cohorts for Air Pollution Effects (Stafoggia et al. 2014). Briefly, we selected hospital discharges with a diagnosis of stroke (ICD-9-CM: 431, 433.x1, 434, 436) and deaths from cerebrovascular disease (ICD-9: 431–436) of subjects who were not discharged from the hospital in the 2 d before death for any cause, or in the month before death for cerebrovascular causes. Subjects who had a stroke from 1 January 1996 to 21 October 2001 were considered prevalent cases and excluded from the incidence analysis (9,696).

Covariates

Covariates included age, sex, and a set of variables reported at inclusion: educational level (university, high school, junior high school, primary school, or less), marital status (married, single, separated/divorced, widowed), occupational status [highly qualified nonmanual (i.e., managers, professors, researchers), other nonmanual, manual workers, other workers (i.e., armed forces, retail sales), housewives, unemployed, retired, other condition (i.e., all subjects who cannot be classified into the previous categories)], birthplace (Rome or other).

To compensate for the lack of information on lifestyles, we adjusted for an area-level indicator of SEP because they have been shown to be predictors of behavioral factors independent of individual indicators (Diez Roux et al. 2003). The index is based on data from the 2001 census in Rome, derived from a factor analysis that included education, occupation, home ownership, family composition, crowding, and immigration at the census-block level (500 inhabitants per block, on average) and categorized into five levels, according to the quintiles of the distribution (Cesaroni et al. 2010).

Residential Greenness

To assess exposure to greenness, we used the leaf area index (LAI) and the normalized difference vegetation index (NDVI).

NDVI is a remote sense indicator based on land surface reflectance, which has been widely used in epidemiological studies to assess the relationship between greenness and health (Crouse et al. 2017; James et al. 2016; Vienneau et al. 2017). It ranges from to 1, with negative values indicating the presence of water bodies, values near zero for bare ground or buildings, and values above zero indicating increasing levels of greenness.

LAI has not been widely used in epidemiological studies. It is a dimensionless quantity defined as the one-sided leaf area per horizontal ground surface area (both in , namely leaf area in per 1 of surface area), and provides structural information on intensity of vegetation, like tree canopy density. LAI measurements range typically from about 0.2 for crop vegetation to more than 10 in dense conifer forests (Asner et al. 2003). Unlike NDVI, LAI is a measure that carries a biological significance, since it represents the amount of leaf surface per of soil. LAI is an important parameter expressing the main ecophysiological processes and vegetation growth, and is used in various ecosystem models (Sprintsin et al. 2009). For example, it was used to calculate air pollution removal from vegetation, making this index important in terms of ecosystem services (Fusaro et al. 2017; Manes et al. 2016; Marando et al. 2016). Its inclusion in an epidemiological study could bridge the gap between different scientific fields, allowing comparisons to be made with the literature on the ecological role of vegetation in air pollution mitigation. Furthermore, it appears easier for stakeholders and citizens to understand and to express in a practical way compared to NDVI. For this reason, we decided to use both the metrics in the present study. Firstly, from the U.S. Geological Survey database (USGS EROS 2015), we downloaded a Landsat 8-Operational Land Imager/Thermal Infrared Sensor remote sensing image with a spatial resolution of , acquired on 17 July 2015. We selected a summer day to combine the period of maximum leaf development with a totally cloud-free image. Through a semiautomatic classification plugin implemented in QGIS (QGIS Development Team 2016), we preprocessed the image and performed the dark object subtraction 1 atmospheric correction. Then, we calculated the NDVI (Rouse et al. 1974), according to the red (R) and near-infrared (NIR) relationship: . We then obtained LAI at a spatial resolution of , following the LAI/NDVI relation () derived from the integration of LAI ground measurements and satellite remotely sensed data in a study carried out in various ecosystems around Lazio coastal areas (Allegrini et al. 2006):

Finally, based on the map obtained (Figure 1), we calculated average LAI in 300- and buffers around each home address. Similarly, we calculated NDVI in 300- and buffers.

We decided to use two different radius measures to evaluate the exposure to greenness in both the proximity of residence () and in the neighborhood ().

Air Pollution and Noise

We attributed air pollution exposure to the coordinates of the baseline home address for each subject.

We estimated annual exposure through a dispersion model called FARM (Flexible Air quality Regional Model), with a spatial resolution of . This is a three-dimensional Eulerian model based on the integration of emission and atmospheric data and includes transportation and multiphase chemistry of pollutants in the atmosphere (Gariazzo et al. 2007, 2011). Data were available for the year 2005.

For annual concentrations, we used a land use regression model, described in detail elsewhere (Cesaroni et al. 2012). In brief, we used Ogawa passive samplers (Ogawa & Co. USA Inc.) to measure concentrations at 78 sites for three 1-wk-long periods in February, May, and October 2007. Then, after attributing the mean of the three values to each sampling location, we fit multivariate linear regression models to predict levels, accounting for several land use, GIS (Geographic information system), and traffic variables. The determination coefficient () from the best-fitting model was 0.704. We validated the model through the leave-one-out cross-validation technique: comparing observed and estimated concentrations of the analysis, the and adjusted was 0.61, and the root mean square error was . Details on the dispersion model, its validation, and the comparison of the results from the and models with actual measurements are reported elsewhere (Cesaroni et al. 2013).

Finally, for all residential addresses, we estimated average traffic noise levels (Lden, day–evening–night equivalent level measured over the 24-h period) by using the acoustic model SoundPLAN 7.4 (SoundPLAN GmbH, Backnang, Germany), based on 2009 traffic flow data. The noise assessment is described in detail elsewhere (Licitra et al. 2016; Ascari 2017).

Statistical Analysis

We evaluated the correlation between different exposure metrics using Pearson’s correlation coefficient.

We used Cox proportional hazard regression models to estimate associations between residential greenness and cause-specific mortality or incident stroke, with age as the time scale and stratifying by sex, to take account of sex differences in baseline rates. We used LAI and NDVI in 300- and buffers around home addresses as a continuous variable per interquartile range (IQR) increments, adjusting for educational level, marital status, occupational status, place of birth, and SEP. We also tested the associations using as exposure a categorical variable based on the quintiles of the distribution, and computed p-trend values using the Wald test for statistical significance across categories.

We explored the shape of relationships between exposures and outcomes by replacing the linear term in the base model with a natural spline with 2 degrees of freedom (Eisen et al. 2004), and compared the two models (the model with the linear term vs. the model with the spline term) using the likelihood ratio test (LRT).

As it has been hypothesized that the benefits of exposure to greenness may be mediated by a mitigation of the exposure to air pollution and noise, we performed a mediation analysis, following a methodology previously described (Imai et al. 2011). The mediation analysis relies on the hypothesis that part of the relationship between the exposure of interest and the outcome in study could be explained by a third actor, the mediator. In order to test this hypothesis, the effect of a variation in the value of exposure (from to ) on the outcome is decomposed in two parts, the causal mediation effect (or indirect effect) and the direct effect, the first one being defined, for each subject i, as

measures the difference that would occur in the observed outcome for subject i () if the value of mediator related to was replaced by the one related to , holding the value of exposure (e) constant. If there is no relationship between exposure and mediator, .

The direct effect can be expressed as

After the average value of these measures over all subjects is estimated, it is possible to express the proportion mediated, which measures the rate between the average causal mediation effect and the average total effect:

It is possible to estimate confidence intervals (CIs) for these measures by running quasi-Bayesian Monte Carlo simulations. In practice, the analysis requires to fit different models: one for the association between exposure and outcome adjusted for the potential mediator (in our case, a Cox model), and one for the association between exposure and mediator (in our case, a linear regression model). Then, making calculations on the coefficients obtained, one can compute the total effect, the natural direct effect, and the natural indirect effect, and then obtain the proportion mediated, which represents the extent of the association mediated by that specific mediator. The results of mediation analyses are valid only if strong and unverifiable assumptions hold, including no unmeasured exposure–outcome confounding, no unmeasured mediator–outcome confounding, no unmeasured exposure–mediator confounding, and no mediator–outcome confounding affected by exposure. We performed separate analyses to assess potential mediation by annual average concentrations of and , and daily average traffic noise (Lden), respectively, with adjustment for the covariates included in the main model. We used the publicly available Mediation package in R version 3.5 (R Development Core Team) (Tingley et al. 2014).

Given the high correlation between the different measures of exposure to greenness and the consistency of results using the different metrics, we explored the possible interaction between exposure and sociodemographic conditions using only LAI in a buffer. We tested for potential effect modification by including interaction terms between the exposure and one potential effect modifier at a time [sex, age class (30–49, 50–69, ), educational level, occupational status (we merged all workers in one category, obtaining a five-level variable: employed, unemployed, housewife, retired, other), and census-block SEP], and assessed heterogeneity performing LRTs between models with and without interaction terms.

Supplementary Analyses

Subjects share similar environmental and social characteristics with residents in the same area. We performed as sensitivity analysis a Cox model with shared frailty to investigate the role of neighborhood (94 urban areas defined by the Municipal Administration) for confounding and clustering in the exposure/outcome association, assuming that subjects residing in the same neighborhood have the same frailty. Across groups, the frailties are assumed to be distributed latent random effects that affect the hazard multiplicatively.

As additional sensitivity analyses, we selected a subset of 7,845 subjects included in our study population, with available information on smoking status as participants to a previous survey, the Italian Studies on Respiratory Disorders in Childhood and Environment (SIDRIA) (Cesaroni et al. 2008). We performed logistic regressions, both with and without adjusting for the other covariates, to investigate the possible association between smoking habits (never vs. ever smokers) and exposure to greenness (LAI in a buffer). We also replicated the survival analysis restricted to SIDRIA participants, adjusting for smoking status. We evaluated only the association with mortality for nonaccidental causes due to the small number of cases.

We assumed that the exposure to residential greenness estimated from a single remote sensing image, acquired on 17 July 2015, could be representative of the overall follow-up period exposure. To test this assumption, we retrieved 14 NDVI summer images of the Rome municipality territory from MODIS [11 or 12 July; MOD13Q1 product (Didan 2015)], for the period from 2001 to 2013 plus 2015, and performed a correlation matrix between NDVI values for each pixel of the images.

We used STATA (version 13; StataCorp) for all statistical analyses except for frailty models, spline models, and mediation analysis, for which we used R software.

Results

We included 1,263,721 residents aged 30 y or more at baseline. We considered lost to follow-up the subjects who migrated from Rome during the study period (9.4%).

Of the participants, 54.5% were females and the mean (standard deviation, SD) age at baseline was 55.1 (15.4) y (Table 1). Mean (SD) exposure to greenness (as LAI in around home address) was 2.7 (0.6), to was , to was , and to traffic noise was 59.8 (8.6) dB.

Table 1

Baseline characteristics of participants according to quintiles of leaf area index (LAI) in around residential address (21 October 2001).

Characteristics	Quintiles of LAI (min–max)					Overall cohort
	1 (1.1–2.1)	2 (2.1–2.5)	3 (2.5–2.8)	4 (2.8–3.2)	5 (3.2–5.8)	Total
n	252,756	252,726	252,849	252,714	252,676	1,263,721
Vital status (%)
Alive at 31 December 2013	71.6	72.8	74.1	75.3	76.8	74.1
Emigrated	9.9	9.5	9.2	9.2	9.3	9.4
Dead	18.6	17.7	16.7	15.5	13.8	16.5
Sex
Females (%)	55.2	55.0	54.4	54.0	53.7	54.5
Age at enrollment
Mean (SD)	56.2±15.9	55.8±15.8	54.9±15.5	54.4±15.1	53.8±14.6	55.1±15.4
Marital status (%)
Married	62.2	64.9	66.6	68.3	69.6	66.3
Single	17.7	15.9	14.9	14.2	13.7	15.3
Separated/divorced	7.1	6.7	6.9	7.0	7.5	7.0
Widowed	13.0	12.5	11.6	10.6	9.3	11.4
Place of birth
Rome	50.3	51.7	52.3	52.4	52.2	51.8
Other	49.7	48.4	47.7	47.6	47.8	48.2
Educational level (%)
University	15.3	15.5	14.6	15.7	20.0	16.2
High school	31.8	32.2	31.5	33.0	36.2	32.9
Secondary school	25.8	26.2	27.2	26.8	23.8	25.9
Primary school or less	27.1	26.1	26.7	24.6	20.0	24.9
Occupational status
Employed NM I	13.0	12.7	12.3	13.3	16.1	13.5
Employed NM II	15.4	15.4	15.4	16.3	17.1	15.9
Employed M	9.0	9.7	10.6	10.5	9.8	9.9
Employed, other	6.5	6.4	6.7	6.6	6.3	6.5
Housewives	21.3	21.5	21.6	20.8	19.9	21.0
Unemployed	4.7	4.8	5.4	5.2	4.9	5.0
Retired	25.0	24.7	23.2	22.7	21.6	23.4
Other condition	5.2	4.8	4.9	4.6	4.4	4.8
Area-level SEP
Very high	14.6	20.4	18.0	20.6	25.5	19.8
High	22.6	19.2	18.3	19.8	22.1	20.4
Medium	26.4	19.0	17.3	18.9	18.7	20.1
Low	22.8	22.8	20.8	19.6	16.0	20.4
Very low	13.6	18.7	25.6	21.2	17.9	19.4
Residential green (mean±SD)
NDVI (300 m)	0.26±0.02	0.31±0.01	0.35±0.01	0.39±0.01	0.46±0.04	0.35±0.07
LAI (1,000 m)	2.43±0.33	2.77±0.34	2.95±0.31	3.14±0.32	3.51±0.42	2.96±0.50
NDVI (1,000 m)	0.32±0.05	0.36±0.05	0.39±0.04	0.41±0.04	0.45±0.05	0.39±0.06
Exposure to air pollution (mean±SD)
PM_2.5 (μg/m^3)	25.0±4.5	23.7±4.1	22.9±4	22.5±4.1	21.0±4.4	23.0±4.4
NO_2 (μg/m^3)	48.5±8.5	45.5±7.8	43.5±7.4	42.0±7.3	38.8±7.7	43.6±8.4
Noise (dB)	62.9±8.2	61.5±7.9	60±8.1	58.4±8.3	56.1±8.7	59.8±8.6

Note: Reported data are complete for all study participants. SD, standard deviation; employed NM I, highly qualified nonmanual workers (i.e. managers, professors, researchers); employed NM II, other nonmanual workers; employed M, manual workers; employed, other: other workers (i.e., armed forces, retail sales); other condition, all subjects who cannot be classified into the previous categories; area-level SEP (socioeconomic position), composite indicator derived from education, occupation, home ownership, family composition, crowding, and immigration at the census-block level; LAI, leaf area index, of leaf surface over of soil area (dimensionless quantity); NDVI, normalized difference vegetation index, remote sensing indicator of green based on spectral reflectance (dimensionless quantity); , averaged values in a buffer around residence; , averaged values in a buffer around residence; , fine particulate matter with aerodynamic diameter ; noise, average traffic noise levels at residence according to the Lden indicator (day–evening–night equivalent level measured over the 24-h period).

Table 1 shows descriptive statistics of the study population by quintiles of LAI in the around the residence. Subjects exposed to higher levels of LAI were younger and better educated, were more likely to be employed in nonmanual work, lived in areas with better socioeconomic conditions, lived farther from high-traffic roads, were exposed to lower levels of air pollution and noise, and showed higher survival rates at the end of follow-up compared to those exposed to lower levels. Table S1 shows descriptive statistics of the study population by quintiles of NDVI in the around the residence. The population characteristics were similar when categorized by quintiles of NDVI vs. LAI.

Over approximately person-years of follow-up, we observed a total of 198,704 nonaccidental deaths, with a crude mortality rate (CMR) of 147.6 deaths for 10,000 person-years. Of these, 81,269 were due to CVDs (CMR: 60.3), 29,654 to IHD (CMR: 22), 18,090 to cerebrovascular diseases (CMR: 13.4), and 12,433 to respiratory diseases (CMR: 9.2). We also observed 29,033 incident cases of stroke, with an incidence rate of 17.7 for 10,000 person-years.

According to Pearson’s coefficients, the metrics for greenness were highly correlated with each other and moderately or poorly correlated with air pollution and noise (Table 2).

Table 2

Pearson’s correlation coefficients of the variables for residential greenness (LAI 300; NDVI 300; LAI 1,000; NDVI 1,000), air pollution, and noise.

—	LAI (300 m)	NDVI (300 m)	LAI (1,000 m)	NDVI (1,000 m)	PM_2.5	NO_2	Noise
LAI (300 m)	1.00	—	—	—	—	—	—
NDVI (300 m)	0.99	1.00	—	—	—	—	—
LAI (1,000 m)	0.76	0.76	1.00	—	—	—	—
NDVI (1,000 m)	0.74	0.75	0.96	1.00	—	—	—
PM_2.5	−0.31	−0.30	−0.40	−0.31	1.00	—	—
NO_2	−0.40	−0.39	−0.43	−0.35	0.78	1.00
Noise	−0.29	−0.29	−0.27	−0.23	0.40	0.50	1.00

Note: —, data not available; LAI, leaf area index, of leaf surface over of soil area (dimensionless quantity); NDVI, normalized difference vegetation index, remote sensing indicator of green based on spectral reflectance (dimensionless quantity); , averaged values in a buffer around residence; , averaged values in a buffer around residence; , average 2005 levels of fine particulate matter with aerodynamic diameter from a Eulerian dispersion model; , averaged 2007 levels of nitrogen dioxide from a land use regression model; noise, average traffic noise levels at residence according to the Lden indicator (day–evening–night equivalent level measured over the 24-h period).

Table 3 shows the results of Cox regression for exposure to residential greenness, expressed as IQR increments in both LAI and NDVI in 300 or around home address, the associations by quintiles of exposure, and p-trend values. When adjusting for individual characteristics and area-level SEP, we observed an inverse association between residential greenness and stroke incidence, and mortality due to nonaccidental, cardiovascular, IHD, and cerebrovascular causes. The associations were strongest in the highest quintiles, and the trends across quintiles were statistically significant. The evidence of an association between exposure to residential greenness and mortality due to respiratory diseases was uncertain. As expected, the associations were very similar both when comparing LAI to NDVI and the two different buffers.

Table 3

Association between residential greenness [expressed as leaf area index (LAI) and normalized difference vegetation index (NDVI) in 300 and around home addresses] and mortality/stroke incidence, Rome Longitudinal Study, 2001–2013.

Outcome	n	Number of events	LAI (300-m)		NDVI (300 m)		LAI (1,000 m)		NDVI (1,000 m)
			HR	(95% CI)	HR	(95% CI)	HR	(95% CI)	HR	(95% CI)
Nonaccidental mortality
Continuous	1,263,721	198,704	0.988	(0.981, 0.994)	0.987	(0.981, 0.994)	0.984	(0.978, 0.990)	0.983	(0.977, 0.988)
Quintile 1	252,756	44,767	1.00	—	1.00	—	1.00	—	1.00	—
Quintile 2	252,726	42,794	0.992	(0.979, 1.005)	0.989	(0.976, 1.002)	1.003	(0.990, 1.017)	1.004	(0.991, 1.018)
Quintile 3	252,849	40,356	0.996	(0.982, 1.009)	0.992	(0.979, 1.006)	0.996	(0.983, 1.009)	0.999	(0.986, 1.013)
Quintile 4	252,719	37,387	0.990	(0.977, 1.004)	0.996	(0.982, 1.009)	0.987	(0.973, 1.000)	0.985	(0.971, 0.999)
Quintile 5	252,676	33,400	0.977	(0.963, 0.991)	0.973	(0.959, 0.987)	0.987	(0.973, 1.001)	0.968	(0.955, 0.982)
p-trend	—	—	0.004		0.005		<0.001		<0.001
Cardiovascular mortality
Continuous	1,263,721	81,269	0.984	(0.974, 0.994)	0.984	(0.975, 0.994)	0.975	(0.966, 0.985)	0.977	(0.969, 0.986)
Quintile 1	252,756	18,768	1.00	—	1.00	—	1.00	—	1.00	—
Quintile 2	252,726	17,919	1.003	(0.983, 1.024)	1.003	(0.983, 1.024)	1.007	(0.986, 1.028)	1.009	(0.988, 1.030)
Quintile 3	252,849	16,587	1.005	(0.984, 1.026)	1.000	(0.979, 1.021)	0.993	(0.972, 1.014)	0.996	(0.976, 1.018)
Quintile 4	252,719	15,032	0.994	(0.973, 1.016)	0.998	(0.977, 1.020)	0.979	(0.958, 1.001)	0.987	(0.966, 1.009)
Quintile 5	252,676	12,963	0.969	(0.947, 0.991)	0.971	(0.949, 0.993)	0.956	(0.935, 0.977)	0.955	(0.934, 0.977)
p-trend	—	—	0.010		0.022		<0.001		<0.001
IHD mortality
Continuous	1,263,721	29,654	0.984	(0.967, 1.000)	0.985	(0.968, 1.001)	0.975	(0.960, 0.990)	0.978	(0.964, 0.993)
Quintile 1	252,756	6,776	1.00	—	1.00	—	1.00	—	1.00	—
Quintile 2	252,726	6,513	1.004	(0.970, 1.039)	0.997	(0.963, 1.032)	1.018	(0.983, 1.054)	1.028	(0.993, 1.064)
Quintile 3	252,849	6,043	0.999	(0.964, 1.034)	0.991	(0.957, 1.026)	1.011	(0.977, 1.047)	1.016	(0.982, 1.052)
Quintile 4	252,719	5,503	0.983	(0.948, 1.019)	0.989	(0.954, 1.025)	1.003	(0.968, 1.040)	1.006	(0.970, 1.043)
Quintile 5	252,676	4,819	0.963	(0.928, 1.000)	0.962	(0.927, 0.999)	0.936	(0.902, 0.972)	0.940	(0.906, 0.976)
p-trend	—	—	0.029		0.056		0.002		0.002
Cerebrovascular mortality
Continuous	1,263,721	18,090	0.964	(0.943, 0.985)	0.965	(0.944, 0.986)	0.956	(0.938, 0.975)	0.962	(0.944, 0.980)
Quintile 1	252,756	4,281	1.00	—	1.00	—	1.00	—	1.00	—
Quintile 2	252,726	4,056	0.996	(0.954, 1.040)	1.002	(0.960, 1.046)	1.010	(0.967, 1.056)	1.021	(0.978, 1.067)
Quintile 3	252,849	3,657	0.973	(0.931, 1.018)	0.977	(0.934, 1.022)	0.951	(0.910, 0.994)	0.941	(0.900, 0.984)
Quintile 4	252,719	3,293	0.958	(0.915, 1.003)	0.956	(0.913, 1.001)	0.942	(0.899, 0.986)	0.964	(0.921, 1.010)
Quintile 5	252,676	2,803	0.930	(0.886, 0.976)	0.939	(0.894, 0.985)	0.916	(0.874, 0.960)	0.916	(0.873, 0.960)
p-trend	—	—	0.001		0.002		<0.001		<0.001
Respiratory mortality
Continuous	1,263,721	12,433	1.014	(0.988, 1.041)	1.011	(0.986, 1.038)	1.004	(0.981, 1.028)	0.998	(0.976, 1.021)
Quintile 1	252,756	2,869	1.00	—	1.00	—	1.00	—	1.00	—
Quintile 2	252,726	2,682	0.980	(0.930, 1.033)	0.955	(0.906, 1.007)	0.977	(0.926, 1.031)	0.964	(0.914, 1.018)
Quintile 3	252,849	2,500	0.987	(0.935, 1.042)	0.973	(0.922, 1.027)	0.982	(0.931, 1.036)	0.979	(0.928, 1.032)
Quintile 4	252,719	2,281	0.987	(0.934, 1.043)	0.982	(0.929, 1.037)	0.960	(0.907, 1.015)	0.948	(0.896, 1.002)
Quintile 5	252,676	2,101	1.034	(0.977, 1.095)	1.012	(0.956, 1.072)	1.004	(0.949, 1.062)	1.008	(0.953, 1.066)
p-trend	—	—	0.320		0.529		0.794		0.888
Stroke incidence
Continuous	1,254,030	29,033	0.977	(0.961, 0.994)	0.976	(0.960, 0.993)	0.982	(0.966, 0.997)	0.981	(0.967, 0.996)
Quintile 1	252,756	6,469	1.00	—	1.00	—	1.00	—	1.00	—
Quintile 2	252,726	6,185	0.984	(0.950, 1.019)	0.985	(0.951, 1.020)	0.984	(0.950, 1.019)	0.987	(0.953, 1.022)
Quintile 3	252,849	5,997	1.000	(0.965, 1.036)	0.988	(0.954, 1.024)	0.984	(0.950, 1.019)	0.988	(0.954, 1.024)
Quintile 4	252,719	5,497	0.974	(0.939, 1.010)	0.977	(0.943, 1.013)	0.949	(0.915, 0.985)	0.951	(0.916, 0.986)
Quintile 5	252,676	4,885	0.951	(0.916, 0.988)	0.946	(0.911, 0.982)	0.967	(0.931, 1.003)	0.967	(0.932, 1.004)
p-trend	—	—	0.013		0.007		0.014		0.015

Note: —, data not available; Models adjusted for educational level, marital status, occupational status, place of birth, area-level socioeconomic position. Stratified by sex, age as time scale. Continuous variables per interquartile range (IQR) increases: ; ; ; . LAI, of leaf surface over of soil area (dimensionless quantity); NDVI, remote sensing indicator of green based on spectral reflectance (dimensionless quantity); , averaged values in a buffer around residence; , averaged values in a buffer around residence; HR, hazard ratio; CI, confidence interval.

When comparing the spline models with the linear term models, the results of LRTs indicated that there was no statistically significant deviation from linearity (LRT p-values ranged from 0.10 to 0.92).

Results from mediation analysis, given that the underlying assumptions hold, suggested that the associations could be partly mediated by a reduction of the exposure to air pollution and noise (Table 4). The proportion mediated varied largely across mediators and outcomes, going from 0 to 92%. It was higher for and and lower for noise, and higher for natural and cardiovascular mortality than for cerebrovascular mortality. To note, no mediation was found for stroke incidence.

Table 4

Mediation analysis: proportion of the associations mediated by environmental hazards (, , and noise).

Outcome	PM_2.5		NO_2		Noise
	Prop. mediated	p-Value	Prop. mediated	p-Value	Prop. mediated	p-Value
Nonaccidental mortality	0.52	<0.001	0.92	0.020	0.27	0.020
Cardiovascular mortality	0.81	<0.001	0.61	<0.001	0.41	0.020
Cerebrovascular mortality	0.39	<0.001	−0.02	0.900	0.26	<0.001
Stroke incidence	−0.09	0.680	−0.01	0.940	0.08	0.460

Note: Analysis made using the R Mediation package (Tingley et al. 2014). Models adjusted for education, marital status, occupation, place of birth, and socioeconomic position. Stratified by sex, age as time scale. , average 2005 levels of fine particulate matter with aerodynamic diameter from a Eulerian dispersion model; , average 2007 levels of nitrogen dioxide from a land use regression model; noise, average traffic noise levels at residence according to the Lden indicator (day–evening–night equivalent level measured over the 24-h period).

Results of effect modification analyses are shown in Figure 2 and reported in full in Table S2. The association between greenness and nonaccidental-causes mortality was significantly different across subgroups by age at inclusion, educational level, census-block SEP, and occupational status, and was stronger among young individuals, those with high or intermediate educational levels and high census-block SEP and those who were employed, retired, or housewives. Occupational status was an effect modifier of the relationship between LAI and mortality for cerebrovascular causes, with stronger inverse associations among housewives and retired people and a positive association for individuals classified as “other.” The association between LAI and respiratory mortality was positive in subjects with a university degree, inverse in those with a high school or junior high school degree, and weakly positive in those with a primary school degree or less (interaction ). Different from the association with cerebrovascular mortality, for which there was not a statistically significant difference between males and females, the association between LAI and stroke incidence was inverse for women and null for men (interaction ). There was also significant effect modification by census-block SEP (interaction ), with inverse associations among those living in census blocks with intermediate, high, and very high SEP, and a positive association for residents of census blocks with very low SEP. It is notable that, with the exception of cerebrovascular mortality, LAI was positively associated with all outcomes among subjects living in census blocks with very low SEP.

Figure 2.

Modification of the association between leaf area index (LAI) exposure in around residence and mortality/stroke incidence by participants’ characteristics, Rome Longitudinal Study 2001–2013. Note: Models adjusted for education, marital status, occupation, place of birth, and socioeconomic position. Stratified by sex, age as time scale. Stratified estimates were derived introducing an interaction term between LAI and the effect modifier and calculating the estimate for each category. Error bars represent the 95% confidence interval (CI). p-Values were derived from likelihood ratio tests. LAI: of leaf surface over of soil area (dimensionless quantity); area-level SEP (socioeconomic position): composite indicator derived from education, occupation, home ownership, family composition, crowding, and immigration at the census-block level; nonaccidental: mortality due to all nonaccidental causes; cardiovascular: mortality for all cardiovascular causes; IHD: mortality for ischemic heart attack; respiratory: mortality for all respiratory causes; stroke: incidence of stroke. See Table S2 for the corresponding numeric data.

Associations estimated with frailty models using neighborhoods as random variance component did not show substantial change compared to those of the standard model (see Table S3).

Results of the analysis restricted to SIDRIA participants are reported in Table S4. The association between exposure to greenness and smoking (ever vs. never smoker) was close to the null and nonsignificant for the buffer [e.g., odds ratio (OR) 0.99 (95% CI: 0.90, 1.08) for LAI 300] and slightly positive and nonsignificant in the buffer [e.g., OR 1.03 (95% CI: 0.95, 1.13) for LAI 1,000]. In the survival analysis, the association followed the same direction observed in the main analysis (i.e., inverse relationship), although it was not statistically significant due to the small number of cases; moreover, results were similar when we furtherly adjusted for smoking status.

Correlations between NDVI values for each pixel on 11 or 12 July for each year in the period 2001–2013 plus 2015 () ranged from 0.75–0.90, which suggests that the 2015 image was representative of values throughout the follow-up period (Table S5).

Discussion

In this longitudinal study of more than a million individuals in Rome, we observed an inverse association between higher amounts of greenness around home addresses and stroke incidence and mortality for nonaccidental, cardiovascular, and cerebrovascular causes. We found weaker associations with mortality due to ischemic heart disease, and no inverse association with respiratory disease mortality. As expected, results were similar when using both LAI and NDVI, and similar or slightly stronger when analyzing the effect of greenness in the larger buffer () around residence compared to the smaller one (). The mediation analysis suggested that the associations with nonaccidental, cardiovascular, and cerebrovascular mortality could be partly mediated by a reduced exposure to air pollution and road traffic noise, while no significant mediation was found for stroke incidence.

Our estimates are consistent with those reported by the currently available meta-analysis (Gascon et al. 2015b) and by most previous studies that explored the association between urban green and mortality (Hu et al. 2008; Lachowycz and Jones 2014; Mitchell et al. 2011; Mitchell and Popham 2008; Wilker et al. 2014), but weaker than those reported by the most recent longitudinal studies (Crouse et al. 2017; James et al. 2016; Vienneau et al. 2017; Villeneuve et al. 2012). Such difference may be due to the broader range of greenness exposure observed in these studies, probably because they were at the national or state level.

In our study, the association with mortality due to respiratory diseases was positive, although not statistically significant. A possible explanation is that more greenness or proximity to green spaces could increase the risk of asthma attacks, as suggested by Dadvand and colleagues in their cross-sectional study on children (Dadvand et al. 2014). However, evidence on the association between exposure to greenness and respiratory health is still mixed (Fuertes et al. 2014).

The negative association we found between residential greenness and stroke incidence, although small, is a novel and major finding of our research. Previous studies have already reported inverse associations between measures of greenness and stroke-related mortality. An ecological study conducted in northwestern Florida found higher risk of stroke mortality in areas with low socioeconomic level, high levels of air pollution, and low exposure to greenness, using dasymetric mapping procedures and Bayesian hierarchical modeling (Hu et al. 2008). A prospective study published in 2014 reported lower all-cause mortality rates for subjects living in areas characterized by higher levels of exposure to greenness in a cohort of stroke survivors (Wilker et al. 2014). Another cohort study from Ontario, Canada, found an association between NDVI and stroke mortality similar to our findings on stroke incidence (Villeneuve et al. 2012).

Several mediators have been proposed to explain associations between exposure to greenness and beneficial health outcomes. Proxy measures of exposure to green spaces have been associated with increased physical activity and reduced obesity and overweight, although the evidence in this field is still mixed, as reviewed in Lachowycz and Jones (2011). A national cross-sectional study of residential proximity to publicly available natural areas in New Zealand and self-reported health concluded that higher levels of physical activity did not fully explain associations between green space and improved cardiovascular and mental health, based on small differences in the association with and without adjustment for physical activity (Richardson et al. 2013). A study conducted in Gothenburg, Sweden, found that green areas were associated with reduced exposure to air pollution and noise; the concentrations of pollutants and noise levels were significantly different between vegetated and nonvegetated sites and before or after leaf emergence (Klingberg et al. 2017). Other studies conducted nationwide in the United States and Italy also supported the association between urban green areas and improved air quality (Manes et al. 2016; Nowak et al. 2006). Other authors hypothesized that exposure to natural environments could exert beneficial effects on mental health. In his review about the restorative benefits of nature, Kaplan proposed that exposure to natural environments may act both by helping mitigate stress and by restoring mental resources from fatigue (Kaplan 1995). A panel study on 1,064 people from the United Kingdom who changed their residence showed that those who moved to greener areas reported a significant improvement in mental health, whereas the contrary happened to those who moved to less green areas (Alcock et al. 2014). In a cross-sectional study conducted in London, the urban street tree density was inversely associated with antidepressant prescription rates (Taylor et al. 2015). Another cross-sectional study from the Netherlands involving 10,069 respondents to the Dutch National Survey of General Practice found that feelings of loneliness and perceived shortage of social support partly mediated the association between percentage of green space and perceived health (Maas et al. 2009). In a 2015 systematic review by Gascon and colleagues, the authors concluded that there was some evidence of a causal relationship between exposure to surrounding greenness and mental health in adults, but it was limited by the heterogeneity of exposure assessment and by the limited number of studies (Gascon et al. 2015a). As a whole, it is probable that all these mediators are involved contemporarily and that they act interdependently (Hartig et al. 2014), making it difficult to isolate the single effect (e.g., benefits on mental health of physical activity in a natural environment with less air pollution). In our study, we estimated the proportion of the association explained by air pollution and noise. Such proportion was larger for air pollutants, arriving up to 92%, a very high value compared to those observed in previous studies, where it was around 4–10% (James et al. 2016; Vienneau et al. 2017). To note, no mediation was found for stroke incidence, suggesting that other mediators not explored in our analysis could play a major role for such disease. However, as noted previously (see “Statistical Analysis” subsection of “Methods”), the results of mediation analysis are valid only if several and unverifiable assumptions hold. The extent to which such assumptions hold in different settings and populations may also contribute to the differences observed among studies.

We found evidence of effect modification in the relationship between greenness and health outcomes by many individual characteristics. We observed stronger associations among females (for cardiovascular mortality and stroke incidence), those with better educational attainment (for nonaccidental and cardiovascular mortality), those who lived in census blocks with better socioeconomic conditions (for all outcomes except cerebrovascular and respiratory mortality), those who were employed (for nonaccidental mortality), housewives, or retired (for nonaccidental and cerebrovascular mortality). The associations were significantly different by age class for nonaccidental and cardiovascular mortality, with the oldest age group closer to the null compared to the middle and youngest age groups. Previous studies reported different associations by age group and socioeconomic factors. The association between greenness and nonaccidental mortality reported by Villeneuve and colleagues were stronger among individuals in the 35–50 age group, and progressively reduced in the 50–65 and age groups (Villeneuve et al. 2012). This may be due to a different use of green spaces by people of different ages, e.g., younger people may use these areas for physical activity more than older people. An observational population study conducted in England found that the association between income deprivation and all-cause mortality was different by categories of exposure to green space: the incidence rate ratio for the most income-deprived individuals compared with the least deprived was 1.93 (95% CI: 1.86, 2.01) in the least green areas, whereas it was 1.43 (95% CI: 1.34, 1.53) in the most green (Mitchell and Popham 2008). In the study by Maas and colleagues we already described, individuals with low income and educational levels showed stronger associations between lack of green spaces and loneliness (Maas et al. 2009). The stronger associations observed among individuals with better socioeconomic conditions could be linked to a different quality and accessibility of the available green areas and parks in more affluent neighborhoods. For cerebrovascular mortality, we also found stronger negative associations among housewives and retired people compared to other categories. These subjects tend to spend more time at home than others, and this could reduce exposure misclassification.

The main limitation of this study is the lack of information on many relevant individual risk factors, such as smoking habits, body mass index (BMI), physical activity, diet, alcohol intake, and social engagement. Some of these factors (like BMI/physical activity and social engagement) are supposed to be on the causal pathway from green to health, and so they should be considered mediators rather than confounders. However, the subgroup analysis on subjects with available information on smoking status showed no evidence of an association between exposure to greenness and smoking habits. Another limitation is that we had only objective measures of the amount of greenness around residences, and no information was available about exposure to greenness in other places frequented in daily life (e.g., the workplace). Moreover, we had no information about actual habits of visiting or using the green spaces in leisure time. Thus, we were not able to explore in detail how different activities connected with green spaces might mediate associations between exposure to residential greenness and beneficial health outcomes. We derived the exposure levels from maps with a spatial resolution, obtaining rather precise measures of the amount of greenness surrounding each address at different buffers (300 and ). Nevertheless, even a resolution of is subjected to limitations, and may misrepresent the available greenness in the study area. This measure was based on a remote sensing image taken at a given moment; thus, it did not allow us to address temporal variations of the exposure. However, we do not expect the amount of greenness to have changed substantially over the follow-up period (this is supported by the high coefficients of the correlation matrix reported in Table S5). Although our study area was limited to the urban and suburban areas of Rome, with the majority of the population being concentrated within the Main Ring Road, Rome has more than 45 million of green surface, and residential greenness based on LAI varied from 1.1 to 5.7 among different census blocks within the study area.

Conclusions

Residential surrounding greenness, as indicated by two different exposure metrics (LAI and NDVI) within 300- and buffer areas, was inversely associated with stroke incidence and with nonaccidental, cardiovascular, and cerebrovascular mortality in our large study cohort of adult residents of Rome. Findings from a mediation analysis suggest that these associations may be partly mediated through a reduced exposure to air pollution and noise. Our results add to current evidence supporting a role of exposure to greenness in promoting better health. Future research should combine objective measures about the amount of greenness around each participant’s home and workplace with more detailed information about the frequency and quality of their interaction with green places. Our findings provide further support for efforts by urban and transportation planners, policy makers, and others to improve health and the quality of life by ensuring that urban residents have access to areas that promote contact with nature and social interactions.

Click here for additional data file.

Click here for additional data file.