Study Areas

This research was conducted from 2009 to 2011 within selected portions of four active US Department of Agriculture (USDA) Forest Service (USFS) cattle grazing allotments located in western Idaho. Pastures scheduled to be occupied by cattle during the 15 June to 15 August time period of each year, according to USFS grazing management plans, formed the study area within each allotment. These study areas (n = 4) ranged in size from 48 to 112 km² and were selected to be representative of the variability in topography, parent materials, soil types, vegetation cover types, hydrology, climate, and livestock management (e.g., herd composition, breeding, calf age at entry), which typically exists in extensive, mountainous cattle grazing lands of the NRM region. Selection was also based on gray wolf presence. Well-established wolf populations resided in and around all four study areas before the study (Nadeau et al. 2008), and monitoring during the study indicated wolf presence reached moderate or high levels in each study area during each study year (see definitions later).

Study Area 1 (112 km²) was selected to typify a management situation where cattle enter the grazing area on 15 June with relatively young calves (2–2.5 mo of age). The dominant landform of this study area was a dissected plateau with steep-walled canyons draining down the bounding slopes. Elevation ranged from 561 to 1 834 m. Slope ranged from 0 to 57 degrees with a mean of 17 ± 11 standard deviation (SD) degrees. Cattle entered this study area on the lower bounding slopes, climbed up the slopes as forage plant phenology progressed, arrived atop the plateau in early July, and then spent the remainder of the study period (until 15 August) grazing among the hillslopes and shallow drainages on the plateau.

Riparian vegetation along drainages at the lowest elevations of Study Area 1 was dominated by willow (Salix sp. L), sedges (Carex sp. L.), and rushes (Juncus sp. L) with Kentucky bluegrass (Poa pratensis L.) and cheatgrass (Bromus tectorum L.) occurring on the stream terraces. Steep walls of the plateau were vegetated by bluebunch wheatgrass (Pseudoroegneria spicata [Pursh] A. Love) and Idaho fescue (Festuca idahoensis [Elmer]) associations with perennial forbs such as arrowleaf balsamroot (Balsamorhiza sagittata [Pursh] Nutt.), parsnipflower buckwheat (Eriogonum heracleoides Nutt.), Cusick's milkvetch (Astragalus cusickii A. Gray), and Snake River phlox (Phlox colubrine Wherry & Constance) occasionally occurring as codominants with the bunchgrasses (Johnson and Simon 1987). The top of the plateau was vegetated by pine savanna or open pine woodlands. Ponderosa pine (Pinus ponderosa Lawson & C. Lawson) and bunchgrasses (e.g., Idaho fescue) dominated the savannas. The open woodlands included a ponderosa pine overstory, a shrub layer of common snowberry (Symphoricarpos albus [L.] S.F. Blake) and/or white spirea (Spiraea betulifolia Pall.), and an herb layer of pinegrass (Calamagrostis rubescens Buckley) and Geyers sedge (Carex geyeri Fernald) or Idaho fescue. Ridge-tops between drainages atop the plateau often lacked forest cover and were thus vegetated by perennial grasslands of bluebunch wheatgrass and Idaho fescue associations. Forests on more mesic slopes and higher elevations contained both Ponderosa pine and Douglas-fir (Pseudotsuga menziesii [Mirb.] Franco).

Soils in Study Area 1 had not yet been mapped. However, soil surveys in nearby areas suggested the bunchgrasses areas in the lower elevations were underlain by loamy-skeletal, mixed, superactive, mesic vitrandic argixerolls, and pachic palexerolls. Open forests atop the plateau were likely supported by loamy-skeletal, isotic, frigid, vitrandic argixerolls. Soils in nonforested areas were likely loamy-skeletal, mixed, superactive, mesic lithic argixerolls.

The Snake River Remote Automated Weather Station (RAWS ID = SRFI1) located west of Cuprum, Idaho at 1 333-m elevation was the most relevant climate station near Study Area 1. Long-term (1998–2016) mean water-yr precipitation at this station was 546 mm (MesoWest 2018). Total precipitation values for the 2009, 2010, and 2011 water yrs were 441, 484, and 537 mm, respectively. Long-term (1998–2016) mean daily air temperatures for the months of June, July, and August were 16.5°C, 23.0°C, and 22.3°C, respectively.

Study Area 2 (48 km²) was selected to represent management situations where calf age at entry was greater (3–3.5 mo) and base elevation (1 011 m) of the grazing area was higher than in Study Area 1. The dominant landform here was also a dissected plateau, but in contrast to Study Area 1, most of the study area extent was situated atop the plateau where the maximum elevation was 1 865 m. As such, although a wide range of slopes (0–62 degrees) were present, this study area was generally flatter ( = 13 ± 8.7 SD degrees) than Study Area 1. In addition, a larger extent of this study area was forested than Study Area 1 and this was likely a consequence of the higher base elevation and associated moisture regime. Cattle entered Study Area 2 at midelevation and dispersed, but, given the limited topography relief in this study area, seasonal movement upslope was generally much less pronounced than at Study Area 1.

The lower elevations of Study Area 2 (up to 1 200 m) occurred along steep canyon walls vegetated by bluebunch wheatgrass and Idaho fescue associations or sagebrush-grasslands dominated by mountain big sagebrush (Artemisia tridentate Nutt. subsp. vaseyana [Rydb.] Beetle) associations (Johnson and Simon 1987). Mid-elevations were rolling mountain slopes vegetated by ponderosa pine savanna on drier exposures and Douglas-fir forest on more mesic slopes. Some ridgetops were nonforested being instead vegetated by bunchgrass or sagebrush associations. The highest elevations (above 1 600 m) were mountain slopes and ridges vegetated by mixed conifer forests of Douglas-fir and grand fir associations. Dry meadows on stream terraces and other flat areas were vegetated by Kentucky bluegrass (Poa pratensis [L.]) and California oatgrass (Danthonia californica Bol.). Tufted hairgrass (Deschampsia cespitosa [L.] P. Beauv.), Hood's sedge (Carex hoodia Boott), and thick-head sedge (Carex pachystachya Cham. Ex Steud.) dominated moist upland meadows. Aspen (Populus tremuloides Michx.) occasionally occurred near springs and other moist areas. Stream riparian areas at lower elevations were dominated by black cottonwood (Populus balsamifera [L.] spp. Trichocarpa [Torr. & A. Gray ex Hook.] Brayshaw), willow (Salix spp. [L.]), and Kentucky bluegrass while willow and sedges (Carex spp. [L.]) dominated at higher elevations.

Soils in the lower elevations of Study Area 2 were loamy-skeletal, mixed, superactive, mesic vitrandic argixerolls and pachic palexerolls (NRCS 2017d). At midelevations, loamy-skeletal, mixed, superactive, mesic lithic argixerolls occurred in the nonforested areas while open forests were underlain by loamy-skeletal, isotic, frigid, vitrandic argixerolls. Soils at higher elevations were not yet mapped. Climatic data from the Snake River RAWS (see earlier) were assumed to also be representative of the climate in Study Area 2.

Study Area 3 (73 km²) was selected to typify situations similar to Study Area 2 with older calves (3–3.5 mo at entry) and higher base elevations than Study Area 1, but in this case, the dominant landform was a mountain and the extents encompassed a large elevational gradient from toeslopes (1 082 m) to the mountain summit (2 478 m). Cattle entered the study area at the very lowest elevations and, as in Study Area 1, progressively followed forage plant phenology upslope reaching the mountain summit area in late July–early August and remaining there until the end of the study period.

Stream riparian areas at the lowest elevations of Study Area 3 dominated by black cottonwood, willow, Kentucky bluegrass, and sedges. Upland slopes at these low elevations (< 1 500 m) were vegetated by bluebunch wheatgrass associations on the drier exposures, while more mesic slopes contained open to moderately closed woodlands dominated by ponderosa pine associations. Douglas-fir and grand fir associations occupied the mesic exposures at midelevations, and pine-bunchgrass savanna vegetated the drier slopes. Open ridgetops and divides between drainages at midelevations were vegetated by bunchgrass grasslands. Vegetation at the highest elevations (> 2 000 m) was dominated by subalpine fir (Abies lasiocarpa [Hook.] Nutt.) associations in forested drainages and fescues (Festuca spp.[L.]) and upland sedges in open areas.

