Impacts of forest biomass operations on forest hydrologic and soil erosion processes ☆ Trees, Forests and People

The use of biomass from forests by humans has been practiced since the dawn of time. Such practices have resulted in severe watershed degradation historically and in today ’ s world. As modern society seeks alternative energy sources, they are once again turning to forests, but now use mechanized harvesting systems to remove the biomass and transport it to processing or utilization sites. This study investigated the impacts of different mechanized biomass operations on sediment displacement, soil hydraulic conductivity, and canopy interception to determine if severe watershed impacts associated with biomass utilization in the northwest U.S.A. may occur as seen elsewhere. Three field sites were selected that used different biomass operation methods: the Payette National Forest (NF) in Idaho using a forwarder system, the Colville NF in Washington State using a tractor-skidder system, and the Flathead NF in Montana using a skyline system. A total of 47 silt fences (5 m × 10 m) were installed to monitor soil movement from biomass operations: 28 silt fences on disturbed plots and 19 on control plots. No soil loss was observed except on the forwarder trails in the Payette NF where animal distur- bances were responsible for some soil movement into the silt fences. Infiltration rates measured using a constant head permeameter showed that the saturated hydraulic conductivity on the disturbed plots for the ground-based operations was lower than the control plots, but unlikely to cause significant changes in runoff on any of the plots except the Payette forwarder trails, where the saturated hydraulic conductivity values averaged only 5.3 mm h (cid:0) 1 , compared to greater than 20 mm h (cid:0) 1 on all other observations. Ground cover increased on the disturbed plots to levels observed on undisturbed plots within two years. The only significant increase in bulk density was observed on the forwarder trails at the Payette NF site where the density averaged 1.41 on the trails, compared to the densities of 1.0 elsewhere. The loss of canopy observed on the three biomass studies led to a complementary study on the effects of canopy loss on snow accumulation and melt on the Priest River Experimental Forest in Northern Idaho. The study found that the canopy only intercepted an average of 7 mm of snow water equivalent during three winter months, but average annual interception of precipitation was 114 mm, indicating that interception of rainfall is greater than snow. Overall, these results suggest that there were no direct detrimental hydrologic impacts due to biomass operations except on the compacted, and sometimes pulverized forwarder trail conditions observed in the Payette NF. The impact of the canopy on rainfall interception is less than the water yield differences due to clear cutting on nearby forested watersheds that were attributed to differences in evapotranspiration.


Introduction
The world's forests have long been sources of biomass for renewable energy and construction materials (McGranahan, 1991;Elliot et al., 2010). Temperate and boreal forests also remove biomass to improve forest health and reduce the risk of wildfire . In recent decades, when world supplies of fossil fuels were compromised, forests were turned to for energy not only in the third world (Badege, 2001) but in developed countries as well (USDA Forest Service, 2003; U.S. Department of Energy (USDOE), 2017). Many developing countries relied on forest biomass for energy and construction, especially in rural areas (Badege, 2001;Hosier and Bernstein 1992;Sattar, 1996). Some African and Southeast Asian countries still consider woody biomass a viable energy source (Naughton-Treves et al., 2007;Sasaki et al., 2009;Owen et al., 2013). Some countries with abundant forest biomass resource, such as Indonesia, considered converting forest biomass into bioenergy (Freppaz et al., 2004;Simangunsonga et al., 2017). Most of the third world studies were conducted on bioenergy potentials (Freppaz et al., 2004;Naughton-Treves et al., 2007;Sasaki et al., 2009;Owen et al., 2013;Simangunsonga et al., 2017) while some studies evaluated soil productivity (Thiffault et al., 2011;de Varies et al., 2021).
Accessing forests to extract biomass for energy has caused damage to other forest resources, such as soil productivity as observed in some developing countries (Badege, 2001;Hosier and Bernstein 1992) and clean water (Elliot, 2010; U.S. Department of Energy (USDOE), 2017). The U.S. Department of Energy (USDOE) (2017) also evaluated the potential impacts of biomass utilization on greenhouse gasses, air pollution and biodiversity. Kort et al. (1998) and Elliot (2010) suggested that removing forest biomass for energy could decrease ground cover and increase compaction, runoff, and soil erosion. Sun et al. (2017) evaluated the impacts of forest biomass removal such as clearcutting and thinning on water yield in the U.S. and estimated water yield increases ranging from 10 mm per year in dry forests to 150 mm in wetter forests in the Northwest and Eastern United States (U.S.) Caputo et al. (2016) assessed the effects of forest biomass harvesting on water and climate regulation in eastern North America and determined that increases in nutrient delivery from watersheds where biomass has been removed would be elevated for about ten years. The U.S. Department of Energy (USDOE) (2017) report estimated small increases in nutrient delivery as well for several years following biomass removal.
In contrast to the third world where biomass harvesting is done manually, in the U.S. and Europe biomass is harvested and transported from the site with heavy machinery (Rummer, 2010). Third world biomass harvesting removes almost all biomass, leaving little residue behind, and in the warmer climates, decomposition of remaining residue was greater than in cooler northern forests (Prescott, 2005). It is not known how much ground cover remains following mechanized harvesting systems removing not only tree boles, but also the slash that is normally left onsite in most mechanized harvesting operations (Elliot, 2010). There have been no studies on the upland impacts of biomass removal using mechanized forest practices on ground cover, soil hydraulic conductivity and soil displacement. The sediment delivery estimated in the U.S. Department of Energy (USDOE) (2017) report was based on watershed studies that would include channel as well as upland erosion.
Forest canopy has played an important role in many hydrologic processes, including evapotranspiration, interception, sublimation, albedo, radiation, and snowmelt (Brooks et al., 1991). The loss of canopy associated with biomass removal may perturbate hydrologic processes increasing runoff, snowmelt rates, and soil erosion, adversely impacting the watershed (Srivastava et al., 2020). The effects of the loss of canopy through timber harvest and thinning on water yield are well documented. Most, but not all studies, reported that water yield increased following timber harvest (Hubbart et al., 2007;Srivistava et al., 2018 and2020;Sun et al., 2017;Troendle et al., 2010). In climates where canopy interacts with snow through interception, sublimation, and melt rates, the effects of the loss of canopy on snow accumulation and melt are less clear (Dobre et al., 2011;Hubbart et al., 2007;Molotch et al., 2007;Troendle et al., 2010). There have been no published studies that measured the impacts of the biomass removal and associated canopy interception on snow processes.
This study assessed the upland impacts of biomass removal using different mechanized forest practices on soil displacement, hydraulic conductivity, ground cover, and rain and snow interception. The findings may aid forest managers in developing best management practices to reduce impacts associated with biomass utilization. To address these biomass impacts, a series of research studies were carried out at three different locations in the Northwestern U.S. measuring changes in ground cover, soil bulk density, soil hydraulic conductivity and soil displacement associated with different biomass harvesting systems. On a fourth site, the effects of canopy loss on snow and rain interception, and snowpack dynamics were measured. This article reports the findings of those studies.

