Snowmobile impacts on snowpack physical and mechanical properties

Snowmobile use is a popular form of winter recreation in Colorado, particularly on public lands. To examine the effects of differing levels of use on snowpack properties, experiments were performed at two different areas, Rabbit Ears Pass near Steamboat Springs and at Fraser Experimental Forest near Fraser, Colorado USA. Differences between no use and varying degrees of snowmobile use (low, medium and high) on shallow (the operational standard of 30 cm) and deeper snowpacks (120 cm) were quantified and statistically assessed using measurements of snow density, temperature, stratigraphy, hardness, and ram resistance from snow pit profiles. A simple model was explored that estimated snow density changes from snowmobile use based on experimental results. Snowpack property changes were more pronounced for thinner snow accumulations. When snowmobile use started in deeper snow conditions, there was less difference in density, hardness, and ram resistance compared to the control case of no snowmobile use. These results have implications for the management of snowmobile use in times and places of shallower snow conditions where underlying natural resources could be affected by denser and harder snowpacks.


Introduction
In the United States snowmobiling accounts for between USD 7 billion (American Council of Snowmobile Associations, 2014) and USD 26 billion (International Snowmobile Manufacturers Association, 2016) in annual revenue, and much of the snowmobile use occurs on public land.The United States National Forest System records about 6 million snowmobile visits annually, accessing about 327 000 km 2 of land (US Forest Service, 2010, 2013a).With continued increases in the number of people participating in winter recreation (Cook and Borrie, 1995;Winter Wildlands Alliance, 2006;US Forest Service, 2010, 2013a;Nagler et al., 2012;Colorado Off-Highway Vehicle Coalition, 2016;Osterberg, 2016), activities including increased snowmobile use may influence snowpack properties in these seasonally snowcovered environments.Of additional concern, is that climate change will result in reduced land available for snowmobiling (Tercek and Rodman, 2016), likely increasing the impact of snowmobile traffic.
There have been limited studies regarding the influence of snowmobile use on snowpack properties (Keddy et al., 1979;Thumlert et al., 2013;Thumlert and Jamieson, 2015).Studies have, however, examined how the snowpack changes due to snow grooming at ski resorts (Fahey et al., 1999;Keller et al., 2004;Spandre et al., 2016a), or to traction and mobility of wheeled vehicles across a snowpack (Abele and Gow, 1975;Shoop et al., 2006;Pytka, 2010).One of the few studies on snowmobile use examined effects on very shallow snow (10-20 cm deep) (Keddy et al., 1979).The authors found a doubling of fresh snow density and a compression of the natural vegetation below the snow (Keddy et al., 1979).Examining deeper snow cover (> 20 cm deep), Thumlert et al. (2013) and Thumlert and Jamieson (2015) examined the distribution of stresses through the snowpack due to type of loading, depth and snowpack stratigraphy (Thumlert et al., 2013).
Changing snowpack conditions from snowmobile use will have other impacts.Aside from the work done by Keddy et al. (1979), there is limited research on how snowmobile activity influences underlying vegetation.The addition of snow due to snowmaking provides an indication of possible changes.Changes from snowmaking include a greater occurrence of soil frost, ice layers may form at the base of the snowpack, and there is often a delay in vegetative growth due to extended snow cover (Rixen et al., 2003).Snowmelt can occur later due to compaction and there is greater heat loss from the densified snowpack and underlying soil, keeping soil temperatures colder longer (Fassnacht and Soulis, 2002;Rixen et al., 2003).
In our research, we specifically examined the effect of snowmobile use on the physical and material properties of the snowpack.The objectives were as follows: (1) to quantify changes to physical snowpack properties due to compaction by snowmobiles; (2) to evaluate these changes based on the amount of use, depth of snow when snowmobile use begins, and the snowfall environment where snowmobiles operate; and (3) to create a simple model to estimate the change in snowpack density due to snowmobile use.This work examines not only changes to the basal snowpack layer, but also to the entire snowpack.The positive economic impact of snowmobiling and increasing winter recreation use from non-motorized activities (such as backcountry skiers, snowshoers, and those on fat bikes) dictates a need to better understand impacts to snow and underlying natural resources in multi-use areas, especially when the information may be used by managers to reduce conflict among recreationists and protect the resource.The Cryosphere, 12, 1121Cryosphere, 12, -1135Cryosphere, 12, , 2018 www.the-cryosphere.net/12/1121/2018/winter recreationists due to the ease of access to backcountry terrain from Colorado Highway 40.Due to heavy use and conflict among users during the winter season, the Forest Service manages Rabbit Ears Pass for both non-motorized and motorized uses.The west side of the pass is designated for non-motorized use and prohibits motorized winter recreation while the east side of the pass is a mixed-use area and is open to motorized use (Fig. 1).This study area was selected to determine if differences in snowpack properties will be observed between the non-motorized and motorized use areas (e.g., Walton Creek versus Dumont Lakes and Muddy Pass in Fig. 1).
Two REP experimental snow compaction study plots were located adjacent to one another within an open meadow north of Colorado Highway 40 at an elevation of approximately 3059 m (Fig. 1).The snow compaction sites were established within an area that prohibits motorized use to protect the study sites from unintended impacts of snowmobilers.Data from the Columbine snow telemetry (SNO-TEL) station, located at an elevation of 2792 m, was used to show how the 2009-2010 winter compared to other winters at REP.The SNOTEL network was established in the late 1970s across the western United States by the Natural Resources Conservation Service to monitor snowpack properties.Initially snow water equivalent and precipitation were monitored, however, temperature and snow depth were added in the 1990s-2000s to aid in operational runoff volume forecasting (see https://www.wcc.nrcs.usda.gov/).Three sites were not experimentally manipulated, i.e., the specific amount of snowmobile use was unknown, and were identified as operational sites along Colorado Highway 40 on REP (Fig. 1 left inset).The "natural" control site was Walton Creek, located west of Rabbit Ears Pass in an open meadow at an elevation of 2895 m within a managed area that prohibits motorized use.Snowshoers, skiers, and snowboarders primarily use this area in the winter to access backcountry terrain.Two sites, Dumont Lakes and Muddy Creek, were located east of REP at an elevation of about 2900 m within an area managed for motorized and mixed uses; the sites were located in open meadows near their respective trailheads (Fig. 1).These trailheads provide backcountry access to snowmobilers.Snowmobile use in the meadows near the trailheads is medium to high, especially on weekends and over holidays (Robert Skorkowsky, personal communication, 2010).The meadow near the Muddy Creek trailhead is more heavily used by snowmobiles than the meadow near the Dumont Lakes trailhead.
Another experimental snow compaction plot was established during the same winter snow season of 2009-2010 at the Fraser Experimental Forest (FEF) near the town of Fraser, Colorado in the Rocky Mountains of central Colorado (Fig. 1).The 93 km 2 experimental forest is a research unit of the United States Forest Service (USFS) Rocky Mountain Research Station (RMRS) located within the Arapaho NF.The FEF snow compaction site was located in a small meadow at an elevation of 2851 m surrounded by lodgepole pine (Pinus contorta) forest.The Fraser Experimental Forest is closed to snowmobile use, but is used to access backcountry terrain by snowshoers, skiers, and snowboarders.The Middle Fork Camp SNOTEL station, located at an elevation of 2725 m, was used to characterize the 2009-2010 winter at FEF.

Experimental snow compaction plots
Snow compaction study plots were established in undisturbed areas at the REP and FEF study areas.Each plot was 22 m wide and 15 m long (Fig. 2a and b).Plots were divided into equal width transects (2 m) and treated with low, medium (FEF only), or high snowmobile use, including a no treatment control transect representing an undisturbed snowpack.Two control transects were used at FEF to represent the undisturbed snowpack (Fig. 2b).Integrating two controls in the FEF study plot allowed for replication and determination of variability.The location of control and treatment plots across each study site were randomly selected.Each transect was separated by a 3 m buffer to eliminate the influence of compaction treatments on adjacent transects (Fig. 2a and b).
Transects were treated by driving a Ski-Doo brand snowmobile weighing about 300 kg including the rider (Fig. 2d) at 10 km h −1 over the length of each transect 5, 25 (FEF only), or 50 times, representing low, medium (FEF only), and high snowmobile use, respectively.Treatments began (Fig. 2c) when non-compacted snow depths were approximately 30 cm (12 inches) for both locations, and when unpacked snow depths equaled approximately 120 cm (48 inches) for REP only (Fig. 2a).Treatments were implemented (Fig. 2e) monthly thereafter, until peak accumulation (Fig. 3).Snowpack sampling was usually performed within a week after each treatment (Figs. 2 and 3).At FEF, snowpack sampling was performed prior to the first treatment to illustrate range of spatial variability across the plots (first set of points in Fig. 4b).

Snow pit analyses and data collection
Snow pit profiles were used to examine the physical properties of the snowpack at both the experimental and operational sites.A vertical snow face was excavated by digging a pit from the snow surface to the ground.Measurements of snow density, temperature, stratigraphy, hardness, and ram resistance were taken vertically along the snowpack profile.Total snow depth was measured from the ground up, and combined with density to yield snow water equivalent (SWE).Physical snowpack properties were compared between non-snowmobile (control) and varying degrees (low, medium (FEF), and high) of snowmobile use (treatment).
www.the-cryosphere.net/12/1121/2018/The Cryosphere, 12, 1121-1135, 2018 Density was measured at 10 cm intervals, from the surface of the snowpack to the ground, by extracting a 250 or 1000 mL snow sample using a stainless-steel wedge cutter (http://snowmetrics.com/) and measuring the mass on an electronic scale with a resolution of 1 g.At least two samples were taken per 10 cm interval.With the 1000 mL wedge cutter, the density of the snow (ρ s in kg m −3 ) was read directly from the scale as the volume of the cutter is 1/1000 of a cubic meter and a gram is 1/1000 of a kilogram.For the 250 mL cutter, the mass measurement results were multiplied by four to obtain density.Snowpack density profiles were created from samples extracted at discrete 10 cm intervals vertically along the working face of the snow pit.The bulk snowpack density was determined by averaging density measurements over the entire depth of the snowpack.A mean of the density measurements for the bottom 10 cm of the snowpack was used to evaluate changes near the snow and ground interface (basal layer).
Temperature measurements were obtained at 5 cm intervals from the top to the bottom of the snowpack using a dial stem thermometer with ±1 • C accuracy.Temperature gradients are well represented by this instrument, and the repeatability of temperature measurements are better than ±1 • C (Elder et al., 2009;American Avalanche Association, 2016).Snowpack temperature profiles and the corresponding bulk temperature gradient were compared.The temperature gradient (T G in • C m −1 ) was calculated as the ratio of the change in temperature ( T in • C) with the distance (d in m) over which the change in temperature occurred.The snowpack temperature gradient was approximated as linear from an upper boundary that was 25-30 cm below the surface to the lower boundary at 0 cm.For this study, the depth below the snow surface where temperature did not fluctuate diurnally was used as the upper boundary to remove bias from diurnal fluctuations (Pomeroy and Brun, 2001).Basal layer temperatures taken at 0 cm were used to compare temperature changes near the snow-ground interface.
Stratigraphic measurements were used to illustrate the evolution of the snowpack over time through characterization of the shape, size, and layering of snow crystals within the snowpack.Classification of grain morphology was based on The International Classification for Seasonal Snow on the Ground (Fierz et al., 2009) and mean grain size was measured and recorded to the nearest 0.5 mm using a hand lens and a crystal card.The crystal forms were identified as precipitation particles, rounded grains, faceted grains, and ice layers.
Hardness is the penetration resistance of the snowpack (Fierz et al., 2009), and is reported as the force per unit area required to penetrate the structure of the snowpack ( McClung and Schaerer, 2006); it is affected by snowpack microstructure and bonding characteristics of the snow grains (Shapiro et al., 1997).Hardness measurements were taken horizontally with a force gauge in each stratigraphic layer using a Wagner Instruments Force Dial gauge (http: //wagnerinstruments.com) with maximum force measurements of 25 and 100 N, and fabricated circular metal plate attachments of 20 cm 2 in area.For each measurement, the circular metal plate was pushed into the snow and the force required to penetrate the snow was recorded.The snow hardness (h i in N m −2 ) for each stratigraphic layer was calculated as the force required to penetrate the snow (F in N) per unit area of the circular metal plate (A in m 2 ).All layers thicker than 5 cm were identified using the 5 cm diameter of the plate.The bulk snowpack hardness (H B in N m −2 ) was determined by weighting each stratigraphic layer hardness measurement by the stratigraphic layer thickness.The hardness associated with the bottom stratigraphic layer for each transect was used to describe hardness changes in the basal layer of the snowpack.The standard ram penetrometer is an instrument with a cone on the end of a rod onto which a hammer of defined weight is dropped from a given height and the depth of penetration is recorded; it was used here to vertically measure the resistance of snow layers to assess the change in ram resistance due to compaction (American Avalanche Association, 2016).A ram profile measurement was taken 0.5 m from the edge of the snow pit wall subsequent to snow pit profile measurements.The mean ram resistance (S B in N) was determined by weighting each ram resistance value obtained from the standard ram penetrometer measurement with the depth sampled.The ram resistance value associated with the bottom layer was measured to describe changes in ram resistance in the basal layer of the snowpack.

Statistical analyses
Data were analyzed using the Mann-Whitney-Wilcoxon rank sum test (Wilcoxon, 1945;Mann and Whitney, 1947).This statistical test is non-parametric and determines whether two independent samples were selected from populations having the same distribution.For this work, the sets of samples compared were density, temperature, hardness, and ram resistance profiles for the five different monthly measurements and the controls (Table 1).A statistical significance was determined for the 95 % (significant) and 99 % (highly significant) confidence interval (p < 0.05, and p < 0.01) and noted with an asterisk in Table 1.

Bulk snowpack density change model
A multi-variate non-linear model was created to estimate the change in bulk snowpack density for various treatments compared to the control (no use) using the following snowpack properties: depth, bulk density, and the number of passes (Fig. 8).The cross-correlation between variables was considered to reduce model over-fitting.The model was calibrated with the experimental data from REP and FEF, and evaluated using data from the operational sites with Walton Creek as the control, Dumont Lakes as medium use, and Muddy Creek as high use.The Nash-Sutcliffe coefficient of efficiency (NSCE, Nash and Sutcliffe, 1970) was used to evaluate the fit of the model.

The measurement winter
The 2009-2010 winter at REP had slightly below average snow depth compared to the 15-year mean, based on the Columbine SNOTEL data averaged from 2003-2017 (Fig. 3a).A peak SWE value of 556 mm on 9 April was 93 % of the historical average.Maximum snow depth measured at the REP snow compaction study plot was approximately 1.5 m and represents a deeper snow cover environment for Colorado.From the Middle Fork SNOTEL data, the 2009-2010 winter snow depth at FEF was also below the 15-year historical average (Fig. 3b).The measured snow depth at the FEF snow compaction study plot never exceeded 1 m, similar to the Middle Fork Camp, and therefore was used to represent a shallower snow cover environment.

Density
The natural variability in density was small at the test sites (Fig. 4).At REP, deep snow (120 cm) compaction treatments were not initiated until after the second sampling date (Fig. 3a), so density for the deep snow low and high use treatments could then be compared to the control; these show minimal difference (Fig. 4).At FEF, there were two sets of control snow pits, and sampling occurred before treatment at all plots (Fig. 3b).These differences in density were greater than those at REP but were still small (Fig. 4).
The mean density values at the FEF plots were almost the same at the end of the sampling period in April (Fig. 5a, column ii).The mean snowpack density increased over the snow season (Fig. 5a), with the exception of the FEF control and at the high use site on 12 February 2010 due to fresh snow are compared to the control on the first two sampling dates (pretreatment, Fig. 3a).At the Fraser Experiment Forest (FEF shown with triangles) the two sets of control snow pits were compared, and all plots were sampled prior to the initial treatment and are compared (see Fig. 5a and b, columns i and ii, respectively).deposition.At the REP snow compaction study site, mean density for high use compaction treatments starting on 30 cm of snow was greater throughout the measurement period than the no use treatment (Fig. 5a, column i, Fig. 6a, columns i, and ii), while the density from low and high use starting on the deeper snowpack of 120 cm was very similar to that measured for no use.The snowpack was more dense for low use on the shallower snowpack (start at 30 cm) than the control for all sampling dates, with the exception of 13 March (Fig. 5a, column i).Density differences are more pronounced for the basal layer (Fig. 5b); for compaction treatments starting at 30 cm, the lowest layers were much more dense than the control (Fig. 6a).Densities for the compaction treatments starting at 30 cm were significantly different than the control and compaction treatments beginning at 120 cm of snow (Table 1a).The density differences between the treatments on the deep snow (120 cm) and the control were not significantly different (Table 1a).
Density increases due to snowmobile use were much greater at Fraser (Fig. 5a, column i and Fig. 5b, column ii) than Rabbit Ears.All treatments at FEF were significantly different than the control, but the difference among treatments was not significant (Table 1a).The density differences among treatments are highlighted in the 10 cm individual density measurements (Fig. 6a) and in the basal layer (Fig. 5b, column ii).

Temperature
Low and high use compaction treatments at the REP snow compaction study site that began on both a shallow snow-pack of 30 cm and on a deep snowpack of 120 cm did not result in significant changes in temperature gradient.The maximum temperature gradients were observed on the earliest sampling date (12 December, Fig. 5c), while they were almost the same for the control, low use, and high use compaction treatments that began on a deep snowpack.Temperature gradients for all treatments decreased throughout the winter season, and were isothermal at 0 • C m −1 by mid to late April (Fig. , 5c, columns i and ii), since the snow had stared to melt (Fig. 3).Overall, temperature gradients were not very different (Fig. 5c) and the variations among treatments were not found to be significant (Table 1b).At FEF, gradients in the high use were greatest after the first treatment and the temperature gradients were essentially the same by March (Fig. 5c, column ii).

Hardness
The snowpack was harder for snowmobile use starting on 30 cm than the control (no use) for both sites (Fig. 5d and e).Mean snowpack hardness did not change much over time (Fig. 5d), except once high use treatments started (6 February) on a deeper snowpack.However, basal layer hardness did decline at REP for both high and low use starting on 30 cm (Fig. 5e, column i).With treatments at FEF, the hardness was always much higher than the control (Fig. 5d, column ii).Hardness initially increased at the REP snow compaction study site following low and high use compaction treatments that began on 30 cm of snow (Fig. 5d, column i), but these were about the same as the control by 17 April, when melt had started.Significant increases in hardness were observed between treatments that began on 30 cm of snow and the control (Table 1c).There was also a significant difference in hardness for deep and shallow initiation depths (Table 1c).In contrast, mean snowpack hardness was not significantly impacted by snow compaction treatments that began on 120 cm of snow (Table 1c).Mean snowpack hardness increased following the initial snow compaction treatments for low starting on 30 cm and high use for starting on both 30 and 120 cm (Fig. 5d, column i).Subsequent compaction treatments did not appear to have a large effect (Fig. 5d, column ii).There were minimal differences by the last sampling date (Fig. 5e, column i).
Snow compaction treatments that began on 30 cm of snow increased basal layer hardness (Fig. 5e, column i), but treatments that began on 120 cm of snow did not impact basal layer hardness (Fig. 5e, column i).For both controls and all treatments that began on 120 cm of snow (Fig. 5e, column i), the maximum basal layer hardness was about 6 kPa.Increased hardness due to snowmobile use showed similar temporal patterns to densification (Fig. 5a and d).At REP, snowmobile use compacted the second layer below the surface, and high use (50 passes) made that layer about 10 times harder than the low use (five passes) snowpack (Fig. 6b, columns i and ii).There was more spatial variability in snowpack hardness (NSCE of 0.50; results not shown graphically) than differences in density (NSCE of 0.93 in Fig. 4) for low and high use compaction treatments versus the control on the first two sampling dates at REP and for the control snow pits at FEF on the pre-treatment date.These larger differences are both attributed to spatial variability, but mostly to the low range of non-treatment hardness values from 0.4 to 5.8 kPa compared to the range of treatment hardness values from 30 to 1157 kPa (Fig. 5d and e).

Ram resistance
Low and high use compaction treatments at REP caused an increase in mean snowpack ram resistance, but the difference was not significant for treatments that began on deep snow (120 cm; Table 1d).After the initial snow compaction treatments mean snowpack ram resistance for low and high use was greater than the control for the entire study period, but by the end of the study period minimal differences were observed between treatments.Basal layer ram resistance increased as a result of low and high use compaction treatments that began on both 30 and 120 cm of snow.Snow compaction treatments at the FEF snow compaction study site caused a significant increase in mean snowpack ram resistance (Table 1d, e.g., Fig. 6c, column iii for the February sampling dates).Basal layer ram resistance increased following the initial snow compaction treatments and continued to increase throughout the duration of the winter season.

Grain size
Smaller crystals in the basal layer were observed for snowmobile use starting on a shallow snowpack compared to the control or starting on a deeper snowpack (Fig. 5f).Rounded grains were observed during the first sampling at REP shallow depth snowmobile start, with faceted grains for the following three sampling dates (Fig. 5f, column i).Rounding facets were observed on the last sampling day at both sites.At FEF, there were 3 to 4 mm faceted crystals prior to the treatments; fragmentation was noted in the faceted crystals found in the basal layer of the treated plots, which began rounding by the last sampling date (Fig. 5f, column ii).The  2009), with the exception of faceted crystals that were fragmented.Note that the snowpack at the low and high use start at 30 cm could not be adequately tested for hardness on the first sampling date at the REP treatment plots.shallower snow at FEF enabled large faceted crystals to grow in the basal layer, up to 9 mm in size (Fig. 5f, column ii).

Operational sites
As illustrated by SWE (Fig. 7d) and snow depth (Fig. 7e), the amount of snow was comparable for the snow pits dug at the three operational sites, even though they were located up to 6 km apart (Fig. 1).Also since these were operational sites, the amount of treatment was not controlled and was based solely on permitted snowmobile use.Generally, patterns of increased density (Fig. 7a), hardness (Fig. 7b), and ram resistance (Fig. 7c) seen at the REP operational sites were similar to the overall patterns seen in the previously presented experiments from REP and FEF (Figs. 5, and 6) with the nonsnowmobile impacted snow pits being less dense (Fig. 7a) and having layers that were less hard (Fig. 7b).From visual inspection of the sites and the measurement results, Muddy Creek had the most snowmobile use and thus exhibited the highest density throughout the winter, and the hardest snowpack for mid-winter (Fig. 7b), but at times the results for Dumont Lakes were similar.

Bulk snowpack density change model
A non-linear bulk snowpack density change model was created using data from the experiments prior to onset of melt conditions (Fassnacht et al., 2010); before the last sampling date (Fig. 3) and prior to when the difference in density between the control and treatments was small (Fig. 5a).Additionally, treatments starting on a deep snowpack at REP were not significantly different than the control (Fig. 5a, Table 1) and were not used in fitting the model.The variables of number of passes per treatment, depth, and bulk density were tested for correlations that might result in model over-fitting.Cross-correlation results were small (R 2 < 0.04), so these variables were used to create the model.Difference in bulk density compared to the control due to snowmobile use is a www.the-cryosphere.net/12/1121/2018/The Cryosphere, 12, 1121-1135, 2018 function of the number of passes per treatment and bulk density, but it is inversely related to snow depth (Fig. 8a).The optimal model had a NSCE of 0.81 (Fig. 8a), which is considered very good (Moriasi et al., 2007).The model was calibrated on the experimental data (Fig. 8a) and applied to the operational sites (Fig. 8b), with no passes occurring equivalent to a density change of 0 kg m −3 .The evaluation results were less optimal, with a NSCE of −0.79 for the four dates tested in December through March (Fig. 8b).The poorer performance of the model at the operational sites is due to an unknown number of snowmobile passes at each site and from limited snowmobile use early in the season (December), resulting in minimal differences between compaction levels at that time (Figs.7 and 8b).Removing the December data points and using only the January through March dates improved the model fit to a NSCE of 0.34 (Fig. 8b).

Observed changes to snowpack properties
Snowpack changes were observed for varying snowmobile use beginning with two different snow depths (REP only in Fig. 5 or 6, columns i and ii) and for two different snowcovered environments (Figs. 5 and 6).A total of 101 snow pits (50 at REP, 15 at the operational sites, and 36 at FEF) were dug and sampled for this work.The increase in density and hardness from snowmobile use is greatest compared to an untreated snowpack in early to mid-season (January) for a deeper snowpack at REP, with density increases of 7-33 % and hardness 4 to 13 times greater than the control (Fig. 5a  and d, column i).For a shallower snowpack at FEF, density increased by 64-76 % and hardness was 500-2000 times greater than the control (Fig. 5a and d, column i).Similar differences were found from ski run grooming in an Australia snowpack with a 400 % increase in hardness early in the snow season but only about a 40 % increase later in the winter (Fahey et al., 1999).Snow grooming increased the average density by up to 36 % compared to non-groomed ski slopes (Fahey et al., 1999;Rixen et al., 2001).
At rest, a snowmobile and its rider exert 4 to 10 kPa of pressure on the underlying snowpack.This assumes a track length from 0.9 to 1.4 m, width of 0.50 m, a snowmobile weight of 200 to 350 kg, and a rider weight of about 100 kg (data from http://www.polaris.com/).There is an increase of less than an order of magnitude due to snowmobile movement.Thumlert et al. (2013), measured stresses of about 10 to 20 kPa at a depth of 30 cm below the surface of a deep snowpack.At 20 cm below the snow surface, the induced stress from a snowmobile is already much less than 10 cm below the surface (Thumlert et al., 2013).Grooming vehicles add a force similar to snowmobiles (Pytka, 2010) based on mass and track size.The snowpack property changes we observed could, therefore, also be representative of impacts from both types of vehicles.Snowpack loading by wheeled vehicles on a shallow snowpack was much greater than that of a snowmobile, peaking at about 350 kPa (Pytka, 2010).In comparison, fresh snow with a density of 100 kg m −3 exerts a pressure of 0.003 kPa on the underlying snowpack (Moynier, 2006).
Compaction due to snowmobile use increased density of the snowpack which influenced snow hardness (Fig. 5d  and e) and ram resistance (Fig. 6c).Compaction altered snow characteristics (Figs.5-7), fragmented faceted grains (Fig. 5f, column ii), and reduced the growth of faceted grains (Fig. 5f).Density measurements for fresh snow (Fassnacht and Soulis, 2002) and/or uncompacted snow (López-Moreno et al., 2013) vary spatially and temporally (Fig. 4), these values can double with just one pass of a snowmobile on a very shallow snowpack (Keddy et al., 1979).The snowpack properties of a shallow snow environment can be more greatly affected by compaction from snowmobile use than those for an area that receives more snow (e.g., Fig. 3b vs. Fig.3a).With more snow accumulation, density also increases, but high levels of snowmobile use will tend to increase the density above what is observed with non-snowmobile impacted snow (Figs.5-7).

Limitations of the measurements
Although snowpack variability over space was limited (Fig. 4), the properties of the snowpack change from site to site and through time.For example, the mean snowpack density was less in February (Fig. 6) than January at FEF (Fig. 5  column ii).From the operational sites, specific hard layers and high values of ram resistance were measured that did not persist until the next monthly sampling (observed in the experimental treatments -not shown graphically).These variations were possibly a combination of naturally occurring spatiotemporal snowpack variability and sampling errors; it can be difficult to obtain reliable hardness measurements in snow disturbed by snowmobiles.Future investigations could focus on specific aspects of this study, such as using a finer temporal resolution, but with fewer treatments.Another source of variability or bias is the type of equipment used for sampling.Density and temperature were measured at 10 cm intervals using the Snowmetrics wedge cutter and dial gauge thermometers.A different sampler could be used to measure the density over each layer and other types of thermometers could be used.Snow-hardness gauges and circular metal plates of known area were used for hardness testing (McClung and Schaerer, 2006), rather than the more simplistic in situ hand hardness test (American Avalanche Association, 2016).However, the hardness of thin layers could not be measured as the circular metal plate used for measurements had a diameter of 5 cm, omitting the possible measurement of these thin layers.Thus, bulk hardness was possibly underestimated.Also, due to the compaction of the snow grains by the high use 30 cm start treatment at REP the hardness could not be measured (Fig. 5, column ii).Different equipment may resolve this issue.

Significance of the changes to snowpack properties from snowmobile use
Snowmobile use was found to have a highly significant effect upon natural vegetation below the snow (Keddy et al., 1979), and by extension from snowmaking as well (Rixen et al., 2003).Ski grooming has been shown to delay the blooming of alpine plants (Rixen et al., 2001) due to later snowmelt and significantly cooler soil temperatures (Fassnacht and Soulis, 2002).Deeper snowpacks were found to not have cooler soil temperatures under the snowpack (Keller et al., 2004), but melted out four weeks later than thinner snowpacks (Keller et al., 2004).Since the changes due to snowmobile traffic on a shallow snowpack were significant (Table 1), the effects of snowmobile use on the soil and vegetation underlying a shallow snowpack should be further investigated.Snowmobile use, starting with a shallow or thin (30 cm) snowpack, resulted in a denser and harder snowpack with a decrease in grain size throughout the season, and rounded crystals or facets observed with the last measurements (Fig. 5f).If compaction penetrates deep enough into the snowpack, it could affect weak layers that cause avalanches (Saly et al., 2016), which are typically composed of soft layers consisting of large facetted grains (e.g., Schweizer and Jamieson, 2003;van Herwijnen and Jamieson, 2007).While this may be useful in very limited and small areas, such as that performed in boot packing programs (e.g., Sahn, 2010) to strengthen snowpacks likely to fail on basal facets, it is very difficult to properly align and reproduce the intensity of repetitive tracks, as done experimentally here (Fig. 2).The effects of snowmobile use for avalanche hazard reduction The Cryosphere, 12, 1121-1135, 2018 through changing snow stability properties requires more investigation.
Other factors acting in concert with snowmobile traffic to affect snowpack properties include wind, snowmaking/grooming, and a changing climate.Without the effects of wind, snow depth will generally be lower for areas with snowmobile traffic (Figs. 2d, e, and 7;Rixen et al., 2001;Spandre et al., 2016a).However, wind is often present in open areas where snowmobiling occurs.Local terrain features and position and extent of canopy cover influence how the wind interacts with the snowpack (Pomeroy and Brun, 2001).In an Australian case study, SWE increased by 45 % in groomed areas (Fahey et al., 1999); at the Rabbit Ears Pass recreational use areas, SWE also increased through time (Fig. 7d) likely due to snow blowing into the depressions created by snowmobile tracks (Fig. 2d); this increased load could further impact the underlying snowpack properties.Further, snowmaking (Spandre et al., 2016a) to supplement natural snow conditions and /or grooming (Fahey et al., 1999;Rixen et al., 2001;Spandre et al., 2016a) compacts the snowpack below it, and alters the underlying snowpack properties (Howard and Stull, 2014;Spandre et al., 2016a, b).Also, a changing climate will likely reduce the extent of snowcovered terrain and decrease the length of the winter recreation season (Lazar and Williams, 2008;Steiger, 2010;Dawson and Scott, 2013;Marke et al., 2015;Schmucki et al., 2015;Tercek and Rodman, 2016;Marty et al., 2017).In addition to possible effects from a changing climate, inter-annual variability of snowpack patterns can be large in Colorado (Fassnacht and Hultstrand, 2015;Fassnacht and Records, 2015;Fassnacht et al., 2017).The effects of this variability should be included in long term motorized use land management considerations.
The significant change to snowpack properties by snowmobiles, except when treatments/use were initiated on a deep snowpack (Table 1), could impact land management decisions for multi-use public lands.The measured depth of influence for a snowmobile is about 90 cm according to work done by Thumlert et al. (2013), but additional work could test starting depths such as 30, 60, and 90 cm in differing snow conditions to identify the depth when snowmobile use has no significant impact.Most ski resorts in the French Alps required a minimum snow depth of 40 cm to offer skiing, with a range from 60 cm in February to 40 cm in April (Spandre et al., 2016b).The US Forest Service (2013b) recommends a minimum of 30 cm before the use of snowmobiles.Increasing the minimum snow depth before allowing snowmobile traffic will reduce changes to the snowpack due to snowmobile use (Table 1).Additionally, the non-linear bulk density change model developed here and applied to operational sites could be used predictively for management needs.This model may be useful in terms of estimating when to limit snowmobile use given changes in specific snow depth and density conditions.
Where the experiments for this study were undertaken, on public lands in Colorado, there are 1.1 to 1.6 million annual snowmobile visits, with an increase from 580 thousand to 690 thousand between 2010 to 2013 in northern Colorado (Routt NF and Arapaho-Roosevelt NF) and southern Wyoming (Medicine Bow NF) (US Forest Service, 2010, 2013a) alone.The annual economic impact of snowmobile use is more than USD 125 million to each state (Nagler et al., 2012;Colorado Off-Highway Vehicle Coalition, 2016).Snowmobile use is likely to continue to increase, and economic gains need to be balanced with potential impacts to the landscape, particularly in those times and places where snowpacks are shallow.

Conclusion
Snowmobiling is a multimillion dollar industry that impacts local and regional economies and public recreation lands.There have been limited studies regarding the influence of snowmobile use on snowpack properties.We examined the effect of snowmobile use on the physical and material properties of the snowpack at sites with varying snowmobile use and seasonal snow conditions.Low, medium, and high snowmobile use was simulated on experimental transects and snowpack sampling results from the treated sites were compared to the snowpack properties observed at undisturbed control sites, and at operational sites with varying levels of use.A non-linear bulk snowpack density change model was developed relating changes in bulk density to snowmobile use as a function of the number of passes, snow depth (inverse relation), and bulk density.The largest differences in snowpack properties occur with snowmobile use beginning on a shallow snowpack (30 cm) compared to no use, which increases snowpack density, hardness, and ram resistance.These increases are directly related to increasing snowmobile use (from low to medium to high).Conversely, snowmobile use that begins on a deep snowpack (120 cm) has a limited effect on the snowpack properties of density, temperature, hardness, and ram resistance as compared to an undisturbed snowpack.These results suggest that from a management standpoint, it may be desirable to limit snowmobile use in shallower snow conditions to avoid increases in density, hardness, and ram resistance that could possibly impact land resources below the snowpack.
written by JTH, SRF, and KJE.The final version of the manuscript was written by SRF and NBHV.SRF generated the figures and created the density model.
Competing interests.The authors declare that they have no conflict of interest.
Acknowledgements.Appreciation goes to Robert Skorkowsky, Kent Foster and Becky Jones of the Hahns Peak/Bears Ears Ranger District of the US Forest Service for their help and support with compaction treatments at the Rabbit Ears Pass study site.Additional thanks goes to James zum Brunnen of the Colorado State University Statistics Department for his assistance with statistical interpretation.Jared Heath would also like to recognize the Colorado Mountain Club for their help supporting this project with a generous grant.Jim Halfpenny, Ned Bair, and two anonymous reviewers provided insight into clarifying this paper.One TC Discussion anonymous reviewer provided very thorough and thoughtful comments that greatly improved this paper, and resulted in the creation of new figures.While the comments from this reviewer provided a challenge, they were appreciated after they had been addressed.TC editor Guillaume Chambon provided additional comments and an important citation that helped reformulate the discussion.
This open-access publication was funded by the University of Göttingen.
Edited by: Guillaume Chambon Reviewed by: Edward Bair and one anonymous referee

Figure 1 .
Figure 1.The snow compaction study plots are located in north-central Colorado.The Rabbit Ears Pass (REP) site is within the Routt National Forest near the town of Steamboat Springs, as are the three operational (non-experimentally manipulated) sites (Walton Creek with no use, Dumont Lakes with low to medium use, and Muddy Pass with high use based on field observations).The Columbine snow telemetry (SNOTEL) station was used to identify the amount of annual snowfall in 2009-2010 compared to the long-term average.The Fraser Experimental Forest (FEF) site is within the Arapaho-Roosevelt National Forest near the town of Fraser.The Middle Fork Camp SNOTEL site was used to represent the year's snowfall.
During the 2009-2010 snow season a set of snow compaction plots were located near Rabbit Ears Pass (REP) in the Rocky Mountains of northern Colorado to the southeast of the town of Steamboat Springs.REP is within the Medicine Bow-Routt National Forest (NF; Fig. 1) along the Continental Divide encompassing over 9400 km 2 of land in Colorado and Wyoming.Rabbit Ears Pass is especially popular during the winter season and is heavily used by snowmobilers and other

Figure 2 .
Figure 2. The sampling design for the snow compaction plots at (a) Rabbit Ears Pass, (b) Fraser Experimental Forest, and photographs of the study plots (c) pre-treatment, (d) during treatment, and (e) after treatment.The colors used for the control and treatment plots are used in Figs.5-8.

Figure 3 .
Figure 3. Mean snow depth from 2003-2017, and for the 2010 water year (WY2010) measured at (a) the Columbine SNOTEL site near Rabbit Ears Pass (REP), Colorado and (b) the Middle Fork Camp SNOTEL near Fraser Experimental Forest (FEF), Colorado, illustrating the dates of treatment and dates of sampling.Data were obtained online from the Natural Resource Conservation Service (NRCS) National Water and Climate Center (http://www.wcc.nrcs.usda.gov/).

Figure 4 .
Figure 4. Spatial variability of mean (yellow) and basal (blue) snowpack density.At the Rabbit Ears Pass (REP shown with circles) the deep snow (120 cm) compaction treatments (low and high use)are compared to the control on the first two sampling dates (pretreatment, Fig.3a).At the Fraser Experiment Forest (FEF shown with triangles) the two sets of control snow pits were compared, and all plots were sampled prior to the initial treatment and are compared (see Fig.5a and b, columns i and ii, respectively).

Figure 5 .
Figure 5.Time series for (i) Rabbit Ear Pass (REP) and (ii) Fraser Experimental Forest (FEF) at the different sampling dates of (a) mean snowpack density, (b) basal snowpack density, (c) snowpack temperature gradient, (d) mean snowpack hardness, (e) basal layer hardness, and (f) mean basal crystal size and shape.The crystal shape is included as per Fierz et al. (2009), with the exception of faceted crystals that were fragmented.Note that the snowpack at the low and high use start at 30 cm could not be adequately tested for hardness on the first sampling date at the REP treatment plots.

Figure 6 .
Figure 6.(a) Density, (b) hardness, and (c) ram resistance profiles for the February sampling dates (6 February at REP and 12 February at FEF) measured at the REP snow compaction study plot for no (control), low, and high use treatments beginning on (i) 30 cm and (ii) 120 cm of snow, and (iii) the FEF snow compaction study plot for no (control), low, medium, and high use treatments beginning on 30 cm of snow.Note that free floating measurements represent overlapping density measurements.The ground is at zero snow depth.

Figure 7 .
Figure 7. Snow-pit data for Walton Creek (no snowmobile use), Dumont Lakes (moderate snowmobile use), and Muddy Creek (high snowmobile use) in the Rabbit Ears Pass recreational use areas illustrating (a) density, (b) hardness, (c) ram resistance, (d) SWE, and (e) snow depth.For (a)-(c), column (i) is the mean snowpack value and column (ii) is the basal layer value.

Figure 8 .
Figure 8. Bulk snowpack density change model for different amounts of use compared to the control of no use (a) calibrated for the two experiment sites (Rabbit Ears Pass, REP and Fraser Experimental Forest, FEF), and (b) applied to the operational sites (Dumont Lakes and Muddy Creek), compared to the no use Walton Creek site.The calibrated model is presented in (a) with the Nash-Sutcliffe coefficient of efficiency (NSCE).The NSCE is presented in (b) for two different time periods: the four pre-melt dates (December-March -4 dates) and the later three pre-melt dates (January-March -JFM).

Table 1 .
Statistical differences (p values) between no snowmobile use (control) and varying snow compaction treatments on snowpack properties at the study plots located at Rabbit Ears Pass (REP) and Fraser Experimental Forest (FEF), Colorado during the 2009-2010 winter season for (a) density, (b) temperature, (c) hardness, and (d) ram resistance.
* * < 0.01 * * < 0.01 * * Statistically significant differences at the p < 0.05 confident level are denoted with an asterisk, and highly significant (p < 0.01) difference are denoted with two asterisks.