Soils at the lower elevations of Study Area 3 were fine-loamy, mixed, superactive, mesic ultic argixerolls and frigid pachic ultic argixerolls (NRCS 2018d). Forest soils at midelevations were loamy-skeletal, mixed, superactive pachic argicryolls. Soils on open, midelevation slopes were clayey-skeletal, smectitic lithic argicryolls. Soils at higher elevations had not yet been mapped.

The NRCS Van Wyck SNOTEL station (ID = 979) located at 1 500 m elevation south of Indian Valley, Idaho was the most relevant climate station near Study Area 3. Long-term (2002–2017) mean water-year precipitation at this station was 663 mm (NRCS 2018c). Total precipitation values for the 2009, 2010, and 2011 water yrs were 650, 742, and 782 mm, respectively. Long-term (2002–2017) mean daily air temperatures for June, July, and August were 15.6°C, 22.4°C, and 21.3°C, respectively.

Study Area 4 (83 km²) had the highest base elevation (1 248 m) of the 4 study areas, and calf age at entry was about 3.5 mo. The dominant landform was a mountain with two parallel ridgelines, a V-shaped stream valley between, and maximum elevation of 2 582 m. Most of the study area extent was forested. Cattle entered at the lowest elevations and moved upslope as the season progressed similar to Study Areas 1 and 3.

The lowest elevations (< 1 600 m) of this study area were vegetated by open ponderosa pine woodlands surrounding occasional open hillslopes dominated by bunchgrasses. Pine associations gave way to Douglas-fir and grand fir associations at midelevations. Shelter-wood silviculture treatments had also created sparsely wooded patches of 10–60 ha distributed occasionally across the midelevations. Stream riparian areas at midelevations were generally forested with grand fir or Engelmann spruce (Picea engelmannii Parry ex Engel.) overstory and contained a shrub layer of willows and herb layer of rushes and sedges. Wet meadows of willow and sedges also occurred at midelevations, but these were very small (< 5 ha). At the highest elevations (> 1 950 m), subalpine fir (Abies lasiocarpa [Hook.] Nutt.) associations dominated the forested areas and fescues and upland sedges occurred in nonforested areas scattered among extensive granite outcrops.

Pine woodlands at the lowest elevations occurred on fine-loamy and loamy-skeletal, mixed superactive pachic argicryolls (NRCS 2018d). Soils on open hillslopes were clayey-skeletal, smectitic lithic argicryolls. Loamy-skeletal, mixed superactive pachic argicryolls and coarse-loamy, mixed, superactive lamellic haplocryepts underlied the Douglas-fir and grand fir associations at midelevations. Soils at higher elevations had not yet been mapped but were coarsely textured and likely derived from granitic parent materials.

The Brundage Reservoir SNOTEL station (ID = 370) located at 1 905 m elevation south of McCall, Idaho was the closest climate station to Study Area 4. Long-term mean annual precipitation (1987–2015) at the site was 1 271 mm (NRCS 2018a). Precipitation totals for the 2009, 2010, 2011 water yrs were 1 270, 1 209, and 1 516 mm, respectively. Long-term (1987–2015) mean daily air temperatures for June, July, and August were 10.5°C, 16.0°C, and 14.9°C, respectively.

Cattle Data Collection

Beef cattle herds in the study areas ranged in size from 195 to 460 cow-calf pairs depending on study area size and stocking rates of about 12–14 ha AUM^–1. These cattle were primarily of British breeds or crosses. Each spring (2009–2011) before entry, 10 mature cows (4–10 yr of age) were randomly selected from the herd associated with each study area. Generally, these mature cows had multiple yrs of experience with the study area landscapes, climates, and herd management actions. Each of the selected cows was fitted with a custom GPS tracking collar (Clark et al. 2006), which recorded the date, time, spatial position, speed, and positional accuracy parameters every 5 min throughout the grazing season (June–October). We assumed, on the basis of our random selection process, that the variability in spatial behaviors expressed by the GPS-collared cows were representative of the respective herds. The GPS collars were retrieved, and data were downloaded as cattle exited the grazing allotments in the fall. Resultant samples sizes, however, were unequal among study areas and yrs due to collar malfunctions and other contingencies. Consequently, to provide a more equal sampling, 3 collared cows from the 10 cows potentially available for each study area–yr combination were selected, on the basis of completeness of GPS data record, for analysis under this study. As the focus of this study was about predicting spatial risk of wolf-cattle encounters during the 15 June–15 August period when cattle were occupying summer ranges and reproductive wolf packs were using rendezvous sites and associated areas, these 36 cattle GPS data sets were truncated to this 2-mo focal period (hereafter referred to as the “rendezvous period”).

Data from these 36 cattle collars were then processed to remove gross GPS positioning error by using a Geographic Information System (GIS) to exclude positions located outside the bounding fences and other limiting perimeter features (e.g., impassable rivers) of each study area. Next, an objective process was applied to flag potentially erroneous positions for further scrutiny. Positions having a Positional Dilution of Precision (PDOP) parameter value of > 10 and/or an instantaneous speed value > 9.3 km/h (5 knots) were flagged as potentially inaccurate. Positions indicating the cow had traveled > 500 m during the 5-min sampling interval (i.e., a sustained velocity of at least 6 km/h) were also flagged. Our previous experience conducting other studies within the region, where continuous visual observation was used to evaluate cattle behavior, indicated that it was rare for range cattle to sustain a velocity of ≥ 6 km/h for 5-min duration (Clark et al. 2017a, 2017c). A custom software package, Kinetic Resource and Environmental Spatial System (KRESS) v. 4, was then used to visualize and individually evaluate the flagged positions in context. Each flagged position was displayed on a high-resolution digital elevation model (DEM) and/ or digital orthophotograph quadrangle (DOQ) background along with ≥ 10 of the preceding and following positions. A line connecting all these positions was overlain to illustrate the movement trajectory. Flagged positions were subjectively evaluated to determine if they substantially departed from the general magnitude and direction of movement along the GPS-based trajectory. Positions that did not deviate substantially from the trajectory were accepted and retained in the data set. Positions exhibiting substantial deviation were further evaluated using the DEM and DOQ to assess whether terrain or other landscape features could logically explain the departure from the trajectory. Positions failing this secondary screening were rejected and removed from the data set. The complete series of error-screening processes excluded an average of 3.4% of the originally recorded positions from each collar data set.

As a final processing step, the custom software program, ASSOC1 (Weber et al. 2001), was used to determine if the three cows selected under each study area–yr combination were behaving spatiotemporally independent of each other (i.e., were not associated). Associated behavior between or among collared cows would violate the independence assumptions of resource selection analyses conducted with these data (Hilbe 2008; Nielson and Sawyer 2013). In this study, dyads that spent > 75% of their time separated by > 75 m from each other were considered nonassociated. These thresholds were based on those used in previous range cattle behavior studies on mountainous grazing lands in the region (Clark et al. 2014, 2016). In all cases, selected cattle in this study were determined to be nonassociated. If an associated dyad had been detected, one member would have been randomly selected for replacement by one of the remaining collared cows in the herd and the association test rerun.

A random sample of 25% of the cattle GPS positions from each study area–yr combinations was extracted and reserved for model validation. The remaining 75% of these data was retained for use as model development data sets.

Wolf Data Collection