Methods
There were two aspects to this study. One aspect explored onsite impacts of biomass removal, while the second aspect evaluated the effects of forest canopy on snow accumulation and melt as there was a notable lack of canopy on all three of the sites where biomass had been removed. Most statistical analyses were carried out in R (R Core Team, 2021), with Excel (Microsoft Corporation, 2018) used to prepare and preview the data, graph results and complete some simple t-tests. For analyses of variance the significance will be presented as "(P <0.0xx)" where "0.0xx" will be the probability that the independent variable or interaction did not significantly cause variability in the dependent variable. For t-tests the significance will be presented as "(α<0.0xx)" where "0.0xx" is the probability that the mean values being compared are not significantly different.

Biomass removal sites
The three biomass removal sites were in the Northwestern U.S. where biomass operations had just occurred ( Fig. 1; Table 1). The harvesting machines used were described in detail in Rummer (2010). In the New Meadows site located in the Payette National Forest in Idaho, a feller-buncher cut trees for saw logs and laid them in piles, and a forwarder transported the merchantable logs to a nearby landing. The remaining treetops and slash were skidded to a landing where the tops were chipped (Fig. 2). The chipper filled dump trucks that shuttled the chips to a staging area where the chips were transferred to larger semitrailer chip vans for transport to a nearby mill for energy generation. Following the forwarder operation, the forwarder trails were obliterated by scarifying and covering with mulch.
A tractor skidder was used in the Summit Pierre site located in the Colville National Forest in Washington State. Normally, the biomass for this site would have been hauled to a biomass cogeneration plant located in Kettle Falls, WA. However, due to the poor biomass market at that time, the chips were not delivered to the cogeneration plant, but were burned at the stockpiled location. The burning of this biomass offsite did not impact the analysis for this site however, as the biomass that would have been used for fuel utilization was removed from this unit prior to this study.
On the third site, a skyline system removed the biomass from the generally steeper Cyclone Reed sale in the Kootenai National Forest in Montana (Table 1). The slash and merchantable logs were pulled to the landings where part of the slash was hauled to a cogeneration plant in Columbia Falls, MT and the merchantable logs were delivered to a post Fig. 1. Locations of three biomass removal and canopy interception study sites the NW U.S. and rail mill in Libby, MT. The remaining slash was left in piles near the road or throughout the unit and later burned. The trails left by dragging the trees up the hill were scarified when the operation was complete.
The study plan was to measure soil displacement, ground cover, bulk density, and saturated hydraulic conductivity on four treatments at each location: control plots with no disturbance, treated plots where biomass was removed but the plot was not compacted by harvesting equipment, access trail plots where biomass was towed or hauled from the site, and off-trail observations adjacent to the trail observations. Plots were installed during the summer following biomass removal.
At each location, between 12 and 21 replicated plots (Table 1) were defined within areas that had been harvested, in adjacent undisturbed areas, and on access trails for completing the four sets of measurements. Displacement plots were nominally 5 m wide and 10 m long, with slope steepness on individual plots ranging from 13 to 43 percent (Table 1 ;  Fig. 3). The aspect and steepness of each plot was measured with a smart phone with the "Spyglass" app ( Fig. 3).

