Informing Sustainable Forest Management: Remote Sensing Strategies for Assessing Soil Disturbance after Wildﬁre and Salvage Logging

: Wildﬁres have nearly become a guaranteed annual event in most western National Forests. Severe ﬁre effects can be mitigated with a goal of minimizing the hydrologic response and promoting soil and vegetation recovery towards the pre-disturbance condition. Sometimes, post-ﬁre actions include salvage logging to recover timber value and to remove excess fuels. Salvage logging was conducted after three large wildﬁres on the Lolo National Forest in Montana, USA, between 2017 and 2019. We evaluated detrimental soil disturbance (DSD) on seven units that were burned at low, moderate, and high soil burn severity in 2022, three to ﬁve years after the logging occurred. We found a range of exposed soil of 5%–25% and DSD from 3% to 20%, and these values were signiﬁcantly correlated at r = 0.88. Very-high-resolution WorldView-2 imagery that coincided with the ﬁeld campaign was used to calculate Normal Differenced Vegetation Index (NDVI) across the salvaged areas; we found that NDVI values were signiﬁcantly correlated to DSD at r = 0.87. We were able to further examine this relationship and determined NDVI threshold values that corresponded to high-DSD areas, as well as develop a model to estimate the contributions of equipment type, seasonality, topography, and burn severity to DSD. A decision-making tool which combines these factors and NDVI is presented to support land managers in planning, evaluating, and monitoring disturbance from post-ﬁre salvage logging.