General wolf presence and wolf-cattle encounter frequency within the four study areas were monitored during 2009–2011 using a combination of telemetry tracking (radio and GPS), scat surveys, camera traps, den/rendezvous site surveys, direct observation, and depredation investigations. Clark et al. (2017b) provided specific descriptions of the wolf-monitoring methods used in the present study. Data from all monitoring sources were used to classify wolf presence and encounter frequency to three levels: low, moderate, and high. Presence/encounter classifications were summarized for each of the 3 summer mo (June, July, and August). Months when no wolf presence was detected, despite rigorous sampling effort, were classified to the low-presence/encounter class. The moderate class was assigned to months when wolf presence was detected but no evidence of cattle pursuit events or depredations were recorded. The high presence/encounter class was reserved for months when pursuit events and/or depredations were documented within the study area.

Table 1

Predictor variables used to develop the a priori set of candidate models evaluated at Study Area 3 for predicting resource selection patterns of mature beef cows on mountainous grazing lands.

Wolf presence and wolf-cattle encounter frequency varied with time and among study areas with Study Area 1 tending to have higher levels during summer mos (June, July, and August) than Study Areas 2 and 4, but all three study areas experienced at least moderate wolf presence during ≥ 1 mo (e.g., July) during all 3 study yrs (Clark et al. 2017b). Although wolf den and rendezvous sites were located about 9 km north of Study Area 3, wolf presence here was generally lower than at the three other study areas.

Monitoring of wolf presence, encounters, and movement patterns was intensified at Study Area 1 during the 2009 study yr. An adult male wolf from the Snake River pack (n = 11 individuals, 5 adults, and 6 pups of the yr) was captured by USDA APHIS Wildlife Services on 22 May 2009 and instrumented with a custom GPS collar (Clark et al. 2006) as a wolf-management response to depredations occurring near Study Area 1. The GPS collar recorded the position of this wolf every 15 min throughout the summer, fall, and winter of 2009. These GPS data were processed to identify and remove gross positioning error in nearly the same fashion as those for the collared cattle. The exception was that these wolf data were not clipped to boundary features as was done for the cattle data. For the purposes of this study, the wolf GPS data set was then truncated to the 15 June–15 August rendezvous period.

This GPS-collared wolf was directly observed, from aircraft and ground, on at least six recorded occasions over the course of the 2009 rendezvous period. During each occasion, the wolf was observed in close company with one to five members of the Snake River pack. Comparison of GPS positions for this wolf with biweekly telemetry flight positions (GPS-based) recorded for the VHF radio-collared, alpha female of the Snake River pack indicated these two wolves were < 250 m of each other during each observation. Consequently, there is strong evidence the movements of the GPS collared adult male wolf were effectively representative of the movements of the Snake River pack in general.

Cattle Resource-Selection Analyses

Resource selection by GPS-collared cattle in the four study areas was modeled using a negative-binomial (NB) regression approach described in detail by Nielson and Sawyer (2013) and previously applied on other mountainous grazing lands in the region (Clark et al. 2014, 2016). These regression models were RSFs as defined by Manley et al. (2002) but differed from logistic regression–based RSFs, which are more typically applied. Standard methods using logistic regression to estimate the exponential RSF (Manly et al. 2002) with samples for used and available locations do not assess the intensity of use but rather treat a habitat unit as used, regardless of whether that area was visited once or multiple times. The NB regression used here goes further and actually does assess intensity of use. Negative-binomial regression only slightly increases the complexity of the sampling process, and model fitting routines are widely available in current statistical software. Model coefficients were estimated using the following equations (1 and 2) (Nielson and Sawyer 2013):

Table 2

Top 12 cattle resource selection function (RSF) models, based on negative binomial regression, for each of 3 study yr at Study Area 3, selected from an a priori set of 184 candidate models based on Akaike's information criterion (AIC) fit scores. Model in bold font was selected as the final cattle (RSF) model, and the model in italics was the runner-up.

The general modeling approach used here involved 5 steps. First, a GIS was used to create virtual, circular plots which were randomly distributed throughout the study area and attributed with values from predictor variable data layers (e.g., terrain slope, vegetation greenness, distance from perennial streams). Second, for population models, GPS data for all collared animals were pooled and then counts were made of any GPS positions that were located within each virtual plot. Third, a candidate set of negative-binomial regression models was fitted to determine which combination of predictor variables best predicted the relative probability of cattle use within the study area. Fourth, a bootstrapping procedure was applied using the individual cattle GPS data sets to derive standard errors and confidence intervals for population model coefficients. Finally, the best negative-binomial regression model was used to create a raster map of predicted probability-of-use classes for the study area.

Specifically, initial sets of circular plots were digitally generated at random locations within each study area. Counts of GPS positions located within these plots were then tallied for each study area–yr combination using a custom script written in the R programming language. Next, zero inflation (i.e., where too many plots contain no GPS positions) was considered because this sampling problem can cause issues in modeling resource selection. In some cases, one could consider using zero-inflated or hurdle models (Nielson et al. 2013), but such necessities are rare. Rather zero inflation cases are often correctable because they usually stem from wrongly identifying which parts of the landscape are truly available to GPS-tracked animals. In the present study, cattle were contained in extensive but fenced pastures, so the availability space was well defined by the pasture boundaries. Nonetheless, checks for zero inflation were advisable and thus conducted by comparing the number of sampled zeros (i.e., plots containing no GPS positions) to the number expected in a standard NB distribution given the mean and variance of the sampling data, along with the sample size. Several choices of plot size and number were iteratively evaluated to find a suitable balance that provided position counts, which well approximated a negative-binomial distribution while also allowing effective detection of variability in animal movement and intensity of use among plots (Nielson and Sawyer 2013). Choices where plots were too small in size and/or too few in number tended to under-sample the GPS positions present while choices with plots too large and/or too many failed to detect variability in animal use. A good compromise was found using 1 500 plots of 500 m in diameter for Study Areas 2, 3, and 4 during all study yrs. Given its larger extent, 2 500 plots of 500-m diameter were used for Study Area 1 during all yrs.

Table 3

Fitted coefficients and statistics for the final cattle resource selection function (RSF) model applied at each of the four study areas for each of 3 study yr. P values significant at the 0.05 alpha level are highlighted in bold font.

A GIS was then used to attribute each plot of a study area with values for nine predictor variables (Table 1). Topographic predictors were derived from 10-m digital elevation models (DEM) sourced from the US Geological Survey. Plots were attributed for the elevation predictor using the mean of all 10-m DEM elevation cells, which were intersected by a plot. Aspect was initially derived as the azimuth bearings for all 10-m cells intersected by a plot. This continuous variable was then converted to a categorical variable with five classes to thus include aspects generalized to four cardinal directions (± 45 azimuth bearing degrees) and flat areas (< 10-degree slope), which served as the reference class. Normalized surface or topographic roughness index values were derived from the 10-m DEM using a GIS to calculate the mean values for 3 × 3 cell neighborhoods (ESRI 2018c). These neighborhood means or roughness values from the resulting grid were then normalized using the following equation:

Table 4

Spearman rank correlation scores (r_s) quantifying the prediction success of the final cattle resource selection function (RSF) model for each study area–yr combination and for application of the model at all four study areas using Global Positioning System data pooled across all study yrs.

A set of 184 candidate cattle RSF models, composed of up to 5 of these predictor variables, was identified a priori (Burnham and Anderson 2002). Model ranking and final selection were then conducted at Study Area 3, where wolf presence and associated predation threat were consistently lower and thus least likely to have had a confounding influence on cattle resource selection patterns compared with other study areas. The intent here was to minimize complex and hidden interactions with confounding factors and select models whose predictive performance was more directly and exclusively dependent on our short list of simple, easily measured environment and habitat predictors (see Table 1). Conducting model selection at Study Area 3, consequently, was expected to yield cattle RSF models, which were robust and generally applicable to other regional sites regardless of wolf presence levels at those sites. All models in the candidate set were fitted to GPS position data from the model development data sets for Study Area 3. Model fitting was conducted separately for each study yr. Fitted models were then ranked on the basis of the Akaike's information criterion (AIC) score (Akaike 1973; Burnham and Anderson 2002). A short list of the 12 best-fitting models was developed for each yr. From these short lists, the model giving the best, most robust AIC fits, when considered across all 3 yr, was selected as the final cattle RSF model.

