Skyglow extends into the world's Key Biodiversity Areas

The proportion of the Earth's surface that experiences a naturally dark environment at night is rapidly declining with the introduction of artificial light. Biological impacts of this change have been documented from genes to ecosystems, and for a wide diversity of environments and organisms. The likely severity of these impacts depends heavily on the relationship between the distribution of artificial night‐time lighting and biodiversity. Here, we carry out a global assessment of the overlap between areas of conservation priority and the most recent atlas of artificial skyglow. We show that of the world's Key Biodiversity Areas (KBAs), less than a third have completely pristine night‐time skies, about a half lie entirely under artificially bright skies and only about a fifth contain no area in which night‐time skies are not polluted to the zenith. The extent of light pollution of KBAs varies by region, affecting the greatest proportion of KBAs in Europe and the Middle East. Statistical modelling revealed associations between light pollution within KBAs and associated levels of both gross domestic product and human population density. This suggests that these patterns will worsen with continued economic development and growth in the human population.


Introduction
The erosion of the night-time through the introduction of artificial lighting, from street-lighting and other sources, has pervasive environmental impacts. These span changes in the physiology and behaviour of individual organisms, in the abundance and distribution of species, in the structure and functioning of ecological communities and in the provision of ecosystem services (Gaston et al., , 2014. A wide diversity of terrestrial, freshwater and marine organisms is influenced, including microbes, plants and many groups of animals (e.g. crustaceans, molluscs, insects, fish, amphibians, reptiles, birds, mammals; Gaston et al., 2013;Bennie et al., 2016).
To date, attention has focussed foremost on the environmental impacts of the direct emissions from sources of night-time lighting. This is, however, only relatively narrowly spatially distributed compared with the skyglow that is caused by upwardly emitted or reflected artificial light being scattered in the atmosphere by water, dust and gas molecules. The latter has been estimated already to extend over 23% of the global land area (Falchi et al., 2016). It can increase background sky brightness to levels comparable to those of late twilight and moonlight, and can obscure the visibility to humans of individual stars and the Milky Way (83% of the human population lives under light polluted skies; Falchi et al., 2016). It is likely to interfere with multiple biological processes including activity patterns of diurnal, crepuscular and nocturnal species (which are variously sensitive to timing of twilight and night-time levels of moonlight; e.g. Moore et al., 2000;Bachleitner et al., 2007), and nighttime orientation and navigation (which often involves the use of stars and other celestial objects; e.g. Foster et al., 2017). Some free-living organisms are able to detect and respond to extraordinarily low levels of night-time light (e.g. Warrant et al., 2004).
Particularly because of the contribution of light that is emitted, or reflected, at relatively shallow angles to the horizontal, skyglow can extend substantial distances (up to hundreds of kilometres) from urban sources of night-time lighting (Luginbuhl, Boley & Davis, 2014). This raises the potential for impacts of artificial night-time lighting to reach many globally important biodiversity areas, even when these are reasonably remote from many other anthropogenic pressures. In this paper, we estimate the extent of this overlap, based on the most recent global modelling of skyglow and the distribution of Key Biodiversity Areas (KBAs).

Data
Estimates of global variation in skyglow were obtained from Falchi et al. (2016). This surface was produced by the modelling of measured upward radiance from artificial sources from satellite imagery (from the VIIRS DNB sensor on the Suomi National Polar-orbiting Partnership satellite), and ground measurements. The data are presented for the entire area (terrestrial and marine) between approximately 85°N and 60°S, at a spatial resolution of 30-arcseconds (~1 km) as an artificial brightness level (mcd/m 2 ). The authors define the level of artificial brightness under which a sky can be considered 'pristine' as up to 1% above the natural background level (ratio of 0.01; artificial brightness = 0.0017 mcd/m 2 ). At a level of 8% or more above natural conditions (ratio of 0.08; artificial brightness = 0.014 mcd/m 2 ) light pollution extends from the horizon to the zenith and the entire sky can be considered polluted. We use these two thresholds in our analysis.
Key Biodiversity Areas are sites that contribute significantly to the global persistence of biodiversity, and use quantitative criteria to identify places that support viable populations of species for which site-scale conservation is appropriate (IUCN, 2016). The current inventory includes Important Bird and Biodiversity Areas, KBAs identified during Biodiversity Hotspot profiles supported by the Critical Ecosystem Partnership Fund and Alliance for Zero Extinction sites (http://www.keybiodiversityareas.org/home). We used boundary data (BirdLife International, 2016) for all 14 979 KBAs for which polygons were available in December 2016 (comprising >95% of all KBAs recognized in March 2018). Only those that lie within the extent of the skyglow data were considered further (14 765). These had a combined coverage of 34 785 579 km 2 , with a median area of 209 km 2 (calculated in ArcGIS 10.3.1, WGS84 datum; ESRI Inc.). The KBA dataset includes characterization by region, and this includes a 'Marine' category identifying those which are predominantly marine based.
Data were also obtained on three further variables that are of potential significance as determinants of the levels of skyglow experienced by KBAs (excluding those on the 'High seas'): gross domestic product (GDP; which provides a measure of economic activity), human population density and protected area coverage: (1) Median GDP per capita, adjusted for purchasing power parity (PPP), was extracted from the gridded (spatial resolution of 5 arc-min) product by Kummu, Taka & Guillaume (2018). GDP grid cells were included in the median calculation if the centre was within the KBA boundary (missing n = 792; 95% of KBAs retained for further analysis). This product uses subnational datasets where possible in combination with national datasets (including the World Bank Development Indicators database and the CIA World Factbook for missing data). (2) Median human population density for each KBA was extracted from the gridded population density data (resolution of~1 km) of the world projected for the year 2015 from census data (CIESIN, 2016) where data were included when the centre of the grid cell was located inside the KBA boundary (missing n = 646). And, (3) protected area coverage for each KBA was calculated using the World Database of Protected Areas (IUCN and UNEP-WCMC, 2016). Protected areas were only included where they had associated spatial boundary data.

Data analysis
Analyses were carried out in the software package R (version 3.5.3) (R Development Core Team, 2014) using the packages 'raster' (version 2.8-19) and 'rgdal' (version 1.3-6). Values of artificial sky brightness were extracted from the global dataset for each KBA. Sky brightness values were stored as a grid (raster) with a spatial resolution of 30-arcseconds (~1 km). Values were included when the centre of the corresponding cell was located within the KBA boundary. Two measures were calculated for each KBA: (1) percentage coverage by pristine night-time skies (artificial brightness of ≤0.0017 mcd/m 2 ratio values ≤ 0.01) and (2) percentage coverage by skies not polluted from the horizon to zenith (artificial brightness of 0.014 mcd/m 2 ratio values < 0.08). In addition, we calculated the percentage coverage of the total area of KBAs by these respective levels of sky brightness.
The percentiles (0-100) were calculated for both the median KBA GDP (per capita, PPP, thousand $) and median KBA population density (people per km 2 ). The median proportion pristine for the KBAs with each combination of GDP and population density percentile was then calculated and plotted. The use of percentiles maximizes the evenness of the sample sizes for each combination.
We used a generalized linear model, with a binomial error structure and logit link function to model whether a KBA had entirely pristine skies or not as a function of GDP per capita, population density, the interaction between the two and the proportion of the KBA that falls within a protected area. To account for strong right skew in population density we applied a log10 transform to that data. Proportion protected had a bimodal distribution with high proportions of KBAs either fully or not protected at all. This variable was dichotomized as either fully protected (100% protected) or not (<100% protected). Model fit is given by McFadden's pseudo-R 2 . Although there was significant pair-wise correlation between the predictor variables, the variance inflation factors were all <3 (Supporting Information Tables S1 and S2).

Results
Over two-thirds (68.6%) of the total number of KBAs assessed contained no area with pristine night-time skies, while less than one fifth (17.8%) had completely pristine night-time skies (ratio values ≤ 0.01 up to 1% above natural conditions; Fig. 1a). Europe had the greatest percentage (94%) of KBAs containing no area of pristine skies, followed by the Middle East (88%) and the Caribbean (77%; Fig. 2a). The only region in which all KBAs had completely pristine skies was the marine region. In Antarctica, 93% of the KBAs were entirely pristine (Fig. 2a).
Nearly one half (46.9%) of all KBAs consisted entirely of area in which night skies were polluted to the zenith (Fig. 1b). However, nearly one third of KBAs (30.7%) were completely free of skies polluted to the zenith (ratio values < 0.08 up to 8% above natural conditions; Fig. 1b). Europe was the region where the greatest percentage (75.6%) of KBAs had night-time skies entirely polluted to the zenith, followed by the Middle East (72.6%) and the Caribbean (59.9%; Fig. 2). In Antarctica, no KBAs were entirely polluted to the zenith.
Of the summed global area of KBAs, over a quarter (26.8%) was not pristine and 14.4% was polluted to the zenith (Fig. 2b). The Middle East was the region with the greatest percentage (80.1%) of the summed KBA area not having pristine skies, followed by the Caribbean (67.9%) and Europe (66.6%; Fig. 2b). These were also regions with the largest percentage of the overall KBA area that is polluted to the zenith (Middle East -49.9%, Europe -49.1%, Caribbean -28.2%; Fig. 2b).
The likelihood of a KBA having pristine skies decreased with increasing GDP and population density, and the interaction between the two, and increased with proportional coverage by protected areas (Table 1). For each ten-fold increase in population density, being pristine was 0.7 x as likely (Odds ratio (OR) = 0.71, 95% Confidence Interval (CI) = 0.66 -0.76). With every increase in GDP by $1000 there was an associated decrease in the odds of being pristine by 7% (OR = 0.93, 95% CI = 0.92 -0.93). KBAs which were fully protected were 1.39 as likely to be pristine compared to those not fully protected (OR = 1.29, 95% CI = 1.04, 1.60). There was also a significant interaction between population density and GDP (Table 1; Fig. 3).

Discussion
Skyglow is often envisaged as an exclusively urban issue. However, both modelling and ground measurements have shown it to be very widespread, often being propagated over  long distances from sources (Kyba et al., 2015;Falchi et al., 2016). Nevertheless, it is perhaps surprising that less than one fifth of KBAs had completely pristine night-time skies, more than two-thirds contained no area with pristine skies and that over a quarter of the total area of KBAs was light polluted. This is especially so considering that the global data on skyglow are likely, if anything, to be conservative estimates of its extent. Such data are primarily estimated using satellite measurements and, as such, are for 'open sky' conditions (cloud-free). They, therefore, do not take into account the amplification of skyglow that can occur by cloud cover (Jechow et al., 2017), and likely under-represent the occurrence of artificial brightness of horizons (which may be important for many organisms, e.g. in influencing predatorprey interactions). These results are especially troubling because it has become increasingly apparent that organisms can respond even to absolutely small (which may nonetheless be relatively large) changes in natural night-time light conditions (Warrant et al., 2004), in perceived day lengths (Gaston et al., 2017), and in artificial night-time lighting (Gaston et al., 2014). The breadth and number of species whose behaviour is influenced by skyglow seems likely to be large.
Unsurprisingly, the likelihood of the skies of a KBA experiencing skyglow tends to increase in countries with higher GDP, and in areas with higher human population  density, consistent with global results (Gallaway, Olsen & Mitchell, 2010). It is possible to have KBAs in regions with relatively pristine skies in areas with high human densities when these populations are economically poor, and likewise in areas with high GDP when these populations are at low density (Fig. 4). However, under virtually all other circumstances KBAs exist under light polluted skies. This strongly suggests that globally the proportion of KBAs experiencing skyglow will almost certainly increase in parallel with developing country economies.
The relatively small areal extent of many KBAs means that most commonly their skies are either entirely light polluted or entirely unpolluted. Whether they have coverage by protected areas is also associated with whether they experience skyglow, consistent with results for upwardly emitted light . However, this seems in most cases unlikely to be a consequence of protection per se. More likely, it is a result of the tendency for protected areas to be distributed away from urban centres (and hence sources of artificial light), and often in regions with reduced competition over land use (Gaston et al., 2008).
Here, we have identified the regions where KBAs are most affected by skyglow. This effectively prioritizes areas for further detailed assessment of the potential risk to species of conservation priority inhabiting these. To assess the risk, the extent of overlap between skyglow and species ranges and occupancy could be calculated. A similar technique has been applied to cacti and mammal species ranges, where the overlap with light pollution was calculated for upwardly emitted light radiation measured directly from satellites Correa-Cano et al., 2018). Further research could also identify areas currently minimally affected by light pollution but at risk of increases in future due to rapidly developing economies or with increasing populations, facilitating the development of mitigation measures in these areas.
Dramatic reductions in anthropogenic pressures on biodiversity are often costly to achieve, and there are commonly substantial lag times between such reductions and biodiversity responses. By contrast, marked reductions in skyglow could be achieved by limiting outdoor artificial lighting to levels and places where it is required by people, which would result in considerable cost savings without undermining the benefits that it brings. Indeed, one might argue that environmental benefits and financial savings from considered lighting policies are closely aligned (Gaston, 2013).

Supporting information
Additional supporting information may be found online in the Supporting Information section at the end of the article. Table S1. Pair-wise Kendall's correlations. Table S2. Generalized variance inflation factor [R, 'vif' function, 'car' package (version 2.1-3)].
Corrections added on 24 February 2020, after first online publication: • The fifth sentence has been revised in the Abstract section on the first page • The third, fourth and fifth sentences have been revised in the first paragraph and the fourth sentence in the third paragraph of the 'Data' section under the 'Materials and methods' section on the second page • The first and fifth sentences have been revised in the first paragraph and the second and sixth sentences in the third paragraph of the 'Data analysis' section under the 'Materials and methods' section on the second page • Figure 1 and legend have been revised on the third page • Paragraphs one, two, three and four have been revised in the 'Results' section under the 'Materials and methods' section on the second, third, fourth and fifth pages • Figure 2 and legend have been revised on the fourth page • Table 1 and legend have been revised on the fourth page • Figure 3 and legend have been revised on the fifth page • Figure 4 and legend have been revised on the fifth page • The third sentence has been revised in the first paragraph in the 'Discussion' section on the fifth and sixth pages