Edited by Jonathan J. Cole, Cary Institute of Ecosystem Studies, Avon, NC, and approved February 12, 2021 (received for review October 13, 2020)
Author contributions: C.K. and D.B. designed research; C.K. performed research; C.K. contributed new reagents/analytic tools; C.K. analyzed data; C.K. wrote the paper; and D.B. provided scientific advising.
- Altmetric
Arctic and boreal regions are undergoing dramatic warming and also possess the world’s highest concentration of lakes. However, ecological changes in lakes are poorly understood. We present a continental-scale trend analysis of satellite lake color in the green wavelengths, which shows declining greenness from 1984 to 2019 in Arctic-boreal lakes across western North America. Annual 30-m Landsat composites indicate lake greenness has decreased by 15%. Our findings show a relationship between lake color, rising air temperatures, and increasing precipitation, supporting the theory that warming may be increasing connectivity between lakes and surrounding landscapes. Overall, our results bring a powerful set of observations in support of the hypothesis that lakes are sentinels for global change in rapidly warming Arctic-boreal ecosystems.
The highest concentration of the world’s lakes are found in Arctic-boreal regions [C. Verpoorter, T. Kutser, D. A. Seekell, L. J. Tranvik, Geophys. Res. Lett. 41, 6396–6402 (2014)], and consequently are undergoing the most rapid warming [J. E. Overland et al., Arctic Report Card (2018)]. However, the ecological response of Arctic-boreal lakes to warming remains highly uncertain. Historical trends in lake color from remote sensing observations can provide insights into changing lake ecology, yet have not been examined at the pan-Arctic scale. Here, we analyze time series of 30-m Landsat growing season composites to quantify trends in lake greenness for >4 × 105 waterbodies in boreal and Arctic western North America. We find lake greenness declined overall by 15% from the first to the last decade of analysis within the 6.3 × 106-km2 study region but with significant spatial variability. Greening declines were more likely to be found in areas also undergoing increases in air temperature and precipitation. These findings support the hypothesis that warming has increased connectivity between lakes and the land surface [A. Bring et al., J. Geophys. Res. Biogeosciences 121, 621–649 (2016)], with implications for lake carbon cycling and energy budgets. Our study provides spatially explicit information linking climate to pan-Arctic lake color changes, a finding that will help target future ecological monitoring in remote yet rapidly changing regions.
Recent widespread changes in Arctic-boreal primary productivity linked to global climate change have been documented in terrestrial (456–7) and oceanic (8, 9) ecosystems. Satellite records have enabled large-scale analysis, providing evidence of uneven changes in vegetation growth in tundra (5, 6, 101112–13) and boreal ecosystems (5, 14, 15), as inferred from vegetation indices. Thus, satellite records have provided critical insight into long-term global change in remote Arctic and boreal regions.
However, less attention has been paid to the large-scale changes in satellite-derived observations of lake color, owing in part to the prevalence of small lakes, which were previously not detectable by earlier generations of satellites (16). Lakes are abundant in permafrost landscapes (1, 17), supporting crucial ecosystem services such as wildlife habitat, biodiversity, nutrient and carbon cycling, and recreation. The lack of pan-Arctic studies of lake color is surprising because half of the world’s total lake surface area is concentrated in Arctic-boreal regions (18), where permafrost thaw has been predicted to dramatically alter lake physical and optical properties (19).
Evidence from field studies has suggested satellite lake color is linked to primary productivity (20) and dissolved organic carbon concentrations, which in turn may be influenced by terrestrial processes (21). For example, shifts in terrestrial vegetation cover are hypothesized to influence lake ecological structure and function (22). As the frequency and timing of spring thaw change, along with total annual precipitation, the magnitude and pathways of water moving within the landscape will impact the routing of terrestrial carbon into aquatic ecosystems. However, the responses of Arctic-boreal lakes to changes in terrestrial ecosystems remains highly uncertain.
Lake color, as observed from space, offers a powerful approach for monitoring at large scales. Satellites such as Landsat observe lake color by quantifying surface reflectance, the fraction of incoming solar radiation emerging from and reflected off a lake’s surface in the visible wavelengths (∼450 to 700 nm). This surface reflectance, referred hereto after as satellite lake color, has been used to infer properties relevant to carbon cycling, including lake clarity (23, 24), phytoplankton stocks (25), and primary production (20, 26). Despite the fact that lake color has been designated an essential climate variable by the Global Climate Observing System, analysis of pan-Arctic trends in satellite lake color remains limited. The majority of satellite lake color studies have been conducted outside Arctic and boreal regions (27), with most high northern latitude work focused on dissolved organic carbon dynamics (282930–31), single regions (32), or physical processes such as changes in inundation extent (33343536–37). This study builds off that work by using 35 y of satellite remote sensing data to systematically investigate changes in lake color, in particular greenness, at the continental scale.
In contrast to satellite remote sensing, limnologists typically measure lake color not as reflectance, but rather as absorbance. Dissolved organic matter contributes color to lakes by absorbing light in shorter wavelengths (∼250 to 440 nm) (38). This absorption can give lakes a “brown” appearance. Decades of field studies in northern lakes have led to the hypothesis that dissolved organic matter concentrations are increasing due to shifts in hydrology and climate (3940–41), resulting in large-scale changes in lake color manifested as increasing light absorption. This increase in absorption, termed browning, has been shown to influence primary productivity (4243–44) because dissolved organic matter absorbs light in the same wavelengths required for photosynthesis (45). Historic Landsat sensors did not measure reliably over water in these shorter wavelengths (46), but they do measure reflectance in green wavelengths (∼560 nm).
Here we explore the hypothesis that satellite lake color, and in particular lake greenness, is changing at the continental scale. We selected reflectance in the green wavelengths, referred to here as greenness, for two reasons. First, green reflectance has ecological significance. Photosynthetic pigments preferentially reflect green light (27). For this reason, satellite greenness has been used for decades to derive primary productivity in the ocean (47, 48) and, more recently, in lakes (26, 4950–51). Second, previous studies have demonstrated the green band is the least influenced by intersensor calibration errors (52, 53), making it a robust choice for time-series analysis. To date, we are aware of no study that has applied this approach to evaluate continental-scale trajectories of pan-Arctic lake greenness, which has been shown to closely track rates of lake primary productivity in diverse northern lake ecosystems (20, 54).
To document trends in satellite lake color, we analyze time series of lake reflectance in green wavelengths (∼560 nm) for over 400,000 lakes in western North America’s Arctic-boreal zone. Using over 54 million observations from Landsat 5, 7, and 8, we generate seasonal composites of lake greenness from 1984 to 2019 for 542,934 lakes larger than 0.1 km2 as a complement to existing studies on growing season greenness in terrestrial ecosystems. At each lake site, we first extracted a time series of lake greenness, defined as growing season surface reflectance (Rs, unitless) in the green wavelengths (∼560 nm), and then calculated the overall trend in lake greenness (∆green Rs, in y−1; Materials and Methods). We discuss our findings in the context of physiochemical, biological, and hydrologic controls. Ultimately, patterns of greenness were compared to historical changes in the landscape and climatic parameters known to influence lake productivity.
Results
Widespread, but Spatiotemporally Variable, Declines in Lake Greenness.
Overall, lake greenness has declined in the past 35 y in North America’s western Arctic and boreal regions by 15% (Fig. 1). Decadal trends and spatial patterns in annual growing season lake color (1984 to 2019) are presented in Fig. 1. Taking all the lakes together, we found a significant negative trend in greenness (Mann–Whitney trend test, Z = −2.88, P = 3.8 × 10−3, N = 472,889), calculated by creating a time series that reflects the average greenness of all lakes for each year (Fig. 1B). Lake greening is rarer, contributing to less than 3% of observed trends, although, in places where greening is observed, the magnitude of the trend is higher (Fig. 2). Greening lakes also had higher green reflectance (Mann–Whitney U test, W = 1 × 108, P < 0.05, n = 469,990).