The final model was then applied to the remaining three study areas to further evaluate its robustness and general applicability. The model was fit using the GPS position data sets (model development) from each study area–yr combination. Prediction performance of the fitted model from each combination was then evaluated with a Spearman rank correlation analysis using the GPS position data sets previously reserved for model validation (Boyce et al. 2002; Sawyer et al. 2009). This analysis examines the relationship between the ranking of predicted probability of cattle use classes (n = 20 ranked classes) and the counts of GPS positions occurring within plots classified to each of these classes. Cases where counts progressively increased with increasing class rank indicated strong correlation and successful model predictive performance. Spearman rank scores were compiled for each study area–yr combination.

The final model was also fitted and evaluated at each study area using GPS data pooled across all 3 study yr to thus provide a set of four population RSF models for use in additional analyses. Standard errors and confidence intervals (90%) were calculated for the coefficients of each these four population RSF models using a bootstrapping routine drawing 1 000 replicate samples from the model development data sets.

Cattle resource selection patterns were mapped for each study area by applying the fitted population model coefficients within a GIS. For each study area extent, spatial data layers for each of the predictors in the final cattle RSF model were loaded in the GIS. The model coefficients were then applied to their respective data layers, using the modeler function within the GIS, to output a resultant raster grid where each 10-m cell was assigned a proportional value representing the predicted probability of cattle use (ESRI 2018d). This raster was then classified into 10 ranked classes using an equal-area classification scheme from which a color-coded cattle RSF map was produced for each study area (ESRI 2018e).

Wolf Resource Selection Analyses

Ausband et al. (2010) applied a simple, three-variable logistic regression model to predict the relative suitability of different habitat types for wolf rendezvous sites at four study areas in mountainous central Idaho. This wolf RSF model included NDVI, normalized surface roughness, and profile curvature as predictors. They used 178 documented rendezvous site locations to validate the predictive performance of their RSF model. In the present study, habitat types within a rectangular region (762 km²) encompassing all four of our western Idaho study areas were classified, using the Ausband et al. (2010) model and their coefficient values, according to relative suitability as wolf rendezvous sites. The same GIS procedures described earlier for cattle RSF mapping were used to create 10-class wolf rendezvous site habitat suitability map with 10-m grid cells for this rectangular region.

Predictive performance of the wolf RSF map was evaluated using a small data set of documented rendezvous sites (n = 14 sites) inside or nearby our four study areas. At each study area, counts of rendezvous sites in the different mapped habitat suitability classes were made. For rendezvous sites not occurring in the very highest suitability class (class 10), the distance between the site and the nearest area of that class was determined. Spearman rank correlation analysis was used to examine the relationship between counts of sites in suitability classes and the ranking of those classes.

Table 5

Fitted coefficients and statistics for the final cattle resource selection function (RSF) model applied for each of the four study areas using Global Positioning System (GPS) data pooled across all study yrs. While the coefficient estimates are based on pooled GPS data, the 90% confidence level (CL) were calculated from bootstrapped samples from individual cattle data sets. P values significant at the 0.05 alpha level are highlighted in bold font.

An additional case-study analysis was conducted at Study Area 1 to evaluate whether habitat types classified as having relatively high suitability for wolf rendezvous sites also tended to have higher wolf presence during the rendezvous period than habitat types of lower suitability. The GPS tracking data set acquired during 2009 for an adult male wolf of the Snake River pack was used in this analysis. GPS positions were acquired at 15-min intervals, yielding 6 591 positions within the study area boundaries during the rendezvous period. Counts of positions located in the different mapped habitat suitability classes were made. Another Spearman rank correlation analysis was conducted for counts and class ranking.

Wolf-Cattle Encounter Risk Mapping

Spatial relationships between cattle resource selection patterns and relative habitat suitability for wolf rendezvous sites were evaluated as a means of predicting the risk of wolf-cattle encounters within the four study areas. For each study area, the raster grid map of classified relative probability of cattle use was overlain on the wolf rendezvous site habitat suitability map in a GIS and coregistration was confirmed. Because both the cattle and wolf RSF maps used the same ranked 10-class value scheme, there were 100 unique combinations of class values possible between coincident raster cells from each map. Counts of occurrences for each of these 100 combinations were made in the GIS (ESRI 2018a). Similarity in the spatial distributions of class values in the cattle and wolf RSF maps was evaluated using a weighted Spearman rank correlation analysis in which the occurrence counts were used as weights. Findings of similarity between the two maps would allow spatial prediction of wolf-cattle encounter risk.

Figure 1.

Map of predicted cattle relative probability of use (10 classes) during the wolf rendezvous period (15 June–15 August), derived using the final cattle resource selection function (RSF) model and cattle Global Positioning System (GPS) data pooled across 3 study yr (2009–2011), for Study Area 1 (112 km²) in mountainous grazing lands of western Idaho. Positions (n = 86 470) from nine GPS-collared mature beef cows (3 cows per study yr) are also displayed relative to roads, streams, and hill-shaded terrain shape.

Wolf-cattle encounter risk maps were created for each study area using the coregistered cattle and wolf RSF maps. Ten risk classes were derived by applying the following equation within a GIS:

Although these 10-class encounter risk maps were logical products of combining the 10-class cattle and wolf RSF maps and were useful for some analyses requiring detailed spatial risk classification, practical experience quickly demonstrated that simpler, more readily interpretable 5-class risk maps would be more effective for management applications. Five-class maps were created by simple combination of class pairs from the 10-class maps.

Risk Map Validation

Rigorously validating conflict-mitigation tools like these wolf-cattle encounter risk maps would typically be challenging as direct observation data for wolf-cattle encounters are rare and difficult to acquire. Fortunately, the GPS tracking data sets collected for wolves and cattle at Study Area 1 provided the opportunity to conduct a case-study validation of the 5-class risk map. Concurrent GPS tracking data acquired in 2009 for 10 collared beef cows (5-min GPS samples) and one collared adult male wolf (15-min samples) representing the movement of the Snake River pack (n = 11 wolves) were used to identify GPS-based, wolf-cattle encounters during the rendezvous period. Concurrent wolf and cattle GPS positions (± 450 s) located within 500 m of each other were considered wolf-cattle encounters. The assumption behind the 500-m threshold distance was that either wolves or cattle or both, through sensory perception, would become aware of the other species presence at this distance, even if direct line of sight was occluded. Researchers investigating wolf-elk relations suggested 1 km was a reasonable threshold for wolf-elk encounters (Middleton et al. 2013; Cusack et al. 2019). Cattle are typically expected to be less vigilant than elk, so the 500-m threshold used here was considered reasonable and appropriate. About 716 occurrences met these encounter criteria; however, not all could be considered independent encounter events. In some cases, multiple collared cows encountered the wolf during the same time frame. If these cows were separated by > 75 m from each other, then each wolf-cow encounter was considered an independent event. The assumption being cattle dyads separated by > 75 m were not behaviorally associated. Otherwise, only a single-encounter event was recorded and this event was randomly assigned to just one of the cows involved. In these and other cases, an encounter might also involve multiple consecutive 15-min wolf GPS positions, thus potentially lasting up to several hours in duration. Only the initial position of an extended encounter was considered an independent event. Encounter events ended with an interruption of at least 15 min in duration. On the basis of these criteria, 200 of the 716 occurrences were identified as independent wolf-cattle encounters. These 200 GPS-based encounter positions were then overlain on the wolf-cattle encounter risk map for Study Area 1, and counts were tallied for encounters occurring in each of the five risk classes. A Spearman rank correlation analysis was used to assess the relationship between counts and class ranking. Distances between encounter locations and the nearest area of the very high risk class (class 5) were compared with those for random points (n = 200). As the distributions for both distance samples were highly left-skewed, a Mann-Whitney U-test was used for the comparison.

Figure 2.

