Effects of large herbivore grazing on relics of the presumed mammoth steppe in the extreme climate of NE-Siberia

The Siberian mammoth steppe ecosystem changed dramatically with the disappearance of large grazers in the Holocene. The concept of Pleistocene rewilding is based on the idea that large herbivore grazing significantly alters plant communities and can be employed to recreate lost ecosystems. On the other hand, modern rangeland ecology emphasizes the often overriding importance of harsh climates. We visited two rewilding projects and three rangeland regions, sampling a total of 210 vegetation relevés in steppe and surrounding vegetation (grasslands, shrublands and forests) along an extensive climatic gradient across Yakutia, Russia. We analyzed species composition, plant traits, diversity indices and vegetation productivity, using partial canonical correspondence and redundancy analysis. Macroclimate was most important for vegetation composition, and microclimate for the occurrence of extrazonal steppes. Macroclimate and soil conditions mainly determined productivity of vegetation. Bison grazing was responsible for small-scale changes in vegetation through trampling, wallowing and debarking, thus creating more open and disturbed plant communities, soil compaction and xerophytization. However, the magnitude of effects depended on density and type of grazers as well as on interactions with climate and site conditions. Effects of bison grazing were strongest in the continental climate of Central Yakutia, and steppes were generally less affected than meadows. We conclude that contemporary grazing overall has rather limited effects on vegetation in northeastern Siberia. Current rewilding practices are still far from recreating a mammoth steppe, although large herbivores like bison can create more open and drier vegetation and increase nutrient availability in particular in the more continental Central Yakutian Plain.


Study area
Climate in Yakutia is highly continental, characterized by very low winter temperatures and relatively high summer temperatures ( Figure S1_1). Temperatures are less extreme in Chersky, which is closer to the Arctic Sea. Here, summer precipitation of July and August is well above the temperature curve in the standard Walter-Lieth climate diagram, while in Yakutsk and Verkhoyansk conditions are closer to summer droughts. Figure S1_1. Climate diagrams of Yakutsk, Verkhoyansk and Chersky, Sakha, Russia. Modified from climate-data.org.
In Yakutia, typical steppe vegetation (Cleistogenetea squarrosae) is only found at sites with special microclimatic conditions (Reinecke et al., 2017): on more or less steep, SW-exposed slopes. Tundra steppes (Carici rupestris -Kobresietea bellardii) were restricted to soil-disturbed hilltops north of the tree line (near Pokhodsk). The number of available sites to study extrazonal steppe vegetation was especially restricted in Chersky, and all of these were out of reach of the large grazers of the Pleistocene Park nearby.
Productivity of vegetation as approximated by above-ground standing biomass (dry weight) in our study regions (see Figure S1_2) ranges between 5 and 62 g/ 40x40 cm²; thus 3.1 to 38.8 g/ m² or 30-400 kg/ ha, with lowest values in the most continental Yana region. Productivity also depends much on habitat type, with meadows and wetlands being most productive and a wide range of productivity values being found in steppes, depending on whether they belong to densely vegetated meadow steppes or sparsely vegetated typical steppes on steep slopes. Figure S1_2. Productivity of vegetation across the study area; given as mean weight of harvested plant biomass extrapolated from subplot of 40x40 cm² per region and vegetation type.

Field sampling
Cattle and horses roamed freely around villages and were not fenced in or herded, thus mimicking natural grazing conditions. Droppings, trails and resting places indicated regular use of study sites by grazers, but at least during our study period floodplain meadows seemed to be the preferred pastures in all regions. Scrub and forest were mostly used as resting places or were frequented during roaming between pastures and had droppings along the way. Small clearings (from cutting, fire or tree fall) with higher herb cover were also occasionally grazed.
Quantification of grazing intensities in this situation was difficult, which is a common problem in rangeland ecology whenever it comes to open range or even mobile grazing systems. Herbivore density of free roaming livestock would be somewhere between that of the Bisonary (where it is highest with 17,7 t/ km²) and the large enclosure of the Pleistocene Park (where it is lowest with 1,0 t/ km²) in the larger vicinity of villages.
Thus, we decided to base our analysis on dung cover. Dung was found to be spread out evenly across the grazed areas and seemed to best reflect the intensity of grazing of different sites ( Figure S1_3).

