Four Decades of Land-Cover Change on the Kenai Peninsula, Alaska: Detecting Disturbance-Influenced Vegetation Shifts Using Landsat Legacy Data

Across Alaska’s Kenai Peninsula, disturbance events have removed large areas of forest over the last half century. Simultaneously, succession and landscape evolution have facilitated forest regrowth and expansion. Detecting forest loss within known pulse disturbance events is often straightforward given that reduction in tree cover is a readily detectable and measurable land-cover change. Land-cover change is more difficult to quantify when disturbance events are unknown, remote, or environmental response is slow in relation to human observation. While disturbance events and related land-cover change are relatively instant, assessing patterns of post-disturbance succession requires long term monitoring. Here, we describe a method for classifying land cover and quantifying land-cover change over time, using Landsat legacy imagery for three historical eras on the western Kenai Peninsula: 1973–2002, 2002–2017, and 1973–2017. Scenes from numerous Landsat sensors, including summer and winter seasons, were acquired between 1973 and 2017 and used to classify vegetation cover using a random forest classifier. Land-cover type was summarized by era and combined to produce a dataset capturing spatially explicit land-cover change at a moderate 30-m resolution. Our results document large-scale forest loss across the study area that can be attributed to known disturbance events including beetle kill and wildfire. Despite numerous and extensive disturbances resulting in forest loss, we estimate that the study area has experienced net forest gain over the duration of our study period due to reforestation within large fire events that predate this study. Transition between forest and graminoid non-forest land cover including wetlands and herbaceous uplands is the most common land-cover change—representing recruitment of a graminoid dominated understory following forest loss and the return of forest canopy given sufficient time post-disturbance.