Introduction
Fire and fuels management has a long history on public lands in western US forests.The occurrence of large wildfires has increased in recent years due to numerous interrelated factors, including a changing climate, decades of fire suppression, human encroachment into the wildland-urban interface (WUI), insects, and disease [1][2][3][4].There is an equilibrium to be achieved with fire on public forested lands, as fire is a natural phenomenon that provides essential ecosystem services, including forest thinning, disease prevention/eradication, and soil nutrient cycling [5].A balancing act is necessary when surrounding communities, air quality, and the secondary effects of severe wildfire (e.g., hydrologic events, vegetation succession changes, etc.) are considered for the use of fire on the landscape.When the natural requirement for fuel reduction cannot be achieved with fire (natural or prescribed), mechanical fuel removal can be implemented.In some circumstances after wildfire, salvage logging is implemented to recover timber value and to remove excess fuels [6].
The extent and degree of the soil and vegetation disturbance by wildfire may have serious ecologic and economic implications on the forest system for years [7].The severity of the burn combined with the resilience of the ecosystem influences the magnitude of secondary fire effects, as well as the trajectory back towards the pre-fire or pre-disturbed state.Vegetation consumption and damage is often the most immediate and visible effect of wildfire.The impact of wildfire varies by pre-fire environmental conditions (e.g., drought, vegetation health, weather) and the intensity and severity at which the fire burned [8,9].Postwildfire precipitation, particularly on severely burned soils, can rapidly lead to severe impacts on water quality and produce potentially life-threatening debris flows in rural areas and in the wildland-urban interface [10].Fire effects on soils go beyond soil loss (erosion); there is also loss of organic matter, and changes in soil chemical, physical, and biological properties [11,12].Mitigation measures, such as mulching or seeding, can be implemented to minimize short-term hydrologic effects [13], while longer-term solutions may include reforestation treatments [2].
The compound effects of severe fire and salvage disturbance can delay recovery by an additional year or more [14,15].Consequences can reach areas far downslope of the logging unit, depending on the degree, scale, and connectivity of the disturbance [16][17][18].Most salvage logging operations employ some ground-based machinery, either for the primary harvest and yarding, or for the eventual transport of the logs offsite [19].Ground-based operations (e.g., feller-bunchers, skidders, tractors, forwarders) often have more damaging ground effects than operations that have aerial components (e.g., skyline or helicopter).Equipment traffic and the disturbance of soil and vegetation from yarding logs can result in a connected network of trails and landings.The configuration of these trails, called skid trails, may look similar to a stream network, or trails may be parallel across (on contour) or with the slope; additionally, the orientation to the slope and the connectivity to roads and nearby streams influences the potential elevated runoff and erosion risk [13,20].Federal regulations and state best management practices (BMPs) [16] optimally exist to ensure logging is completed with minimal disturbance to the forest floor, existing roads, and surrounding unlogged areas.Natural resource values to be protected include water, vegetation, soils, timber, wildlife, fish, and productivity [21].
Land management planning agencies have protocols in place to manage, minimize, and mitigate potential disturbance from logging activities [22].Contracts for logging are awarded to private timber companies and contractors, and there is an expectation that contractors follow the protocol to minimize detrimental soil impacts.For example, there are Forest Service directives to maintain soil productivity by limiting the amount of detrimental soil disturbance (DSD) to less than 15% of the areal extent of a logging unit [23,24].DSD is defined in [19], and is a combination of the compaction, rutting, and displacement of soils after logging.Factors impacting DSD include burn severity, soil parent material, topographic features, the logging equipment used (e.g., ground or skyline harvest; tracked or rubber tire logging equipment; helicopter), the harvest season (winter or summer), and site/soil resiliency [22,[25][26][27].DSD can exist prior to the logging operation from past forest activities, and effects can be cumulative and long-lasting from previous logging, severe fire, and salvage operations [28].There is a particular need to consider the effects of logging operations on burned slopes, where ground cover is sparse and where there is direct connectivity to roads and streams [13,20].
Since the passage of the National Forest Management Act of 1976 (NFMA) [24], regions and forests in the US Forest Service have implemented standards and guidelines to maintain soil productivity by limiting forest soil disturbance.Minimizing the total area of disturbance is an operational strategy for minimizing the impact to forest soils and restoration expenses [29].The USDA Forest Service monitors soil disturbance as a means to evaluate regulatory compliance.In the Lolo National Forest Monitoring and Evaluation Report [30], the monitoring reported findings where logging units were assessed for detrimental soil disturbance against an areal threshold of 15% DSD.Exceeding 15% DSD indicates potential long-term impairment to soil productivity.The main soil quality indicators that led to DSD were soil compaction and displacement in the 10-30 cm portion of the soil profile.The postharvest monitoring requires a minimum 30-point sampling protocol within each logging unit, which can increase in sample size depending on the desired confidence Forests 2023, 14,2218 3 of 25 interval [31,32].On-the-ground sampling efforts take time and resources, particularly when estimating the areal extent of various disturbances across a harvest unit.
Remote sensing is often the best solution for landscape-scale environmental monitoring.Satellite imagery can be used to evaluate disturbance and to monitor change over time, and there are standardized protocols for mapping and monitoring burned landscapes [8,9].There is a need to monitor the disturbed environment after wildfire and salvage logging over both the spatial scale of the salvage site and the temporal scale of the recovery period.
While the Normalized Difference Vegetation Index (NDVI) is conventionally a greenness index (equation in the Methods) [33], its simplicity and utility has resulted in other uses in burned and disturbed landscapes, and generally, as an inverse of exposed soil [34][35][36].NDVI values above 0.4 are generally live vegetation, with lower values ranging from nonphotosynthetic vegetation (NPV) to a mix of soil and NPV, to nearly bare soil in the 0.35 down to 0.1 range.In another study, Ref. [34] evaluated the use of a post-salvage ground-cover treatment adding logging slash to reduce runoff and erosion, and the ability to evaluate it remotely with satellite imagery.The relationship between the amended ground cover and NDVI was significant, and NDVI was used as a proxy for estimating and mapping ground disturbance as it related to erosion potential.They found an NDVI threshold of 0.32 to be a breakpoint between salvaged areas that had sufficient ground cover (generally >60%) to minimize postdisturbance runoff and erosion.
Salvage logging operations were conducted on several fires in the Lolo National Forest (NF) following wildfires in 2016 and 2017.Because of the environmental and legal requirements for logging on public lands, there is substantial documentation available concerning potential detrimental impacts to the disturbed lands.Public comments on the Environmental Impact Statement (EIS) included recommendations for more research to set realistic and consistent thresholds for determining forestry activities that are detrimental to forest productivity and water quality for drinking or for aquatic habitat.Specifically, they call for (1) interpreting results of after-activity "effectiveness" monitoring, including the severity and areal extent of soil disturbance, and (2) assessing the tradeoffs in risks to soil capacity between activities to reduce fuels and wildfire hazard compared to the consequences of wildfire.
We propose to investigate satellite imagery and the use of NDVI alongside groundcover and ground-disturbance measurements to detect the change in bare soil and woodslash cover multiple years after salvage-logging operations.Significant and/or prolonged changes in ground cover will be considered a proxy for estimating forest-floor disturbance.The objectives of this study are to (1) evaluate post-salvage landscapes for detrimental soil disturbance (DSD) 4-5 years after the logging, (2) relate ground-disturbance conditions to very-high-resolution satellite imagery (<2 m), (3) use the relationship between NDVI and ground cover to classify and map soil-disturbance conditions from the pre-fire to current period, and (4) merge these complimentary image and field-analysis components in an innovative methodology that will deliver a protocol and product to identify the areas at greatest risk for DSD following salvage logging.An overarching goal of this project is to aid land managers in minimizing soil degradation after wildfire and salvage logging while still allowing for sustainable and efficient use of the resources within impacted forests and watersheds.

Fire Information
The Copper King Fire burned 12,000 ha, starting in July 2016 on the Lolo NF in Sanders County, MT (Figure 1).Thirty percent of the burned area was classified as moderate and high soil burn severity (burn-severity information for all three fires at: [37]).Soil burn severity was classified from Landsat-derived (30 m pixels) dNBR (differenced Normalized Burn Ratio) values, which were field-validated to reflect the effects of the fire on the soil [9].The Rice Ridge Fire burned 65,000 ha in July 2017 in the Flathead and Lolo NFs (Missoula, Powell, and Lewis and Clark Counties).More than half (56%) of this fire was classified as moderate or high burn severity.The Sunrise Fire burned 10,500 ha, also in July 2017, on the Lolo NF (Mineral County), 34% of which was classified as moderate and high severity.All three fires are geographically dispersed across mountainous western Montana, separated east-to-west by about 160 km.Vegetation regimes in this region are generally mixed conifer, with warmer west and southerly slopes dominated by ponderosa pine (Pinus ponderosa) and Douglas fir (Pseudotsuga menziesii), and cooler or higher elevation slopes dominated by subalpine fir (Abies lasiocarpa), lodgepole pine (Pinus contortus), and grand fir (Abies grandis).Fire regimes in this region are generally mixed severity with a return interval of 5-25 years for low-intensity fire on the dry ponderosa pine forest type and 100+ years for high-intensity stand-replacing fire in the subalpine forest type [38,39].Soils in this region are derived from metasedimentary deposits and are well-drained; volcanic ash-capped soils dominate, and textures are largely silt loams (USDA NRCS 2020; web soil survey online application).
Sanders County, MT (Figure 1).Thirty percent of the burned area was classified as moderate and high soil burn severity (burn-severity information for all three fires at: [37]).Soil burn severity was classified from Landsat-derived (30 m pixels) dNBR (differenced Normalized Burn Ratio) values, which were field-validated to reflect the effects of the fire on the soil [9].The Rice Ridge Fire burned 65,000 ha in July 2017 in the Flathead and Lolo NFs (Missoula, Powell, and Lewis and Clark Counties).More than half (56%) of this fire was classified as moderate or high burn severity.The Sunrise Fire burned 10,500 ha, also in July 2017, on the Lolo NF (Mineral County), 34% of which was classified as moderate and high severity.All three fires are geographically dispersed across mountainous western Montana, separated east-to-west by about 160 km.Vegetation regimes in this region are generally mixed conifer, with warmer west and southerly slopes dominated by ponderosa pine (Pinus ponderosa) and Douglas fir (Pseudotsuga menziesii), and cooler or higher elevation slopes dominated by subalpine fir (Abies lasiocarpa), lodgepole pine (Pinus contortus), and grand fir (Abies grandis).Fire regimes in this region are generally mixed severity with a return interval of 5-25 years for low-intensity fire on the dry ponderosa pine forest type and 100+ years for high-intensity stand-replacing fire in the subalpine forest type [38,39].Soils in this region are derived from metasedimentary deposits and are well-drained; volcanic ash-capped soils dominate, and textures are largely silt loams (USDA NRCS 2020; web soil survey online application).and surface fuels largely consumed, leaving disturbed soils exposed and vulnerable to erosion.Additionally, the risk to human safety, property, and other resources is elevated.Natural resources at risk for all three sites in the immediate fire aftermath include bull trout habitat, native plant communities, and noxious weed invasion.All three fires had some fire-induced soil water repellency due to soil heating, as well as some steep drainages that were deemed to be at risk for hydrologic impacts [37].

Field Sites
Across the three fires, the salvage units that were selected for field sampling had a range of burn severities, and both tractor and skyline logging (Table 1).Four units with ground-based tractor harvesting were sampled on the Copper King Fire.Salvage logging on these units was completed in May of 2018.Elevation in these sites was about 1200 m with average slopes of 25%.Three ground and skyline harvest units were logged on the Rice Ridge Fire in summer 2018.Elevation in these units was 1600-2000 m with slopes of 16%-27%.Three skyline units were salvage logged in summer 2018 after the Sunrise Fire.These units had a mean elevation of 1500 m and ~25% slopes.The logging system and topographical characteristics were compiled from the Lolo NF and the Forest Service Activity Tracking System (online: https://data.fs.usda.gov/geodata/edw/edw_resources/fc/S_USA.Actv_CommonAttribute_PL.gdb.zip;accessed on 22 March 2023).The year and month of logging completion was confirmed visually from Sentinel-2 images. 1 All skyline units were specified as single-span. 2Rubber-tire tractors were specified on the Copper King units but not on the Rice Ridge unit. 3Harvest descriptions in [40].

Field Sampling
We followed the Forest Soil Disturbance Monitoring Protocol [32] to evaluate DSD in 7 logging units on the Copper King and Rice Ridge Fires in the fall of 2022 (3-4.5 years after salvage logging was completed).Randomly located transects were sampled in 7 monitoring units across the 2 burned and salvage-logged areas.Each random transect had 30 points where the soil-disturbance class was estimated, which was used to calculate the percent of the detrimental soil disturbance (DSD) for the unit.Field crews assessed the degree of soil disturbance by examining the deformation of the soil structure and alteration of the A horizon in the upper soil solum, in addition to the amount of alternation of the forest floor [31].The visual clues served as surrogates for quantitative measures of the bulk density and organic matter, which facilitate soil physical, biological, and chemical functions.To assess deep-soil compaction (10-30 cm), field crews used a shovel inserted into the soil to note the depth of compaction.If needed, a slice of soil was removed to check for platy and massive structures, both of which indicate compaction.Surface visual assessments can also serve as a surrogate for compaction within the profile.For example, if the organic horizons were not as deep in the trafficked areas when compared to adjacent undisturbed soil, then it was likely that compaction occurred.Similarly, loss of the organic horizons or surface mineral horizon could indicate soil displacement that can lead to an alteration of soil properties.
If the point along the transect was classified as either disturbance class 2 or 3 (Table 2), the point was considered detrimentally disturbed.Due to time constraints, we used a shorter walk-through assessment collecting information on unit disturbance, site characteristics, and soil pit observations on 3 units on the Sunrise Fire.At all of the units in Table 1, we also sampled random transects to quantify ground cover.Transects generally ran parallel to the slope gradient and were 50 to 400 m in length depending on unit shape; multiple shorter transects were sometimes used to get at least 10 points per unit.Transect endpoint locations were collected with a Trimble GPS unit (Trimble Inc., Westminster, CO, USA).Each transect was sampled at the endpoints and every 10 m; e.g., a 50 m transect had 6 ground-cover points.Ground-cover points were measured using a quadrat (scaled frame) photo-analysis approach in which 4 photos each contained 0.25 m 2 (Figure 2).These photos were converted to the percent of ground-cover classes (e.g., mineral soil, litter, vegetation, gravel, rock, woody debris) using the Cover Monitoring Assistant program (CMA; available at http://resources.nativerevegetation.org/; accessed on 19 April 2023).

Spatial Data Layers
The soil-burn-severity layers based on the dNBR from GTAC (Geospatial Technology Application Center; https://burnseverity.cr.usgs.gov/baer;accessed on 3 March 2023) were used to calculate the majority soil-burn-severity class (low, moderate, or high) of each field point and salvage unit.LANDFIRE layers EVC 2016 (Existing Vegetation Cover; https: //landfire.gov/;accessed on 3 March 2023) and EVT 2016 (Existing Vegetation Type), which were also derived from 30 m Landsat images, were used to approximate the major landcover attributes in the 1-2 years prior to the fires.The SSURGO (Soil Survey Geographic Database; https://websoilsurvey.nrcs.usda.gov/app/;accessed on 3 March 2023) database was used to pull the majority soil information of the salvage units.The United States Geological Survey (USGS) 10 m digital elevation models (DEMs) were downloaded from the National Map Downloader (https://apps.nationalmap.gov/downloader/;accessed on 20 February 2023) for elevation, slope, and aspect data.

Spatial Data Layers
The soil-burn-severity layers based on the dNBR from GTAC (Geospatial Technology Application Center; https://burnseverity.cr.usgs.gov/baer;accessed on 3 March 2023) were used to calculate the majority soil-burn-severity class (low, moderate, or high) of each field point and salvage unit.LANDFIRE layers EVC 2016 (Existing Vegetation Cover; https://landfire.gov/;accessed on 3 March 2023) and EVT 2016 (Existing Vegetation Type), which were also derived from 30 m Landsat images, were used to approximate the major land-cover attributes in the 1-2 years prior to the fires.The SSURGO (Soil Survey Geographic Database; https://websoilsurvey.nrcs.usda.gov/app/;accessed on 3 March 2023) database was used to pull the majority soil information of the salvage units.The United States Geological Survey (USGS) 10 m digital elevation models (DEMs) were downloaded from the National Map Downloader (https://apps.nationalmap.gov/downloader/;accessed on 20 February 2023) for elevation, slope, and aspect data.

Remotely Sensed Imagery
NDVI was calculated by the following equation, where R is the red band and NIR is the near-infrared band: Index values were extracted at all field plot and salvage unit locations, and the mean NDVI for each was calculated for each site and year.
All commercial satellites utilized in this study are owned by MAXAR Technologies, Longmont, CO, USA (https://www.maxar.com/;accessed on 5 August 2023).The WorldView-2 (WV2) satellite has a pixel size of 1.8 m and collects spectral data across eight visible and near-infrared (NIR) bands.A typical skid trail is 4 m wide; therefore, it is likely to be captured in one or two WorldView pixels.GeoEye-1 (GE1) has the same ground resolution as WV2, and has five spectral bands spanning the visible and NIR range.WorldView-3 (WV3) has a higher ground resolution, with 1.24 m pixels and more spectral bands, including shortwave infrared bands.MAXAR's commercial satellites must either be tasked (requested) to acquire imagery over the area of interest (e.g., a wildfire), or the archive can be searched (https://discover.maxar.com/;accessed on 9 March 2023) for previous collections, which are available for purchase or free for US federal government employees via an existing contract.Images are downloadable or delivered ready for analysis.We used a combination of pre-collected and archived images along with one new collection per fire.For the archived data, we opportunistically selected the best available images with dates that aligned annually and phenologically with dates of ground-cover measurements taken in fall 2022 (Table 3).The goal was to assemble a range of images for each site: pre-fire, post-fire, and post-salvage image(s).The tasked new collection in September 2022 coincided with the field ground measurements.Images were delivered georeferenced, radiometric-and sensor-corrected, and ready for orthorectification.Residual geometric error was corrected by automatic orthorectification in ENVI (ENVI 4.4, ITT Visual Information Solutions, Boulder, CO, USA).The resulting root mean square error (RMSE) of the position for each geometric correction was kept below one pixel.A dark-object subtraction [41] was applied to each image to normalize the atmospheric effects between annual images.
2.6.Statistics and Data Analysis 2.6.1.Overview Statistical software programs R (Ver.4.3.1)[42] and SAS (Ver.9.4; SAS Institute, Cary, NC, USA) were used to perform statistical analyses.Model diagnostics included testing residuals for normal distribution, plotting residuals versus predicted values, plotting residuals versus order, and quantile-quantile plots of residuals versus predicted values [43] to ensure model assumptions were met.A combination of parametric (correlation and mixed model) and nonparametric (classification tree) analyses were run to meet multiple objectives.
The correlation analysis was run as an exploratory exercise to evaluate the relationship between DSD and other field measures and NDVI to identify important variables moving forward.The parametric mixed model provided insight on variables' contribution to the post-salvage NDVI values, in particular the differences between the factor variables (e.g., harvest method and harvest season).The nonparametric classification tree was primarily used to identify breakpoints in NDVI values associated with the levels of DSD and to provide classification thresholds.Some analyses were done with the individually sampled field points (n = 29 or 39), while other time values were aggregated to the salvage unit means (n = 7 or 10); the higher or lower n indicates whether or not DSD values were recorded at the field points; DSD was not explicitly measured at the Sunrise Fire.Other analyses were done with the full dataset of all units that were salvaged after the three fires (n ~300); it is specified in each Results section and in each figure in which the dataset was used.

Model Details
In SAS, a linear mixed-effects model [44] was developed for the response variable post-salvage NDVI using the full dataset of unit means.Fire, equipment, month, year, SBS class, aspect, slope, soil type, EVC and EVT tested as fixed effects, and site as the random effect.As fixed effects were found insignificant (p > 0.05), they were iteratively removed from the model.For each response, an autoregressive correlation structure was applied.
The significance of differences across multiple comparisons among fixed effects were tested using least-squares means.
A classification tree analysis was performed to evaluate the relationship between field-measured variables (field points only) and their contribution to the determination of DSD using the R software (ver.4.3.1)(R Development Core Team, 2022) package 'rpart' [45].Classification-tree methods are nonparametric and function as recursive partitioning regressions.They operate by partitioning variables into increasingly homogeneous response-variable classes (e.g., detrimental soil disturbance) or breakpoints (e.g., NDVI values).The weight of each predictor variable is calculated, and each predictor variable is used in order of significance to divide the data [46].All data at a terminal node are assigned the same class of disturbance.The classification tree was first built using all available variables and narrowed down as variables proved unimportant.Variables that did not meaningfully improve the division of the response variable into classes or breakpoints were eliminated until only significant variables remained.

Ground Cover
Total ground cover included logging slash, litter, and vegetation; total ground cover was inversely correlated with bare mineral soil, which is the primary factor affecting postdisturbance runoff and soil erosion [47,48].
An average of 12% soil (soil + rock + gravel) cover was reported on all sampled plots regardless of prior disturbance, salvage logging, or equipment method, etc.The range of soil was also the smallest compared to the other cover categories, ranging from 0% to 25% (Figure 3, left).There was more nonphotosynthetic vegetation (NPV) in the plots than live vegetation; however, this was likely due to the timing of the field sampling in early fall.Organic ground cover of at least 60% is generally considered to be the threshold for elevated runoff or erosion events [49][50][51].In the seven units sampled for DSD, the percent of the area that was still found to have detrimental disturbance 2-4 years after salvage logging ranged from 1% to 20% (Table 4).This was a significant finding considering the likely higher degree of disturbance immediately after the logging disturbance.The three units with more than 10% DSD remaining were all primarily ground-based tractor-logged.There was a strong correlation between the exposed soil that was measured via the photo method and DSD that was measured using other methods (Figure 3, right).Many other factors besides the percent of soil exposed are taken into consideration when DSD is calculated.Table 4. Unit means of ground cover and detrimental soil disturbance (DSD) from the field-sampled units; DSD was not explicitly measured at the Sunrise Fire.Mean NDVI values were calculated and In the seven units sampled for DSD, the percent of the area that was still found to have detrimental disturbance 2-4 years after salvage logging ranged from 1% to 20% (Table 4).This was a significant finding considering the likely higher degree of disturbance immediately after the logging disturbance.The three units with more than 10% DSD remaining were all primarily ground-based tractor-logged.There was a strong correlation between the exposed soil that was measured via the photo method and DSD that was measured using other methods (Figure 3, right).Many other factors besides the percent of soil exposed are taken into consideration when DSD is calculated.Table 4. Unit means of ground cover and detrimental soil disturbance (DSD) from the field-sampled units; DSD was not explicitly measured at the Sunrise Fire.Mean NDVI values were calculated and extracted from each salvage unit using Sentinel-2 satellite imagery.Nonphotosynthetic vegetation (NPV).

Imagery
Four sets of imagery were evaluated for each fire (Table 3; Figure 4).The dates of the pre-fire images ranged from 2014 to 2017, and one post-fire (but prior to salvage) image was also identified for each fire, and the dates of those images ranged from 2016 to 2018.The post-fire to post-salvage date range is more difficult to identify with certainty due to the timing of salvage logging (i.e., the salvage completion dates are available but not the start dates).On the ground-based units in particular, it was fairly easy to determine whether or not salvage logging had been initiated from the 2 m WorldView imagery (Figure 4).Therefore, we feel confident that we have identified post-salvage imagery for each of the units that we sampled, regardless of the dates in the tables.
NDVI values from the pre-fire images were similar, ranging from 0.6 to 0.7, indicating green vegetation typical of a mixed-conifer forest.Post-fire values dropped significantly to ~0.25 across all fires (Table 4; Figure 5).Interestingly, the post-salvage values had the widest range, both lower and higher than the post-fire values.By 2022, the sites had begun to homogenize again, with less variability in NDVI (~0.4), and values were increasing towards pre-fire levels.
The relationship between the post-salvage NDVI and percent of DSD was significant (Figure 6).As DSD increased, NDVI values decreased steeply.For contrast, the 2022 NDVI values were also plotted against DSD, and this emphasized the degree of vegetation recovery and homogenization of NDVI values in the 2-3 years between the post-salvage and 2022 imagery.
The correlation analysis (Table 5) highlights several significant relationships between the ground-collected and image (NDVI) variables.Soil burn severity was negatively correlated with NPV cover and post-fire and post-salvage NDVI.As the burn severity increased, these other variables decreased.Soil cover was positively correlated with DSD; more soil and rock exposure indicated a higher degree of detrimental soil disturbance.The NPV was only positively correlated with the post-salvage NDVI; there was likely a high volume of timber slash and other woody debris in the post-salvage images.DSD was strongly negatively correlated with both the post-fire and post-salvage NDVI (r = −0.62,−0.73 (respectively); p < 0.0001).Collectively, these results indicate the capability to use NDVI to map DSD as it relates to more exposed soil, as well as the presence of NPV, which may indicate logging activity.
whether or not salvage logging had been initiated from the 2 m WorldView imagery (Figure 4).Therefore, we feel confident that we have identified post-salvage imagery for each of the units that we sampled, regardless of the dates in the tables.
NDVI values from the pre-fire images were similar, ranging from 0.6 to 0.7, indicating green vegetation typical of a mixed-conifer forest.Post-fire values dropped significantly to ~0.25 across all fires (Table 4; Figure 5).Interestingly, the post-salvage values had the widest range, both lower and higher than the post-fire values.By 2022, the sites had begun to homogenize again, with less variability in NDVI (~0.4), and values were increasing towards pre-fire levels.The relationship between the post-salvage NDVI and percent of DSD was significant (Figure 6).As DSD increased, NDVI values decreased steeply.For contrast, the 2022 NDVI values were also plotted against DSD, and this emphasized the degree of vegetation recovery and homogenization of NDVI values in the 2-3 years between the post-salvage and 2022 imagery.The relationship between the post-salvage NDVI and percent of DSD was significant (Figure 6).As DSD increased, NDVI values decreased steeply.For contrast, the 2022 NDVI values were also plotted against DSD, and this emphasized the degree of vegetation recovery and homogenization of NDVI values in the 2-3 years between the post-salvage and 2022 imagery.The correlation analysis (Table 5) highlights several significant relationships between the ground-collected and image (NDVI) variables.Soil burn severity was negatively correlated with NPV cover and post-fire and post-salvage NDVI.As the burn severity increased, these other variables decreased.Soil cover was positively correlated with DSD;

Classification
We used regressive partitioning trees on the individual field-point data to attempt to determine breakpoints in the data that indicated a high/moderate/low degree (or high/low) of disturbance.Again, using the post-salvage (PS_NDVI) values in the partitioning trees, two breakpoints were clearly identified (Figure 7).NDVI values greater than 0.32 had a low proportion of DSD, while values lower than 0.22 had a high proportion of DSD (Figure 8).Moderate DSD was between 0.22 and 0.32 NDVI.We were unable to extrapolate this analysis to the greater dataset (all salvaged units) because we only had DSD values for the points that were field-sampled.
We repeated the partitioning tree analysis with the 2022 NDVI (NDVI_22) data to see how NDVI thresholds might change over time.The breakpoint for the low DSD increased to 0.38, while the breakpoint for the high DSD increased to 0.31.This is intuitive, as all sites have had vegetation recovery in the 2-3 years between the two sets of images.There was no natural separation for a 'moderate' DSD class using the 2022 NDVI values.This also makes sense, as the sites are likely becoming more similar as the vegetation recovers.
DSD values for the points that were field-sampled.
We repeated the partitioning tree analysis with the 2022 NDVI (NDVI_22) data to see how NDVI thresholds might change over time.The breakpoint for the low DSD increased to 0.38, while the breakpoint for the high DSD increased to 0.31.This is intuitive, as all sites have had vegetation recovery in the 2-3 years between the two sets of images.There was no natural separation for a 'moderate' DSD class using the 2022 NDVI values.This also makes sense, as the sites are likely becoming more similar as the vegetation recovers.

Mixed-Model Results
A mixed model was used to evaluate the influence of several variables on the postsalvage NDVI values.All available topographical, vegetation, soil, logging equipment, soil burn severity, and logging completion date variables were initially considered.Nonsignificant variables were iteratively eliminated until significant and meaningful variables remained (Table 6).

Mixed-Model Results
A mixed model was used to evaluate the influence of several variables on the postsalvage NDVI values.All available topographical, vegetation, soil, logging equipment, soil burn severity, and logging completion date variables were initially considered.Nonsignificant variables were iteratively eliminated until significant and meaningful variables remained (Table 6).Factor variables (fire, equipment, season, SBS) were evaluated for their relative differences by class.Tractor logging had a significantly lower NDVI (0.21) compared to skyline logging (0.25).The salvaged units on the Sunrise Fire had the lowest post-salvage NDVI (0.12), followed by Copper King (0.22) and Rice Ridge (0.36).The units classified as low SBS had the highest post-salvage NDVI (0.27) compared to the moderate-and high-SBS units (0.22, 0.21; respectively).Of the continuous variables, both elevation and slope were significant, but compared to the other variables had relatively minor effects on the model.
The season variable was determined by the month that the salvage logging was completed, and there was a significant effect on the post-salvage NDVI values (Figure 8).When salvage logging was completed in winter months (December-March), NDVI was higher.NDVI values were significantly lower in early summer or growing-season months (April-July), and the lowest NDVI values occurred when logging was completed in autumn months (August-November).This indicates less disturbance when logging was done over snow, more disturbance during the growing season, and the greatest disturbance during the fall, when soils may be occasionally wet and susceptible to rutting and compaction.

Imagery
Using the ranges of NDVI values that indicated the post-salvage DSD from the thresholding analysis, we classified NDVI values and overlaid them on true-color imagery (Figure 9).The skid trails are highly visible and show the network of traffic patterns from the salvage logging.Areas of greater disturbance (red/yellow pixels) likely had more traffic passes or a greater concentration of activity.Areas of relatively less disturbance (yet still greater than the surrounding undisturbed land) are green/blue pixels and likely had fewer traffic passes.
These images can be useful for managers, as they highlight the areas that were detrimentally disturbed, as well as provide some framework for the degree of disturbance of the whole unit.The skid-trail networks can also be digitized for further analysis for road-network connectivity and disturbance, as it pertains to greater watershed disturbance after wildfires.In order to quantify the percent of each of these units that was likely detrimentally disturbed, we plotted NDVI values on histograms with bins (Figure 10) closely corresponding to the values in the legends in Figure 9. NDVI values less than  The thresholding analysis (Figure 7) found higher NDVI values to indicate likely DSD in an image that was several years postdisturbance compared to shortly after the event (0.31 and 0.22, respectively).The 2022 Copper King imagery was overlaid with classified NDVI values thresholded at the higher limits to estimate the areal disturbance of the salvage units (Figure 11).Using 2022 imagery, the breakpoint for DSD is ~0.31, giving 10%-25% DSD in these four units (U4 18%; U12 26%; U13 10%; U37 25%).While these values may seem high for 4+ years post-salvage disturbance, we measured 3%-13% exposed soil and 3%-20% DSD at the field plots (Table 4, Figure 2); thus, the field data correspond with the image classifications.The thresholding analysis (Figure 7) found higher NDVI values to indicate likely DSD in an image that was several years postdisturbance compared to shortly after the event (0.31 and 0.22, respectively).The 2022 Copper King imagery was overlaid with classified NDVI values thresholded at the higher limits to estimate the areal disturbance of the salvage units (Figure 11).Using 2022 imagery, the breakpoint for DSD is ~0.31, giving 10%-25% DSD in these four units (U4 18%; U12 26%; U13 10%; U37 25%).While these values may seem high for 4+ years post-salvage disturbance, we measured 3%-13% We were able to observe the spread and degree of disturbance in the salvage units on all three fires (Figures 9, 11 and 12; Rice Ridge, Copper King, and Sunrise, respectively) using the post-salvage classified NDVI values.Somewhat varying NDVI values indicate different degrees of disturbance on all three fires, and that can primarily be attributed to the time between the fire and salvage logging, and the image that was used to evaluate disturbance.The time between salvage-logging completion and the imagery was ~1 month for the Rice Ridge and Sunrise Fires, and NDVI values at or below 0.21 highlighted ground disturbance directly due to logging operations (e.g., skid trails, patches of newly exposed soil) (Figures 9 and 12).While on the Copper King Fire, more than 4 years passed between the salvage and the high-resolution image we acquired; consequently, more revegetation had occurred, and higher NDVI values (~0.3) indicated remaining soil disturbance (Figure 11).exposed soil and 3%-20% DSD at the field plots (Table 4, Figure 2); thus, the field data correspond with the image classifications.We were able to observe the spread and degree of disturbance in the salvage units on all three fires (Figures 9, 11, and 12; Rice Ridge, Copper King, and Sunrise, respectively) using the post-salvage classified NDVI values.Somewhat varying NDVI values indicate different degrees of disturbance on all three fires, and that can primarily be attributed to the time between the fire and salvage logging, and the image that was used to evaluate disturbance.The time between salvage-logging completion and the imagery was ~1 month for the Rice Ridge and Sunrise Fires, and NDVI values at or below 0.21 highlighted ground disturbance directly due to logging operations (e.g., skid trails, patches of newly exposed soil) (Figures 9 and 12).While on the Copper King Fire, more than 4 years passed between the salvage and the high-resolution image we acquired; consequently, more revegetation The results from the analyses are summarized along with management considerations in Table 7.Our goal was to provide functional guidance for post-fire-salvage-logging planning.Mindful planning considering the topography, background soil burn severity, time since the fire, and logging season can all significantly lessen DSD.Although immediate ground-cover measurements were not taken, we often found 0%-30% ground cover after high-severity fires in the immediate post-fire period [14,52,53], with salvage logging potentially decreasing that value even further [20,34,54].On these fires, the wildfire appears to have reduced cover by at least half in our study sites, which we can estimate from NDVI values (Figure 5).Sites were salvage-logged 2-3 years after the fire, and ground-cover values showed high variability in this period; some sites had less cover post-salvage than post-fire due to tree removal and soil exposure on the skid trails, while other post-salvage sites have more cover, likely attributed to vegetation regrowth, logging slash debris, and the general 'recovery' in the postdisturbance period.Moreover, salvage impacts are highly localized even within salvage units; areas which are impacted are often highly impacted, while other areas are untouched by the salvage operation.By 2022, all units were at least 5 years post-fire and at least 3 years post-salvage, and cover had increased to 80%-90% on most units, with the Sunrise units having the least total ground cover at ~75% (Table 4), which is likely due to a site factor, such as soil depth or localized rainfall.The values found on these sites are typical for a mixed-conifer forest in the western US, as Ref. [34] found total vegetation cover on skid-trail plots at a site in eastern Washington at least doubled between the first and second post-fire years.Another study in the same region found vegetation recovery was slow after delayed salvage logging (3-4 years post-fire), particularly on sites burned at high or moderate burn severity [14], with exposed soil still high at ~50-75%.Although the percent of soil exposed is only a proxy for DSD, there was a significant relationship between DSD measurements, the percent of exposed soil, and NDVI values (Table 5; Figures 3 and 6).

Using Imagery to Map DSD
Matching the spatial scale of remotely sensed imagery to the scale of the element or process of interest is an essential part of a successful mapping campaign.Wildfire disturbances are often mapped at the 10-30 m (pixel) scale because of the frequently available Sentinel-2 or Landsat satellites, which have appropriate bands for detecting changes in the vegetation and soil conditions due to the wildfire [55].The on-the-ground heterogeneity of a wildfire disturbance can certainly be finer than 10 or 30 m, but for a full-fire burn severity map, this type of imagery is the standard protocol [8,56,57].For this study, however, a typical logging skid trail is 2-4 m wide, which makes very-high-resolution imagery, such as WorldView, better suited for evaluating fine-scale soil disturbances.Currently, WorldView and other commercially available imagery must be ordered, but it is becoming more widely used and available.To estimate salvage logging completion dates, we used 10 m Sentinel-2 imagery and found that, in highly disturbed-ground salvage units (such as Rice Ridge Unit 26, see Figure 4), we were able to identify skid trails easily with true-color and NDVI images.On units that were less disturbed (e.g., Copper King Units 4 and 37; see Figure 11) and did not have obvious skid trails, it was more difficult to determine when salvage had occurred.The scale of the disturbance was finer than what was captured with 10 m pixels.
We applied NDVI as a disturbance indicator or as proxy for the percentage of DSD in a salvage unit.We were able to estimate the areal percentage of disturbance in each salvage unit using classified NDVI values.Skid trails, heavily trafficked roads and intersections, and log landings are the major contributors of disturbed soil after a salvage operation [13].The patterns of trail networks and landings that are the result of logging activities are visible in imagery and are enhanced by classified NDVI values (Figures 4c and 9b).Best management practice guidelines aim for less than 15% DSD in a unit.We found values ranging from about 5% to 25% in various units, some of which are higher than desirable.However, we feel that our estimates are likely within the range of variability of what would actually be measured on the ground, and feel comfortable recommending that managers use the estimates from classified NDVI values for post-salvage disturbance estimates.Researchers found that both logging methods and timber operators contribute to the degree of disturbance after logging, and there can be a need for quality control [25].An areal estimate of the disturbance from satellite imagery presents an efficient proxy to estimate DSD to supplement soil field assays.Field assessments help to clarify and ground-truth the image estimates according to the soil type and ecological recovery.A complaint in the public record from the Lolo NF salvage sales was that only 7 units (out of 300+) were evaluated for DSD [30], which were different than the 7 units that were evaluated as part of this study.A remote assessment would allow a much greater percentage of sites to have a consistent evaluation.
The degree of disturbance estimated is strongly related to the time since the disturbance has passed, and we found that we were able to map disturbance patterns (skid trails and patches) even 4+ years after logging on the Copper King Fire (Figure 11).However, the rate at which revegetation occurs is often rapid after the first postdisturbance growing season [52,53], and is likely captured by one or two images in the proceeding year(s).Thus, there is no critical need in most cases for weekly or monthly images, like one could get with either Landsat or Sentinel time series, or at even higher temporal resolution with Harmonized Landsat-Sentinel time series [55].For both spatial and temporal considerations, we feel that very-high-resolution imagery (<5 m) is most appropriate for mapping logging skid trails and post-salvage disturbance.With a modest amount of effort, we were able to acquire very-high-resolution imagery at sufficient periodic timesteps to capture the change in the forest environment due to wildfire and salvage logging at the desired spatial and temporal resolution.
On the Lolo NF, presalvage soil monitoring found around 2%-4% DSD in potential logging units, primarily from existing jammer road systems [30].Because it may not be possible to return a disturbed system fully to its pre-disturbance state, the goal is often to rehabilitate along a trajectory that aims to come as close to the pre-disturbance state as possible.The same Lolo NF Monitoring Report [30] found 25% or more DSD after an initial assessment of high-burn-severity units on the Copper King Fire.Other Lolo units (not part of this study) averaged around 7.5% DSD, which is well below the BMP limit of 15%.This range of DSD values across sites is not surprising, as many factors, such as soil burn severity (i.e., prior extreme disturbance), slope, equipment, intensity, and harvest season, all play a role in the degree of disturbance from logging activities [16,22,27].The range of values emphasizes the need for monitoring at a broad scale, which is much more feasible with remote sensing techniques rather than relying on field visits alone.
The factors that we found to be influential for DSD after salvage logging broadly fell under three categories: topographical, prior disturbance, and logging activity.To-pographically, higher slopes and higher elevations both lead to an increased chance of DSD (Tables 6 and 7).There are BMPs for slope, but not for elevation.A review by [28] indicated that the greatest disturbance on machine-operated skid trails occurred on slopes greater than 20%.The authors of [26] found that bulk density and rut depth increased on steeper slopes due to unbalanced load distribution.These results agreed with what we found (Table 7), i.e., that there is an increased risk of DSD at slopes greater than 25%, and even more so on slopes greater than 35%.We found that, as elevations increased, so did DSD, and this is likely the result of soil types and vegetation density differing by elevation.Lower vegetation density, lower vegetative ground cover, and rocky soils at higher elevations could all be elements signaling the potential for greater disturbance after salvage logging.Similarly, lower elevations than the study area could prove to have slower recovery due to less vegetation vigor and overall production.
Prior disturbance for salvage logging is primarily soil burn severity, which is specific to the effects of the fire on the ground surface rather than the canopy [8,9], but could also include land use such as grazing or prior logging activity [58].High soil burn severity in particular indicates severe effects of the fire on the soils, which can mimic and compound the disturbance from logging: the removal of surface organics, reduced infiltration, increased risk of runoff, and erosion [9,22].Other studies have found logging traffic increased the soil bulk density by 25% up to a soil depth of 10 cm [20,34,59,60].In addition to soil disturbance leading to a potential hydrologic response, long-term tree recovery and stand performance may also be compromised [21,61].We calculated the dNDVI (pre-NDVI-post-NDVI) for both post-fire and post-salvage NDVI from the values in Table 4.The mean dNDVI for the post-fire timeframe was 0.39, while the post-salvage dNDVI was 0.40.This was interesting because it hinted at the possibility of an additional impact from the secondary disturbance of salvage logging, or at least showed little or no recovery between fire and salvage.Postsalvage dNDVI was more strongly correlated (r = 0.71) to DSD than the post-fire dNDVI (0.61) on the plots where DSD was measured.While it is difficult to determine from these values alone the proportion of disturbance that can be attributed to either fire or salvage, or the accumulation of both, it is clear that NDVI is capable of identifying areas of varying degrees of disturbance.A decrease in NDVI values can be attributed to the destruction of vegetation (increase in red values) as well as more exposed soil (decrease in NIR values).Fire and logging can impact both vegetation and soils in the perturbed areas.
We use the term 'logging activity' broadly to incorporate the seasonal timing of the logging operation, equipment type, and the intensity or impact of the logging on the unit.Seasonally, we found that winter logging (December-March; Figure 8) had the least impact on soils.These results are generally agreed upon in the literature, as frozen soils are more resistant to compaction and rutting [61,62].Reeves et al. [25] also proposed that seasonality was a main factor, and included the interaction of the harvest season × land type as a key variable in their model, along with slope and aspect.We evaluated two equipment types in this study, skyline (cable suspension of logs) and tractor (log skidders) [19].The result was that the tractor-skidder units were significantly lower in NDVI values and higher in DSD values (Table 6).This was an intuitive result, but also advantageous when supported by data and statistics, as well as other studies [22].Reeves et al. [25] discussed equipment operator skill and sale administrator knowledge on local conditions having an impact on the degree of soil disturbance during the logging activities.Equipment operators' awareness and understanding of the need for adhering to BMPs likely plays a role in the disturbance outcome.

Management Considerations
Identifying areas more susceptible to high impact and greater disturbance from harvest activities allows land managers to develop alternative strategies to meet management objectives and can help prioritize the allocation of monitoring resources.Decision tools, such as the one in Table 7, that incorporate geospatial layers like NDVI, which indicates levels of disturbance, can play a crucial role in project planning and are a key component of adaptive management strategies.All of the factors that we have identified as influential for predicting detrimental soil disturbance are somewhat ameliorated with planning.Naghdi [26] suggested that the adverse effects of salvage logging traffic can be partially mitigated with adequate preparation, including spatially pre-planning machine-tracked trails.Particular attention should be paid to minimizing traffic passes and eliminating skidding on steep (>25%-35%) slopes.When possible, logging on wet or severely disturbed soils should be avoided.If one or more of the conditions in Table 7 cannot be avoided, then minimizing the number of confounding disturbance factors should be considered.
To include a real-life example from the Lolo NF, the timber industry asked for modifications from the Lolo Forest Plan in order to harvest more lands; for example, harvesting on slopes >35%-50% [30].The leadership team responded in a way that allowed some variances to the BMPs with appropriate monitoring to allow for additional timber harvest to meet other objectives (V.Archer, personal communication).We are not suggesting rigid rules for harvest, but rather a set of factors to consider when planning in order to meet a wide range of land-management objectives.Our methods presented here provide an operational, repeatable process for monitoring soil disturbance after salvage logging over many years, particularly in the most sensitive 1-3 years following the salvage disturbance.Because NDVI is sensitive to changes in both green vegetation and soil exposure, it would be appropriate to apply these methods to prescribed fires as well as fuel treatments (e.g., thinning), which are less disruptive to the impacted environment.NDVI is capable of detecting minor disturbances [36] and should be an efficient method for other types of disturbance mapping.Limitations to the application might include situations where the impact to the soil is much higher than that of the canopy (occlusion), or a situation where vegetation is affected by disease or pest, which might have zero associated soil disturbance.Prior knowledge or historical images of the study area would help inform disturbance history and increase the likelihood of success from the methods presented here.The decision table and GIS procedure presented here move towards making both planning and monitoring decisions easier and more accessible for land managers.There is a need for the standardization of protocols and thresholds to apply technology across multiple landscapes, and the methods presented here should be replicated on other fire and logging environments for evaluation.

Conclusions
Our primary goal was to develop an operational methodology for remotely evaluating detrimental soil disturbance (DSD) after post-fire salvage logging.We found strong correlations between exposed mineral soil and DSD measured in the field, and subsequently found strong correlations between post-salvage NDVI values and DSD.Three to four years post-salvage, we found a range of exposed soil of 5%-25% and DSD of 3%-20% across low-, moderate-, and high-soil-burn-severity units.NDVI values from the 2 m WorldView-2 imagery <0.22 and <0.31 soon after salvage logging and 3 years post-salvage, respectively, indicated high DSD.We modeled DSD based on NDVI and other logistical logging factors (topography, season, soil burn severity, equipment), and made management recommendations for minimizing DSD from salvage logging.While our general management recommendations are not exceedingly novel, the additional geospatial analysis and threshold NDVI values provide a repeatable and scalable method for mapping and monitoring logging disturbance over large areas and many years.Our management considerations table is formatted as a decision-making tool to deliver functional guidance for post-fire salvage-logging planning.Mindful planning considering the topography, soil burn severity, time since fire, and logging season can all significantly lessen DSD.The classified NDVI maps can be used to evaluate salvage disturbance, delineate existing skid trails, and eventually incorporate skid trails and soil-disturbance information into watershed-modeling tools.

Figure 1 .
Figure 1.Locations of the three fires sampled.Soil-burn-severity maps from classified and fieldvalidated dNBR maps (Landsat, 30 m pixels), the salvage sale units, and the field units are shown for each fire.About one-third of the Copper King and Sunrise Fires, and just over half the Rice Ridge Fire, were classified as moderate and high soil burn severity by post-fire assessment teams [9].Most trees are killed under high-burn-severity conditions, with understory

Figure 2 .
Figure 2. Ground cover photos from Copper King Unit 12 (top row), Sunrise Unit 22 (bottom right), and Rice Ridge Unit 26 (bottom left) in September 2022.These photos represent typical groundcover conditions from the most disturbed unit we sampled on each fire; bare soil ranges from 6% to 60%.The time since salvage logging was more than four years for the Copper King unit, and three years for the Rice Ridge and Sunrise units.

Figure 2 .
Figure 2. Ground cover photos from Copper King Unit 12 (top row), Sunrise Unit 22 (bottom right), and Rice Ridge Unit 26 (bottom left) in September 2022.These photos represent typical ground-cover conditions from the most disturbed unit we sampled on each fire; bare soil ranges from 6% to 60%.The time since salvage logging was more than four years for the Copper King unit, and three years for the Rice Ridge and Sunrise units.

2. 5 .
Remotely Sensed Imagery NDVI was calculated by the following equation, where R is the red band and NIR is the near-infrared band: NDVI = (NIR − R)/(NIR + R)

Figure 3 .
Figure 3. Box and whisker plot of the unit values of the ground-cover categories from the 2022 field photo analysis (left; n = 10), and a scatterplot comparing the percent of exposed soil and the percen of DSD measured in fall 2022 (right; n = 7).On the boxplots, minimum and maximum values are represented by the whiskers; the first and third quartiles are shaded; the median is in the shaded box.

Figure 3 .
Figure 3. Box and whisker plot of the unit values of the ground-cover categories from the 2022 field photo analysis (left; n = 10), and a scatterplot comparing the percent of exposed soil and the percent of DSD measured in fall 2022 (right; n = 7).On the boxplots, minimum and maximum values are represented by the whiskers; the first and third quartiles are shaded; the median is in the shaded box.

Figure 4 .
Figure 4. Rice Ridge Fire Unit 26; ground-cover field-plot locations are overlaid in the figure.Image years are (a) 2014 (pre-fire), (b) 2018 (post-fire), (c) 2019 (post-salvage), and (d) 2022.The reported completion date of salvage logging was August 2019; the 2022 image shows about 3 years of recovery.

Figure 4 .
Figure 4. Rice Ridge Fire Unit 26; ground-cover field-plot locations are overlaid in the figure.Image years are (a) 2014 (pre-fire), (b) 2018 (post-fire), (c) 2019 (post-salvage), and (d) 2022.The reported completion date of salvage logging was August 2019; the 2022 image shows about 3 years of recovery.

Forests 2023 , 27 Figure 5 .
Figure 5. Box and whisker plot showing field unit means (n = 10) of NDVI values over time and across the fire and salvage logging disturbances.Pre-fire NDVI (Pre_NDVI), post-fire NDVI (PF_NDVI), post-salvage (PS_NDVI), and 2022 NDVI (NDVI_22) image dates are in Table 3.On the boxplots, minimum and maximum values are represented by the whiskers; the first and third quartiles are shaded; the median is in the shaded box.

Figure 5 .
Figure 5. Box and whisker plot showing field unit means (n = 10) of NDVI values over time and across the fire and salvage logging disturbances.Pre-fire NDVI (Pre_NDVI), post-fire NDVI (PF_NDVI), post-salvage (PS_NDVI), and 2022 NDVI (NDVI_22) image dates are in Table 3.On the boxplots, minimum and maximum values are represented by the whiskers; the first and third quartiles are shaded; the median is in the shaded box.
quartiles are shaded; the median is in the shaded box.

Figure 6 .
Figure 6.Scatterplot of NDVI values compared to ground-measured percent of detrimental-soildisturbance (DSD) values at the field units (n = 7).Post-salvage and 2022 NDVI values are plotted.

Figure 7 .
Figure 7. Regressive partitioning trees from the rPart package in R. The tree on the left is thresholding the field points (n = 29) into DSD classes using the post-salvage NDVI data (PS_NDVI), while the tree on the right is using 2022 NDVI data (NDVI_22).

Forests 2023 , 27 Figure 7 .
Figure 7. Regressive partitioning trees from the rPart package in R. The tree on the left is thresholding the field points (n = 29) into DSD classes using the post-salvage NDVI data (PS_NDVI), while the tree on the right is using 2022 NDVI data (NDVI_22).

Figure 8 .
Figure 8. Mixed-model results showing the difference in NDVI values by the season salvage logging was completed.Completion dates for each field unit are listed in Table 1 and are generally in the fall; however, this figure was created from the full dataset (n ~ 300), where salvage unit completion dates were more varied.Winter (December-March); early Summer (April-July); Autumn (August-November).The dashed horizontal line at NDVI = 0.22 is a breakpoint below which indicates detrimental soil disturbance from the classification tree analysis in Figure 7. Minimum and maximum values are represented by the whiskers; the first and third quartiles are shaded; the median is in the shaded box.Significant differences (p < 0.05) between NDVI values by season are indicated by capital letters A and B.

Figure 8 .
Figure 8. Mixed-model results showing the difference in NDVI values by the season salvage logging was completed.Completion dates for each field unit are listed in Table 1 and are generally in the fall; however, this figure was created from the full dataset (n ~300), where salvage unit completion dates were more varied.Winter (December-March); early Summer (April-July); Autumn (August-November).The dashed horizontal line at NDVI = 0.22 is a breakpoint below which indicates detrimental soil disturbance from the classification tree analysis in Figure 7. Minimum and maximum values are represented by the whiskers; the first and third quartiles are shaded; the median is in the shaded box.Significant differences (p < 0.05) between NDVI values by season are indicated by capital letters A and B.

0. 22 27 Figure 9 .
Figure 9. Classified NDVI values from Rice Ridge Fire imagery in the identified range of high to low DSD, which highlight the skid trails from salvage logging.Panels are Rice Ridge (a) Unit 20, (b) Unit 26, and (c) Unit 61.

Figure 9 .
Figure 9. Classified NDVI values from Rice Ridge Fire imagery in the identified range of high to low DSD, which highlight the skid trails from salvage logging.Panels are Rice Ridge (a) Unit 20, (b) Unit 26, and (c) Unit 61.

Figure 10 .
Figure 10.Histograms showing the distribution of pixels from the 2019 NDVI-classified images from the Rice Ridge (RR) Fire in Figure 9; (a) Unit 20, (b) Unit 26, and (c) Unit 61.The frequency × 100% gives the percentage of each of the classes in the units, with NDVI values below 0.22 indicating potential detrimental soil disturbance (DSD).

Figure 10 .
Figure 10.Histograms showing the distribution of pixels from the 2019 NDVI-classified images from the Rice Ridge (RR) Fire in Figure 9; (a) Unit 20, (b) Unit 26, and (c) Unit 61.The frequency × 100% gives the percentage of each of the classes in the units, with NDVI values below 0.22 indicating potential detrimental soil disturbance (DSD).

Figure 11 .
Figure 11.Classified NDVI images from the Copper King Fire using 2022 imagery.Field units: (a) Unit 4; (b) Units 12, 13; (c) Unit 37.These images represent approximately 4 years of recovery after the salvage operations.

Figure 11 .
Figure 11.Classified NDVI images from the Copper King Fire using 2022 imagery.Field units: (a) Unit 4; (b) Units 12, 13; (c) Unit 37.These images represent approximately 4 years of recovery after the salvage operations.

Figure 12 .
Figure 12.True color imagery from the Sunrise Fire.(a,b) are 2018 images approximately one year after the fire and only one-month post-salvage; classified NDVI values are overlaid on panel (b), showing widespread post-salvage disturbance; (c) shows revegetation in the 2022 image.

Figure 12 .
Figure 12.True color imagery from the Sunrise Fire.(a,b) are 2018 images approximately one year after the fire and only one-month post-salvage; classified NDVI values are overlaid on panel (b), showing widespread post-salvage disturbance; (c) shows revegetation in the 2022 image.

Table 3 .
Satellite 1 Some post-fire images are potentially midsalvage harvest depending on the timing.

Table 6 .
Results from the mixed model of the form: Post-salvage NDVI = Fire • Equipment • Month • SBS • Slope • Elevation.SBS: soil burn severity.The model includes all tractor and skyline units (n = 318) with complete datasets (no missing values) for all three fires.

Table 6 .
Results from the mixed model of the form: Post-salvage NDVI = Fire • Equipment • Month • SBS • Slope • Elevation.SBS: soil burn severity.The model includes all tractor and skyline units (n = 318) with complete datasets (no missing values) for all three fires.

Table 7 .
Decision tool for evaluating risk or likelihood of detrimental soil disturbance (DSD).