Spatiotemporal patterns in Arctic and boreal lake color change from 1984 to 2019. (A–D) Changes in lake greenness vary spatiotemporally, and overall lake greenness declined from 1984 to 2019. (A) Map showing the location of study lakes (red) and distribution of global circumpolar permafrost landscapes (dark gray). (B) Time series of median greenness (Rs, unitless) from Harmonized Landsat Collection 1 Surface Reflectance dataset (1984 to 2019) of study lakes (N = 472,889). Lake trend lines were calculated using a robust Theil-Sen estimator. Trends in lake greenness (∆green Rs y−1) vary spatiotemporally. Lakes with significant (P < 0.05) negative (C) and positive (D) trends below 75° N.


Latitudinal profiles of lake counts and color. (A) Distribution of study lakes by latitude and (B) greenness trend for all lakes (lighter colors) and significant lakes (darker colors). Plots were made by binning statistics by every 1° of latitude.
Calculating individual significance per lake, 26% of lakes showed statistically significant (P < 0.05) shifts in color, and spatiotemporal differences are revealed in Fig. 1 C and D. Of the lakes with significant trends, 97% show declining greenness. Lakes within the zone of continuous permafrost showed stronger negative trends (Fig. 3). Significant differences existed between biomes (H = 16,230, P < 0.05, n = 472,889; SI Appendix, Figs. S5 and S6 and Tables S1 and S2). Tundra, temperate grasslands, and shrub biomes showed similar and greater declines in greenness than forested biomes (SI Appendix, Table S1 and S2). Areas of particularly apparent rapid change are located in the Yukon–Kuskokwim delta (Fig. 1C), where studies have linked terrestrial browning to saltkill from coastal flooding (55). Shrubification of previously barren tundra regions indicates that energy and nutrients are present enough to support more complex vegetation, which in turn can alter local energy budgets, impacting the rate of permafrost thaw (56) and eventual transfer to aquatic ecosystems (19). The smaller changes in lake greenness within forested areas of Arctic boreal North America could be explained by net losses of forests and increasing drought (13), which could decrease hydrologic connectivity and terrestrial dissolved organic carbon (DOC) loads available for transfer (21) to aquatic systems.