Introduction
The boreal forest is a global resource in terms of biodiversity, silviculture production [1], and global carbon stocks [2]. Boreal forests, which encompass about 30% of Earth's forested area [3], are a massive carbon reservoir, composing one third or between 367.3 and 1715.8 Pg of Earth's terrestrial carbon [4]; Although humans have inhabited the region for millennia, until recently the human impacts to land cover have been dwarfed by natural disturbances. The Dena'ina and Alutiiq people settled the Kenai Peninsula around 1000 A.D. and engaged in vibrant subsistence economies [34]. Settlers from Russia and the United States arrived in the late 17th and 18th centuries to support the fur trade,  Although humans have inhabited the region for millennia, until recently the human impacts to land cover have been dwarfed by natural disturbances. The Dena'ina and Alutiiq people settled the Kenai Peninsula around 1000 A.D. and engaged in vibrant subsistence economies [34]. Settlers from Russia and the United States arrived in the late 17th and 18th centuries to support the fur trade, salmon canneries, and fox farming. The human population increased after gold was discovered in 1895; in the 1940s, when homesteading opened; and again in the 1950s, when oil and gas were discovered [35]. These activities began and continue to result in larger human-induced land-cover change including forest removal for infrastructure and agricultural resources and incidental human-caused wildfires. The Chugach National Forest was established in 1907 and the Kenai Moose Range (now a National Wildlife Refuge) was established in 1941. Highway construction to link the western peninsula to Anchorage began in 1946. Wildfires and spruce bark beetle outbreaks are significant disturbance agents within forests of the Kenai Peninsula [23,24,36]. Reconstruction of historical mean fire intervals (MFI) for the Kenai Peninsula indicate relatively short MFI of approximately 80 years in black spruce forests of the northern half of the peninsula [37] and much longer, 400-600-year MFI in white and Lutz spruce forests [24]. Major wildfires that have occurred on the peninsula within the period of detailed fire history records include the 1947 Kenai fire (128,727 ha), 1969 Swanson River fire (34,522 ha), 2015 Funny River fire (79,565 ha), and the 2019 Swan Lake fire (67,656 ha) as reported in the Alaska Large Fire Database (https://www.frames.gov/catalog/10465) [38]. Spruce bark beetles have also been an historically significant disturbance agent within white, Lutz, and Sitka spruce (Picea sitchensis (Bong.) Carrière) forests, with large outbreaks resulting in significant tree mortality approximately every 50 years [24]. Notably, spruce beetle outbreaks that began in the 1980s and continue today have affected the majority of Lutz spruce forests in the southern portion of the Kenai Peninsula [23,39]. Other natural disturbance agents (wind, fire, insect, fungus, and geomorphic processes) have affected all forest patches since post-glacial forest establishment [29].

Approach
We classified the entire Kenai Peninsula Lowland study area using temporally constrained satellite imagery, a random forest classifier, and training data derived from manual interpretation of high-resolution aerial imagery. Available data limited our approach to producing land-cover classification for three eras: circa 1973, circa 2002, and circa 2017. In each of these eras, two conditions were met. First, high-resolution aerial imagery was available for ground-truthing and model development within parts of the study area. Second, satellite imagery of sufficient quality for both summer (June, July, August) and winter (February, March, April) seasons was available for landscape classification. After historical land-cover classifications were developed for each era, land-cover change for three periods of assessment were tabulated. These periods include the early 1973-2002 period, the later 2002-2017 period, and the net change experienced over the entire 1973-2017 period.

Vegetation Classification
Vegetation type was both interpreted and modeled as one of nine general classes: needleleaf forest, broadleaf forest, mixed forest, alpine, barren, shrub, herbaceous, wetland, or water. These lifeform classifications are ecologically significant to the landscape and are also common to other land-cover classifications such as the National Land Cover Database (NLCD) [40]. We did not classify urban or developed land-cover types and instead used the extent of developed land covers from the 2011 NLCD (https://www.mrlc.gov/data/nlcd-2011-land-cover-alaska-0) to mask built environments from each era so that only consistently naturally vegetated areas were considered.

Image Acquisition
We used Landsat imagery acquired between 1973 and 2017 to classify vegetation into nine general classes based on vegetation height, density, and photosynthetic activity. We inspected the entire Landsat Collection-1 Level-1 archive housed within the U.S. Geological Survey (USGS) Center for Earth Resources Observation and Science (EROS) EarthExplorer web application for our study area and seasons of interest. From all Landsat platforms, we ultimately used scenes from Landsat 1 Multispectral Scanner (MSS), Landsat 7 Enhanced Thematic Mapper Plus (ETM+), and Landsat 8 Operational Land Imager (OLI) defining three eras: circa 1973, circa 2002, and circa 2017. These eras were solely the result of cloud-free scene availability during both summer and winter seasons. These three eras are the only eras where such conditions were met within our study area. Landsat Collection 1 data products are deemed suitable for time-series pixel level analysis and are organized into two tiers. Tier 1 contains the highest quality scenes considered suitable for time-series analysis, while Tier 2 consists of scenes potentially suitable for time-series analysis even though they do not meet Tier 1 criteria (https://lta.cr.usgs.gov/). Our analysis used both Tier 1 and Tier 2 scenes (Table 1). To facilitate the primary delineation between forest and non-forest areas, we used the blue band (Landsat 8 B2, Landsat 7 B1) from late winter scenes as a measure of vegetation density and general structure. We chose the blue band based on the snow's highest reflectance in the blue portion of the electromagnetic spectrum. For Landsat 1, which was not designed to acquire data in the blue band, we used the infrared band (Landsat 1 B7), which is only slightly less sensitive to reflectance of light from snow. Detection of forested area is based on the difference in albedo between forest canopy and open ground, where dense forest canopies mask the reflective ground snow, and non-forested areas are moderately to highly reflective in relation to decreased height and density of vegetation overtopping snow. Using this method, we found that land-cover types that are often misclassified, such as graminoid wetlands and low stature black spruce, were much more readily discernible when the winter scene was considered.  Two overlapping winter scenes were used for each season for both the 2017 and 2002 eras. The 1973  era was prone to winter cloud cover, so five scenes from 1974-1975 were mosaicked together (Table 1). We used level 1 digital number (DN) values for all winter scenes except for those used in the c1973 era. For this era, each scene was manually clipped to remove cloud and shadow contamination; the DN values were standardized by using dark object subtraction based on ice-free offshore water and then rescaled so that large snow-covered lowland lakes common to overlapping scenes had a DN value of 140. These scenes were then ordered and mosaicked in ArcMap [41] using the FIRST operator, a method that determines the pixel value from the first raster dataset encountered in a mosaic list.
We used the Level 2 Enhanced Vegetation Index (EVI) spectral index derived from clear-sky summer scenes for the 2017 and 2002 eras available from USGS EROS Science Processing Architecture (ESPA, Equation (1)). This index is a ratio of red (R) and near infrared (NIR) signatures that incorporates the blue (B) band to reduce atmospheric noise, and results in a single-band product that scales the "greenness" of vegetation and permits classification of lifeforms based on the amount of photosynthetically active vegetation: As an example, EVI can be used to differentiate broadleaf (highly photosynthetically active) and evergreen (less photosynthetically active) vegetation [42,43]. As ESPA does not provide Level 2 spectral indices for Landsat 1 or Landsat 2 from the 1973 era, we used the Transformed Vegetation Index (TVI) based on the red and infrared bands (Equation (2)). This index also functions as a measure of greenness by detecting differences in photosynthetic activity [44]. Scenes were selected to capture peak greenness (July-August) and minimal cloud contamination. Two overlapping scenes were also used for each summer season for each era.

Model Development
The study area was divided into 30 × 30-m cells (0.09 ha) based on the Landsat 7 and Landsat 8 pixel resolution. Each cell was assigned a feature list containing pixel values from the winter and summer Landsat mosaics and values for slope, aspect, and elevation, which were derived from the Alaska 2-Arc-second digital elevation model (DEM) with a 60-m resolution [45]. A random forest classifier was developed for each era based on 3600 training points (400 points for each land-cover type). Training points were randomly distributed within training areas, and training areas were opportunistically distributed to capture the regional and geomorphic extent of each land-cover type to the extent possible given availability of aerial imagery. Training areas for the 2017 era were guided by plot data used in the development of the 2017 Kenai Vegetation Mapping Project [46]. Training areas were drawn around homogenous stands of vegetation identified by these points; these totaled 15,929 hectares and were distributed across the entire study area. Land cover was verified using high-resolution imagery acquired in 2016 by the Kenai Peninsula Borough (http://maps.kpb.us/gis/rest/services/2016_Imagery/MapServer).
For the 2002 era, we used previous interpretations of aerial imagery completed circa 1999 when the majority of the Kenai Peninsula's land cover was manually classified. The shapefiles delineating land-cover type were acquired from the Kenai Peninsula Borough's interactive Land Cover Viewer (http://maps.kpb.us/gis/rest/services/Vegetation/MapServer). This approach gave our 2002 era training area a net area of 844,920 hectares. For the 1973 era, training areas totaled 34,275 hectares and were restricted to the northern half of the study area due to coverage of high-resolution air photos acquired in 1975 and provided upon request by the Kenai National Wildlife Refuge. No reference data were available, and training areas were placed where homogenous vegetation patches could be confidently identified and manually delineated. Each training point was assigned feature list values and a land-cover type based on the corresponding cell values from the winter and summer Landsat mosaics, values for slope, aspect, and elevation derived from the same 60-m DEM, and manual interpretation of high-resolution aerial imagery. These training data were read into the randomForest package in R [47] and used to create a predictive model (classifier). The random forest classifier was then run for every cell within the entire study area. The classifier looked at the feature list for every pixel within the study area to determine the most probable land-cover type ( Table 2). An independent model was developed and executed for each era [48]. Each model used 500 trees and tested two variables at each node. All other parameters remained in the default settings offered in the R package [47]. In addition to classifying vegetation type, we derived the classification probability for each pixel that represents the proportion of decision trees that "voted" for the majority classification [48]. In our classifications, 11.2% probability represented the minimum proportion of votes (56 out of 500) needed to achieve a majority vote. To produce a conservative estimate of land-cover transition rates, we summarized our estimates of land-cover area and transition rates for pixels with a probability greater than 50%. We bound these estimates with error bars representing probabilities >25% (upper) and probabilities >75% (lower).

Land-Cover Assessment and Change Detection
Our region of analysis encompasses a total area of 9300 km 2 . We excluded urban and developed areas mapped in the 2011 NLCD, which removed approximately 160 km 2 from within the study area for a final area of 9142 km 2 . The proportion of this final area that is classified depends on the probability level used in our summary statistics. All 9142 km 2 are classified with a minimum probability of 11.2%. When we only consider model predictions where the probability of the majority predicted land cover was greater than 25% (over twice the minimum required for a majority vote), 2-3% of our final study area becomes indeterminable. In other words, 0.4-0.5% (37 km 2 -46 km 2 ) of our study area becomes indeterminable depending on the era. If we restrict our analysis to pixels where the probability of the majority predicted land cover is greater than 75%, our confidence in model predictions is maximized but 64-81% (5864 km 2 -7403 km 2 ) of our study area becomes indeterminable depending on the era. When using the intermediate threshold of a >50% probability, 27-38% (2505 km 2 -3503 km 2 ) of our study area becomes indeterminable. This 50% probability threshold provides a conservative estimate and allows us to more confidently approximate the spatial extent of various land-cover types and provide margins of error associated with vegetation transition types and transition localities ( Figure 2).

Land Cover by Era
For the 1973 era, 71% (6457 km 2 ) of the study area was classified. Forest covered approximately 39% of the study area (2494 km 2 ) (Figure 3), while non-forest land-cover types excluding water covered 61% of the study area and totaled 3429 km 2 (Table 5, Figure 3). Of the forest cover, roughly equal parts were needleleaf (16% or 1024 km 2 ) and mixed species (16% or 1064 km 2 ) forest ( Figure 2). Broadleaf forest made up roughly 6% of the total forest area (407 km 2 ) ( Figure 2). Approximately 43% and 20% of the total non-forested area was made up of graminoid and shrub cover types, respectively, To produce a spatially explicit land-change product, a feature list was produced for each grid cell within the study area that contained the predicted land-cover type for each era. Using the ArcGIS environment, this feature list could be queried to group specific sequences of land-cover classes into interpreted ecological processes (e.g., deforestation). For example, if a pixel was classified as needleleaf, mixed, or broadleaf forest in one era and then transitioned to a non-forested class in a subsequent era, it was termed "deforestation". Grid cells that were not forested in an era and transitioned to a forest class in a later era were considered "reforestation." We considered pixels that remained in wetland classification for the 1973 and 2002 era but transitioned to a woody classification (forest or shrub) in 2017 to represent afforestation of wetlands. We considered pixels that remained in shrub classification for the 1973 and 2002 era but transitioned to forest in 2017 to represent afforestation of shrub lands and a possible elevational increase in tree line. The total areas of each transition type were summarized within the ArcGIS environment by counting the number of grid cells that shared specific transition types. Lastly, raster datasets were produced to graphically display the spatial extent of forest loss and gain among the eras.

Classification Accuracy
When using random forest classification routines, a subset (approximately 36%) of the training data is automatically reserved and employed in model testing. The test statistic is reported as the out of bag estimate of error rate. The error rate for 1973 was 14.5%. Out of bag estimates of error rate for 2002 and 2017 were 30.67% and 10.06%, respectively. General classification accuracy was good across all three eras with classification accuracies of >80% for forest (needleleaf, broadleaf, or mixed) vs. non-forest (Table 3). User accuracy for classifying forest was highest for 2017 at 97.7% and lowest for 2002 at 85.4%. The 1973 era attained user accuracy of 94.5 %. General classification accuracy appeared to be influenced by training area as the larger the training area, the lower the user's and producer's accuracy. Table 3. Confusion matrices summarizing out of bag error estimation rates for forest and non-forest land-cover types.

1973
Referenced Classification accuracy was >50% for all vegetation types across all eras ( Table 4). The 2002 era had the lowest user's accuracy (UA) compared to the 2017 and 1973 eras. This is likely because the 2002 era has a broader training area with more scattered training locations and greater variation in the attributes associated with the training data points. The 2017 era had the highest UA across all vegetation types, with all values above 80%. This may represent greater classification potential from the increased radiometric resolution of the OLI sensor found in Landsat 8 as compared with ETM+ sensor used in Landsat 7 and MSS sensor used in Landsat 1 and 2.

1973-2002
Between 1973 and 2002, 18% or 695 km 2 of non-forest pixels transitioned to forest. The majority (54% or 373 km 2 ) of these transitions were in pixels classified as graminoid (Table 6). Transition to forest from shrub (46% or 321 km 2 ) was the second most common transition. The majority (50% or 351 km 2 ) of the forest gained during this period was needleleaf forest. Mixed species forest was the next most common (29% or 200 km 2 ) followed by broadleaf forests (21% or 145 km 2 ).

1973-2017
Between 1973 and 2017 (the net change within our study area), 23% of the area (886 km 2 ) classed as non-forest in 1973 had transitioned to forest by 2017. The majority of the net forest gain between 1973 and 2017 was needleleaf (52%) totaling 459 km 2 . Mixed species forest was the next most common (25% or 225 km 2 ), followed by broadleaf forests (23%, 202 km 2 ). Reforestation occurred in nearly equal proportions of area classified as shrub (47% or 416 km 2 ) and graminoid (52% or 463 km 2 ) land cover ( Table 6). A portion (~0.3%) of this forest gain came from pixels previously classified as alpine.
Afforestation, the intrusion of shrubs and trees into historically non-forested areas, accounts for a small proportion of the observed forest gains in this time frame. Conversion of persistent wetland into forest made up 0.2% or 0.6 km 2 ( Table 6). Of this area, 76% converted to needleleaf, 8% converted to mixed forest, and 17% converted to broadleaf forest. Wetland transition into shrubs occurred over an additional 2.6 km 2 . Transition from shrub to forest accounted for 0.3% or 0.9 km 2 of forest gain ( Table 6). The majority of shrub to forest transition was to needleleaf forest (47%) followed by transition to broadleaf forest (29%) and mixed species forest (24%).

Four Decades of Land-Cover change on the Kenai Peninsula
Our analysis indicates that substantial land-cover change has occurred across the Kenai Peninsula over the last 44 years with the exception of water extent, which remained stable between 1973 and 2017 with a slight peak in 2002 (Figure 3). Although deforestation events related to spruce bark beetle outbreaks, wildfires, and timber harvest have occurred on the Kenai Peninsula over the past four decades [24,36,49], our classifications suggest that the peninsula has actually become more forested over the period of analysis, because significant and prolonged forest loss in the southern half of the peninsula ( Figure 5) has been offset by forest regrowth in the northern half. Forest regrowth is characterized as increased cover of needleleaf and broadleaf forest between 1973 and 2017. These transitions were primarily within the large fire perimeters of the Kenai fire of 1947 and the Swanson River fire of 1969.
Additional reforestation occurred as the infilling and maturation of forest cover around and within existing forest patches. Prominent reforestation occurred within the 1994 Windy Creek fire, located in the central portion of the study area. We also noted what appeared to be reforested logging patches in the southern portion of the study area. Some reforestation occurred within areas previously classified as alpine. It is possible that these pixels represent an encroachment of trees into higher elevations; however, we suspect these are erroneous classifications due to the improved ability of Landsat 8 over earlier platforms to detect small and scattered trees characteristic of elevational tree lines. Transition areas from wetland to forest and/or shrub are scattered throughout the study area and predominantly clustered within the large wetland complexes north and south of the communities of Kenai and Soldotna. Additional transitions are located within the wetland complex north of Skilak Lake. The transition of forest to shrub occurred along mountain slopes where forest stands blend into shrub thickets. This pattern aligns with the elevational increase in tree line that has been documented elsewhere on the peninsula [50]. bark beetle outbreaks, wildfires, and timber harvest have occurred on the Kenai Peninsula over the past four decades [24,36,49], our classifications suggest that the peninsula has actually become more forested over the period of analysis, because significant and prolonged forest loss in the southern half of the peninsula ( Figure 5) has been offset by forest regrowth in the northern half. Forest regrowth is characterized as increased cover of needleleaf and broadleaf forest between 1973 and 2017. These transitions were primarily within the large fire perimeters of the Kenai fire of 1947 and the Swanson River fire of 1969. Additional reforestation occurred as the infilling and maturation of forest cover around and within existing forest patches. Prominent reforestation occurred within the 1994 Windy Creek fire, located in the central portion of the study area. We also noted what appeared to be reforested logging patches in the southern portion of the study area. Some reforestation occurred within areas previously classified as alpine. It is possible that these pixels represent an encroachment of trees into higher elevations; however, we suspect these are erroneous classifications due to the improved ability of Landsat 8 over earlier platforms to detect small and scattered trees characteristic of elevational tree lines. Transition areas from wetland to forest and/or shrub are scattered throughout the study area and predominantly clustered within the large wetland complexes north and south of the communities of Kenai and Soldotna. Additional transitions are located within the wetland complex Interestingly, mixed-forest cover decreased between 1973 and 2017. We hypothesize that this was due to transitions from mixed forest to needleleaf forest via maturation of coniferous tree canopies in mixed-forest patches, which sways classification towards a needleleaf designation. Mixed forest to broadleaf transitions are coincident with large forest patches in the northern half of the peninsula and may represent disagreement in classification between the eras. The increased spatial and spectral resolution in Landsat 7 and Landsat 8 offer greater ability to differentiate broadleaf forests from the more common mixed-forest types.
Known disturbances from insect outbreaks, wildfire, and timber harvests align well with our detection of land-cover change, giving us confidence that the methods employed in our study capture true change. Land-cover changes are primarily attributed to pulse disturbances that are readily detected by both Landsat sensors and human observers. Between 1973 and 2002, forest loss was concentrated within the southern half of the peninsula (Figure 4), primarily associated with known spruce bark beetle mortality, resultant salvage logging operations, and the 1996 Crooked Creek wildfire [38].  [38]. Smaller additional areas of forest loss are distributed across the southern half of the peninsula and likely represent human-and insect-caused forest loss. Some deforestation within the northern half of the study area can be attributed to prescribed fire and forest mastication efforts completed within the Kenai National Wildlife Refuge. These disturbance events were often perceived as negative impacts on the biophysical and socioeconomic systems by people living in the area [51]. The loss of forest and conversion of forest to graminoid-dominated cover has received considerable attention from both the scientific community and general public [36]. Less perceptible to human observation, but well-captured in our long-term analysis, are the slow shifts in vegetation primarily due to forest succession. This land-cover change has not received the same attention as that of forest loss, likely due to the slower process of natural succession and the fact that the majority of forest succession is occurring away from populated areas. The increase in forest area is largely a reflection of transitions from shrubs to forest during the 1974-2002 period. Within this time frame, reforestation occurred primarily as recovery from the 1947 Kenai fire. The greater extent of reforestation in the north as compared with the southern portion of the peninsula may be partially explained by disturbance history and species responses. Wildfires have historically been more frequent and extensive in the northern portion of the peninsula [24,38], whereas bark beetle-caused forest mortality has been more extensive in the southern portion ( Figure 6). Both disturbance agents (fires and beetle outbreaks) have occurred throughout the study area, suggesting that the region is broadly vulnerable to disturbance; however, wetter and cooler conditions in the south relative to the north [33] may explain the reduced frequency and extent of fires, while drier and warmer conditions in the northern peninsula have facilitated wildfire spread. Low-elevation forests, in which the principal needleleaf species are white and/or Lutz spruce, are the most vulnerable to beetle attack [52] and were historically most common in (or in the case of Lutz spruce, restricted to) the southern portion of the Kenai Peninsula. Black spruce forests, most extensive in the northern portion of the peninsula, were found to be least vulnerable to beetle damage [52] but are especially flammable and therefore regularly disturbed by wildfire [37,53]. For many of the woody species found on the Kenai Peninsula, exposure of mineral soil improved rates of seedling establishment and survival [54][55][56], which is more likely to occur with forest fires that create tree throws (subsurface depressions left by tree roots) and combust the surficial organic horizon [57][58][59]. In contrast, bark beetle-caused tree mortality can lead to a release of Calamagrostis grass, which is negatively correlated with seedling establishment [60]. This may be slowing reforestation dynamics in the southern portion of the peninsula [21] and highlights unique species, disturbance, and biogeographic patterns of vulnerability and resilience. in the southern portion of the peninsula [21] and highlights unique species, disturbance, and biogeographic patterns of vulnerability and resilience. The vast majority of the reforestation we report in this study is natural succession of forest patches where forest cover is the historical land-cover type ( Figure 7A,A') and instances of non-forest status only occur following pulse disturbances such as wildfire ( Figure 7B,B'). A subset of this reforestation, however, represents true afforestation-the establishment of forest communities in locations that are historically devoid of forest. This is the case for locations that are traditionally wetland or upland non-forest that have clearly transitioned to a forested status ( Figure 7C,C'). We found it necessary to restrict our analysis to our first two eras to identify pixels that were persistently non-forested and that transitioned to forest between 2002 and 2017. We further emphasize that site specific investigations are required to document true encroachment of forest cover into areas that have no known disturbance or appear to have been historically persistent non-forested areas (Figure The vast majority of the reforestation we report in this study is natural succession of forest patches where forest cover is the historical land-cover type ( Figure 7A,A') and instances of non-forest status only occur following pulse disturbances such as wildfire ( Figure 7B,B'). A subset of this reforestation, however, represents true afforestation-the establishment of forest communities in locations that are historically devoid of forest. This is the case for locations that are traditionally wetland or upland non-forest that have clearly transitioned to a forested status ( Figure 7C,C'). We found it necessary to restrict our analysis to our first two eras to identify pixels that were persistently non-forested and that transitioned to forest between 2002 and 2017. We further emphasize that site specific investigations are required to document true encroachment of forest cover into areas that have no known disturbance or appear to have been historically persistent non-forested areas ( Figure 7C,C').
Large and significant shifts in vegetation continue to occur across the study area. During the summer of 2019, numerous wildfires resulted in further natural forest loss. The largest of these fires, the +67,000-hectare Swan Lake fire, largely burned at stand-replacing severity in black spruce forests that had reforested since the 1947 Kenai fire but also burned climax stage spruce and birch forests and subalpine hemlock forests with no known fire history. Additionally, the fire burned above the tree line along the northern Kenai Mountains.
The coming decades are likely to bring further shifts in land cover associated with changing climate and amplification of natural disturbances. Alaskan ecosystems are highly vulnerable to climate change [61], because northern latitudes are warming at twice the rate of the global average [62]. Warming temperatures that increase spruce susceptibility to beetle attack, shorten beetle life cycle, and increase the number of individuals that mature in a given year can synchronize beetle attacks at landscape or regional scales [23,39,63], putting remaining or regenerating spruce forests at risk of sustained insect outbreaks. Above-average temperatures and prolonged drought are contributing factors to increased incidence of wildfires in the boreal forest [64], a pattern that may be amplified at local scales where widespread flushes of Calamagrostis grass and increased overall surface fuel loads following spruce beetle outbreaks have increased fuel flammability [65]. Recent disturbances and future projections mean that existing estimates of forest and non-forest cover, both in our study and others, are no longer representative of rapidly changing, on-the-ground conditions. This necessitates the continued systematic mapping of prominent vegetation for the purposes of inventory and monitoring, land and resource management, and an understanding of the resilience of boreal forest ecosystems to global change. The coming decades are likely to bring further shifts in land cover associated with changing climate and amplification of natural disturbances. Alaskan ecosystems are highly vulnerable to climate change [61], because northern latitudes are warming at twice the rate of the global average [62]. Warming temperatures that increase spruce susceptibility to beetle attack, shorten beetle life cycle, and increase the number of individuals that mature in a given year can synchronize beetle attacks at landscape or regional scales [23,39,63], putting remaining or regenerating spruce forests at

Assumptions and Limitations
Two primary components of our methods influence the accuracy of our model outputs. First, we assume that a change in the spectral parameters of pixels from one era to the next represents true change in vegetation type. High magnitude shifts in spectral conditions are most likely to represent changes in land cover, but subtle variation in spectral reflectance can occur for a location despite no actual land-cover change. This variation originates from the natural seasonal variation in vegetation (e.g., the extent of leaf out, foliar senescence, vegetation health, and physically wet vegetation), the physical variability of atmospheric conditions, and the configuration of the observing remote sensing platform. We have attempted to minimize the effects by limiting our data to cloud-free scenes and scenes that capture peak greenness (mid-summer) and maximal snow extent and depth (late winter/early-spring). For example, our methods depend on the use of winter clear-sky scenes to deduce forest from non-forest and infer forest canopy structure. Under ideal conditions, this offers great power in defining land-cover types given the significant differences in reflected spectra. In dense evergreen forests, snow-covered, closed-canopy trees could be interpreted due to increased reflectivity [66]. While we strived to use the best available scenes from the Landsat campaign, it is possible that snow-covered trees introduced error and inflated our reforestation or deforestation estimates. We attempt to reduce the detection of false transitions by restricting our summary to cells that had high probabilities (over 50%) of being correctly classified in each era, maximizing our confidence in the assigned land-cover type and identified land-cover transitions. However, the exclusion of less probable pixels may mask useful but more subtle patterns.
Second, the location and extent of training areas influences model performance when classifying land-cover types. This was detected around the maximal and minimal elevational limits of the alpine land-cover training areas. As a result, restricted distribution of training areas could artificially restrict the distribution of predicted land-cover types. The historical nature of our analysis also prohibited site visits to verify field conditions as the vegetation has likely changed over time. We attempted to minimize this effect by opportunistically establishing training areas throughout the study area to the best of our ability to confidently interpret land-cover types using high-resolution imagery. Furthermore, we rely on pre-existing datasets for ground-truthed reference data in the 2002 and 2017 eras. Any error associated with these data at the time of collection has been passed along into our classification accuracy. While overall agreement between the respective eras of classifications, known disturbance events, and land-cover patterns was good, an obvious inconsistency can be found with the distribution of alpine regions in the 2002 era. While the alpine and subalpine shrub in the Caribou Hills and north-western Kenai Mountains are generally consistent between 1973 and 2017, there are spatial differences in the 2002 era. This is due to alpine areas having been classified most broadly as "alpine" and more specifically as "alpine-shrub" in the Kenai Peninsula Land Cover dataset. When deriving our reference data for this era, we incorporated these areas as alpine assuming the broadest classification would be the most accurate. However, upon comparison, classifying these reference areas as shrub may have improved the era-to-era agreement of vegetation at these specific locations.