Greenness, mortality and mental health prescription rates in urban Scotland-a population level, observational study

BACKGROUND: Recent studies have shown an association between vegetation around dwellings and mortality, with mental health as a possible mediator. OBJECTIVES: Examine whether there is an association between greenness and mortality or greenness and the proportion of the population being prescribed drugs for anxiety, depression or psychosis in urban areas of Scotland. METHODS: Two greenness maps were prepared based on Landsat 8 Normalised Difference Vegetation Index data from 2013 to 2016, one for summer and one for winter. Greenness was sampled from these maps around each of 91,357 urban postcodes. The greenness data was averaged by 4,883 urban Data Zones covering 71% of the Scottish population and compared with mortality and prescription rate data from the Scottish Index of Multiple Deprivation. RESULTS: The areas least green in the summer were found to have higher mortality rates but no association was found between mortality and winter-greenness. The largest relative differences of mortality were around 9%. High levels of summer-greenness were associated with an increase in mental health prescription rates but areas with the highest differences between summer and winter-greenness had lower prescription rates than other ‡ © Hyam R. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. areas. The largest relative difference in prescription rate was 17%. All models controlled for overall deprivation. It is hypothesised that the year round greenness of mown grass is associated with increased mental health prescriptions and obscures the benefits of other kinds of vegetation on both mortality and mental health. DISCUSSION: There is an association between greenness and mortality and greenness and mental health. The association is both statistically significant and large enough to be of importance for policy making. Higher levels of non mown grass vegetation may be preferable for human wellbeing but more detailed understanding of the diversity of plant life in urban areas and how people related to it is required to make more specific recommendations.


Introduction
There is considerable body of work showing the positive effect greenness and green space has on human health in urban areas, recent reviews of evidence include, Amano et al. 2018, Berto 2014, Fong et al. 2018, Gascon et al. 2016, Hartig et al. 2014, Maas et al. 2009, Rojas-Rueda et al. 2019, WHO 2016, Zhang et al. 2017. The current state of knowledge is well rehearsed in these comprehensive reviews and will not be repeated here.
Two recent papers from North America have shown a decrease in mortality with increased greenness around subjects' home addresses, as measured by Normalised Difference Vegetation Index (NDVI) , Crouse et al. 2017. These studies were large in scale, both involved tracking over one hundred thousand individuals for eight to ten years. A third study (Xu et al. 2017) took and ecological approach to mortality at the neighbourhood level in Hong Kong and found significant associations between NDVI and cardiovascular and diabetes mortality.
Taking these papers as an inspiration, this ecological study examines whether this effect is visible at the population level in urban areas of Scotland. Three broad research questions are addressed.

1.
Is there a detectable relationship between mortality rates and greenness near peoples' homes in urban Scotland? 2.
Is there a relationship between mental health prescription rate and greenness that may indicate a mediating effect on mortality? 3.
Are urban areas of multiple deprivation also deprived of greenness?

Methods
The study is primarily concerned with four variables: Greenness (satellite measured NDVI), Deprivation (the Scottish Index of Multiple Deprivation, SIMD), Mortality (Standardised Mortality Ratio) and the proportion of population being prescribed drugs for anxiety, depression or psychosis. Mental health prescription rate is included because mediation analysis by  suggested that mental health was important in explaining the association between greenness and mortality. Other recent studies have also shown links between green space (Helbich et al. 2018) and street trees (Taylor et al. 2015) and mental health prescription rates and it has already been established that urban areas have higher mental health prescription rates compared to rural areas in Scotland (McKenzie et al. 2013).