Greening trends binned by bioclimatic and physical controls on lake color. Box plots showing lake greenness trend binned by ecoregion (Top), lake depth, permafrost extent, and ground ice content (Bottom). Numbers on the left indicate the number of lakes in each class.
Taken together, these results show the greatest changes in lake greenness are concentrated in the shrub and tundra ecosystems underlain by continuous permafrost regions (Figs. 1C and 2). Negative trends in lake color observed here are supported by evidence of increased organic carbon export to lakes (57, 58), potentially due to thermokarst-activated flow paths between aquatic and terrestrial systems (3). Other factors shown to drive increased light absorption include warming and precipitation changes (42, 596061–62), recovery from acidification (40), altered biogeochemical interactions (62), and changes in land use (41, 63). Stronger absorption of light by dissolved materials can impact the amount of light available for photosynthesis (64), reducing productivity and weakening the strength of the Arctic carbon sink (65).
However, changes in lake color are complex and may be driven by a combination of water column properties, bottom substrate, and lake depth (66). Our analysis indicates that shallow lakes are more vulnerable to changes in color. To assess the impact of depth, we compared trends between shallow (<10 m depth) and deeper (>10 m depth) lakes. The majority of lakes in the study are shallow as indicated by a mean depth of 3.5 m. Shallower lakes exhibited 34% greater declines in greenness than deeper lakes (Fig. 3), and the difference was significant (Mann–Whitney U test, W = 560,284,235, P = 1.43 × 10−30, n = 469,990) compared to deeper lakes (>10 m, n = 2,899). It is important to note, however, that lake color derived from satellites is restricted to the near surface (67), so this analysis is restricted to the surface layer for deeper lakes. Here, measures of lake depth have been estimated using a static geostatistical model conditioned on topography (68). Assessing the impact of changes in lake depth on color would require multitemporal lake depth estimates. Regional lake levels have been characterized using AirSWOT, an experimental airborne Ka‐band radar interferometer (69). Most pan-Arctic studies have focused on changes in lake area and surface water extent (34, 35, 38, 70, 71), not on lake depth, precluding our ability to capture the historical change in lake water level or depth. Nonetheless, these results imply that shallower lakes may be undergoing more dramatic transformations.
We implemented an extensive screening method to produce high-quality time series. The main sources of uncertainty in this analysis are discrepancies arising from using multiple generations of Landsat sensors and the correct identification of water pixels using Landsat’s Quality Assessment (QA) Bands. To minimize uncertainty from using multiple Landsat sensors, we selected the green band, which has the lowest calibration errors across sensors (53). Another possible bias is the fewer number of observations in the Landsat-5 archive due to downlink and data quality issues. Sparse Landsat-5 observations in the 1980s and 1990s can bias the interpretation of trends. Our analysis showed that 93% of lakes had >20 y of data (SI Appendix, Fig. S4), and lakes with fewer than 10 y were excluded (Materials and Methods). To assess bias introduced by using composites derived from multiple Landsat sensors, we calculated trends with and without harmonization coefficients (72) and with and without the 1991 to 1993 Pinatubo anomaly and found the overall trend remained negative and significant in both cases (SI Appendix, Fig. S7). A series of conservative filters (SI Appendix, Fig. S1), including Landsat’s QA band and an independent water mask derived from the Global Surface Water Dataset (73), were used to ensure analyzed pixels were not clouds, vegetation, or snow/ice (Materials and Methods). These extra efforts were undertaken to provide evidence that the results were robust and free from potential systematic errors (74).
Lake Greenness Declines Most in Warming and Wetting Regions.
The continental-scale shifts in lake color observed here appear to track regional changes in climatic conditions. To evaluate lake color trends in the context of climate conditions, we analyzed air temperature and precipitation trends using the ERA5 climate reanalysis dataset (75) (Fig. 4 A–F). Consistent with other studies, the past three decades (1984 to 2019) are characterized by an overall rise in western North American Arctic and boreal air temperatures (2) by 0.031 °C (SD ± 0.028 °C) per year (Fig. 4A), which is lower but on the same order of magnitude as estimates made using a process model for the time period of 1982 to 1998 (76).