Weather on biomass removal sites
A weather station was installed at each biomass removal site recording precipitation, wind speed and direction, temperature, and relative humidity. On the New Meadows, ID site, an additional weather station was installed 1200 m from the weather station near the treated area in a clearing in the adjacent undisturbed forest. The station in the treated site had more canopy in the vicinity of the station than the station in the nearby clearing.
On the Montana site, there was a lapse in data in the first year when the weather station was severely damaged by gun fire and the data logger stopped recording in the third winter. Missing data were estimated from the Ashley Divide SNOTEL station located 22.1 km southeast of the research site at an elevation of 1469 m (Snow Telemetry; USDA NRCS, 2017), similar to the research site at 1402 m.

Soil displacement
Soil displacement was measured by installing silt fences on the downhill side of plots assumed to be the width of the silt fence (5 m) and 10 m below a natural boundary, usually the crest of a small hill, to collect sediment displaced by surface erosion processes (Fig. 3). The installation of a geotextile silt fence was a low-cost technique to retain displaced sediment to estimate water erosion on relatively small plots (Robichaud and Brown, 2002). On the Idaho site, 21 plots were installed (7 control, 7 on the trail, and 7 off the trail), 14 on the Montana site (7 control and 7 on the trail), and 12 on the Washington site (6 control and 6 on the trail) (Table 1). In May or June of each year, accumulated sediment, if any, was collected from the uphill side of each silt fence and weighed on site. A sediment sample from each fence was collected to Table 1 Site characteristics of biomass utilization study sites in the NW U.S. (Fig. 1   determine the water content and dry mass by oven drying at 105 • Celsius in the Moscow, Idaho Forestry Sciences Laboratory (MFSL).

Mineral soil exposure and duff amount
Ground cover was measured at three points within each soil displacement plot following plot installation, and during the summer for three years following the operation. The measurement locations for the ground cover within the displacement plots were determined by using a random number generator on a scientific calculator. For each location, a random number between zero and 1.0 was multiplied by length and width of the plot, which were nominally 10 m long and 5 m wide. The datum (0, 0) was set at the SE corner of the silt fence plot. At each location, a 1-m square frame with a grid of wires at 0.1 m was placed on the ground, and the presence of cover (vegetation, branches, rock, duff) and exposed mineral soil were recorded to determine the fraction of each cover category. The three cover observations were averaged for each plot. The exposed bare mineral soil was for detailed statistical analysis as that measurement is frequently cited in forest plans to evaluate disturbance impacts. An analysis of variance including all the interactions was applied to the mineral soil exposure observations to evaluate differences among sites and years since treatment.
On the snow interception study, we also measured duff depth and mass (Mg ha − 1 ) at four points around each of the 12 rain gages, and used the average of the four for each rain gage location during the first year of the study. A t-test was used to determine if there were differences in depth and mass between treatments (clearing vs. forest).

Soil bulk density
Soil bulk density was measured using an undisturbed soil core when the displacement plots were installed the summer following the harvest. A cylindrical metal sampler with a 4.8 cm diameter by 10.1 cm long (total volume 183 cm 3 ) was driven into the soil with a sliding hammer. In each soil displacement plot, a core was collected from 0 to 10 cm, and a second core from the same hole from 11 to 21 cm. The sampler was carefully removed to preserve the soil core (Blake and Hartge, 1986;Page-Dumroese et al., 1999). A single set of samples was collected from each displacement plot. Samples were placed in cans and taken to the MFSL where they were weighed, dried at 105 • Celsius, and weighed again to determine both soil water content and dry bulk density. Soil texture was determined by wet sieving and pipetting methods (Guy, 1969).

Hydraulic conductivity
Adjacent to each of the displacement plots, and on treated areas and off trail locations for all sites, hydraulic conductivity was measured at six locations using a constant head permeameter ( Fig. 4; NRCS, 2010). The initial measurements were made when the plots were installed and a second set the following year during the summer. Six permeameters were fabricated in the workshop at the MFSL but modified from the NRCS (2010) design by replacing the top glued cap with a cam locking cap to allow ease of refilling (Fig. 4). Two telescoping legs were designed to support each permeameter so that one person could operate multiple permeameters at the same time. The scale for the permeameter was calibrated in inches of water infiltrated (NRCS, 2010). Following the NRCS (2010) procedure in the field, a small 28-mm diameter hole about 110 mm deep was excavated. A supporting collar for the permeameter was positioned at ground level above the hole. The rate of flow into the soil from the permeameter was determined by the volume change of water in the permeameter reservoir for a given time increment. When the flow rate reached equilibrium for two consecutive readings, the measurements were stopped, and the final flow rates were converted to mm h − 1 . This value was assumed to be a measurement of the saturated hydraulic conductivity (NRCS 2010).
The results were averaged across sites for each year (year 1 and year 2) and an analysis of variance carried out to determine differences between years and interactions among the sites within years. The values were then averaged for both years for each site, and a second analysis of variance carried out to determine if there were differences between sites or interactions within the treatments by site.

Canopy effects on snowpack
One feature of the three biomass removal sites was the absence of forest canopy following the biomass removal operation. To measure the hydrologic effects of loss of canopy, an ideal site with a climate similar to the three biomass sites was the Priest River Experimental Forest (PREF) (Fig. 1). The proximity of the site to the MFSL made it feasible for monthly trips to the site during winter weather to carry out a snow survey. Within the PREF, a site was identified that had been subjected to a simulated wildfire in October 2006 followed by salvage logging in the summer of 2007, so there was minimal canopy (Elliot and Glaza, 2009). There was also some historic snow course data from an earlier study on that site (Dobre et al., 2011 and2012). Adjacent to the cleared site was a western red cedar (thuja plicata)/ponderosa pine (pinus ponderosa) forest that had not been disturbed for more than 80 years (Fig. 5). To potentially link results of the study to satellite imagery, the center point  of each 30-m Landsat pixel was overlaid on the site (https://www.usgs. gov/landsat-missions/landsat-8). Six rain gauges were located on the center points of Landsat pixels in the area where the canopy had been removed, and six gauges were installed beneath the adjacent forest canopy area in the center of Landsat pixels (Fig. 5). Each year in November, the rain gages were "winterized" by the addition of a reservoir of environmentally safe antifreeze above the rain gage, so the rain gage measured the overflow from the antifreeze when either rain or snow fell into the reservoir. The antifreeze reservoirs were removed in April of each year. A previous study had suggested that snow melt rate may be related to soil water content (Dobre et al., 2011), so in addition to measuring precipitation with 12 winterized rain gages, we installed two soil water content sensors at each rain gage, placed approximately 10 cm deep. In previous studies, the soil water content gauges had not been reliable, hence the justification for installing two such sensors at each rain gage site. In addition to the 12 rain gauges and 24 soil water probes, an additional 32 snow course points were established on the 30-m grid (Fig. 5). Two weather stations monitoring temperature, humidity and wind speed and direction were installed, one in the center of the clearing and another near the road on the edge of the forested site (Fig. 5). To measure incoming and outgoing long wave and short-wave radiation, two radiation stations were established near the PREF headquarters where a constant electrical supply was available to keep the radiation sensors warm so they would not get covered with snow. One station was set in a clearing, and the other was set beneath a mature canopy similar to the cedar/pine canopy on the snow course site. Each station had four radiation sensors, long-wave and short-wave sensors pointing up measuring incoming radiation, and long-and short-wave sensors pointing down measuring outgoing radiation. Net radiation was determined by subtracting the total outgoing radiation from the incoming radiation. Radiation data were collected for the water years of 2014, 2015, and 2016.
During each of the winters, snow depth and density were measured at monthly intervals starting in January if there was measurable snow at some of the snow course points (Natural Resource Conservation Service (NRCS) 2017). Rain and snow canopy interception were assumed to be the difference between the rain gauges in the clearing, and those beneath the canopy.
In October 2014, the ground cover was measured at each rain gage at three nearby locations using the same method as the other studies. In addition, the depth of ground vegetative residue (duff) was measured at four points around each rain gage, and all the organic material removed from within a 30-cm ring to be dried and then weighed in the laboratory. The mass of material collected within the rings was converted Mg ha − 1 .
A least square means analysis was applied to the precipitation gage data that included an analysis of variance determining whether there were differences between the years or the treatments and the interaction between years and treatments. An analysis of variance on the effects of year and treatment was carried out for the snow water equivalent and soil water content observations. T-test analyses were carried out to compare the temperatures and net radiation values between the clearing and beneath the forest canopy. A correlation analysis was conducted to determine if there was a significant relationship between snow water equivalent and the soil water content averaged from the preceding ten days.

Biomass removal sites
On the Idaho Biomass site, skid trails were easily identified. On the Washington site, the scattering of slash following biomass removal masked the locations of the skidding trails. The skidder tended to use the unit's total area for multiple passes, making it difficult to discern the offtrail locations when collecting bulk density samples within the treated area. In Montana, few trails were visible that could be attributed to the skyline operation for plot installation. Thus, only the "Treated" samples for ground cover, bulk density and soil displacement were collected in the Washington and Montana sites to compare to the control sites. We were however able to identify enough locations for "Trail" samples for measuring hydraulic conductivity on the Washington and Montana sites.

Precipitation
The observed weather data for the biomass removal sites are summarized in Table 2. The water year runs from October 1 of the year preceding the water year until September 30 of the water year. The properties of the largest daily precipitation amounts, and the storm with the greatest intensity each year were determined, as these events are generally associated with the largest erosion events in each year (Robichaud et al., 2016). The Idaho site was the wettest, with the average annual precipitation for the four years of record of 915 mm on the treated site and 983 mm with the clearing station. The Washington (540 mm) and Montana (500 mm) sites were drier. The highest daily precipitation amounts were 115 mm in Idaho, 93 mm in Montana and 33 mm in WA. The Idaho storm was an approximately a 100-year return period event, as was the Montana storm. The maximum 10-and 30 min precipitation depths did not occur the same days as the maximum depths that were associated with long duration precipitation events.

Soil displacement
The only site where soil loss was observed was on the Idaho Site (Fig. 6). The soil loss rates on the Idaho plots were low with an average yearly rate on the control sites of 24 kg ha − 1 . Average soil displacement rates of more than 40 kg ha − 1 were measured on the trail plots and 4 kg ha − 1 on the off-trail plots in the first year of the Idaho study. There was no soil displacement observed on any of the control or treatment plots on the other two sites.

Mineral soil exposure
The Idaho Biomass site was our first installation and served as a pilot study for the other two sites. After the first year, we determined that there was no significant difference between the control and the off-trail plots, with minimal soil displacement from the control plots and no soil displacement on the off trail plots. The bare mineral soil exposure averaged 0.03 for both treatments. In subsequent years, we measured more mineral soil exposure on the control plots than on the off-trail plots on the Idaho site (Tables 3 and 4). Table 3 shows the detailed results from the Idaho site for all the cover categories. We never recorded any ash or char for any treatment, and less than 0.01 rock cover. The control plots had more branches (0.122) than the off-trail plots (0.095), but less other organic material such as needles, twigs, and duff (0.841 vs, 0.897, respectively). We observed similar ground cover conditions at the Montana and Washington sites, with considerable ground cover for the off-trail and the control, so we only measured the control plots and the "treated" plots on those sites. Table 4 presents a summary of the mineral soil fraction measured for all three sites averaged over four years and for the average bare mineral soil fraction for each year for each site. An analysis of variance showed that there were significant differences among the sites, with the Idaho site having the greatest soil exposure, among the years, with the bare soil declining each year, and due to treatment with the trail having more soil exposure than the control or off-trail (P<0.001 for all comparisons). There was also a significant interaction between site and treatment, suggesting that the difference between the trail and the control in Idaho was much greater than at the other two sites; between treatment and year, suggesting that the rate of recovery was not the same for all treatments; and between site and year, suggesting that the Idaho likely recovered faster than the other two sites (Table 4).
On the snow canopy study, the average ground cover was 0.957 in the clearing, and 1.00 in the forest. Since all 12 of the forest ground cover measurements were zero, the bare soil in the clearing was significantly greater than the bare soil in the forest (α<0.01). The average depth of duff in the clearing was 1.63 cm and in the forest was 4.83 cm. The average duff mass was 31.45 Mg ha − 1 in the clearing and 52.29 Mg ha − 1 in the forest treatment. Both the depth and the mass were significantly different (α<0.01).

Bulk density
Observed soil bulk densities on the biomass utilization sites ranged from 0.68 g cm − 3 on the control treatment for the Washington site to 1.41 g cm − 3 on the treated plots for the Montana site (Table 5). There were significant differences due treatment, between sites, and between depths (Table 5). There were no significant interactions. The trail/ treated observations had higher bulk densities than the control or off trail values, whereas the MT site had the highest overall bulk densities (1.24).

Table 2
Weather Observations on the three biomass removal sites for the biomass utilization study in the NW U.S. The Idaho site had two weather stations, one in the treated area (Idaho T) and one in a nearby clearing (Idaho C).    trails were more dispersed and covered with slash, so samples were taken from "Treated" areas that may have included hidden skid trails. trails were more dispersed and covered with slash, so samples were taken from "Treated" areas that may have included hidden skid trails.

Hydraulic conductivity
Individual saturated hydraulic conductivity values ranged from zero to nearly 160 mm h − 1 and the average observations for each treatment, site and year are shown in Table 6. The data set was highly skewed with the Shapiro-Wilks normality test resulting in a probability of normality less than 6 × 10 − 10 . The data were transformed by adding 1 mm h − 1 to the conductivity values to remove the zero values, and then taking the cube root to achieve normality (P = 0.16) before applying an analysis of variance to the data set. The trail/treatment observations were the lowest at all sites, with a mean of 6.3 mm h − 1 on the Idaho site, and the greatest mean value for trails of 35.9 mm h − 1 on the Washington Site ( Table 6). The Idaho control site had an unexpectedly low mean value of only 12.7 mm h − 1 , less than the treated and off-trail values of 25.7 mm h − 1 . On the Washington site, the treated mean value (43 mm h − 1 ) was also greater than the control mean value (42.1 mm h − 1 , Table 6). The hydraulic conductivity values measured in the second year after the biomass removal were significantly greater than the first year (P < 0.001). There were also significant differences among the sites (P < 0.001). However, the high variability within the observations meant that there were no significant differences among the treatments nor any significant interactions.

Rain and snow interception
In the snow and rain interception study at the PREF during the winter months there was little net interception (Table 7, Fig. 7), although there was variation from day to day. Snowfall on cold days resulted in positive interception, whereas shortly after a snowfall event, melting snow sometimes caused a "negative" interception because rain gauges under the canopy measured precipitation whereas rain gages in the clearing did not (Fig. 7). In the winter of 2015, there was more precipitation measured in the forest (Table 7). During wet periods in the spring and fall when temperatures were above freezing, there was considerable canopy interception observed from some rainfall events (Fig. 7). During the three winter months, the observed daily precipitation collected by the rain gages beneath the canopy averaged 7 mm less than the observed values in the clearing, but the difference was not statistically significant (Table 7). Table 7 shows that for the annual totals, the clearing rain gages measured significantly more precipitation than the forest during the three-year study, with an average annual interception of 114 mm (P<0.01). There were also significant differences among years for total precipitation (P < 0.001). In particular, the precipitation depths measured in the clearing in 2016 and 2017 were greater than the clearing in 2015, and greater than all the forest rain gauges (Table 7).

Air temperatures and net radiation extinction
On the PREF canopy study, the mean temperature measured by the forest weather station was 8.00 • C which was significantly less than in the clearing where the mean temperature as 8.14 • C (P < 0.001). Fig. 7 shows the day-to-day variation in temperature from the weather station in the forest. The differences between the two treatments were too small to display in a graph.
The forest canopy had a significant impact on net radiation. A paired t-test determined that the net radiation in the clearing with mean value of 78.35 w m − 2 was significantly greater (P < 0.001) than in the forest where the mean radiation was only 8.38 w m − 2 . Fig. 8 shows the daily net radiation values. Differences in radiation were much greater in the Table 6 Saturated Hydraulic Conductivity (mm h − 1 ) for four treatments on three sites 2 years following the treatment for the biomass utilization study in the NW U.S. There were significant differences between the two years (P<0.001), and among the three sites (P<0.001). The differences among the treatments were not significant, nor were either of the interactions.   Fig. 7. Cumulative and daily interception of snow and rainfall (mm), and average daily temperature ( • C) in the canopy interception study at the Priest River Experimental Forest.

Fig. 8.
Net radiation (w m − 2 (incoming short wave + long wave) -(outgoing short wave and long wave)) from a weather station in a clearing and a weather station in a nearby forest in the canopy interception study at the Priest River Experimental Forest.
summer than in the winter. We experienced some technical problems with the instrumentation including falling branches in the forest, loss of power, and animal damage, so of the 1103 days in the study, there were 252 days with missing data in the clearing, and 481 days of missing data in the forest. Table 8 summarizes the snow water equivalent from the monthly snow course observations. There was no snow on the ground in March 2015. There were significant differences between the SWE values observed in the clearing (8.9 mm) compared to under the forest canopy (7.03 mm). There were also significant differences among the years, with 2015 having the lowest average value (4.17 mm) and 2017 the greatest (11.6 mm). There were no significant interactions among the data.

Snow water equivalent (SWE) and soil water
The mean soil water content was the same for the forest and the clearing (0.19 m 3 m − 3 ). Fig. 9 shows the variation of soil water during the year. There was no consistent trend, with the forest soil generally wetter than the clearing in the winter of 2015, but drier in the winters of 2016 and 2017. The forest soils tended to dry out earlier in the summer but wetted up earlier in the autumn. The Pearson correlation between SWE and soil water was not significant with a value of − 0.131 (P = 0.629). Although not significant, the negative correlation value suggests that the soil was drier when the snowpack was deeper, but as the snow melted and SWE declined, soil water increased.

Discussion
The main finding of this study is that there is minimal soil displacement on trails, and negligible displacement on the treated areas (Fig. 6) despite some very large rainfall events on the Idaho and Washington sites ( Table 2). The low erosion rates are likely because of the high level of ground cover that remained on all three sites following the biomass removal (Figs. 3 and 4; Table 3). Rather than a decrease in ground cover as hypothesized by Elliot (2010) due to biomass removal, the off-trail or treated ground cover was similar to the control treatments (Tables 3 and 4). On a small watershed study near the Idaho site, Covert (2003) reported that observed sediment delivery from a 2-ha watershed ranged from 0 to 90 kg ha − 1 over a four-year period, and much of the sediment was attributed to channel erosion on a skid trail immediately upstream of the inlet to the runoff collector (Covert, 2003). The larger runoff contributing area in the Covert (2003) study may explain why soil displacement was greater in that study than we observed, which was less than 50 kg ha − 1 and frequently zero (Fig. 6). Our observations were similar to those reported by Elliot and Glaza (2009) who measured zero erosion following forest activities, and Elliot and Miller (2004) who measured sediment displacement rates from less than 1 to more than 200 kg ha − 1 , with the greatest displacement attributed to pocket gophers Much of the soil displacement measured on the Idaho site trails was associated with post treatment animal activity, including pocket gophers depositing sediment in front of the silt fences, and ongoing trail disturbance from elk. The Idaho trail became the preferred path for browsing elk, reducing vegetation recovery (Table 4). Elliot and Miller (2004) had attributed soil displacement on a fuel management study in Northern Idaho to gopher activity. In that study, the gophers were a greater problem on the control plots than on the skid trails and were the dominant source of displaced sediment on the control plots. The Idaho site also had the greatest precipitation, which may be why it was the only site with measurable erosion in the control treatment.
On the Idaho site, the harvesting system used "small" dump trucks to shuttle the chips from the site to a staging area where larger chip vans were loaded. This practice resulted in a pulverized soil condition on the access roads, aggravating the erosion potential (Fig. 10). Traffic was not sufficient on the forwarder trails on this site to result in such severe disturbance. This increased risk of erosion on access roads confirmed what Elliot (2010) had hypothesized, that increased traffic associated with fuel management may increase erosion risk. On the Washington and Montana sites, where the dump truck shuttling chips was not a part of the system, there did not appear to be any road impacts.
The higher bulk density observations in the Montana site may be a reflection of the greater content of clay in this glacially-derived soil. It is also possible that the soils may have been wetter during the harvesting operation, or that there may have been compaction from previous logging operations (Froehlich et al., 1985). Freohlich et al. (1985) observed Table 8 Snow water equivalent (mm) for the clearing and the forest for each month (no snow observed in March 2015) and year for the canopy interception study at the PREF. There were significant differences among the years (P < 0.001) and between treatments (P = 0.006), but no differences in SWE among the months, nor for interactions.   similar bulk densities to the Montana treated site on locations that had previously been in skid trails. They reported that on loamy soils, bulk densities of skid trails were slow to recover to pre harvesting levels, so this high value may be due to a previous operation, or to the biomass harvesting activity. The increase in bulk density (Table 5) on the Montana site, however, did not appear to cause a decrease in hydraulic conductivity (Table 6). The two rain gages on the Idaho biomass site were not consistent (Table 2), with the three-year mean precipitation on the treated site, with some canopy recording 848 mm and the nearby control site with less canopy averaging 983 mm. This 135 mm difference in measured precipitation was similar to the more intense canopy study at the PREF site in Northern Idaho, where there was a significant difference of 114 mm in measured precipitation annually between the canopy and the clearing ( Table 7), and that difference was attributed to canopy interception. In the Idaho biomass study, we averaged 80 events per year with greater than 2 mm of precipitation on the biomass site each, suggesting an average single storm interception of 1.7 mm for the partial canopy on that site. In the PREF study, the average number of events with greater than 2 mm of precipitation was 68, leading to an estimated average single storm interception of 1.7 mm, the same as estimated for the Idaho biomass site further south. Typical rainfall interception rates for forests are 3.6 to 3.8 mm per event (Ward and Elliot, 1995), for complete canopies, suggesting that in both the Idaho sites, the canopy cover was incomplete. Amatya et al., 2016 reported that conifer forest interception rates were typically 18 to 25% of the annual precipitation, which translates to 176 to 245 mm for the biomass site, and 142 to 198 mm for the PREF site, also greater than the interception amounts we observed of 135 mm for the biomass site and 114 for the PREF site. The lower estimated interception rates we observed may also be due to less interception during the winter months (Table 7).
The rain gage results from the canopy study (Table 7 and Fig. 7) showed that there was no significant difference between rain gages in the clearing and in the forest during the winter, averaging 7 mm more in the clearing but annually there was a significant difference of 114 mm more precipitation measured in the clearing. This indicates that interception was greater from rainfall events than from snowfall events. The snow course results showed that the snow water equivalent (SWE) in the forest treatment was significantly less than in the clearing, with an average difference of 1.9 mm. The rain gage values cannot be directly compared to the snow course observations as they are not measuring the same variables, with the rain gauges indicating rain or snow input, whereas the snow course measures the net effect of previous SWE and snow input less snow melt and sublimation leaving the snowpack. This means that even though the clearing received more snowfall increasing the SWE, it may also have experienced a higher melt rate due to greater radiation (Fig. 8) or slightly higher air temperatures (8.14 vs 8 • C) reducing the SWE in the clearing.
The interception by the forest is lower than snow interception followed by canopy sublimation as reported by Molotch et al. (2007) who reported a difference of 20 mm during the winter 3 months. Troendl et al. (2010) reported that in a cold snow zone, canopy interception may be 25 to 35% of the winter snowpack. Canopy interception and sublimation are unlikely to be the cause of an average annual increase in runoff of 270 mm (Hubbart et al., 2007) or 284 mm (Srivastava et al., 2020) due to clear cut timber harvest observed at the Mica Creek Experimental Watershed, located 140 km southeast of the study site. Both the Hubbart et al. (2007) and Srivastava et al. (2020) attributed increased stream flow from clear cut plots to decreased evapotranspiration and not snow interception. All the sites in this study were in a maritime/continental climate, whereas the Molotch et al. (2007) study was conducted in Colorado U.S., at a higher elevation with colder winter temperatures and lower winter water vapor pressures in a continental climate.
The net radiation results (Fig. 8) had some gaps in the data due to equipment failure, loss of power, or in the forest instrumentation, damage from wildlife or falling tree limbs. During the summer months, radiation was much greater in the clearing, frequently greater than 150 w m − 2 compared to beneath the forest canopy where net radiation seldom exceeded 20 w m − 2 . In the winter months however, net radiation was often similar under the canopy compared to the clearing, likely due to much shorter days with lower solar radiation and more cloud cover. Some days had a negative net radiation, and on those days, the clearing sometimes lost more energy to radiation than the forest, hence the more negative values (Fig. 8). This may have resulted in a colder snowpack in the clearing, reducing subsequent sublimation or melt rates from radiation in the days that immediately followed. Additional climate factors such as wind speed may impact interception and distribution of rain and snow on along the boundaries of forest and clearing, and further analysis of factors impacting melt rates, such as relative humidity may also be influencing the observed SWE values. Such analyses are beyond the scope of this article.
The soil water probes were near the surface (10 cm deep), intended to measure the effects of soil water on snow melt. There was, however, no significant correlation between SWE and soil water content. The data from the soil water sensors (Fig. 9), however, provide additional insight into the effects of a forest canopy and duff on soil water. Soil water is the net result the addition of precipitation and snowmelt less lateral flow, deep seepage and evapotranspiration. At the start of the study (October 2014), soil water contents were similar for the forest and the clearing (0.2 m 3 m − 3 ), near the overall average value of 0.19 m 3 m − 3 . In the dry summer season, soil water appears to drop more rapidly in the forest than in the clearing, likely reflecting the higher ET values from mature forests Hubbart et al., 2007;Srivastava et al., 2020). However, by the end of the summer, soil water content at 10 cm is less in the clearing, which may reflect the lack of canopy to reduce surface heating from solar radiation and a less deep duff layer reducing surface evaporation in the clearing, resulting in drier soils around the shallow probes. When autumn rains begin to wet the soil, the water content appears to rise more rapidly for the forest than the clearing (Fig. 9), which may be due to higher infiltration in the undisturbed forest soils compared to soils that had previously been disturbed by fire and salvage logging in 2006. Our observations from the biomass sites, however, did not show significantly different hydraulic conductivity values associated with biomass removal (Table 6), but none of the biomass sites had been burned. The constant head permeameter we used (Fig. 4) measures the hydraulic conductivity beneath the duff layer, and not on the duff surface. Other studies have shown that the presence of a mulch layer increases infiltration (Mannering and Meyer, 1963). The significantly greater mass of duff in the forest may also have retained some of the precipitation, leaving the soil beneath drier. Raaflaub and Valeo (2009) measured duff porosity values around 0.8 m 3 m − 3 , suggesting the clearing duff with a depth of 1.6 cm could store up to 1.3 cm of precipitation, whereas the forest duff with a depth of 4.8 cm could store as much as 3.8 cm, reducing the precipitation reaching the soil in the forests, hence the drier forest soils in the winters of 2016 and 2017 (Fig. 9). Fig. 9 shows that the maximum difference between the clearing and the forest in the winter is only about 0.05 m 3 m − 3 or a difference of 0.5 cm of water in the top 10 cm of soil. The high ground cover fraction measured in three biomass sites suggest that even though the overstory is gone, the effects of the duff layer will still be influencing near-surface hydrology. These results suggest that the role of forest duff on hydrology, and the effects of forest management on duff warrant further investigation.
This study showed that the mechanized removal of biomass followed by forest regeneration as practiced in the U.S. appears to cause minimal environmental damage compared to adverse impacts from overutilization reported in some third world countries (Badege, 2001;Hosier and Bernstein 1992;Sattar, 1996). In addition, in the U.S. both public and private forest managers follow best management practices intended to minimize upland soil erosion and downstream water quality deterioration (https://www.stateforesters.org/bmps/), such as leaving sites covered in slash as observed in this study.

Summary and conclusions
A study was carried out on three sites in the Northwest U.S. to measure the impacts of biomass utilization on soil and water resources. Ground cover was minimally impacted from the biomass removal, but cover was reduced on the trails the year following the treatments. Soil compaction and a decrease in hydraulic conductivity were observed on the trails, but not on the treated areas. On two of the three sites, no soil displacement was measured on either the treated forested areas or on the trails. Only on one site was there significant soil loss measured on the trails. This was likely due to a higher level of traffic and wildlife disturbance on the trails following the biomass removal, as well as greater precipitation.
The loss of canopy due to biomass removal was evaluated in more detail. We found that canopy interception was greater for rainfall than for snowmelt, and the combination of interception under the canopy and likely greater melt rates in the clearing led to a net result of 1.9 mm more snow water equivalent measured in the clearing, a small but statistically significant difference.
Future areas for further research are to evaluate edge effects on SWE and interception, and the effects of the duff layer on soil water dynamics. Both factors will be influenced by biomass utilization.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.