Map of predicted cattle relative probability of use (10 classes) during the wolf rendezvous period (15 June–15 August), derived using the final cattle resource selection function (RSF) model and cattle Global Positioning System (GPS) data pooled across 3 study yr (2009–2011), for Study Area 2 (48 km²) in mountainous grazing lands of western Idaho. Positions (n = 151 280) from nine GPS-collared mature beef cows (3 cows per study yr) are also displayed relative to roads, streams, and hill-shaded terrain shape.

Figure 3.

Map of predicted cattle relative probability of use (10 classes) during the wolf rendezvous period (15 June–15 August), derived using the final cattle resource selection function (RSF) model and cattle Global Positioning System (GPS) data pooled across 3 study yr (2009–2011), for Study Area 3 (73 km²) in mountainous grazing lands of western Idaho. Positions (n = 129 925) from nine GPS-collared mature beef cows (3 cows per study yr) are also displayed relative to roads, streams, and hill-shaded terrain shape.

Depredation Risk

A small data set of cattle depredations (n =16), confirmed as wolf caused by USDA APHIS Wildlife Services, was used to evaluate the efficacy the wolf-cattle encounter risk maps for predicting the relative spatial risk of wolf depredation. At three of the four study areas, counts were made of depredations occurring within each of the five mapped risk classes. No depredation data were available for Study Area 4. The very small sample sizes of depredations at individual study areas did not allow specific analyses. Instead, depredation count data for all study areas were pooled and a Spearman rank correlation analysis was conducted on these pooled data.

At the three relevant study areas, distances between depredations locations and the nearest area of the very high risk class were calculated for comparison with distances for random points (n = 200). Distance data were pooled across study areas, and the depredation and random sample means were compared with a Mann-Whitney U-test.

Cattle Resource Selection

Table 2 lists the top 12 cattle resource selection models for each of the 3 yr at Study Area 3. These top models were selected, on the basis of AIC fit score rankings, from the set of 184 a priori candidate models. Just two models, one a simpler derivative of the other, were found to consistently occur among the top models across all 3 study yr. Both of these negative-binomial regression models contained five predictors—slope, distance from roads, distance from perennial streams, ponderosa pine cover, and aspect—but differed regarding inclusion of the quadratic forms for the slope and streams predictors. The more complex model, presented in bold font in Table 2, was selected as the final cattle RSF model for this study as it consistently had better AIC fits than the simpler, runner-up model (presented in italics in Table 2).

The final cattle RSF model was then evaluated for robustness by fitting to GPS data from all four study areas for each of the 3 study yr. Model fit results for each of these 12 study area–yr combinations are presented in Table 3. While the relative importance of the individual predictors within the model varied considerably among these situational cases, the model remained remarkably effective for accurately predicting the relative probability of cattle use (hereafter referred to as “predicted use”) within this broader, more diverse scope. Spearman rank correlation analysis, using the GPS data sets previously reserved for model validation, yielded high prediction success (r_s > 0.93) for 11 of 12 cases (Table 4). Furthermore, in these 11 cases, the model performed on par or even somewhat better in Study Areas 1, 2, and 4 than in Study Area 3, where the model was developed. While even in the exceptional case, Study Area 4 during 2011, the Spearman score was still quite high (r_s = 0.81), the evident reduction in prediction success relative to the other cases suggested additional influential factors may have been at play in this specific situation. Questions about this case are explored later. Meanwhile, it is clear from the Spearman rank results that the predictive performance of this cattle RSF model was quite robust across a broad spatiotemporal scope.

As a side note, for a separate publication, the authors fit all 184 models from the original candidate set to each of the 12 study area by yr combinations. In addition, they also fit a much larger candidate set (1 139 models) involving both the original set plus models, which included other predictors based on terrain indices (e.g., topographic roughness [TRI], position [TPI], and wetness [TWI]) and slope-position classifications. Some models certainly performed better for some study area–yr combinations than others. However, this more extensive testing again confirmed that the original final model was the best in terms of overall robustness in performance across all combinations.

Figure 4.

Map of predicted cattle relative probability of use (10 classes) during the wolf rendezvous period (15 June–15 August), derived using the final cattle resource selection function (RSF) model and cattle Global Positioning System (GPS) data pooled across 3 study yr (2009–2011), for Study Area 4 (83 km²) in mountainous grazing lands of western Idaho. Positions (n = 127 961) from nine GPS-collared mature beef cows (3 cows per study yr) are also displayed relative to roads, streams, and hill-shaded terrain shape.

Figure 5.

Predicted cattle relative probability of use response to the terrain slope (degrees) predictor, in the final cattle resource selection function model, for each of four study areas in mountainous grazing lands of western Idaho.

With a robustly effective cattle RSF model then in hand, a simpler composite set of results was then needed for further analysis and presentation. Therefore, the final model was fitted for each study area using GPS data pooled from across all 3 study yr. Spearman rank analysis was again applied and yielded scores (r_s) of 0.99, 0.95, 0.98, and 0.92 for Study Areas 1, 2, 3, and 4, respectively, thus confirming the final model performed effectively as a population model (see Table 4). Predicted cattle resource selection patterns for each study area are presented in Table 5 and Figures 1-4.