Climate trends: overall trends in temperature (A, C, and E) and precipitation (B, D, and F) for annual (A and B), spring (C and D), and summer (E and F) time periods, and the proportion of lakes whose basins are in each climate category (G). Box plots of distributions of green trend by basin climate trend category calculated by basin for each lake using the GLCP dataset (Materials and Methods; H).
We then linked each lake to the climate reanalysis dataset by matching each lake to its smallest basin using the Global Lake Area, Climate, and Population (GLCP) dataset (77) (Materials and Methods and SI Appendix). We calculated that 60% (n = 224,422) of lakes are located in wetter basins and 40% (n = 152,603) of lakes are conversely in drier basins. We separated temperature and precipitation trends into four categories: warmer and wetter, warmer and drier, colder and wetter, and colder and drier (Fig. 4G). We found a significant difference in greening trends between these four categories (Kruskal–Wallis test, H = 21,326, P < 0.05, n = 377,025; Fig. 4H). Lakes in warmer and wetter basins had a 2.5 times greater mean decline in green reflectance than lakes in colder and drier basins (Fig. 4H). This result is consistent with the current hypothesis that precipitation and permafrost warming are altering terrestrial-to-aquatic hydrologic connectivity. This impacts downstream delivery of terrestrial solutes (78), including organic carbon, to lakes (21, 79), which, in turn, could negatively influence primary productivity and further induce net heterotrophy (80).
The decline in lake greenness under warmer and wetter conditions is likely caused by an overall increase in terrestrial-to-aquatic hydrologic connectivity from thawing, which can result in higher DOC concentrations in lakes (58). While greenness is typically associated with photosynthetic pigments, absorbing organic matter, called chromophoric dissolved organic matter (CDOM), can also dominate absorption in green wavelengths, directly influencing the greenness (45) (SI Appendix). Higher DOC concentrations exert a two-way influence on primary productivity. At lower concentrations, increased terrestrial inputs can supply nutrients such as phosphorus and nitrogen, thereby stimulating phytoplankton growth (81). However, high concentrations of colored DOC can reduce light availability, shading out producers (44, 82) and weakening the strength of CO2 uptake (45).
Seekell et al. (64) established a DOC threshold of 4.6 mg L−1 in Alaskan and Swedish lakes. Above this tipping point, the DOC–primary production relationship switches from positive to negative.
According to new estimates of global lake DOC concentrations, 89.5% of lakes are above that tipping point (83). Individual field studies substantiate this finding; interior Alaskan lakes have average DOC concentrations as high as 29.2 ± 25.7 mg L−1 (21). Other studies showing similarly high values (>10 mg L1) have been performed in the Mackenzie Delta (84) in northern Alberta and the Northwest Territories, Canada (85), and in central Greenland (86). These studies provide limnological context that suggests DOC concentrations could be impairing light penetration and impeding photosynthesis. These findings suggest that, rather than spurring productivity, DOC shading is dominating. Shading could reduce primary productivity of both phytoplankton and benthic producers (87).
Elevated CDOM has been linked to the suppression of benthic primary production in shallow lakes. Benthic producers have been estimated to be responsible for up to 98% of Arctic lake gross primary productivity (GPP) (88) and are limited primarily by light (89). Given this study’s shallow average lake depth (3.5 m), the overwhelming evidence from field surveys that light reaches shallow lake bottoms (20, 9091–92), and the numerous field studies testifying that heightened DOC can suppress GPP (44, 45, 64), it is likely this decline in greenness reflects, in some part, changes in benthic reflectance. Declining greenness could be a symptom of a larger shift from clear-water systems dominated by benthic production to brown-water ecosystems that are net-heterotrophic (82, 93).
Further Implications.
The continental-scale shifts in lake greenness, and their links to large-scale climatic changes, are not necessarily surprising, since temperature and precipitation patterns have been shown to shape other lake properties, such as dissolved organic carbon concentrations (94) and patterns of greenhouse gas (GHG) emissions (95). Yet our observations have several important implications for water quality, primary productivity, and GHG emissions.
First, changes in color can impact lake surface temperatures (96), which alters mixing regimes and lake stratification (97). Greater stratification can restrict nutrient upwelling and lower bottom-water oxygen concentrations, reducing both water quality and primary production (98). Second, changes in lake greenness may have major implications for lake C storage and GHG emissions. Lake ecosystems are both sources and sinks of carbon (99100101102–103), varying significantly between and within regions. This balance depends to varying degrees on rates of primary production, which have been linked to greenness (20, 54). Lakes that are less green may store less C in sediments because burial rates correlate with primary productivity (104). Changes in greenness may also translate to changes in methane fluxes because productivity and CH4 emissions are correlated (105106107–108). Positive correlations between methane fluxes and chlorophyll-a have been observed in northern temperate lakes (109). In the Mackenzie Delta, methane has been demonstrated to be strongly related to both macrophyte biomass and dissolved organic carbon (110), two variables that can also strongly influence lake color. The importance of substrates in controlling methane emissions has been confirmed in laboratory and field studies of northern-latitude lakes (109, 111). As a result of this relationship, surface reflectance has even been used as a proxy for methane concentration in Siberian rivers (112). However, relationships to carbon dioxide fluxes are more complicated because of hydrologic loading of external CO2. Also, stronger stratification in deeper lakes may starve pelagic phytoplankton of nutrients, causing reduced primary production (113). However, large lakes comprise only ∼1% of the lakes in this study and are likely to already be contributing less methane due to greater methane oxidation in the water column and lower sediment organic substrate availability (95).
Another implication of this decline in greenness is a decrease in GPP. Previous work has established a relationship between airborne and satellite observations of Arctic-boreal lake greenness and GPP (20). Using published relationships between Landsat greenness and in situ GPP, declines in greenness between the first and last decades of this analysis translate to a drop in carbon uptake of 58 Tg C during the growing season (Materials and Methods and SI Appendix, Table S4). Strong evidence exists suggesting the growing season is lengthening in Arctic lakes (114), perhaps by as much as 8.1 d earlier per 100 y (115) across the Northern Hemisphere. Earlier ice-out has been linked to earlier and longer phytoplankton blooms in Finnish lakes (116). However, Landsat’s 16-d return period limits the data available to pinpoint shifts in growing season length. These results, based on an anniversary window approach (117), cannot inform on net annual carbon fluxes because this analysis is constrained to a narrow window during the ice-free growing season and therefore does not include shifts in phenology, under-ice, or nongrowing season dynamics. Field measurements across the entire year should be prioritized to fully constrain net annual carbon fluxes. However, our results indicate a potential 16% reduction in carbon uptake from declining gross primary productivity during the growing season. This finding implies that the mechanism by which lakes take up carbon during the growing season is weakening. Thus, increased DOC accelerates emissions both directly by fueling CO2 and CH4 emissions (118) and indirectly by suppressing photosynthetic carbon uptake (64).
Conclusion
In this study, we present a comprehensive analysis of changing greenness for Arctic and boreal lakes in western North America. We find that lakes with declines in greenness overwhelmingly outnumbered greening lakes, and we show this trend occurs over a much larger region than previous large-scale studies of lake color (40, 119).
As rising temperatures accelerate permafrost thaw, more soil organic carbon and nutrients are expected to be released across terrestrial–aquatic interfaces (57). Evidence increasingly indicates widespread, but spatially complex, changes in Arctic and boreal terrestrial productivity (120). This, in combination with vegetation transitions (121), is predicted to have nonlinear consequences for lake primary productivity (122) and aquatic ecosystems more broadly. The decline in greenness documented here is consistent with the growing number of paleolimnological and biogeochemical field studies that show changes in light absorption resulting from altered terrestrial–aquatic connectivity (22, 42, 45, 58, 64, 119, 123). These findings indicate changes in ecosystem structure and function that in turn impact lake carbon cycling dynamics, which are known to be globally significant (124).
Many field studies, constricted by time and funding, trade off between broad geographic coverage versus higher temporal resolution sampling. In contrast, satellite remote sensing used here encompasses a broad range of diverse geographies while simultaneously preserving high spatial resolution (<100 m). This balance ensures the trends detected are comprehensive while also preserving the uneven change that can occur within and between regions. Our results show that declines in greenness are most evident in warmer and wetter ecosystems, providing additional support to the hypothesis that, in fact, hydrologic change is altering links between terrestrial and aquatic systems (21). These findings warrant continued monitoring of high northern latitude lakes to characterize the significance of lake color change within the broader carbon cycling dynamics of Arctic and boreal landscapes. The boreal Arctic is changing unevenly, and we show that lakes are part of this uneven change.
Materials and Methods
The overall workflow for calculating growing-season composites from 1989 to 2019 for each lake in this study is presented in SI Appendix, Fig. S1. Satellite remote sensing analyses were performed in Google Earth Engine (125); statistics were calculated in Python (126) using a suite of packages and spatial joins were conducted in QGIS (127). Data visualizations were created using all three platforms.
Satellite Data.
The Landsat program offers the longest continuously operating Earth observation satellite record. Landsat Collection 1 Tier 1 surface reflectance from Landsat 5 (n = 66,825), 7 (n = 65,902), and 8 (n = 29,993) was used. Scenes acquired over the study domain (total n = 278,284 scenes; SI Appendix, Fig. S2) were accessed through Google Earth Engine (73), a parallelized geospatial analysis platform featuring colocated storage and compute for large geospatial datasets, including the Landsat archive. The global coverage and analysis-ready Level 1 Surface Reflectance product (46) has been shown to be suitable in previous studies for inland water color analysis including for lake clarity (24), chlorophyll-a and turbidity in rivers (128), and suspended sediment (129). Uncertainties stemming from orbital decay and radiometric intercalibration are noted, and the green band was selected as less influenced by remote sensing uncertainties relative to other bands (52, 53).
Satellite Imagery Masking.
Clouds, cloud shadows, snow, and ice were identified using the CFMask algorithm flags (130131–132) included in the Level 2 Pixel QA band (46). Additional thresholds were imposed to further mask clouds and cloud shadows using Landsat 8’s aerosol band, and haziness was filtered out using the atmospheric opacity bands (threshold of 0.3; L5, L7). Saturated pixels were identified using the radsat_qa band (133), and scenes with low solar zenith angles (<60°) were excluded to minimize glint (134).
Growing-Season Composites.
Growing-season greenness was calculated per pixel as the median value during the peak growing season, which was defined as June and July. Past work confirmed that this is the window of time in which maximum photosynthesis (135, 136) occurs with a minimum of interference from clouds and ice (20, 137). This approach avoids uncertainties introduced by changes in the growing season length or seasonal phenology and was conservatively selected to restrict observations to the ice-free period (after June 1) and to avoid fall cloudiness (SI Appendix, Fig. S3). After filtering lakes to cloud-free observations, lakes were observed an average of 2.8 times each growing season, resulting in over 53 million observations total from 1984 to 2019. These observations were used to generate seasonal composites for each year for each lake (n = 584,596).
Data Sampling.
Growing-season greenness for each lake was identified by sampling a 3-by-3 pixel area centered on each lake centroid calculated from HydroLAKES polygons within the study domain. The median, SD, and pixel count of the green band inside each 3-by-3 pixel box were calculated for each lake at the native scale of the Landsat green band (30 m). Samples without a majority of valid pixels in the 3 × 3 pixel box (n < 5 of 9) were discarded to further reduce influence from adjacent clouds or land (128). Growing season greenness for each lake for each year was then exported for visualization and statistical analysis. The dataset was further screened to exclude negative pixels (<1% of the dataset) and lakes with less than 10 y of data (2% of the dataset). To conservatively ensure that each observation was cloud-free and over water, a final filter using the Pixel QA (pixel_qa) band was imposed to identify only lakes identified as cloud-free and water. Lake centroids were further intersected with the Global Surface Water Dataset (73), and only those identified as having only permanent surface water (transition_class = 1) inside the 3 × 3 grid were maintained in the dataset.
Time Series Analysis.
To assess changes in lake surface reflectance over time, we examined time series of summer (June/July) composites from 1984 to 2019. June and July were selected based on precedent in the literature (135, 136) and to avoid high cloud cover in August (SI Appendix, Fig. S3). This method was required because Landsat’s 16-d return period is insufficient to detect ice-out reliably. The majority (70%) of lakes had over 30 y of observations available (SI Appendix, Fig. S4). Long-term trends in lake color over the four decades were analyzed as change in greenness per year (Δgreen y−1) from 1984 to 2019. Trends were calculated using Theil-Sen’s Slope Estimator from the SciPy package (138), and slopes were tested for significance using a Mann–Kendall test, which is designed to identify monotonic trends and is been widely used to identify terrestrial greening and browning trends (15). This nonparametric approach accounts for gaps in observation years.
Matching Lakes to Basins.
In order to calculate basin-level climate for each lake, lake basins were assigned using the harmonized GLCP dataset (77). The GLCP offers basin-level climate and population data paired to the HydroLAKES dataset. Lakes were matched to basins using the HydroBASINS dataset (version 1.c format 1), which are basin polygons derived from the 15–arc-second–resolution HydroSHEDS dataset. The HydroLAKES lake polygons used in this study, identified by a unique ID (hylak_id), and their corresponding basin IDs were selected from the GLCP dataset. Level 12 basins only were selected as they hold the majority (83%) of lakes (SI Appendix, Table S3) and to exclude artifacts in higher-level basin assignments that are known problems in the GLCP (77). This step resulted in 30,902 unique basins being identified and used for climate analysis.
Ecoregion and Permafrost Extent.
Ecoregions and biomes were assigned to each lake using the geographic boundaries given in the World Wildlife Fund Terrestrial Ecoregions (139) dataset (SI Appendix, Fig. S5). Ecoregions reflect distinct assemblages of flora and fauna that are in turn shaped by climate and landform.
The 825 ecoregions and 14 biomes featured in the WWF dataset were accessed from https://www.sciencebase.gov/catalog/item/508fece8e4b0a1b43c29ca22 as a polygon shapefile (WGS84‒EPSG:4326). Ecoregion and biome classifications were assigned using a spatial join between the WWF layer and HydroLAKES centroids. The majority (76%) of lakes were located in the tundra or boreal forests/taiga biomes, with the remainder found in temperature conifer forests, temperate shrublands, or broadleaf and mixed forests. Shapefiles of circumpolar permafrost extent and ground ice content (140) projected in Lambert Azimuthal Equal Area were downloaded from the National Snow and Ice Data Center (https://nsidc.org/data/ggd318) and reprojected to WGS-84. Lakes were assigned an ecoregion and permafrost extent by spatial-joining the HydroLAKES centroids to the ecoregion and permafrost extent reprojected layers. Fig. 3 shows results for the 472,889 lakes stratified by permafrost extent class, ground ice content, lake depth, and biome classification.
Climate Reanalysis Data.
Trends in temperature and precipitation were calculated from European Centre for Medium-Range Weather Forecasts fifth-generation global atmospheric reanalysis data (ERA5) with a native spatial resolution of 0.25° × 0.25° dating back to 1979 (75). ERA5 data were released in 2019 and are freely available in Google Earth Engine and the Copernicus Climate Data Store. Studies have shown ERA5 to have satisfactory performance in accurately modeling storm surges (141), Canadian prairie climate (142), and global discharge reanalysis (143) compared to previous generations and other models. ERA5 monthly aggregated precipitation and temperature values were summarized annually and seasonally over the study domain. Summer (June, July, and August) and spring (April, May, and June) trends were identified by summarizing data using a Sen’s slope estimator within those time periods by Level 12 basin assigned using the GLCP.
Gross Primary Productivity.
Declines in GPP were calculated from the average greenness for the first (1984 to 1994) and last (2009 to 2019) decade using the relationship between Landsat growing-season greenness and field measurements of GPP identified in Kuhn et al. (20) (SI Appendix, Table S4). Average growing-season GPP rates (in g O2 m−2 d−1) were calculated for the first and last decade and multiplied by the length of growing season (60 d) used in the remote sensing analysis. The resulting areal GPP (in g O2 m−2) was then multiplied by the total lake area within the domain (417,152.40 km2) to estimate the total change in O2 (Tg) uptake between the two decades. Total lake area was calculated by summing the individual HydroLAKES lake area for lakes passing quality control. Photosynthetic uptake of carbon was calculated using the carbon-oxygen photosynthetic quotient of 1.2 where, for every 1 mol of oxygen produced, 1.2 mol of carbon is assimilated (144).
Acknowledgements
This work was supported by the NASA Earth and Space Science Fellowship (80NSSC18K1336), the University of Washington College of the Environment Integral Environmental Big Data Award, and the NASA Terrestrial Ecology Program through Arctic and Boreal Vulnerability Experiment Grants (80NSSC18K1336 and NNX15AU04A). We thank the eScience Institute, the developers at Google Earth Engine, and NASA and the US Geological Survey for keeping the Landsat archive free and open. We thank Dr. Uma Bhatt for her assistance with ERA5 data, Drs. Matthew Bogard and Ben Miller for helpful comments, Drs. Jill Deines and Jon Wang for their help with the Landsat cloud masks, and M.M.C. and anonymous peer reviewers for their contributions to this work.
Data Availability
Data for this research are presented in the manuscript and Supplementary Information. In addition, tables of the Landsat time series and trends have been deposited in the publically available NASA ABoVE Science Cloud (https://daac.ornl.gov/cgi-bin/dsviewer.pl?ds_id=1866) (145).
References
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
Declining greenness in Arctic-boreal lakes
