Evaluation of Light Pollution in Global Protected Areas from 1992 to 2018

: Light pollution, a phenomenon in which artiﬁcial nighttime light (NTL) changes the form of brightness and darkness in natural areas such as protected areas (PAs), has become a global concern due to its threat to global biodiversity. With ongoing global urbanization and climate change, the light pollution status in global PAs deserves attention for mitigation and adaptation. In this study, we developed a framework to evaluate the light pollution status in global PAs, using the global NTL time series data. First, we classiﬁed global PAs (30,624) into three pollution categories: non-polluted (5974), continuously polluted (8141), and discontinuously polluted (16,509), according to the time of occurrence of lit pixels in/around PAs from 1992 to 2018. Then, we explored the NTL intensity (e.g., digital numbers) and its trend in those polluted PAs and identiﬁed those hotspots of PAs at the global scale with consideration of global urbanization. Our study shows that global light pollution is mainly distributed within the range of 30 ◦ N and 60 ◦ N, including Europe, north America, and East Asia. Although the temporal trend of NTL intensity in global PAs is increasing, Japan and the United States of America (USA) have opposite trends due to the implementation of well-planned ecological conservation policies and declining population growth. For most polluted PAs, the lit pixels are close to their boundaries (i.e


Introduction
Artificial light pollution is a general environmental concern by changing the natural form of brightness and darkness in the ecosystem, which further affects the phenology of insects and biodiversity in ecosystems [1], especially in protected areas (PA). The nighttime light (NTL) reflects human activities and can be monitored by satellites [2]. Given the rapid global urbanization [3,4], NTL has imposed noticeable pressures on natural systems and wildlife, causing ecological problems such as artificial light pollution [5]. PA plays a crucial role in global biodiversity conservation [6] and has been extensively used in international programs such as the Aichi biodiversity targets [7] and Sustainable Development Goals (SDGs) [8,9]. In these programs, PAs were regarded as critical elements to maintain the Remote Sens. 2021, 13, 1849 3 of 16

The Harmonized Global Nighttime Light Data
We collected the harmonized global NTL data from 1992 to 2018 as our primary data source of light pollution [2,47]. The harmonized NTL data integrate satellite observations from the Defense Meteorological Satellite Program Operational Line-scan System (DMSP/OLS) and the Visible Infrared Imaging Radiometer Suite Day Night Band (VI-IRS/DNB) data. The temporal spans of DMSP and VIIRS NTL data are from 1992 to 2013 and from 2012 to 2018, respectively, with a near-global coverage (i.e., 65 • S to 75 • N) [48]. The raw records in DMSP NTL data are digital number (DN) values (i.e., ranging from 0 to 63), notably different from the radiance in the VIIRS/DNB data. Additionally, their spatial resolutions are different, i.e., 30-arc and 15-arc in DMSP/OLI and VIIRS/DNB, respectively. To extend the interval of DMSP NTL observations, the harmonized global NTL data were generated using the kernel density approach [26,49] and the sigmoid function, which aggregates the 15-arc resolution to 30-arc seconds and converts annual VIIRS observations to the DMSP-like NTL data [2], respectively. Given that there are more uncertainties when the DNs in NTL data are lower, we screened out pixels with DN values less than 7 in our study as suggested in [2].

Methodology
We evaluated the light pollution of global PAs from 1992 to 2018 using the harmonized NTL and the WDPA data ( Figure 1). First, we categorized all PAs as three major types of light pollution: continuously polluted PAs, discontinuously polluted PAs, and nonpolluted PAs (Figure 1a). Then, we investigated the temporal trend of NTL in polluted PAs (i.e., continuously and discontinuously), where the NTL change in different regions around the PA was explored (Figure 1b). Details of the definition and temporal trends of NTL are presented in the following sections. more than 30,000 PAs were used in our study.

The Harmonized Global Nighttime Light Data
We collected the harmonized global NTL data from 1992 to 2018 as our primary data source of light pollution [2,47]. The harmonized NTL data integrate satellite observations from the Defense Meteorological Satellite Program Operational Line-scan System (DMSP/OLS) and the Visible Infrared Imaging Radiometer Suite Day Night Band (VIIRS/DNB) data. The temporal spans of DMSP and VIIRS NTL data are from 1992 to 2013 and from 2012 to 2018, respectively, with a near-global coverage (i.e., 65°S to 75°N) [48]. The raw records in DMSP NTL data are digital number (DN) values (i.e., ranging from 0 to 63), notably different from the radiance in the VIIRS/DNB data. Additionally, their spatial resolutions are different, i.e., 30-arc and 15-arc in DMSP/OLI and VIIRS/DNB, respectively. To extend the interval of DMSP NTL observations, the harmonized global NTL data were generated using the kernel density approach [26,49] and the sigmoid function, which aggregates the 15-arc resolution to 30-arc seconds and converts annual VIIRS observations to the DMSP-like NTL data [2], respectively. Given that there are more uncertainties when the DNs in NTL data are lower, we screened out pixels with DN values less than 7 in our study as suggested in [2].

Methodology
We evaluated the light pollution of global PAs from 1992 to 2018 using the harmonized NTL and the WDPA data ( Figure 1). First, we categorized all PAs as three major types of light pollution: continuously polluted PAs, discontinuously polluted PAs, and non-polluted PAs ( Figure 1a). Then, we investigated the temporal trend of NTL in polluted PAs (i.e., continuously and discontinuously), where the NTL change in different regions around the PA was explored (Figure 1b). Details of the definition and temporal trends of NTL are presented in the following sections.

Definition of Light Pollution Categories
We classified all PAs as three major types according to the spatial extent and year of lit pixels in/around them ( Figure 2). The definition of light pollution in our study is the occurrence of NTL light within the 50 km buffer zone of PA, as suggested in [36]. The first Remote Sens. 2021, 13, 1849 4 of 16 is PAs without any lit pixels (non-polluted) inside the PA during the whole period (1992 to 2018) (Figure 2a). By contrast, the second is continuously polluted PAs, which have had lit pixels in and around their buffer zones since 1992 (Figure 2b). A representative of this pollution category can be referred to as the case in Paris ( Figure S1), where the mean NTL around these PAs has been continuously increasing since the 1990s. Unlike those continuously polluted PAs, the discontinuously polluted PAs have lit pixels around during the middle of period 1992 to 2018, as illustrated in the case of Zhangjiajie (China) ( Figure S2) due to the rapid urban sprawl process. Considering that most non-polluted PAs are distributed in high-latitude regions, we focused on these polluted types (i.e., the continuously and discontinuously polluted).

Definition of Light Pollution Categories
We classified all PAs as three major types according to the spatial extent and year o lit pixels in/around them ( Figure 2). The definition of light pollution in our study is the occurrence of NTL light within the 50 km buffer zone of PA, as suggested in [36]. The firs is PAs without any lit pixels (non-polluted) inside the PA during the whole period (1992 to 2018) (Figure 2a). By contrast, the second is continuously polluted PAs, which have had lit pixels in and around their buffer zones since 1992 (Figure 2b). A representative of this pollution category can be referred to as the case in Paris ( Figure S1), where the mean NTL around these PAs has been continuously increasing since the 1990s. Unlike those contin uously polluted PAs, the discontinuously polluted PAs have lit pixels around during the middle of period 1992 to 2018, as illustrated in the case of Zhangjiajie (China) ( Figure S2 due to the rapid urban sprawl process. Considering that most non-polluted PAs are dis tributed in high-latitude regions, we focused on these polluted types (i.e., the continu ously and discontinuously polluted).

The Temporal Trends of Nighttime Light
We explored the temporal trend of NTL change using the Theil-Sen slope and Mann-Kendall test in polluted PAs (i.e., continuously and discontinuously polluted). We meas ured the slope of NTL change using the Theil-Sen trend approach, a robust and simple linear regression method [50]. We also calculated the significance level using the p-value provided by the Mann-Kendall approach [51]. Combining these two approaches for trend analyses is extensively used in studies such as ecological and environmental change [52,53].
To identify those buffers with lit pixels that have different NTL intensities around PA, we proposed two intervals of buffers (i.e., the first polluted and the high-intensity intervals) that can well-characterize the spatial extent of those sensitive regions. The firs buffer in a specific year was defined as the nearest buffer with lit pixels. Thus, the firs polluted interval (Figure 3a) was defined as between the furthest (e.g., the 30th buffer in 1992) and the nearest (e.g., the buffer in 2018) during the whole period. This interval re flects the NTL impact on PAs from a spatial distance, and those PAs having closer buffers with lit pixels are easier to be affected by the light pollution. Differently, we also defined the interval using the intensity of NTL (Equation (1)), in which both the DN values of li

The Temporal Trends of Nighttime Light
We explored the temporal trend of NTL change using the Theil-Sen slope and Mann-Kendall test in polluted PAs (i.e., continuously and discontinuously polluted). We measured the slope of NTL change using the Theil-Sen trend approach, a robust and simple linear regression method [50]. We also calculated the significance level using the p-value provided by the Mann-Kendall approach [51]. Combining these two approaches for trend analyses is extensively used in studies such as ecological and environmental change [52,53].
To identify those buffers with lit pixels that have different NTL intensities around PA, we proposed two intervals of buffers (i.e., the first polluted and the high-intensity intervals) that can well-characterize the spatial extent of those sensitive regions. The first buffer in a specific year was defined as the nearest buffer with lit pixels. Thus, the first polluted interval (Figure 3a) was defined as between the furthest (e.g., the 30th buffer in 1992) and the nearest (e.g., the buffer in 2018) during the whole period. This interval reflects the NTL impact on PAs from a spatial distance, and those PAs having closer buffers with lit pixels are easier to be affected by the light pollution. Differently, we also defined the interval using the intensity of NTL (Equation (1)), in which both the DN values of lit pixels and their distances to PAs were considered (Figure 3b). Similarly, we calculated the minimum (e.g., the 50th buffer in 1992) and maximum (e.g., the 20th buffer in 2018) buffers from 1992 to 2018 and defined their interval as the high-intensity interval (e.g., the 20th to 50th buffer in Figure 3b).
where I i is the intensity of i th buffer, B i is the i th buffer, and NTL sum is the sum of NTL in the buffer ring i (i.e., B i − B i−1 ). pixels and their distances to PAs were considered (Figure 3b). Similarly, we calculated the minimum (e.g., the 50th buffer in 1992) and maximum (e.g., the 20th buffer in 2018) buffers from 1992 to 2018 and defined their interval as the high-intensity interval (e.g., the 20th to 50th buffer in Figure 3b). (1) where Ii is the intensity of i th buffer, Bi is the i th buffer, and NTLsum is the sum of NTL in the buffer ring i (i.e., Bi − Bi-1).

Figure 3.
Conceptual illustration of buffer intervals around PA revealed by the first polluted buffer (a) and the high-intensity buffer (b). Note: yellow and purple were used to indicate different pollution intervals (i.e., the first pollution and high-intensity interval).
We analyzed the temporal trend of NTL in different regions for continuously and discontinuously polluted PAs. Due to urban sprawl and human activities, PAs are likely to be affected by NTL over the past decades, i.e., the first polluted buffer becomes close to PAs with increasing intensity of NTL over the years. Therefore, for continuously polluted PAs, we analyzed the temporal trend of NTL inside the PA and the high-intensity interval because the first polluted buffer was merged with PA since 1992. However, for discontinuously polluted PAs, we investigated the temporal trend of NTL in the first polluted interval and high-intensity interval (a detailed case of the high-intensity interval can be referred to Seattle in Figure S3).

Spatially Explicit Distribution of Light Pollution Categories
Those PAs in Europe, the USA, and East Asia have had more severe light pollution than other regions over the past decades ( Figure 4). Globally, numbers of PAs defined as non-polluted, continuously polluted, and discontinuously polluted are 24,650, 8141, and 16,509, respectively (Figure 4a), suggesting that almost 50% of global PAs have been affected by NTL. To better visualize their distributions, we generated kernel density maps using the centroid of PA with different light pollution categories and then composited them at different channels ( Figure 4b). Most light pollution was observed in regions with latitudes from 30°N to 60°N, exhibiting three critical regions of Europe, the USA, and East Asia. The continuously polluted PAs (red color in Figure 4b) are widely distributed in Europe, suggesting these PAs have been affected by NTLs since the 1990s. This is largely We analyzed the temporal trend of NTL in different regions for continuously and discontinuously polluted PAs. Due to urban sprawl and human activities, PAs are likely to be affected by NTL over the past decades, i.e., the first polluted buffer becomes close to PAs with increasing intensity of NTL over the years. Therefore, for continuously polluted PAs, we analyzed the temporal trend of NTL inside the PA and the high-intensity interval because the first polluted buffer was merged with PA since 1992. However, for discontinuously polluted PAs, we investigated the temporal trend of NTL in the first polluted interval and high-intensity interval (a detailed case of the high-intensity interval can be referred to Seattle in Figure S3).

Spatially Explicit Distribution of Light Pollution Categories
Those PAs in Europe, the USA, and East Asia have had more severe light pollution than other regions over the past decades ( Figure 4). Globally, numbers of PAs defined as non-polluted, continuously polluted, and discontinuously polluted are 24,650, 8141, and 16,509, respectively (Figure 4a), suggesting that almost 50% of global PAs have been affected by NTL. To better visualize their distributions, we generated kernel density maps using the centroid of PA with different light pollution categories and then composited them at different channels ( Figure 4b). Most light pollution was observed in regions with latitudes from 30 • N to 60 • N, exhibiting three critical regions of Europe, the USA, and East Asia. The continuously polluted PAs (red color in Figure 4b) are widely distributed in Europe, suggesting these PAs have been affected by NTLs since the 1990s. This is largely due to the relatively close distance between PAs and cities in Europe. The numbers of continuously and discontinuously polluted PAs are 4096 and 4601, respectively, with a total area of 1.69 million km 2 . There is a distinct separation of continuously polluted (i.e., east) and discontinuously polluted (i.e., west) PAs in the USA. The total area of polluted PAs in the USA is 0.95 million km 2 , and the areas of 1674 continuously and 4028 discontinuously polluted areas are 0.3 and 0.65 million km 2 , respectively.
The widely distributed discontinuously polluted PAs indicate that urbanization with increased NTL use has posed noticeable challenges to the protection of PAs. For example, there are many global urban clusters greater than 100 km 2 in the USA [54]. The category of continuously polluted is the dominant type of light pollution in East Asia, particularly in South Korea and Japan. It is worth noting that the WDPA contains limited PAs in China, making it challenging to evaluate China's light pollution. Moreover, five hotspots deserve attention worldwide, including Caracas (the capital of Venezuela) and Rio de Janeiro (Brazil) in South America, West Africa, South Asia, and Australia ( Figure 4b). The urbanization of South America began in the early years and the deforestation mainly occurred in Amazon [55], which caused many continuously polluted PAs during the whole period. Besides, with rapid urbanization [56], the PAs with non-polluted were gradually translated to discontinuously polluted, and most of them are located in West Africa, South Asia, and Australia. The widely distributed discontinuously polluted PAs indicate that urbanization with increased NTL use has posed noticeable challenges to the protection of PAs. For example, there are many global urban clusters greater than 100 km 2 in the USA [54]. The category of continuously polluted is the dominant type of light pollution in East Asia, particularly in South Korea and Japan. It is worth noting that the WDPA contains limited PAs in China, making it challenging to evaluate China's light pollution. Moreover, five hotspots deserve attention worldwide, including Caracas (the capital of Venezuela) and Rio de Janeiro (Brazil) in South America, West Africa, South Asia, and Australia (Figure 4b). The urbanization of South America began in the early years and the deforestation mainly occurred in Amazon [55], which caused many continuously polluted PAs during the whole period. Besides, with rapid urbanization [56], the PAs with non-polluted were gradually translated to discontinuously polluted, and most of them are located in West Africa, South Asia, and Australia.

Temporal Trends of NTL in Different Light Pollution Categories
There is an increasing temporal trend in polluted PA worldwide over the past decades, especially in Europe and South America ( Europe account for 51% of the world, dominated by increasing NTL (Figure 5a). Among these PAs, 24% of them show a significant increased trend with p-values below 0.05. Most of the discontinuously polluted PAs with a significant increase are distributed in South America, West Africa, India, and South Asia (Figure 5b) (Table 1). For example, the proportion of PAs with a significant increase in South America is larger than 50%, suggesting a notable growth of NTL luminance. The difference between Europe and South America in different pollution categories indicates the rapid growth of city light in South America's developing areas due to urban sprawl. At the continental scale, the polluted PAs (i.e., continuously and discontinuously polluted) show a distinct spatial heterogeneity ( Figure 5), especially in South America, West Africa, South Asia, and Oceania. The dominant trend of NTL change in South America is increasing in continuously and discontinuously polluted PAs, suggesting a considerable intervention of human activities related to the NTL luminance (e.g., urban sprawl). The discontinuously polluted PAs in Africa show an increasing trend (accounting for 50%), especially in West Africa.

The distance of Light Pollution to the Protected Areas
Globally, most polluted PAs are close (i.e., less than 10 km) to the first polluted buffer and the high-intensity buffer (Table 2). Here, the averaged first polluted buffer and the high-intensity buffer were defined as their averages in different years. Then, we divided these PAs into three categories according to the distribution of the averaged first polluted buffer and the high-intensity buffer of all polluted PAs. We found most PAs with the av-  Global polluted PAs with a decreasing temporal trend are mainly distributed in the USA and Japan, and such a pattern is consistent in both the continuously and discontinuously polluted categories ( Figure 5). The proportion of polluted PAs with a decreasing trend in Japan is about 85%, of which 35% of PAs are statistically significant. The proportion of discontinuously polluted PAs with a decreasing trend reaches 90% in Japan, and around 50% of these PAs pass the significance test. This phenomenon is mainly attributed to the well-planned conservation policies of light pollution in these two regions. For example, Japan is the first country to develop guidelines for mitigating the effect of light pollution [57]. The continuously polluted PAs with decreasing trends are mainly (~65%) distributed in the Great Lakes region of the USA. The discontinuously polluted PAs in the USA are also dominated (75%) by a decreasing trend. This phenomenon is related to the phase of urbanization and population loss; e.g., there is a noticeable outflow of the population in the Great Lakes region [58].

The Distance of Light Pollution to the Protected Areas
Globally, most polluted PAs are close (i.e., less than 10 km) to the first polluted buffer and the high-intensity buffer (Table 2). Here, the averaged first polluted buffer and the high-intensity buffer were defined as their averages in different years. Then, we divided these PAs into three categories according to the distribution of the averaged first polluted buffer and the high-intensity buffer of all polluted PAs. We found most PAs with the averaged first polluted buffer of less than 10 km are in Europe, the USA, South America, and Asia (Figure 6a). In the USA, there is a noticeable spatial heterogeneity of the distance of light pollution to the PAs; i.e., PAs in the east of the USA are more likely to be affected due to the widely distributed cities. According to the global urban boundaries [54], cities with over 100 km 2 in the USA are mainly distributed in the east. The averaged first polluted buffer of PAs in East Asia (e.g., Japan) is closer (less than 10 km) than those in South Asia (e.g., India) (10 km~25 km). While in South America, the first polluted buffer of most PAs in the Amazon rainforest is above 25 km. The PAs that are far from NTL (beyond 25 km) are mainly distributed in Australia, New Zealand, and the west USA. Additionally, PA types need attention due to their functions ( Table 2). For example, PA types of IV (habitat/species management area) and V (protected landscape/seascape) account for 64% of all PAs, with the averaged first polluted buffer less than 10 km. The light pollution is limited to PAs with Ia (strict nature reserve) due to its high priority of conservation. Those PAs close (less than 10 km) to the averaged high-intensity polluted buffer are mainly in Europe and Brazil (Figure 6b), accounting for 38% of the global polluted PA. Compared to the first polluted buffer, the high-intensity interval represents more substantial impacts of human activity on PAs. Similarly, most (76%) of these PAs belong to the IV (habitat/species management area) and V (protected landscape/seascape), suggesting these types are more easily affected by human activities.

The Temporal Trends of NTL in High-Intensity Intervals
Hotspots of global PAs can be revealed from their spatial maps and temporal trends of light pollution in high-intensity intervals, while showing that the continuously polluted PAs are highly correlated with urbanization ( Figure 7). Spatially, PAs with a high value of NTL in the high-intensity polluted intervals from 1992 to 2018 are mainly distributed in the USA, Europe, East Asia, Venezuela, and Brazil (Figure 7a). However, when regarding their temporal trends of NTL over the past decades, the USA and Japan have decreas- Figure 6. The distance of PAs to the averaged first (a) and high-intensity (b) polluted buffers. The visualization of maps is composited by channels of distance to buffers (buffer <= 10, R; 10 < buffer < 25, G; buffer >= 25, B). Those red points with distances beyond 40 km are mainly caused by unstable light sources (e.g., ports, logging, and shipping).

The Temporal Trends of NTL in High-Intensity Intervals
Hotspots of global PAs can be revealed from their spatial maps and temporal trends of light pollution in high-intensity intervals, while showing that the continuously polluted PAs are highly correlated with urbanization ( Figure 7). Spatially, PAs with a high value of NTL in the high-intensity polluted intervals from 1992 to 2018 are mainly distributed in the USA, Europe, East Asia, Venezuela, and Brazil (Figure 7a). However, when regarding their temporal trends of NTL over the past decades, the USA and Japan have decreasing NTL (Figure 7b). Meanwhile, we found that the temporal trends of light pollution are increasing in the Middle East, Indian, and South Asia. In the USA, the temporal trends of NTL in polluted PAs translated from decreasing in the north to increasing in the south (Figure 7c). The main polluted category in Europe is continuously polluted (Figure 7d). The temporal trend of NTL is increasing, except northern Europe, where the population has been slightly decreasing over the past decades [59]. In Japan, the continuously polluted PAs are close to the urban, whereas the time trend decreases (Figure 7e). This phenomenon mainly contributes to well-planned ecological conservation policies; e.g., Japan is the first country to develop guidelines for light pollution [57]. The trend of light pollution in South America is increasing, especially in Venezuela ( Figure 7f) and Brazil (Figure 7g). What these places have in common is that continuously polluted PAs are near urban areas larger than 100 km 2 .
Remote Sens. 2021, 13, x FOR PEER REVIEW 11 of 17 is the first country to develop guidelines for light pollution [57]. The trend of light pollution in South America is increasing, especially in Venezuela ( Figure 7f) and Brazil ( Figure  7g). What these places have in common is that continuously polluted PAs are near urban areas larger than 100 km 2 . can be found in enlarged views. Note: the mean value of NTL and temporal trend of annual NTL sum in the high-intensity interval were visualized by the kernel density approach to identify these five global hotspots.

The Influence of Policies on Light Pollution
The dominant decreasing trend of NTL in PAs in Japan and the USA is closely related to protection policies and organizations. The first governmental guideline about light pollution was published in March 1998 in Japan [57]. The Ministry of Environment in Japan has published light pollution control guidelines to help local authorities take efficient control measures and enforce special environmental regulations. For instance, the Okayama Figure 7. Distribution of global polluted PAs with high value of NTL in the high-intensity interval (a) and their temporal trends of NTL over the past decades (b). Detailed cases in (c-g) in (a) can be found in enlarged views. Note: the mean value of NTL and temporal trend of annual NTL sum in the high-intensity interval were visualized by the kernel density approach to identify these five global hotspots.

The Influence of Policies on Light Pollution
The dominant decreasing trend of NTL in PAs in Japan and the USA is closely related to protection policies and organizations. The first governmental guideline about light pollution was published in March 1998 in Japan [57]. The Ministry of Environment in Japan has published light pollution control guidelines to help local authorities take efficient control measures and enforce special environmental regulations. For instance, the Okayama prefecture enacted the first regulation in Japan to deal with light pollution to protect spacebased radio astronomy observations in 1989. Moreover, in the USA, some organizations have been founded to mitigate the negative impacts of light pollution [60], such as the New England Light Pollution Advisory Group (NELPAG), the International Dark-Sky Association (IDA), and the Dark Skies Advisory Group (DSAG) [61]. With years of hard work, the status of light pollution in Japan and the USA has been mitigated considerably, which is helpful to realize the sustainable development of ecology and environment in/around urban domains.

The Ecological Impact of Light Pollution
Many studies have demonstrated that light pollution could impact wildlife on foraging, reproduction, migration, and communication in natural systems. Artificial light can extend the period of foraging behavior of animals, which could introduce disturbances in the ecosystem and communities [1]. Additionally, artificial light may confuse organisms accustomed to navigating in a dark environment [62]. The sea turtles hatched from nests would be disorientated by nighttime light, which would affect their reproduction [63]. Birds can be disoriented and collided with others by lights at night [64]. Many insects, such as moths, are likely to be attracted by nighttime light [18].
Compared to continuously polluted PAs, the short-term nighttime lights revealed from the discontinuously polluted PA (e.g., near Zhangjiajie in Figure S2) have a greater impact on the ecosystem, e.g., the wild animals need to adapt to such emerging nighttime lights. This is also an important reason for humans to prevent light pollution from spreading to PA. To identify the sprawl of light pollution in the future, we provided early warning of global light pollution for biodiversity and wildlife by mapping the light pollution source (i.e., the first polluted buffer and the high-intensity buffer) and their distances to PAs. These new insights into light pollution are not available in previous studies [14,36].
Light pollution continues to expand on Earth. However, for many species and ecosystems, its impacts have not been examined on a global scale. Human pressures on the natural system, commonly referred to as threats to biodiversity, are from diverse human activities [65], such as built environment, road, pasture, and nighttime light, etc. Especially, the impact of nighttime light on wildlife is complex, which can change the habits of some species (e.g., owl, insects, and migratory birds, etc.). Specially, the International Union for Conservation of Nature (IUCN) Red List will be employed to explore the impact of light pollution on threatened species on a continental scale in the future. This is helpful for biodiversity conservation, i.e., the INCU Red List was widely used to evaluate the vulnerability of PA in the face of climate, human footprint, and agricultural expansion [66,67].
Additionally, blue light is a new emerging source of light pollution due to the increased use of LED lighting [60,68]. Unfortunately, such kinds of lights cannot be well-detected by DMSP/OLI and VIIRS/NDB [69], and its impact was not considered in our study. Some new remote sensing platforms have recently been used to monitor blue light pollution using advanced sensors [69] and new indices [70], which will be considered in future works of global light pollution.

The Relationship between Light Pollution and Urbanization
Only 26% of 24,650 polluted PAs were found to have significantly decreasing trends, whereas 56% of them were found to have increasing trends, and these results are related to the rapid urbanization. We found that the areas of polluted PAs with the increasing temporal trend are relatively large. The continental polluted PAs (9.90 million km 2 ) with increasing trends (56%) account for at least 6.58 million km 2 in extent. Here, we used the global artificial impervious area (GAIA) [71] to study potential implications for the planet's urbanization and biodiversity conservation (Table 3). We found the ratios of PA with the increasing trend are consistent with the increasing rates of impervious surface area (ISA) (Figure 8). The increasing rates of ISA from 1992 to 2018 in Africa and Oceania are approaching 200%, but the ratios of PAs with the increasing trend are significantly different (e.g., 74% in Africa and 45% in Oceania). This indicates a lack of attention to light pollution in PAs during the urbanization process in Africa and South America. Japan handled the conflict between urbanization and conservation of PA well, illustrating the importance of positive protection policies. Moreover, given that the rapid global urbanization in China [72,73] and the overall trends (Figure 8), PAs will be facing a great threat of light pollution in the near future. 3, x FOR PEER REVIEW 13 of 17 area (ISA) (Figure 8). The increasing rates of ISA from 1992 to 2018 in Africa and Oceania are approaching 200%, but the ratios of PAs with the increasing trend are significantly different (e.g., 74% in Africa and 45% in Oceania). This indicates a lack of attention to light pollution in PAs during the urbanization process in Africa and South America. Japan handled the conflict between urbanization and conservation of PA well, illustrating the importance of positive protection policies. Moreover, given that the rapid global urbanization in China [72,73] and the overall trends (Figure 8), PAs will be facing a great threat of light pollution in the near future.

Uncertainty
Here, we discussed the uncertainty of the evaluation of light pollution using the harmonized global NTL dataset [2]. The uncertainty of this dataset mainly lies in two aspects, including the employed data and model. On the one hand, the DMSP and VIIRS NTL, as two different kinds of data, have different observation records (i.e., digital number and radiation), transit time (i.e., 9:00 pm and 1:30 am), and spatial resolution (i.e., 30-arc and 15-arc). On the other hand, a sigmoid function was employed to convert VIIRS radiance data to the DMSP-like data. When the DMSP-like digital number values from VIIRS are low, the uncertainly of the model is relatively large. Therefore, a technical validation was performed, and it was discovered that excluding low luminance regions using the thresh- Figure 8. The relationship between light pollution (i.e., the ratio of polluted PA with an increasing trend) and urbanization (i.e., the increasing rate of impervious surface areas).

Uncertainty
Here, we discussed the uncertainty of the evaluation of light pollution using the harmonized global NTL dataset [2]. The uncertainty of this dataset mainly lies in two aspects, including the employed data and model. On the one hand, the DMSP and VIIRS NTL, as two different kinds of data, have different observation records (i.e., digital number and radiation), transit time (i.e., 9:00 pm and 1:30 am), and spatial resolution (i.e., 30-arc and 15-arc). On the other hand, a sigmoid function was employed to convert VIIRS radiance data to the DMSP-like data. When the DMSP-like digital number values from VIIRS are low, the uncertainly of the model is relatively large. Therefore, a technical validation was performed, and it was discovered that excluding low luminance regions using the threshold of 7, the derived NTL result from VIIRS is closer to the DMSP data [2]. Some nighttime light with DN values less than 7 in small towns and settlements might be neglected in our study. Overall, the approach in this dataset outperforms other approaches in terms of temporal consistency. For example, we compared the results with other studies on global light pollution [23]; a similar temporal trend distribution appeared (e.g., the east of the USA and Japan have a clear decline trend in NTL from the DMSP).
There are additional limitations between satellite-based global light pollution evaluation and ground light pollution measurement. First, the outdoor lighting facilities and technologies rapidly change in the spectrum, direction, and total flux of lamps [60]. As a result, the night sky will change, and the traditional NTL observations such as DMSP/OLI and VIIRS/DNB may fail to detect new light sources such as blue light. Second, due to the viewing angle, polarization, and wavelength, the ground light pollution cannot be directly observed by satellites [38,69]. Future satellite mega-constellations are expected to significantly increase light pollution [74], which brings great uncertainly to future research on light pollution.

Conclusions
In this study, we developed a framework to evaluate the light pollution in global PAs using WDPA and globally harmonized NTL data over the past decades (1992 to 2018). We identified three categories of light pollution (i.e., 5974 non-polluted, 8141 continuously polluted, and 16,509 discontinuously polluted) in PA and investigated the global distribution of these PAs. Then, we analyzed the temporal trend of NTL in polluted (i.e., the continuously and discontinuously polluted) PAs. Specifically, NTL trends in the first polluted interval and high-intensity intervals were analyzed with consideration of protection policies and urbanization.
Those PAs affected by NTL are mainly distributed in regions with latitudes from 30 • N to 60 • N, representing three key regions (i.e., Europe, the USA, and East Asia). Overall, because of global urbanization, the temporal trend of NTL in global PAs is increasing, except for the USA and Japan, showing a slightly decreasing trend. This is probably related to socioeconomic and demographic factors. Most PAs are close to lit pixels when regarding their distances to the first or the high-intensity polluted buffers (i.e., less than 10 km), especially for the PA types of habitat/species management area (IV) and protected landscape/seascape (V). These PAs are more likely to be affected by light pollution than other regions under future global urbanization scenarios.
Characterizing the spatial relationship between PAs and NTL and the trend of NTL intensity in/around PAs is of high importance for biodiversity conservation. Our study provides a global perspective on light pollution in worldwide PAs and identified different types and hotspots of PAs using temporally consistent NTL time series data. Global urbanization is the main driver raising the risk of light pollution in/around PAs, of which their distance to NTL has become close over the past decades. As a result, more attention is required to mitigate these risks to support sustainable development. It is also worthy to note that the number of PAs in China is not so many in WDPA, although China has experienced prevalent urban sprawl over the past decades [72].
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/rs13091849/s1. Figure S1: An illustration of continuously polluted category. The case area in Paris (France) (a) and change of nightlights (b). Figure. S2: An illustration of the discontinuously polluted category in Zhangjiajie (China) (a) with NTL dynamics (b) and temporal trend analysis (c). Figure.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.