Because the relationship between predicted cattle use patterns and suitable wolf rendezvous site habitat types was a principal interest of this study, a performance comparison between the final cattle RSF model and the simple, three-variable model described by Ausband et al. (2010) was conducted. The latter model was fitted for each study area using the pooled cattle GPS data and Spearman rank analyses were conducted for model performance. While this simpler model, with an entirely different set of predictors, did not perform as well as the final model, it did in some cases successfully predict cattle use patterns (r_s range 0. 46–0.97, see online  supplementary Table S1, available online at  https://doi.org/10.1016/j.rama.2019.08.012).

Specific Responses to Predictors

Predicted cattle responses to the individual predictors within the final model were then evaluated across the four study areas based on fits using GPS data pooled across yrs. Generally, predicted cattle use decreased in a curvilinear fashion with increasing slope (Fig. 5). Cattle at Study Area 2 were predicted to be the most sensitive to slope, exhibiting a sharp decline as slope increased from 2.2 degrees to about 10 degrees and then use levels flattened out with minimum predicted use occurring at about 22 degrees. Cattle at the other three study areas were predicted to exhibit a more sigmoidal decline in use with increasing slope. Peak use occurred on slopes of < 10 degrees, and slopes of > 25 degrees were largely avoided. Only the quadratic form of slope was significant at Study Areas 1 and 3. It is notable that predicted use at Study Area 4 initially increased with slope before peaking at about 7 degrees and then declining.

Figure 6.

Predicted cattle relative probability of use response to the distance to nearest road (m) predictor, in the final cattle resource selection function (RSF) model, for each of four study areas in mountainous grazing lands of western Idaho.

In most cases, predicted cattle use increased with distance from roads (Fig. 6). Predicted use at Study Area 1 exhibited a slight initial decline with distance to a minimum at 445 m and then remained low out to about 1 500 m, where use then inclined sharply to a maximum at 1 866 m. A gentler curvilinear increase in predicted use with distance occurred at Study Area 2, where the rate of increasing use was initially moderate and steepened at about 750 m. Maximum use occurred at 1 274 m. Although predicted cattle use at Study Area 3 did initially increase with distance from roads, maximum use occurred at about 1 665 m and then use declined with minimum use occurring at about 4 011 m. Predicted use, however, was relatively high (> 0.6) along the entire distance range for Study Area 3. Distance to roads was not a significant predictor of cattle resource selection at Study Area 4.

Predicted cattle use generally increased with distance from perennial streams, but there was some variability among study areas in the details of this response (Fig. 7). Like the response to distance from roads, predicted cattle use at Study Area 1 initially declined slightly, reaching a minimum at 628 m, remained quite low out to about 1 500 m from streams, and then increased abruptly reaching a maximum at 2 004 m. These responses to roads and streams seem to run counter to the mapped predicted cattle resource selection patterns for Study Area 1 (see Fig. 1), which suggest cattle prefer areas near roads and streams. This apparent discrepancy is further explored in the Discussion section later. Predicted cattle use exhibited nearly linear increases, at least initially, with distance from streams at both Study Areas 2 and 3. However, while maximum predicted use was reached at 1 212 m for Study Area 3, predicted use at Study Area 2 peaked at 1 602 but then declined with further distance from streams. Study Area 4 was again a somewhat exceptional case, where predicted use started relatively high (> 0.8) near streams and initially increased with distance but peaked at 334 m and then declined with increasing distance to a minimum level at 1 562 m.

Figure 7.

Predicted cattle relative probability of use response to the distance to nearest stream (m) predictor, in the final cattle resource selection function model, for each of four study areas in mountainous grazing lands of western Idaho.

Predicted cattle response to terrain aspect varied among study areas. Cattle at Study Area 1 were 0.587 times (0.368, 1.08 CL 90%) as likely to use areas with northerly aspects compared with flat areas (< 10 degrees slope). Predicted use of other aspects did not differ from that of flat areas. At Study Area 2, cattle were 1.54 times (1.28, 1.93 CL 90%) as likely to use northerly aspects, 0.596 times (0.395, 0.789 CL 90%) to use easterly aspects, and 0.674 times (0.527, 0.838 CL 90%) to use westerly aspects as they were flat areas. No significant effect of aspect on predicted cattle use was detected at Study Area 3. Cattle at Study Area 4 were less likely to use of easterly (0.644 odds with 0.470, 0.846 CL 90%), southerly (0.425 with 0.237, 0.647 CL 90%), and westerly aspects (0.515 with 0.320, 0.758 CL 90%) than flat areas.

Predicted cattle use at Study Area 1 increased by 0.733% (0.379%, 1.22% CL 90%) with each percentage point increase in areal coverage by the ponderosa pine vegetation type. At Study Area 4, predicted cattle use increased by 1.79% (1.02%, 2.45% CL 90%) with each percentage point increase in ponderosa pine. No significant effects of ponderosa pine cover on predicted cattle use were detected at Study Areas 2 and 3.

Figure 8.

Map of predicted relative habitat suitability for wolf rendezvous sites (10 classes), derived using the Ausband et al. (2010) wolf resource selection function (RSF) model and their coefficients, for Study Area 1 (112 km²) in mountainous grazing lands of western Idaho. Locations of documented wolf rendezvous sites (n = 6) are also displayed relative to roads, streams, and hill-shaded terrain shape.

Cattle resource selection patterns at the four study areas can be summarized with just a few lines. Cattle were predicted to strongly select for areas of relatively low slope (< 10 degrees). These areas could be near streams and/or roads but most often occurred as meadows, benches, saddles, or flat ridgetops quite distant from streams and roads. In some cases, areas dominated by ponderosa pine were selected for, particularly if those areas were also relatively flat. Cattle quite strongly avoided steep slopes (< 20 degrees), particularly those that were sparsely vegetated or had relatively low forage productively.

Wolf Resource Selection

Figure 8 is a map of predicted relative habitat suitability for wolf rendezvous sites (10 classes), derived using the Ausband et al. (2010) wolf RSF model and their coefficients, for Study Area 1. Predicted habitat suitability maps for the other three study areas are provided online in the supplementary materials ( Figs. S1-S3, available online at  https://doi.org/10.1016/j.rama.2019.08.012). Although Ausband et al. (2010) had effectively validated this wolf resource selection model for their study areas in central Idaho, further validation was required for the present application in western Idaho. Unfortunately, documented rendezvous site location data were not available within most of the present study areas. Data were available, however, for Study Area 1, where six documented rendezvous sites were occupied by the Snake River wolf pack in 2009 (Fig. 8). These rendezvous sites were well distributed across the extent of the study area and thus served as a small but effective validation set. When overlain on the habitat suitability map for Study Area 1, all six sites occurred in areas classified to either class 10, very highest suitability (5 sites) or class 9, very high rendezvous site suitability (1 site). In this latter case, the documented site location was within 11 m of areas classified to be the very highest suitability class.

While documented wolf rendezvous sites did not occur within the remaining three study areas, there were eight documented sites nearby (4–21 km distant) and being located within the wolf RSF-mapped region encompassing all the study areas were thus suitable for use in model validation. One rendezvous site occurred 4 km northwest of Study Area 2 and was in an area classified to the very high suitability class (9). This site was located only 2 m from an area of the very highest class. Four other rendezvous sites were located 6–9 km north of Study Area 3. Two of these four sites occurred in the very high class (class 9) and were within 10 m and 20 m of areas classified to the very highest suitability class. One of the remaining sites occurred in the high-moderate class (class 7) and was located within 51 m of areas classified to the very highest suitability class. The fourth site occurred in the low (class 3) suitability class but was located < 103 m from areas predicted to have the very highest suitability for wolf rendezvous sites. Three additional rendezvous sites occurred about 14–21 km north or northwest of Study Area 4. Two of these three sites were in areas classified to the very highest suitability class. The third site occurred in high suitability (class 8), was adjacent to (8 m) areas in the very high class, and was within 82 m of areas in the very highest suitability class. A Spearman rank correlation score (r_s = 0.77), based on the total of 14 documented rendezvous site locations, indicated there was a strong correlation between counts of sites within suitability classes and the ranking of these classes.

Figure 9.

Map of predicted wolf-cattle encounter risk (five classes), derived from the combination of maps for predicted cattle relative probability of use and for relative habitat suitability for wolf rendezvous sites, for Study Area 1 (112 km²) in mountainous grazing lands of western Idaho. Also shown are locations of Global Positioning System–detected wolf-cattle encounters (n = 200) and confirmed wolf-caused depredations (n = 11) relative to roads, streams, and hill-shaded terrain shape.

Predicted wolf rendezvous–site habitat suitability for the four study areas can be summarized quite simply. Habitat types with the very highest suitability were areas of flat or gently sloping concave terrain. Contrary to descriptions by Ausband et al. (2010), these habitat types were most often forested. The conifer and, in some cases, aspen or cottonwood tree canopy in the present study resulted in the relatively high NDVI or vegetation greenness values that Ausband et al. associated with highly suitable rendezvous site habitat. While often occurring along perennial stream courses, these habitat types were also frequently located in the heads of small side canyons where stream flow was seasonal or ephemeral. Those very highest suitability habitat types that did occur in open areas were typically located at relative high elevation near the heads of major drainages. Habitat types of very lowest suitability occurred on steep slopes, particularly southerly or westerly aspects, which were relatively sparsely vegetated. Apart from most very highly suitable habitat types being forested, these findings align quite well with those reported by Ausband et al. (2010) for their model application in central Idaho.

These results then beg the question, if habitat types can be effectively classified and mapped with regards to their relative suitability as wolf rendezvous sites, would one then expect to observe higher wolf presence during the rendezvous season (15 June–15 August) in habitat types of the very highest rendezvous-site suitability compared with those of lower suitability classes? This question was evaluated in Study Area 1 using GPS tracking data from the adult male wolf representing the Snake River pack during the rendezvous period of 2009. Of the 6 591 positions acquired in this study area during the rendezvous period, about 26.3% (1 734) of the positions were located in areas classified to the very highest rendezvous site suitability class (class 10). An additional 27.3% (1 798) of the positions occurred in the very high class; thus combined, more than half (53.6%) of the positions from this wolf occurred in the top two suitability classes. A Spearman rank analysis comparing counts of positions within classes to ranking of those classes yielded a very high score of 0.988, indicating spatial presence levels of this wolf, and probably the Snake River pack in general, were strongly correlated with predicted habitat suitability for rendezvous sites.

Predicting Wolf-Cattle Encounters

Coregistered cattle and wolf resource selection maps were evaluated as a means of predicting the spatial risk of wolf-cattle encounters within the four study areas. The weighed Spearman rank correlation analysis indicated strong similarity existed between the cattle and wolf maps for all study areas except Study Area 2, where the similarity was weak (r_s = 0.20) but still significant (P < 0.05). Spearman scores for Study Areas 1, 3, and 4 were 0.70, 0.51, and 0.58, respectively. These results suggested a positive relationship existed, during the 15 June–15 August time period, between predicted cattle resource selection patterns and wolf rendezvous site habitat suitability patterns. On the basis of the final RSF model, cattle were predicted to select for meadows, benches, saddles, and other relatively flat and productive areas. These same areas were classified by the Ausband et al. (2010) model as very highly suitable for wolf rendezvous sites. Cattle were predicted to avoid areas of steep slopes, especially those only sparsely vegetated, and these same areas were also predicted to have very low suitability as wolf rendezvous sites.

Wolf-cattle encounter risk maps were created for each study area using these coregistered cattle and wolf maps. Figures representing these 10-class encounter risk maps are included in the online supplementary materials associated with this paper ( Figs. S4– S7, available online at  https://doi.org/10.1016/j.rama.2019.08.012). Simpler, more easily interpretable 5-class risk maps, however, were found to be more useful for management applications. The 5-class wolf-cattle encounter risk map for Study Area 1 is illustrated in Figure 9.

Case-Study Validation of Risk Map

Concurrent wolf and cattle GPS tracking data acquired in 2009 during the rendezvous period were used in a case-study validation of the wolf-cattle encounter risk map for Study Area 1. About 62% of the 200 independent encounters, identified with these concurrent GPS data, occurred in the very high risk class (class 5; see Fig. 9). The top two risk classes accounted for 84% of all encounters. The moderate-, low-, and very-low-risk classes contained 24, 8, and 0 encounters, respectively. The Spearman rank correlation analysis indicated a strong correlation (r_s = 0.99) between encounter counts and risk class ranking. Wolf-cattle encounters occurred a mean distance of 32.7 ± 46.5 SD m (range 0–253 m) from areas classified to the very high risk class. Mean distance between random points (n = 200) and areas of the very high risk class was considerably longer at 237 ± 309 SD m (range 0–1 650 m, P < 0.0001). Combined, these results clearly demonstrate the risk map was highly successful at predicting where, within the rugged and diverse landscape represented by Study Area 1, wolf-cattle encounters were most likely to occur.

Predicting Depredation Risk

The case study described above revealed the spatial risk of wolf-cattle encounters can be accurately predicted, and although this knowledge is quite important for management, a critical question is still unanswered. Does increased spatial knowledge of encounter risk lead to greater predictability of where depredation events will occur within rugged and remote cattle grazing lands? Fortunately, a small data set of confirmed wolf depredations (n = 16) was available to address this question at three of the four study areas.

At Study Area 1, five of 11 (45.5%) confirmed depredations occurred in areas classified to the very high wolf-cattle encounter risk class (class 5). About 72.7% of the 11 depredations occurred within the top two risk classes. One depredation each occurred in the remaining three risk classes. The mean distance between depredation locations and areas classified to the very high risk class was 25.9 ± 42.6 m SD (range of 0–132 m).

Two cattle depredations were confirmed in Study Area 2, and both occurred in the high-risk class (class 4). These depredations were a mean distance of 21.2 ± 27.4 m SD (range of 12.5–51.2 m) from the nearest areas of the very highest risk class.

Three depredations occurred at Study Area 3, two in areas classified to high risk (class 4) and one in low risk (class 2). The mean distance between depredation locations and areas classified to the very highest risk class was 34.6 ± 51.1 m SD (range of 2.6–93.5 m). No confirmed depredation data were available for Study Area 4.

The Spearman rank correlation analysis using depredation count data pooled across all study areas indicated the wolf-cattle encounter risk map effectively predicted spatial risk of wolf-caused depredation (r_s = 0.77). The Mann-Whitney U-test also using pooled data demonstrated the mean distance between depredation locations and nearest areas of the very high risk class (28.3 ± 40.3 m SD) was shorter than between random points and nearest areas of this class (181 ± 219 m SD; range 1.6–1 650 m; P < 0.0001). Although the sample sizes of confirmed depredation locations used for validating the wolf-cattle encounter maps of these three study areas were quite small, these results do lend some evidence that increased spatial knowledge of cattle-wolf encounter risk will help us better understand and predict where wolf-caused livestock death and injury losses are most likely to occur on mountainous cattle grazing lands.

Previous Spatial Risk Modeling Work

Spatial risk modeling is a common analytical means for predicting and mapping relative risk across landscapes. These analyses are often used to assess potential consequences of environmental contamination (Li et al. 2007; Carlon et al. 2008), natural disasters (Jaiswal et al. 2002; Keef et al. 2009), disease transmission (Tachiiri et al. 2006; Beck-Worner et al. 2007), and crime and human conflict (Kennedy et al. 2011; Rustad et al. 2011). Risk modeling is also an effective means of spatially quantifying the risk of predator-livestock conflicts (Mech et al. 2000; Miller et al. 2015a, 2015b). Considerable work has been conducted in the Upper Midwest to predict spatial risk of wolf-caused depredations (Treves et al. 2004, 2011; Edge et al. 2011; Treves and Rabenhorst 2017), but research in the NRM is limited (e.g., Bradley and Pletscher 2005; Hanley et al. 2018a, 2018b). A matched-pair analysis approach, contrasting sites affected and unaffected by depredation, was developed and applied in Wisconsin and Minnesota (Treves et al. 2004) and also applied in Michigan (Edge et al. 2011) to predict wolf depredation risk at township (92 km²) and farm-vicinity scales (10 km²). Risk mapping was done with a state-wide scope, but only at the township scale, which was certainly useful for broad-scale management planning and policy making but finer-scale mapping, is also critically needed. Working in Idaho and Montana, Bradley and Pletscher (2005) also contrasted depredated and nondepredated sites but at the finer, fenced-pasture scale (~201 ha) and determined predictors such as elk presence, pasture size, cattle herd size, vegetation cover, and distance from human residences affected wolf depredation risk while husbandry predictors did not. However, neither the work in the Midwest nor the NRM effectively incorporated wolf and cattle spatial behavior into modeling approaches except through rather coarse measures (e.g., wolf pack range extents). Rather than relying entirely on data from past depredations, which no doubt carry with them some degree of hidden bias, the resource-selection modeling approach applied in the present study began by first understanding wolf and cattle space-use behaviors, then identified consequent spatial overlaps, and finally predicted wolf-cattle encounter and conditional depredation risks based on these modeled behavioral responses. The limitations, implications, and impacts of this RSF modeling approach to spatial risk prediction are discussed in this section.

Robustness of Cattle RSF Model

The Spearman rank score for model predictive performance was notably lower for Study Area 4 in 2011 (r_s = 0.81) than for any other study area-yr combination (r_s > 0.93). Additional factors were likely affecting cattle selectivity in this case. Plots of the 2011 cattle GPS positions relative to GIS layers for elevation and vegetation cover types revealed cattle largely remained in the lower-elevation portions of the study area and generally occupied open or sparsely wooded areas throughout most of the summer grazing season. Specifically, these positions were clustered on bunchgrass-dominated hillslopes, openings within ponderosa pine forest, areas where shelter-wood silviculture treatments had created sparse tree canopy, and at midelevations, in small meadows surrounded by extensive granite outcrop. These selectivity patterns were a clear departure from those of 2009 and 2010 when cattle were initially clustered under forest canopy at lower elevations, moved progressively up through forested slopes and stream courses as the season advanced, and finally were distributed among subalpine meadows in late summer. Median elevation of cattle positions in 2011 (1 571 m) was about 377 m and 442 m lower than the medians for 2009 and 2010. Combined, these responses to elevation and vegetation cover type suggest climate, specifically an unusually deep and lingering snow pack, played a role in shaping cattle resource selection patterns in 2011. Peak snow water equivalent (SWE) for the 2011 water yr at the Brundage Reservoir SNOTEL station (ID = 370) near Study Area 4 was 980 mm on 3 May 2011, which was 383 mm and 246 mm greater and occurred 20 and 31 d later than peak SWE during 2010 and 2009, respectively (NRCS 2018a). The last SWE value (i.e., indicating the last measurable presence of snow) recorded at this station in 2011 was on 27 June, which was 20 and 28 d later than 2010 and 2009, respectively.

Given the unusual snow conditions in 2011, cattle at Study Area 4 seemed to have selected for openings and other sparsely wooded areas where snow had melted out earlier than areas shaded by forest canopy and/or at higher elevations where snow would have lingered well into the summer, thus delaying cattle access and forage production. Study Area 4 has the highest base elevation of the four study areas and thus was the most susceptible to late snow-cover anomalies. Although annual precipitation amounts were higher than long-term averages at some of the other climate stations in 2011, the unusual snow conditions and resultant cattle selective patterns observed at Study Area 4 did not extend to the other study areas (MesoWest 2018; NRCS 2018c). Therefore, even when confronted with an anomalous situation, the cattle RSF model actually performed quite well at Study Area 4 during 2011, regardless of its comparison with the other study area–yr combinations.

Distance from Roads and Streams

The map of predicted cattle use patterns for Study Area 1, at first glance, seems to indicate cattle would select for areas located near roads and perennial streams (see Fig. 1). Conversely, a look at the fitted response curves for the road and stream distance predictors (see Figs. 6–7) seems to suggest cattle would make little use of areas within 1 500 m of roads and perennial streams. Clearly, there is some complexity here that needs to be worked out. Cattle were generally predicted to select for the eastern portion of this study area, including areas that happened to be near roads and streams, and to largely avoid the western portion. In so doing, cattle would also avoid areas near and up to a considerable distance (e.g., 1 500 m) from the roads and extensive stream network situated in the western portion of the study area. Therefore, while distance to roads and streams were significant components of the cattle RSF model at Study Area 1, the nature of their influence was complex and likely interactive with terrain slope and other predictors. Generally, slope differed considerably between the eastern portion, located on the relatively flat or rolling top of the plateau, and the western portion of the study area, which occurred on the complex boundary slopes of the plateau. Cattle were predicted to select for areas of < 10 degrees slope, which occurred primarily in the eastern portion, and avoid areas of steep slopes (> 20 degrees), which generally occurred in the western portion. Furthermore, among available flat or gently sloping areas, cattle were predicted to select specifically for those areas dominated by ponderosa pine cover. As such, slope and pine cover were the primary and secondary drivers of cattle selectivity and distance to roads and streams seemed to have played tertiary, and less readily interpretable, supporting roles at Study Area 1.

Applicability Scope

Inclusion of ponderosa pine predictor into the final cattle RSF model may seemingly limit the effective scope of this model and the cattle-wolf encounter risk mapping supported by it. Although ponderosa pine occurs widely throughout the western United States (NRCS 2018b), this vegetation type is generally confined to relatively low-moderate elevations. This raises the question “Will the cattle RSF model perform adequately for grazing areas located entirely above the ponderosa pine zone?” Unfortunately, this question cannot be fully resolved without additional testing at grazing areas of this type. However, given the effective performance of the cattle RSF model at Study Areas 3 and 4, which contained extensive areas in the mixed conifer and spruce-fir zones well above the ponderosa pine zone, it is reasonable to expect the model will perform well enough for initial applications. Further evaluations of the model under these applications may then reveal the need for refinement or specificity (e.g., perhaps replacement of the ponderosa pine predictor with mixed conifer). A follow-on application of this research is planned for two study areas in central Idaho, where cattle have been GPS tracked since 2005 and which are located entirely above the ponderosa pine zone (Breck et al. 2012). Five additional study areas in northeastern Oregon with cattle GPS tracking data will also be used to evaluate the robustness of the cattle RSF model (Clark et al. 2017a, 2017b).

This discussion prompts an additional point concerning elevational limits. Ausband et al. (2010) cited unpublished data indicating documented wolf rendezvous sites in Idaho all generally occurred below 2 765 m elevation. In fact, these workers constrained validation of their habitat suitability maps to below this elevational threshold. This wolf selectivity response is probably explained by thermal constraints. Nights are typically quite cold at high elevations, even in summer; thus, rendezvous sites established above this threshold would expose wolf pups to adverse thermal stress (but see Mech 1993). Thermal conditions and other elevational effects may also constrain the distribution of wild ungulate prey, which could, in turn, also influence wolf selectively. If wolf rendezvous sites and thus concentrated wolf presence during the rendezvous period are constrained below this or a similar elevational threshold, then by extension the wolf-cattle encounter risk mapping concept presented here would also be constrained by that threshold.

Limited Validations

Validation of the wolf RSF model for rendezvous site habitat suitability was limited by the rather small sample of documented rendezvous sites (n = 14) associated with the four study areas. However, Ausband et al. (2010) thoroughly validated the model at their central Idaho study areas, about 10–150 km east of the present study areas, using a large sample of 178 rendezvous sites. Results from the present study largely conform to those from Ausband et al. (2010), thus providing some confidence this simple three-variable model is quite robust. Nevertheless, there were some discrepancies to address. Although Ausband et al. described most of their surveyed rendezvous sites as being in wet meadows, most areas classified to the very highest suitability class at the four western Idaho study areas were forested rather than meadow settings. In addition, 11 of the 14 documented rendezvous sites in the present study were located in forested areas. Of the remainder, one was on an open, relatively dry slope and only two were in wet meadows. The wolf biology literature certainly provides evidence wolves use forested areas as rendezvous sites (Ballard and Dau 1983; Theuerkauf et al. 2003; Capitani et al. 2006). Therefore, it seems both wet meadows and forested areas, particularly those on flat or gently sloping terrain near perennial water sources, can provide very highly suitable rendezvous sites. Discrepancies between wolf resource selection observations in central and western Idaho may then simply stem from differing availability of wet meadow areas rather than the robustness of the wolf RSF model.

The wolf-cattle encounter risk maps developed in this study were validated using a large sample of GPS-based encounter observations (n = 200) but only at one of the four study areas. Collection of concurrent and intensive GPS observations of sympatric wolves and cattle was arguably the most effective and least-biased approach for validating the encounter-risk maps. In fact, there were few viable alternatives (e.g., direct visual observation, automated camera systems, backtracking) and these would have suffered from site accessibility limitations, view occlusion, and other sampling biases. Unfortunately, concurrent wolf-cattle GPS data are extraordinarily challenging to acquire given the cooperation, commitment, and consensus required among affected government agencies, cattle producers, and other stakeholder groups to implement politically charged data collections of this sort. In this study, these challenges were only surmounted at one study area and only for one study yr. Therefore, on one hand, the case study conducted at Study Area 1 in 2009 was a rigorous validation of encounter risk mapping within that limited scope. On the other hand, robustness of the encounter risk mapping concept could not be effectively evaluated and thus awaits further validation efforts.

Use of wolf-cattle encounter risk maps for predicting spatial risk of wolf-caused depredation was validated at three of the four study areas, but sample sizes for confirmed depredations within each of these study areas were small. However, when considered in composite, the data set of 16 total depredations represents a moderate, more statistically effective sample size for validation. The three main findings from this data set were 1) 75% of these depredations occurred in the top two risk classes, 2) a strong positive relationship existed between among depredation counts and risk class ranking (r_s = 0.77), and 3) depredations occurred closer to areas of very high encounter risk than would be expected at random; all lend considerable support to conclusions that wolf-cattle encounter risk maps are indeed effective for predicting the spatial risk of wolf depredation. Further efficacy evaluations are certainly needed, particularly in other regions experiencing wolf-cattle conflict, to assess applicability and robustness of the encounter risk mapping concept. These evaluations could be facilitated and greatly enhanced by routine collection of site coordinates (GPS) during formal depredation investigations.