Snow physical properties may be a significant determinant of lemming population dynamics in the high Arctic

Cyclic population fluctuations are common in boreal and Arctic species but the 23 causes of these cycles are still debated today. Among these species, lemmings are Arctic 24 rodents that live and reproduce under the snow and whose large cyclical population 25 fluctuations in the high Arctic impact the whole tundra food web. We explore, using 26 lemming population data and snow modeling, whether the hardness of the basal layer of 27 the snowpack, determined by rain-on-snow events (ROS) and wind storms in autumn, 28 can affect brown lemming population dynamics in the Canadian high Arctic. Using a 7-29 year dataset collected on Bylot Island, Nunavut, Canada over the period 2003-2014, we 30 demonstrate that liquid water input to snow is strongly inversely related with winter 31 population growth (R 2 ≥ 0.62) and to a lesser extent to lemming summer densities and 32 winter nest densities (R 2 = 0.29 to 0.39). ROS in autumn can therefore influence the 33 amplitude of brown lemming population fluctuations. Increase in ROS events with 34 climate warming should strongly impact the populations of lemmings and consequently 35 those of the many predators that depend upon them. Snow conditions may be a key factor 36 influencing the cyclic dynamics of Arctic animal populations. 38 Graphical Abstract: Lemming burrows in an Arctic snowpack

temperature fluctuates around 0°C. This may form melt-freeze layers which, even though 90 they do not resemble ice layers found in wetter climates, might still impede lemming 91 movement and access to food. The occurrence of these events in autumn therefore needs 92 to be considered together with wind effects to evaluate the hardness of the basal snow 93 layer. 94 We evaluated the hypothesis that the hardness of the basal snow layer is a 95 determinant of lemming population dynamic in the High Arctic. We predicted that a hard 96 bottom snow layer due to rain-on-snow events and/or wind storms in autumn negatively  brown lemmings because they are far more abundant than collared at our study site, 112 especially in years of peak population (Gauthier et al., 2013;Fauteux et al., 2015). three pictures a day of mesic habitat. Given these data sources, the number of 123 independent snow start dates that could be obtained ranged from 1 (in 1993) to 6 (in 124 2014). 125 Regarding satellite images, we used data from MODIS, Globsnow-2 and Landsat. 126 Satellite images give more representative data but cloud cover is often a problem Track Scanning Radiometer (A-ATSR) and ATSR-2 sensor at 0.01° spatial resolution for 135 years 1997135 years and 1998135 years (Metsamaki et al., 2015. This dataset was complemented by 30 m 136 spatial resolution true color images from Landsat satellites (5, 7 and 8) for years 1993 to 137 2014. Landsat data were collected from the USGS web site 138 (https://earthexplorer.usgs.gov/).

139
Despite these multiple data sources, the snow start date was sometimes far from 140 certain and we attributed a reliability index for the snowpack onset date of each year 141 based on the coincidence between dates given by satellites, snow gauges and time-lapse 142 photographs. This index was 100% when at least 4 dates were identical and there was no 143 conflicting date. The index was 90% when 3 dates were identical and a 4 th date was most 144 likely the same, and there was no conflicting data. This happened for example when a 145 signal for a snow gauge was low, or when a satellite image indicated significant but 146 incomplete snow cover. The index was 80% when 3 dates were identical with no 147 conflicting dates (typically some methods yielded no date because of persistent cloud 148 cover in satellite data). An index lower than 80% indicated either fewer than 3 reliable we therefore consider this date 100% certain. Snowpack onset dates as well as the dates 157 when snow height reached 5 cm as determined by Crocus are reported in Table 1.

159
In the absence of field measurements in autumn, we evaluated the hardness of the 160 basal snow layer using a snow physics model. Detailed snowpack models initially 161 developed for avalanche hazard forecasting such as Crocus (Vionnet et al., 2012) are 162 appealing as they have a strong focus on snow mechanical properties. Our approach was 163 to quantify the main processes contributing to snow hardness, wind drift and the amount 164 of liquid water in snow, and to use both these variables as proxies for snow hardness. We 165 used Crocus to quantify drifting snow as implemented by Vionnet et al. (Vionnet et al.,166 2012). Essentially, the procedure calculates a unitless driftability index D i that depends 167 on wind speed and type of snow at the surface of the snowpack. If D i >0, then snow is 168 wind-drifted and Crocus simulates fragmentation and compaction of surface snow by 169 wind-induced snow transport. The higher the value of D i , the more the snow is drifted 170 and consequently the denser and harder the snow formed in drifts will be (Kotlyakov, 171 1961;Domine et al., 2011). Regarding the impact of melting, our approach was to 172 quantify the amount of liquid water, W i , in kg m -2 (equivalent to mm of column water) 173 formed in the snowpack, either by melting of deposited snow because of warm air 174 temperature or strong radiation, or brought by rain (Vionnet et al., 2012). Since liquid 175 water is preferentially found at junctions between snow grains (Colbeck, 1973), the more 176 liquid water is present, the stronger the bonds between grains will be when the water 177 refreezes and therefore the harder the snow will be.
D r a f t radiation, precipitation, specific humidity). When wind speed and air temperature were Lemming population data 202 We estimated brown lemming summer densities using a capture-mark-recapture The population growth rate over the winter was obtained by dividing the June 226 lemming density of year t+1 by the August density of year t on the natural log scale (due 227 to some zero values, we added 0.05 to the denominator, as 0.05 roughly corresponds to 228 half the lowest density that could be estimated on our grids). In June 2010, lemming 229 density was likely underestimated due to a persistent snow cover on our live-trapping 230 grid. To correct for that, we used the average of the June/July density estimates for our  Table 2 show the values of snow hardness indices SD i and SW i and lemming 260 population data for years where both data types were judged sufficiently reliable. The 261 quality of the snow data increased over the years because of the increase in sources of 262 data to determine the snow onset date, so that the analysis includes mostly recent years 263 starting in 2004. Likewise, lemming population monitoring methods improved over time.

264
In particular, live-trapping data are the most reliable (Fauteux et al., 2015) and these data 265 were obtained starting in 2004. We therefore concentrated on lemming data obtained over  (Table 1) 293 was accompanied by winds reaching 10 m s -1 which did not produce the expected 294 vigorous drifting and the resulting snow hardening (Fig. 3). It is noteworthy that this 295 intense ROS event resulted in the lowest lemming populations recorded (2013 data in 296 Table 2) even though the wind contributions appears as zero (2012 SD i data in Table 2, 297 and Fig. 3). 298 We therefore propose that it is not possible to evaluate the effect of wind 299 independently from that of ROS, since both sometimes happen at the same time and ROS 300 masks wind effects, while the reverse is not true. The presence of liquid water thus 301 appears to be the dominant factor leading to snow hardening in autumn at our study site.

302
Melting due to warm air or strong radiation could lead to more complex considerations 303 but here all liquid water occurrences were caused by ROS events. Based on these 304 considerations, we now explore the relationship between the occurrence of liquid water in 305 the snowpack in autumn and brown lemming demographic variables.

306
The winter population growth rate is probably the variable we expect to be the        Table 3. R 2 with p values in parentheses for the relationships between the 6 brown 528 lemming population variables of Table 1 and snow variables SDi (wind drift) and SWi 529 (liquid water formation) that promote the formation of hard snow in the basal layer.

530
Values that were significant or nearly so (p  0.07) are in bold.    D r a f t June year t+1) in two summer habitats (wet and mesic) and total nest density sampled at snow-melt in two wintering habitats (mesic and gullies). Standard errors are indicated for lemming and nest counts.