Figure S1_3
: Scatterplots of dung cover in % per grazing herbivore (bison, horse, small mammals, other large herbivores like cattle, musk ox and moose) in all plots, grouped by study region. Figure S1_4: Scatterplots of dung cover in % per grazing herbivore (bison, horse, small mammals, other large herbivores like cattle, musk ox and moose) in all plots versus slope inclination, grouped by study region.
We measured 22 grazing-related plant functional traits; (  Mongolia (Jigjidsuren and Johnson, 2003). In a few cases (<10 values), trait data was missing and we thus had to use nearest neighbor imputation to fill out single missing values.  We collected biomass and soil samples from each vegetation plot. In grassland habitats (meadow, steppe, tundra steppe) three subplots with a size of 40 x 40 cm² were randomly selected across the vegetated area of each plot to account for spatial variability. We took one soil and one plant biomass sample per subplot in these habitats. In forest, scrub and tundra habitats, only one soil sample and no biomass samples were taken per plot. Aboveground biomass of each subplot was cut approximately 1 cm above ground using scissors, excluding dead standing biomass, and then air dried. We sampled the topsoil below the litter layer using a 100 cm 3 core cutter.
Slope inclination was estimated in the field as percent inclination. Slope aspect (N, NE, E, SE, S, SW, W, NW) was measured in the field, using a Garmin GPS (Garmin GPSMAP 64s) by walking straight downhill for a few meters until the direction was reliably given. In addition, we checked google earth maps of the location for correctness. We then derived northernness and easternness from aspect in degrees (360°; with 0°=N, 90°=E) by taking the cosine (northerness) and sine (easterness), thus transforming the degrees to a value between -1 and 1.
Heat load was calculated according to McCune (2007). This variable estimates temperatures on a land surface, based on the amount of potential direct incident radiation (DIR), which depends on latitude, slope aspect and slope inclination, while taking into account the time of the day that surface is subjected to this radiation.
The intensity of grazing by each grazing animal (bison, horse, cattle, small mammals) was estimated based on the density of droppings within each vegetation plot (in %). Small mammals were mostly represented by ground squirrels (Urocitellus parryii). Defecation is considered part of a grazing effect (fertilization), apart from the actual intake of plant biomass. Other proxies, which have proven useful in other studies, for example in Mongolia and Tibet, like distance to town or water well, did not work in our setting, as livestock was not bound to settlements (except cattle to some degree). We also estimated approximate intensity of grazing in the field, but the simple index (high, medium, low) we developed from this information did also not prove useful. Instead, dung density was crucially evaluated already in the field and found to be the best approximation for grazing intensity, even when considering animal movement. Other studies have also shown that dung density is a useful indicator (e.g. Wang et al. 2018).
Macroclimatic variables (Bio 1/ 7/ 10/ 12/ 15/ 18/ 19) were extracted from WorldClim (Hijmans et al., 2005). Annual Mean Temperature (Bio 1) and Annual Mean Precipitation (Bio 12) give basic information on climate; Mean Temperature of Warmest Quarter (Bio 10) and Precipitation of Warmest Quarter (Bio 18) give information on vegetation-relevant summer conditions; and Temperature Annual Range (Bio 7) and Precipitation Seasonality (Bio 15) give information on seasonal differences in climate, thus its continentality. We used GPS coordinates of plots to extract the spatially explicit climate data from the WorldClim model.

Sample processing
Soil samples were initially dried in the lab for 48 hours at 40°C. Samples were then sieved using a 2 mm coarse screen, using the fine material for further analysis. We measured pH (H2O) and electric conductivity (EC) after 1h and 24 h. We measured the C/N ratio through combustion in a CN analyzer (Vario Elementar). To assess the amount of plant available nutrients (Ca, Mg, K, P) we prepared soil extractions following the Olsen P method (Sims, 2000). Nutrient contents in these extractions were measured by spectrometry (ICP-OES, Institute of Soil Science, Hannover University). Rest water was measured after drying of samples at 105°C for 24h. The carbonate content was first assessed with a quick test using 10% HCl, and samples showing a reaction were further analyzed using a calcimeter following Scheibler's method (ON L 1084(ON L -99, 1999. Rest water content was used to calibrate nutrient contents per g soil and carbonate content to correct C/N measurements.
Plant biomass was cut into pieces of 1-3 cm length using ceramic scissors and then separated about 2 (1-3; depending on amount of plant material) times using a dividing cross. A mixed sample of the biomass was then finely ground (Leuphana University of Lüneburg; Umweltanalytisches Labor, IHI Zittau, University of Dresden   Table 1 for abbreviations. and is shown. See Table 1 for abbreviations.