Data and study area<\/p>\n
The foundational dataset for this global analysis of ALAN dynamics was the Black Marble product suite of NASA (Collection 1), derived from the DNB sensor onboard the VIIRS instrument3<\/a>. The DNB observes light in the wavelength range of 500\u2013900\u2009nm, with an equatorial local overpass time of roughly 1:30\u2009a.m. Standardized quality assurance of the Black Marble data ensures consistency across time, geography and sensors (Suomi-NPP and the NOAA-20\/NOAA-21VIIRS DNB), enabling cross-mission compatibility vital for continuous monitoring applications51<\/a>. We acquired two specific daily products (that is, VNP46A2 and VNP46A1) from the NASA Level-1 and Atmosphere Archive and Distribution System Distributed Active Archive Center (LAADS DAAC), covering the period from 1 January 2013 to 31 December 2023.<\/p>\n The primary daily NTL data for detecting changes was the VNP46A2 product, which provides daily ALAN radiance values (nW\u2009cm\u22122\u2009sr\u22121) corrected for atmospheric influences and bidirectional reflectance distribution function (BRDF) effects from lunar illumination geometry and diverse surface reflectance variability3<\/a>. This rigorous design improves radiometric stability and distinguishes Black Marble from earlier NTL products, supporting quantitative, science-quality analysis. Data are provided in geographic coordinates (WGS84) at a nominal spatial resolution of 15 arc-seconds (about 460\u2009m at the equator), finer than the intrinsic sensor resolution of about 750\u2009m. Complementing this, we used the VNP46A1 products to provide the pixel-level quality assessment flags of cloud, snow\/ice and solar\/lunar contamination information21<\/a>. These products also supplied detailed viewing geometry, specifically the sensor viewing zenith angle (VZA), a key parameter for the VZA-stratified COntinuous monitoring of Land Disturbance (VZA-COLD) change detection algorithm used to mitigate angular effects on observed NTL4<\/a>.<\/p>\n To ensure high-quality inputs for the change detection algorithm and improve processing efficiency, two pre-processing steps were conducted: systematic filtering of low-quality daily images and masking of persistently dark areas with no historical artificial light (Supplementary Information Section\u00a01<\/a>). The ALAN change detection results were generated under the linear latitude\/longitude geographic projection, consistent with the input atmospheric- and lunar-BRDF-corrected Black Marble data. This projection was also retained for visualizations of the global maps. However, all quantitative analyses in this study, such as area estimates and accuracy assessments, were conducted following the MODIS\/VIIRS 500-m nominal resolution sinusoidal equal-area projection to ensure accurate area-based calculations.<\/p>\n The geographic scope of this study was defined as the global terrestrial regions located between 70\u00b0\u2009N and 60\u00b0\u2009S. This latitudinal range was chosen because it encompasses the vast majority of the landmasses of Earth and virtually all primary human settlements and areas of substantial ALAN52<\/a>. Polar regions were excluded from the analysis because of data acquisition challenges (for example, polar day, extensive snow and ice). Oceans and large inland water bodies were masked using the Terra Moderate Resolution Imaging Spectroradiometer (MODIS) Land Water Mask product in 2014 (MOD44W Collection 6.1)53<\/a>. The World Bank Official Boundary data (accessed 1 January 2025) were used for regional aggregation (continents,\u00a0countries, and territories). The core period for ALAN change detection was 1 January 2014 to 31 December 2022. Data from 2013 served for model initialization, and data from 2023 for confirming end-of-series changes (see details in Supplementary Information\u00a0Section\u00a01<\/a>).<\/p>\n Definition of ALAN change types<\/p>\n In this study, ALAN changes are categorized into two primary types based on their temporal patterns: abrupt and gradual (Supplementary Fig. 4c<\/a>). These are further classified as either brightening or dimming, depending on the direction of change in radiance. Abrupt ALAN changes are short-term shifts characterized by sudden step-like changes in NTL radiance or a structural break in the time series. These changes typically unfold over a span of weeks to months and often correspond to discrete events such as urban construction or demolition (Fig. 1c<\/a>), natural disasters (Fig. 1h<\/a>), or armed conflicts (Fig. 1d<\/a>). Abrupt changes also include changes that caused sudden redirections of a longer-term trend, such as the onset of economic recession or a surge in foreign investment or immediate policy actions (for example, beginning of rapid lighting installation), leading to a noticeable inflection in the trajectory of ALAN dynamics (Fig. 1j<\/a>). By contrast, gradual ALAN changes represent long-term, continuous trends that unfold over periods exceeding 1\u2009year. These changes exhibit a relatively stable, directional pattern, either brightening or dimming, without abrupt discontinuities. Gradual changes typically reflect sustained socioeconomic or demographic processes such as rural expansion (Fig. 1i<\/a>), economic transformation or the systematic rollout of new lighting technologies (for example, LED retrofitting). Unlike abrupt events, they result in a smooth and persistent evolution in night-time radiance over time.<\/p>\n Each ALAN change is further classified by its direction (that is, brightening and dimming) based on Extended Data Table 5<\/a> (equations (4) and (6)). Brightening corresponds to a\u00a0positive model-estimated change in radiance, whereas dimming corresponds to a negative change. This directional classification is essential for disentangling the complex global patterns of illumination gain and loss, enabling a more complete understanding of the dynamic behaviour of ALAN across space and time.<\/p>\n Although these persistent changes are the focus of our analysis, transient fluctuations that dissipate within days to a few weeks, such as those caused by temporary outages, meteorological anomalies or daily variations in power supply, are not analysed. These ephemeral changes, which typically return to the original status of the NTL intensity within 1\u2009month, fall outside the analytical scope of this study. Our emphasis is on sustained alterations in ALAN radiance rather than day-to-day variability in the signal. The workflow of our method and analysis is shown in Extended Data Fig. 5<\/a>.<\/p>\n ALAN change metrics calculation<\/p>\n The core methodology for detecting changes in ALAN was based on the VZA-COLD algorithm, detailed in ref.\u20094<\/a> and Supplementary Information Section\u00a02<\/a>. VZA-COLD mitigates viewing geometry effects by stratifying daily NTL observations for each pixel into four VZA intervals (0\u00b0\u201320\u00b0, 20\u00b0\u201340\u00b0, 40\u00b0\u201360\u00b0, 0\u00b0\u201360\u00b0). For each stratum, a harmonic time series model is continuously fitted using robust regression to capture intra-annual seasonality and inter-annual trends in ALAN radiance5<\/a>. Abrupt ALAN changes are confirmed if 14 consecutive observations are anomalous (residual\u2009>\u200975% change probability) in any VZA stratum. This ensures that only sustained and statistically robust deviations from the local baseline are identified. This design is well-suited to capturing sudden human-driven ALAN change events, which typically manifest as persistent and significant radiance deviations. Gradual ALAN changes are identified in segments between these breaks if the linear trend coefficient (bi in equation (3) in Supplementary Information Section\u00a02<\/a>) of the fitted 0\u00b0\u201360\u00b0 VZA model is statistically significantly different from 0 (P\u2009<\u20090.05).<\/p>\n To reduce false positives arising from small fluctuations or noise, a post-processing step filters all detected changes to have a minimum magnitude threshold (1.0\u2009nW\u2009cm\u22122\u2009sr\u22121). This conservative threshold excludes minor light variations, such as day-to-day fluctuations in natural airglow. For abrupt changes, detections with change magnitude (measured in DNB radiance) below this threshold were filtered. For gradual changes, only those segments with an absolute change in DNB radiance above this threshold over the duration of the segment were considered. This dual criterion ensures that identified trends are both statistically significant and physically meaningful. Several key adaptations, including a dynamic harmonic period for high-latitude areas to handle polar day and snow impacts and a skipping-update strategy to reduce computation, were implemented to improve model stability and reduce commission errors caused by polar-light effects, enabling robust global-scale processing. An incremental online processing framework using a 4-year moving window was adopted for continuous updates and memory management (for details, see Supplementary Information Section\u00a02<\/a>).<\/p>\n To quantitatively characterize detected ALAN changes, we derived a suite of pixel-based metrics, including \u2018change time\u2019, \u2018change area\u2019, \u2018change intensity\u2019, and \u2018radiance change\u2019 for abrupt, gradual and total ALAN changes, based on outputs from the adapted VZA-COLD algorithm (Extended Data Table 5<\/a>). For all subsequent spatial analyses and to ensure accurate area calculations, all derived ALAN change metric maps were reprojected from their native geographic coordinates (WGS84) to a 500-m nominal resolution sinusoidal equal-area projection. The nearest-neighbour resampling method was used for this reprojection to preserve the original pixel values, particularly for categorical change information and derived intensities, while ensuring that area calculations were not distorted by variations in pixel size across different latitudes.<\/p>\n For each detected abrupt change, the change time was recorded as the date when the breakpoint first appeared in the fitting model. This timestamp marks the transition point between the before and after segments of the harmonic fit, often reflecting a real-world transition in lighting conditions. The change area for all pixel-based metrics was standardized to 0.2147\u2009km2, corresponding to the actual spatial footprint of the\u00a0sinusoidal equal-area grid cell, enabling consistent area-based aggregation over the globe. The change intensity, described Extended Data Table 5<\/a> (equation (4)), quantifies the magnitude of an abrupt change. This was calculated to represent the net, long-term impact after the abrupt event on the stable NTL level. It was defined as the difference between the predicted overall ALAN radiance (derived from the de-seasonalized ai\u2009+\u2009bix terms of the 0\u00b0\u201360\u00b0 VZA model) of the change pixel at the start date of the time segment immediately following the detected breakpoint, and the predicted overall ALAN radiance at the end date of the time segment immediately preceding the breakpoint (Supplementary Fig. 4c<\/a>). A positive value for this metric indicates a step brightening in the stable light level, whereas a negative value signifies a step dimming. To quantify the total effect of an abrupt change in a spatially meaningful way, we derived the radiance change by multiplying the change intensity by the area of the pixel in Extended Data Table 5<\/a> (equation (5)). This metric expresses the sum of radiance change in ALAN and is important for scaling pixel-level transitions to broader spatial summaries.<\/p>\n For time segments identified as having gradual ALAN change (Supplementary Fig. 4c<\/a>), the change time was marked as the dates from when the gradual change was first observed to the end date of the gradual time segment. The change area, again, was fixed at 0.2147\u2009km2. The change intensity was calculated on an annual basis, described in Extended Data Table 5<\/a> (equation (6)). This metric represents the area-averaged total change in NTL radiance attributable to the sustained gradual trend over that specific calendar year. It was calculated using the slope term from the 0\u00b0\u201360\u00b0 VZA interval model for that particular gradual change segment (denoted as MG in equation (6) in Extended Data Table 5<\/a>) multiplied by the number of days within that calendar year that belonged to that gradual change segment. This provides a measure of the intensity and direction of the gradual trend (for example, the speed of ongoing urbanization or the rate of decline in a depopulating area). The corresponding radiance change was obtained by multiplying the annual change intensity by the fixed pixel area (Extended Data Table 5<\/a>, equation (7)), enabling consistent comparison and aggregation of gradual radiance changes across large regions or over multi-year periods.<\/p>\n To capture the cumulative impact of all detected ALAN dynamics at the pixel level, we define a set of metrics for the total change (the sum of all types of changes) as the overall alteration at a pixel resulting from the combined effects of both abrupt and gradual ALAN changes. The total change intensity metric, described in Extended Data Table 5<\/a> (equation (8)), represents the overall ALAN change intensity experienced by a pixel. It was calculated by summing the differences in radiance caused by all detected abrupt changes and all annual gradual changes occurring within a defined period (for example, for each year for time series analysis, or cumulatively over the entire 2014\u20132022 study period for total change mapping). This provides an integrated measure of the net impact of all ALAN dynamics on a pixel during the study period. The total radiance change metric, as described in Extended Data Table 5<\/a> (equation (9)), serves to quantify the total radiance impact of all types of ALAN changes across space.<\/p>\n Accuracy assessment and unbiased area estimation<\/p>\n Accuracy assessment and unbiased area estimates followed established protocols for land change mapping54<\/a>,55<\/a>. An independent validation dataset was generated using a stratified random sampling design to ensure statistically robust estimations of accuracy by appropriately weighting the estimators to account for the large differences in mapped areas across classes. The sampling frame consisted of all terrestrial pixel-years (a 500-m nominal resolution sinusoidal equal-area pixel in a specific calendar year) from 2014 to 2022, totalling more than 635 million potential units. Two separate sampling schemes were implemented to target abrupt and gradual changes, respectively, as both can occur within the same year. For the abrupt change assessment, strata were defined annually as abrupt change detected or non-abrupt change (including stable and gradual change pixels). For the gradual change assessment, strata were defined as gradual change detected or non-gradual change (including stable and abrupt change pixels). Sample sizes were allocated based on mapped area proportions to achieve the accuracies of the target user of about 70% for the smaller change strata and around 95% for the larger non-change strata, with a target standard deviation of 0.005 for overall accuracy and a minimum of 200 samples per change stratum. This resulted in the selection of 2,071 sample units for the abrupt change assessment and 1,902 sample units for the gradual change assessment (Supplementary Fig. 5<\/a>).<\/p>\n Reference data for each validation sample were generated through careful visual interpretation by trained analysts using a custom web-based application developed in Google Earth Engine. This tool facilitated consistent interpretation by showing the full VZA-stratified DNB NTL time series, corresponding high-resolution optical imagery (from Google Earth and PlanetScope) and original daily NTL image chips for each sample unit. Interpreters assigned a reference classification label (that is, Abrupt Change, Non-abrupt Change, Gradual Change, and Non-gradual Change) based on a holistic assessment of all available evidence for that specific pixel-year. During this process, analysts explicitly identified discrepancies between the map and reference labels, recording specific instances of commission errors (false positives) and omission errors (false negatives) to ensure robust accuracy estimation. For samples in which the mapped change agrees with the reference data, interpreters also qualitatively noted the potential direct causal drivers defined in Extended Data Table 4<\/a>, based on the information from the remote sensing images and other open-access resources, such as the VIIRS Nightfire gas flaring data from SkyTruth56<\/a>, conflict data from UCDP (Uppsala Conflict Data Program), earthquake data from CrisisWatch, hurricane records from NOAA (National Oceanic and Atmospheric Administration), social media (for example, X and Weibo), and relevant news and reports. To ensure objectivity in this driver attribution, we implemented a strict quality control protocol: each sample was independently interpreted by two trained analysts. Any discrepancies between their driver labels were flagged and adjudicated by a third, senior analyst to produce the final consensus label.<\/p>\n Based on a comparison of map labels to reference labels, confusion matrices were constructed separately for the abrupt and gradual change maps (Extended Data Table 1<\/a>). From these matrices, unbiased estimators of overall accuracy, user accuracy (a measure of commission error) and producer accuracy (a measure of omission error) were calculated with 95% CIs, appropriately weighted by the area proportions of each stratum to avoid bias from large stable areas54<\/a>,57<\/a>. To obtain statistically robust and unbiased estimates of the cumulative area experiencing different types of ALAN change (globally and regionally), we applied area-weighted estimators based on the reference samples54<\/a>. This standard method leverages the stratified validation sample to correct map-based pixel counts for classification errors, providing unbiased area estimates with associated 95% CIs (Extended Data Table 3<\/a>). Note that apart from these unbiased cumulative area estimates, all other area-related spatial analyses in this study rely on direct mapped pixel counts.<\/p>\n Spatial and temporal analysis<\/p>\n For analysing broad spatial patterns and for visualization purposes, the 15-arc-second pixel-level ALAN change metrics were aggregated into regular linear latitude\/longitude grid cells. Two grid resolutions were primarily used: 0.5\u00b0 and 2\u00b0. Within each of these grid cells, summary statistics such as the sum of change area, the area-averaged change intensity and the sum of radiance change were calculated (see definitions in Extended Data Table 5<\/a>). This gridded aggregation provided a uniform spatial framework that facilitated the visual analysis of broad spatial patterns across the globe (0.5\u00b0) and the calculation of zonal statistics (with 2\u00b0), such as those used for creating the latitudinal and longitudinal profile plots (Fig. 3<\/a> and Extended Data Fig. 2<\/a>). Moreover, for analyses tied to geopolitical or socioeconomic contexts, the equal-area 500-m nominal resolution sinusoidal pixel-level results were aggregated based on administrative boundaries (continents,\u00a0countries, and territories) obtained from the World Bank dataset. Total change area, area-averaged change intensity and total radiance change were calculated for each administrative unit.<\/p>\n To evaluate how the extent and intensity of ALAN change evolved over the 2014\u20132022 study period, we performed temporal trend analysis on the annually aggregated statistics. This was done at global and, where data permitted, country and territory levels. The non-parametric Theil\u2013Sen regression method58<\/a>,59<\/a>,60<\/a> was used to estimate the slope (specifically, the median slope) of the trend in annual ALAN change over the 9-year period. Trends were considered to be statistically significant if the P-value from the Mann\u2013Kendall test was less than 0.05. These trend tests were applied separately to metrics for net ALAN change (sum of brightening and dimming), brightening ALAN change and dimming ALAN change (Fig. 4<\/a>).<\/p>\n Although the\u00a0actual change values (for example, change area in km2, change intensity in nW\u2009cm\u22122\u2009sr\u22121) provide direct measures of alteration, understanding the rate or relative change compared with the initial state of illumination is often crucial for contextualizing the impact of ALAN changes, especially when comparing regions with vastly different baseline light levels. Using simple annual composites from the start year (that is, 2014) as a baseline for calculating relative change can be misleading, as these empirical data can be affected by short-term variability, seasonality or even undetected early changes that occurred before the model initialization was complete. Therefore, to establish a more robust and stable baseline representing the initial NTL conditions, we generated synthetic baseline NTL radiance values and area for 1 January 2014 (Extended Data Table 2<\/a>). This was achieved by using the overall trend terms of the harmonic model that was fitted to the 0\u00b0\u201360\u00b0 VZA interval data, primarily based on the 2013 observations used for model initialization (ai\u2009+\u2009bix in Supplementary Information\u00a0Section\u00a02<\/a>, equation (3)). We defined the global lit area of the\u00a02014\u00a0baseline (or initial global lit area) as all regions in which this synthetic baseline radiance exceeded 1.0\u2009nW\u2009cm\u22122\u2009sr\u22121. This calculation effectively removes the influence of seasonality captured by the harmonic terms and represents the de-seasonalized, stable baseline radiance that was predicted by the initial model fit, thus being free from the influence of any ALAN changes detected after the model was initialized. These synthetic baseline radiance values provide a consistent and unbiased representation of initial lighting conditions across all pixels. Relative changes observed over the study period (for example, percentage radiance change in ALAN) were then calculated by comparing the observed cumulative radiance change to these synthetic baseline radiance values (Fig. 3e<\/a>).<\/p>\n","protected":false},"excerpt":{"rendered":"Data and study area The foundational dataset for this global analysis of ALAN dynamics was the Black Marble…\n","protected":false},"author":2,"featured_media":381866,"comment_status":"","ping_status":"","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[23],"tags":[5563,4068,85,132152,46,4069,55527,141,145],"class_list":{"0":"post-381865","1":"post","2":"type-post","3":"status-publish","4":"format-standard","5":"has-post-thumbnail","7":"category-space","8":"tag-environmental-impact","9":"tag-humanities-and-social-sciences","10":"tag-il","11":"tag-interdisciplinary-studies","12":"tag-israel","13":"tag-multidisciplinary","14":"tag-optical-sensors","15":"tag-science","16":"tag-space"},"_links":{"self":[{"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/posts\/381865","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/users\/2"}],"replies":[{"embeddable":true,"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/comments?post=381865"}],"version-history":[{"count":0,"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/posts\/381865\/revisions"}],"wp:featuredmedia":[{"embeddable":true,"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/media\/381866"}],"wp:attachment":[{"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/media?parent=381865"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/categories?post=381865"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/www.newsbeep.com\/il\/wp-json\/wp\/v2\/tags?post=381865"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}