Sampling rationale
A wide range of population statistics in Scotland are compiled to geographical units called Data Zones (DZ) (Suppl. material 1). One of these statistics is a summary statistic, the Scottish Index of Multiple Deprivation (SIMD). SIMD ranks the 6,976 DZs by overall deprivation based on indicator statistics in seven domains (Employment, Income, Health, Crime, Housing, Education, Access). Two of the indicators in the Health domain are Standardised Mortality Ratio (SMR) and Depression (The proportion of population being prescribed drugs for anxiety, depression or psychosis).
DZs were originally based on primary school catchment areas and then manually adjusted to produce polygons that are a compact shape, take account of physical boundaries, were felt to be homogeneous and have a population of between five hundred and one thousand (Office of the Chief Statistician 2004). DZ are effectively communities just large enough to assure anonymity in openly published statistical data. Because the DZ creation aimed to produce a comprehensive and contiguous division of the country all areas have to be attributed to a zone even if they are not a populated place. This produces problems when sampling for greenness of a DZ. As natural boundaries were followed green spaces were not divided between DZs but attributed to single zones. A park, field or golf course will be included within the DZ for a nearby community simply because it has to be put somewhere. Many DZs have large green spaces tagged onto them that are not necessarily representative of the lived experience of the people in that DZ whilst others appear to have little green space despite being adjacent to large, open areas. An example of this in the centre of Edinburgh is DZ S01008684 which includes not only an area of dense housing and the Scottish Parliament buildings but also the adjacent 260 hectare Holyrood Park. The park is not included in other adjacent DZs which have very little greenness. This is a common problem. The overall greenness of a statistical area may not reflect the exposure to greenness of a resident of that area especially if they live on the edge of it (Xu et al. 2017). Assessing green space by looking at the NDVI level for whole DZs would therefore be inappropriate. It would grossly overestimating some DZs and underestimating others. An approach was therefore devised that examined circular buffer areas around exemplar addresses based on the small user postcodes database. These results were then averaged for each DZ for comparison with the SIMD data that is released at the DZ level. Where greenness of a DZ is discussed in this study it is the average greenness of addresses sampled within that DZ at a particular buffer size not the overall greenness of the geospatial polygon defining the DZ. Geographic coordinates of postcodes was considered a good proxy for locations of a subset of actual homes in urban areas. This was not considered to be a reasonable assumption to make for rural postcodes and so the study was restricted to urban areas (Fig. 1).  looked at 250 m and 1,250 m buffer zones around dwellings. Crouse et al. (2017) examined greenness in buffers of 250 m and 500 m and found stronger associations with the 500 m buffers. They suggested further studies should explore this relationship. Both studies used data from the NASA MODIS (Moderate Resolution Imaging Spectroradiometer) instrument which has a 250 m resolution which the Nyquist-Shannon sampling theorem would suggest is at the limit of resolution at the 500 m buffer size. This study uses data from the Landsat 8 satellite, which has a higher 30 m resolution and examines buffers around addresses of 100 m, 250 m and 500 m (Fig. 2).  took seasonal variability into account when calculating overall exposure to greenness during the study period but not the subjects' experience of seasonality. Crouse et al. (2017) only looked at average, maximum summer NDVI because much of the ground will be covered by snow in the Canadian study area in winter months. This study includes average summer and average winter NDVI as well as the difference between these two calculated on a per address basis. This measure of seasonality was introduced as an initial attempt to assess nature of vegetation making up the greenness around addresses.   21, 205/21, 206/21, 204/20, 205/20, 206/20 WRS-2 (NASA 2018) for the calendar years 2013 to 2016. These scenes cover most of Scotland and include all the major urban areas. A full list of the 333 products is given in Suppl. material 2.
All of Scotland is over 54° North and so for many satellite images the sun is at too low an angle to give reliable surface reflectance data especially in the winter months. Scotland also has an oceanic climate so the ground is often obscured by cloud or mist. To build a detailed, contiguous NDVI map of the whole country therefore requires combining images taken on many satellite passes especially if points are to be sampled multiple times to overcome measurement errors. The images downloaded from USGS were therefore combined. A cloud free version of each NDVI image was created by setting the pixels that corresponded to cloud, snow or water in the Quality Assurance Assessment band to NA. These cloud free images were then combined into a single, mosaic stack of images to cover all of the study area and then averaged down to a single layer as a tiff image. This was done for two seasonal periods, Winter (October, November, December of 2013, 2014, Addresses were considered to be urban if they were scored as less than three on the Scottish Government's eight fold Urban Rural Classification, falling in settlements of greater than ten thousand people. The urban/rural classification was at the granularity of the postcode but analysis was taking place at the level of the DZ. Some DZs contained postcodes with different urbanisation scores. Subsequent analysis was therefore carried out on two sets of DZ, mixed-urban (n=4,883) which contained postcodes of different scores and a subset of these, pure-urban (n=2,246), that only contained postcodes scored as 1 in the Urban Rural Classification.

Sampling NDVI
The geospatial point data for each postcode was extracted from the small users postcode dataset. The gBuffer function from the rgeos package in R was used to create circular polygons around each of these points with radii of 100, 250 and 500 m. The NDVI values of the pixels within these buffer zones was averaged for both the summer and winter maps. A further three values were added for each point by calculating the difference between summer and winter at each buffer size as an absolute value. Each postcode then had nine associated values. The mean and standard deviations for each of these values for postcodes within each DZ were calculated resulting in nine measures (and standard deviations) for each DZ that were taken forward into further analysis. DZs were ranked into quintiles of average greenness both for the full set as well as the pure-urban DZs (Suppl. material 4).

Controlling for Deprivation
Throughout the analyses the model was controlled for deprivation using SIMD which is a summary statistic containing both the SMR and the Depression (prescription rate) data as two of its indicators. This biases the results to some degree. The possibility of constructing two new measures of deprivation (SIMD minus Mortality and SIMD minus Depression) was considered. It is apparent from the initial results that the full SIMD has a very small influence on the associations being measured. Any new measures would have even smaller effects, would not change the conclusions being drawn and would be a potential source of error so only the full SIMD, as published, was used here.

Results
The total number of postcode sample points was 91,357. Because of missing data not all of these had NDVI measurements for every pixel. Samples where the buffer had greater than 5% or more missing data were excluded from further analysis. This resulted in 90,689 for the 100 m buffer size, 90,060 for the 250 m buffer size and 89,537 for the 500 m buffer size going forward into the analysis. The average number of sample points per datazone for each combination of season and buffer size was 19.27 (standard dev 10.4).
The DZs used in the study have a total population of 3,776,237, the pure-urban subset represent 1,795,419 people. These are 71% and 34% of the Scottish population of 5.3 million respectively.
There is a constraining relationship between summer-greenness, winter-greenness and seasonal-difference. An urban area that has no vegetation in the summer is unlikely to develop it in the winter to any great extent and an area that is green in the winter is likely to continue being green in the summer. Despite this a minority of DZs showed small negative seasonality. For the 250 m buffer this was 10% both for mixed and pure-urban samples. For mixed-urban the average negative NDVI was 24% of the average positive NDVI, 20% for pure-urban areas. Aerial photos on Google Maps of sites of negative seasonal NDVI change showed a predominance of deciduous trees over lawns from which it was assumed that a major cause of negative NDVI may be less productive tree foliage obscuring highly productive grass during the summer months. This assumption warrants further investigation but for the purposes of this study the absolute value of NDVI difference between summer and winter was used.
The constraints were also evident in the overlap in membership of the NDVI quintiles. For the 250 m buffer the greenest winter quintile shared 62% of its DZs with the greenest summer quintile because areas, typically of mown grass, that were green in the winter were also likely to be amongst the most green areas in the summer. For the least green quintiles 71% of DZs were shared between winter and summer because areas that lacked vegetation were never going to be green either in winter or summer. Of the 449 most seasonal pure-urban DZs only 2% were also in the most winter-green quintile whilst 43% were in the least winter-green quintile.
Following exploratory Pearson's correlation analyses between all variables an Analysis of Covariance (ANCOVA) was used to model mortality/prescription rate conditional on different levels of greenness whilst controlling for deprivation.
There did not appear to be any meaningful correlations in the exploratory stage of analysis. For summer-greenness r ranged from 0.014 to 0.151 for pure-urban for both mortality and prescription rate. It was lower for mixed-urban. For winter-greenness r ranged from 0.007 to 0.031. Scatter plots of the NDVI against each of the variables didn't reveal any obvious patterns. The only measure that might be considered weakly correlated was prescription rate against winter-greenness with r ranging from 0.287 to 0.230 for pure-urban areas. This suggests that increased greenness near the home in the winter may be associated with increased (not decreased) use of prescriptions for mental health.
James et al. 2016 specifically reported a 12% lower rate of all-cause non-accidental mortality for subjects living in the highest quintile of greenness compared to those in lowest quintile and so comparisons were made between greenest and least green quintiles for summer and winter as well as the most seasonal and lease seasonal quintiles controlling for mental health prescription rate (SIMD) using the two linear models "Mortality G reenness Quintile + SIMD Rank" and "Depression ~ Greenness Quintile + SIMD Rank" in the R function lm. For the pure-urban datasets this gave significant differences between the extreme quintiles for all combinations of greenness measure and buffer size apart from mortality versus seasonality (See Tables 1, 2  The mortality rate (SMR) is centred on a value of 100 so unsurprisingly the mean for all DZs in the SIMD dataset is 99.766 (sd=44.866). The mean for DZs in this study was 104.301 (sd=49.060) and 107.853 (sd=54.336) for the pure-urban areas showing an underlying increase in mortality rates in more urbanised areas but with a great deal of variability. The largest differences shown in Table 1 therefore represent a relative difference in mortality of about 9%. Results for the mixed-urban DZs were similar if slightly stronger and better supported. For the 250 m buffer areas summer Δ was -11.007 (t= -5.761983, p= < 0.001), winter was -8.532 (t= -4.507, p= <0.001) whilst seasonal-difference still failed to produce a significant result.
The mean prescription rate for all DZs in the SIMD data is 17.55% (sd=5.01). The mean for all DZs in this study was 18.07% (sd=5.26) for mixed-urban and 17.67% (sd=5.66) for pure-urban areas. The largest differences shown in Table 2 therefore represents a relative difference in prescription rate in the region of 17%. Results for prescription rate in the mixed-urban DZs were similarly distributed but smaller and less well supported compared to those for pure-urban areas. For the summer 250 m buffer areas summer Δ was 1.34% (t = 9.35790, p < 0.001), winter was 1.46% (t = 10.368, p < 0.001) and seasonal difference was -0.36% (t = -2.532, p = 0.01).
Finding significant differences between the highest and lowest quintiles but not an overall linear relationship suggested a more complex pattern. In an ANCOVA for quintiles of all combinations of greenness measure against mortality and prescription rate for both mixedurban and pure-urban sets controlling for SIMD there were significant differences (p < 0.01) in all combinations apart from two cases. No significant difference between quintiles was found for any buffer sizes in winter-greenness verses mortality in the pure-urban subset of DZs and in seasons verses mortality for the 100 m buffer size in the mixed-urban set.
To understand the size of the statistically significant differences found Cohen's d was calculated for the differences between highest and lowest quintiles for 250 m buffers for pure-urban DZ. Summer-greenness was associated with a small to medium changes in both mortality and prescription rate (Cohen's d 0.265 and 0.296). Seasonal-difference was associated with medium size differences in mortality (Cohen's d 0.467). Winter-greenness and seasonal-difference were associated with medium to large differences in prescription rates, albeit in opposite directions (Cohen's d 0.728 and 0.764).
To help interpret these relationships Tukey's HSD test was carried out for each combination of variables for the five quintiles G1 (least green or least seasonal) to G5 (greenest or most seasonal). Significance in these comparisons is taken as p < 0.01.

Mortality and summer-greenness
For pure-urban DZs at 250 m G1 was significantly higher than all other quintiles (p < 0.01) with the difference ranging from 11.550 to 15.231 SMR. G2 through G5 didn't vary from each other significantly at all. This pattern was stronger at 100 m buffer size but broke down at 500 m buffer size. For the mixed-urban DZs the pattern was similar but didn't break down at the 500 m buffer size.

Mortality and winter-greenness
As per the ANCOVA for pure-urban DZs there were no significant differences between quintiles at any buffer size. For the mixed-urban DZs G5 was significantly lower than other quintiles at 500 m but this broke down at 100 m and 250 m with low levels of significance of any differences between quintiles compared.

Mortality and seasonal difference
For pure-urban DZs G5 was significantly lower than all the other quintiles for 250 m which held for 100 m and 500 m. For mixed-urban DZs G1 and G5 had significant lower mortality rates at larger buffer sizes only but the confidence levels were low and the difference disappeared for the 100 m buffer size.

Mental health prescription rate and summer-greenness
For pure-urban DZs there were significant differences between most of the quintiles with the largest being between G1 and the rest. Notably higher levels of greenness were associated with higher prescription rates. For mixed-urban this did not hold. There were significant differences between some combinations of quintiles but confidence levels were low apart from at 500 m where G5 stood out as having small but significant reduced prescription rates compared to the other quintiles.

Mental health prescription rate and winter-greenness
For pure-urban DZs there were significant differences between most of the quintiles which were maintained at all buffer sizes. Areas of higher winter-greenness had higher prescription rates. These differences are visualised in Fig. 3. There is a suggestion of a bimodal distribution in the density plot. The relationship changed for mixed-urban DZs where G5 no longer had significantly higher rates.

Mental health prescription rate and seasonal difference
The significant differences between quintiles approximated a mirror image of that for pureurban DZs, winter-greenness for 250 m and 500 m buffers. These differences are visualised in Fig. 4 which also indicates a possible bimodal aspect to the distribution. At 100 m only G5 had a significantly lower rate with no significant differences between the other quintiles. The same pattern of variation was found in the mixed-urban DZs.

Discussion
Mechanisms for how the benefits of greenness and greenspace on health may occur fall into three broad areas: Physiological mechanisms involve vegetation ameliorating the hostile aspects of the built environment (Beckett et al. 1998, Dadvand et al. 2012, Nowak et al. 2006); psychological mechanisms center on two main theoretical frameworks Attention Restoration Theory (Kaplan 1995, Basu et al. 2018 and Stress Recovery Theory (Ulrich 1983, Ulrich et al. 1991; Social mechanisms including increased community interactions and reduction of loneliness and isolation (Ward Thompson et al. 2016). There is likely to be a complex relationship between these mechanisms. Having a pleasant, safe place to walk with friends requires both physiological, psychological and social services to be facilitated by urban greenery. Most studies, like this one, have been in developed countries and it is uncertain if the results can be extrapolated to developing regions (Amano et al. 2018).
Urban areas in Scotland fall within a single, oceanic climatic zone at or close to sea level. Although there are many deciduous trees in some areas there is also a large amount of mown grass which grows nearly year around. One would therefore expect a lower level of variance in NDVI than in the North American studies which cover a range of climatic zones and altitudes. The published figures suggest this may be the case. The mean NDVI for 250 m buffer zone of the current study is 0.48 (sd=0.08) while James et al. (2016) had 0.47 (sd=0.12) and Crouse et al. (2017) 0.58 (sd=0.11) although quantitative comparison of NDVI between studies is error prone and no formal conclusions should be drawn from these numbers it appears that even with the low variance in Scotland there was still sufficient information to observe associations between greenness and mortality and greenness and prescription rates.
The association of mortality rates with greenness needs to be interpreted in the context of the constrained pattern of seasonal differences, that areas green in winter are also green in summer but areas not green in summer are never green in winter. No overall relationship between winter-greenness and mortality was found, only a relationship between summergreenness and mortality, yet a large amount of that summer-greenness must be the result of the same vegetation that causes winter-greenness. This is confirmed by the seasonaldifference data where differences in mortality are again found suggesting that the relationship is only with summer-greenness. Only the least summer-green areas have a detectable rise in mortality rate suggesting that there is an absolute level of summergreenness beyond which no benefit to mortality rates occurs. The influence of greenness on mortality has to approach an asymptote at some level and, because of the low variance in Scottish urban greenness, it is plausible the majority of the urban areas are at or above this level. Any planning policy might therefore only look at increasing summer-greenness in a restricted range of areas. Alternatively differences in mortality at higher levels of summergreenness may be obscured by associated winter-greenness and if this could be controlled for the relationship may be more linear. This is supported by the most seasonal areas having a fall in mortality rate compared to the rest.
The findings on mental health prescription rate were initially surprising. It was expected that, as a potential mediator of mortality rate, they would be found to vary either in the same direction as the mortality rate or not at all. The observed relationship could be entirely spurious but the size of the differences and level of statistical support as well as the robustness between pure-urban and mixed-urban datasets make this unlikely. There could be a direct causal link between winter-greenness and prescription rates but it seems unlikely that mown lawns cause mental health problems or attract people with poor mental health. A more likely explanation is a confounding variable that leads to the geospatial association of a this kind of vegetation and mental illness. Deprivation would be the lead candidate but an attempt was made to controlled for this in the analysis. Identification of the confounder will require other studies but the restricted range of climatic characteristics of urban Scotland rule out some of the physiological candidates such as heat stress. Pollution may be a candidate in association with busy roads but psychosocial mechanisms may also play a part (Dzhambov et al. 2018).
Despite varying in the opposite direction to mortality, mental health prescription rate might still be a mediator of mortality rate or share a confounding variable. Prescription rates fall with increased seasonal difference just as mortality rates do. More deciduous areas have lower mortality rates and reduced prescription rates suggesting, again, that the nature of the vegetation is as important as the amount. A recent study from the Netherlands (Helbich et al. 2018) also showed a nonlinear relationship between mental health prescription use and greenness. The study was based on the detailed Dutch land use classification scheme rather than direct measures of NDVI and used a single notion of greenness but another study in New York (Reid et al. 2017) at a similar resolution compared grass with tree cover and found that tree cover had a stronger association with self reported health than grass, even when controlling for closeness to public parks. Other studies have found that tree cover in particular is negatively associated with depression (Taylor et al. 2015) or positively associated with school performance (Donovan et al. 2018, Sivarajah et al. 2018, Li et al. 2019 or vascular health (Donovan et al. 2015, Donovan et al. 2013. A lab based study has recently shown a higher restorative potential for seasonal images even in designed digital landscapes (Kuper 2018). But some studies fail to find associations with natural environments in general and mental health prescribing (Gidlow et al. 2016) .
This more complex relationship between greenness and health has implications for the suggested mechanisms as to how nature effects well-being. Publicly accessible, mown parkland might appear to increase the overall greenness in a community but have little effect on noise and air pollution or provide meaningful social spaces that are visited throughout the year. In contrast street trees and a patchwork of private or shared gardens may be more beneficial, provided they are tended in a way that maintains sufficient volume of seasonally diverse vegetation. Looking at new mothers Feng and Astell-Burt (2018) concluded "... it may be how mothers perceive green space nearby and what those spaces enable them to do, rather than simply how much there is overall, that is important for promoting mental health in the postpartum period."

Conclusions
There is a relationship between greenness around peoples' homes and standardised mortality rates at the DZ level (Research question 1). Least green areas have higher mortality rates. This relationship only holds for levels of greenness in summer, not winter. It only exists for the least green 20% of areas, the other 80% being indistinguishable. Areas with highest levels of seasonal difference have lower mortality rates compared to others suggesting that the kind of vegetation plays a role. The sizes of differences in mortality are small to medium representing a relative difference of around 9%.
Higher levels of greenness around homes is associated with higher mental health prescription rates especially greenness in the winter (Research question 2). There is a stronger relationship closer to the home. The opposite is true of seasonal difference. Higher seasonal difference is associated with lower prescription rates. The differences are medium to large in size with the largest seasonal difference having a relative difference in prescription rate of 17%.
Because relationships between greenness and mortality and mental health prescription rates hold despite controlling for deprivation it is unlikely areas of higher deprivation are also deprived of raw greenness in general (Research question 2). This does not exclude the kind of vegetation playing a role in deprivation or lack of greenness playing a role under certain circumstances.
This is an ecological, observational study and no causation can be concluded from the fact that variables measured vary together however the context, level of support and size of the variation suggests that the variables may form parts of chains of causation rather than their association merely being chance. Further investigation of these relationships requires patient level data on the one hand and more granular information about kinds of vegetation and people's lived experience of that vegetation on the other. Our urban environments are complex cultural constructs which require a truly interdisciplinary approach to understand.

Suppl. material 2: Landsat Products Used
Authors: USGS Data type: Geospatial Brief description: List of products downloaded from USGS and used to create summer and winter maps. Download file (13.34 kb)

Suppl. material 3: Small User Postcodes
Authors: Scotland Registry Data type: Geospatial Brief description: List of small user postcodes as downloaded from Scotland Registry. Download file (7.07 MB)

Suppl. material 4: Results Data Tables
Authors: Hyam, Roger Data type: Geospatial Brief description: Nine data tables containing the NDVI and SIMD data used in the analysis. Download file (1.18 MB)