Prioritizing winter habitat quality for greater sage-grouse in a landscape influenced by energy development

Prioritizing habitats that provide the best options for the persistence of sensitive species in human-modified landscapes is a critical concern for conservation. Linking occurrence and fitness parameters across multiple spatial scales provides an approach to address habitat prioritization for species of concern in disturbed habitats. To demonstrate the usefulness of this approach, we generated resource selection and survival risk models as a framework to quantify habitat value for wintering female greater sage-grouse (Centrocercus urophasianus) inhabiting a 6,093-km study area in northwest Colorado and south-central Wyoming, USA, being developed for oil and natural gas reserves. Our approach allowed us to evaluate the relative influence of anthropogenic development and environmental attributes characterizing a large landscape on habitat selection and habitat-specific survival in winter for female sage-grouse. When combined, these models provided a spatial representation of habitat quality to inform management and conservation of critical wintering habitats. We used 537 locations from 105 radio-marked female grouse obtained from 18 fixed-wing flights across winters 2007–2008, 2008–2009, and 2009–2010. Wintering sage-grouse selected areas with higher wetness potential (0.75-km scale), intermediate (quadratic form) total shrub cover (18.83-km scale), higher variability in shrub height (18.83-km scale), and less heterogeneity in Wyoming big sagebrush (Artemisia tridentata wyomingensis; 4.71-km scale) cover and total shrub cover (18.83-km scale). Anthropogenic surface disturbance (0.75-km scale) was negatively associated with occurrence. Winter survival for female grouse was positively correlated with heterogeneity in big sagebrush cover at the 0.75-km scale, but negatively correlated with heterogeneity in total shrub cover at the 18.83-km scale. We did not detect an association between anthropogenic variables and female winter survival. However, displacement of sage-grouse in the energy extraction area may have masked our ability to identify anthropogenic variables potentially influencing survival. Our winter habitat quality map indicated highly effective winter habitat (high occurrence-low survival risk) was limited, only representing 17.1% of our study area. Consequently, displacement from these limited, high-quality winter habitats could have profound consequences to population persistence.


INTRODUCTION
Identification of habitat requirements for critical life stages in an animal's annual life cycle is important when prioritizing habitat for a particular species.Habitat selection studies provide information regarding habitat characteristics that promote species occurrence.Identification of habitats that facilitate proportionally high use by individuals is particularly important for species of conservation concern.In addition, species may utilize a spectrum of habitat quality within their range (Donovan and Thompson 2001), making inferences about habitat use difficult.Moreover, animal occurrence alone may be a misleading indicator of population fitness (Van Horne 1983, Aldridge and Boyce 2007, Chalfoun and Martin 2007) as it relates to habitat performance (i.e., fitness consequences relative to habitat selection choices).Habitat quality is not simply the relative presence or density of individuals, but a function of an occupied habitat's conduciveness to survival and reproduction (Van Horne 1983, Hall et al. 1997).Thus, linking survival to habitat use is one way to accurately assess habitat quality to prioritize conservation efforts focused on population persistence of declining species.
Range-wide declines in the greater sage-grouse (Centrocercus urophasianus; hereafter, sagegrouse) have led to its recent designation as ''warranted but precluded'' from listing under the U.S. Endangered Species Act of 1973 (U.S. Fish and Wildlife Service 2010).Identification of seasonal habitat needs and prioritizing the best habitats for its conservation are thus of critical importance.Sage-grouse utilize an array of habitat characteristics within the sagebrush (Artemisia spp.) biome for nesting, brood-rearing, and wintering (Schroeder et al. 1999, Connelly et al. 2004).Habitat selected by sage-grouse in winter generally occurs in flatter topography dominated by large, continuous patches of sagebrush in areas away from conifer (Eng and Schladweiler 1972, Doherty et al. 2008, Carpenter et al. 2010).Sage-grouse rely exclusively on sagebrush as a food source in winter (Patterson 1952, Wallestad et al. 1975), and characteristics such as variation in topography and sagebrush height are believed to be essential in providing access to sagebrush forage under varying snow conditions and depths (Beck 1977, Remington and Braun 1985, Hupp and Braun 1989, Homer et al. 1993, Schroeder et al. 1999, Connelly et al. 2000, Frye et al. 2013).Many individual sagegrouse within populations make large movements (.10 km) from summer to winter habitats (Beck et al. 2006, Bruce et al. 2011, Fedy et al. 2012).Thus, it is likely that grouse undergoing long-distance movements from summer to winter range seek habitats that provide specific conditions favorable to winter survival.
Average annual range-wide breeding-age mortality rates for female sage-grouse range from 55 to 75% (Schroeder and Baydack 2001).Fallwinter survival estimates for sage-grouse tend to be higher (78-97%; Beck et al. 2006, Baxter et al. 2013); however, winter survival can be a limiting factor for female sage-grouse and may be depressed during severe winters (Moynahan et al. 2006, Anthony andWillis 2009).Low winter survival is attributed to greater snow accumulation during extreme cold that excludes birds from forage, and thermal and predator concealment cover (Moynahan et al. 2006, Anthony andWillis 2009).In addition, winter survival may also be lower for high elevation populations that experience extreme weather conditions (Anthony and Willis 2009).Because sage-grouse are a relatively long-lived species with moderate reproductive rates, population persistence is particularly sensitive to adult female survival (Johnson and Braun 1999, Schroeder et al.1999, Taylor et al. 2012).A recent meta-analysis indicated that female survival is the most significant demographic rate, followed by chick survival, and finally nest success (Taylor et al. 2012).Therefore, adult female survival in winter has consequential implications to sage-grouse population persistence (Moynahan et al. 2006).
Linking habitat occurrence with demographic parameters is particularly important in humanmodified landscapes where habitat quality may be compromised when habitats are lost or fragmented (Jones 2001, Schlaepfer et al. 2002, Aldridge and Boyce 2007).Animals may avoid previously suitable habitat when anthropogenic features are introduced to landscapes (Dyer et al. 2001, Johnson et al. 2005, Sawyer et al. 2006).Avoidance behavior associated with exposure to energy development is well documented in prairie grouse (Holloran 2005, Pitman et al. 2005, Aldridge and Boyce 2007, Doherty et al. 2008, Pruett et al. 2009, Carpenter et al. 2010, Dzialak et al. 2012).In addition, variation in avoidance response to energy field infrastructure among sage-grouse individuals has been documented, showing avoidance during the day, but not at night, suggesting avoidance of human activity (Dzialak et al. 2012).Sage-grouse exhibit strong fidelity to seasonal habitats, including wintering areas, and show little flexibility to change once these areas have been selected (Connelly et al. 2004).Thus, it is imperative that the quality of critical wintering areas is not compromised by anthropogenic development that could lead to wintering areas becoming population sinks (Delibes et al. 2001, Aldridge andBoyce 2007).
Research on other avian (Chalfoun and Martin 2007, Liebezeit et al. 2009, Gilbert and Chalfoun 2011) and non-avian (Johnson et al. 2004) species has demonstrated variations in fitness rates (survival and reproduction) related to anthropogenic features.Predation is the proximate cause of mortality for sage-grouse through all life stages (Schroeder et al. 1999, Connelly et al. 2011) and can be negatively correlated with sagegrouse fitness parameters (Bui et al. 2010, Hagen 2011, Dinkins 2013).Anthropogenic features and associated activities have been shown to alter predator composition and abundance when introduced to undisturbed habitats (Knight andKawashima 1993, Restani et al. 2001).The influence of anthropogenic features on major predators of sage-grouse, positive or negative, has potential to benefit (creating local refugia from predation) or hinder sage-grouse populations (increase predator abundance).
Several studies have reported avoidance of energy infrastructure and negative effects of development on sage-grouse population parameters during reproductive life stages (Lyon and Anderson 2003, Holloran 2005, Aldridge and Boyce 2007, Walker et al. 2007, Holloran et al. 2010, Kirol 2012).However, no study has explored winter female survival as it relates to habitat-specific survival and energy development.Because female survivorship is a key element of population productivity, our study in northwest Colorado and south-central Wyoming, USA sought to understand the suite of ecological conditions of winter habitats that were most critical in promoting population persistence within a developing oil and natural gas field.The specific objectives of our study were to: (1) develop winter resource selection and habitatspecific survival models for female sage-grouse in an area of ongoing oil and gas extraction, (2) evaluate the relative influence of environmental characteristics and anthropogenic features on winter habitat selection and survival of female sage-grouse, and (3) spatially depict (map) habitat quality to identify areas of high sagegrouse conservation and management importance.Linking winter habitat occurrence with survival provided us a means to prioritize habitat quality for sage-grouse within a landscape undergoing energy development.Prioritization of habitat quality is one way that conservation efforts may best focus limited resources to promote population persistence for species of concern.

Study area
Our study area encompassed 6,093-km 2 in southeastern Sweetwater and southwestern Carbon Counties, Wyoming, and northern Moffat County, Colorado, USA (41817 0 55.874 00 N, 107844 0 6.904 00 W; Fig. 1).The study area was within the Wyoming Basin Sage-Grouse Management Zone (U.S. Fish and Wildlife Service 2010).The study area was characterized by a vast sagebrush steppe with elevations ranging from 1,909 to 2,529 m and low average annual precipitation (24 cm; Natural Resource Conservation Service 2006).The region encompassing the study area was typified by cool temperatures with average daily temperatures ranging between a low of À168C and a high of 0.58C in midwinter and between 138C and 248C in midsummer (U.S. Bureau of Land Management 2006).Temperature extremes ranged from À468C to 388C with the frost-free period generally occurring from mid-May to mid-September.Precipitation was evenly distributed throughout the year with minor peaks in May, July, and October.The snowiest months were in December and January with an average of 98.6 cm of snow falling during the year (U.S. Bureau of Land Management 2006).Because the study area encompassed a large landscape characterized v www.esajournals.orgby wide variation in elevation and topography, site-specific climatic conditions varied across the 3 winters encompassing our study.
The study area was composed of Bureau of Land Management (60%), lands owned by the states of Colorado and Wyoming (10%), and private ownership (30%).The study area encompassed approximately 50 discrete oil and gas fields of varying size with the majority in Wyoming (Wyoming State Geological Survey 2012).Within our study area we verified 2,512 established wells (oil, coal bed natural gas, and conventional natural gas) and approximately 6,918 km of paved and improved gravel roads by the end of the study in 2010.Additional land uses in our study area included livestock grazing and hunting, which typically occurred outside of winter.

Field procedures and monitoring
We captured female sage-grouse in spring (mid-March through late April) and late summer (late-August through September) 2007, 2008, and 2009 using established spot-lighting and hoopnetting protocols (Giesen et al. 1982, Wakkinen et al. 1992).In the spring we captured grouse from 14 leks that were evenly distributed across the study area to ensure equal capture effort and to obtain a random sample of the population (Manly et al. 2002).We secured VHF radio transmitters (Model A4060; Advanced Telemetry Systems Incorporated, Isanti, Minnesota, USA) to female grouse with a PVC-covered wire necklace.Transmitters weighed 22 g (;1.4% of mean female sage-grouse body mass); had a battery life expectancy of 789 days; and were equipped with motion-sensors (radio-transmitter pulse rate increased in response to inactivity after 8 hours).Female sage-grouse were captured and handled according to the University of Wyoming Institutional Animal Care and Use Committee approved protocols (03032009) and Wyoming Game and Fish Department Chapter 33 permits 572 and 699.
Sage-grouse were located with a fixed-wing airplane with aerial telemetry antennas at approximately 3 to 4 week intervals during winter from 1 November to 15 March 2007March -2008March , 2008March -2009March , and 2009March -2010. .Our winter season delineation corresponded to large movements of radiomarked females (see Fedy et al. 2012) and was corroborated by the winter focal period identified by Carpenter et al. (2010).We recorded use locations using Global Positioning System (GPS) receivers.We estimated error in aerial locations by relocating n ¼ 7 randomly placed transmitters at fixed locations in a blind trial.On each flight, we attempted to relocate each fixed location transmitter.We calculated location error as the median linear distance between estimated flight locations and true transmitter locations.

Landscape predictor variables
We considered a suite of predictor variables on the basis of a priori information from previous landscape-scale research (Homer et al. 1993, Aldridge and Boyce 2007, Doherty et al. 2008, Carpenter et al. 2010, Doherty et al. 2010, Kirol 2012; Table 1).These variables encompassed environmental and anthropogenic categories that we evaluated at three spatial scales: 0.490-km radii (0.75 km 2 ), 1.224-km radii (4.71 km 2 ), and 2.448-km radii (18.83 km 2 ) based on estimated median location error, weekly movements, and biweekly movements estimated from aerial relocations, respectively.
We used remotely-sensed sagebrush products to estimate shrub height for all shrub species, and percent canopy cover of sagebrush (all Artemisia species combined), big sagebrush, Wyoming big sagebrush, and all species of shrubs (Homer et al. 2012).In our study area, for example, shrub canopy cover may reflect a combination of shrubs co-occurring such as Wyoming big sagebrush, antelope bitterbrush, and greasewood.We calculated mean and standard deviation (SD) for shrub height and percent cover for each shrub attribute within each spatial scale, with SD serving as a proxy for habitat variability or heterogeneity (Kastdalen et al. 2003, Aldridge and Boyce 2007, Carpenter et al. 2010).Lastly, we assessed quadratic relationships in shrub variables to evaluate potential non-linear relationships.
We used a 10-m digital elevation map (DEM; U.S. Geological Survey 2011) to calculate slope, aspect and elevation and to calculate a Vector Ruggedness Measure (VRM) and a Topographic Wetness Index (TWI).VRM uses the variation in slope and aspect to create a single measure of terrain ruggedness (Sappington et al. 2007).We rescaled VRM values by multiplying the original values by 1000 for ease of interpretation.TWI measured wetness potential based on drainage of the local upslope and local slope (integration of slope and aspect; Sorensen et al. 2006, Theobald 2007).TWI also incorporated solar insulation to identify differences in north-and south-facing aspects and is useful for predicting soil moisture, plant productivity, and species density (Zinko et al. 2005, Sorensen et al. 2006).TWI is similar to the compound topographic index (CTI), but unlike CTI, it incorporates solar radiation.
We obtained energy well data for our study area, including type, location, status, production, and spud date (initiation of drilling) from the Wyoming Oil and Gas Conservation Commission (2010) and the Colorado Oil and Gas Conservation Commission (2010) databases.Because en-ergy development was ongoing during the study we time-stamped variables associated with fossil fuel development based on the spud dates of wells associated with these variables to accurately characterize when they were established.Therefore, spatial analysis for each winter only included infrastructure established prior to that winter, which enabled us to depict temporal additions to human infrastructure.In addition, we used 2006 and 2009 imagery obtained from the National Agriculture Imagery Program (NAIP; U.S. Department of Agriculture 2010) to inspect the analysis area and validate energy infrastructure.We quantified Euclidean distance and exponential decay functions (Fedy and Martin 2011) from grouse use and available locations to nearest improved gravel road and nearest energy well.Decay variables corresponded to each spatial scale in the form e (Àd/a) where d was the distance in km to each feature, and a was set to the radii of each of the three scales (Leu et al. 2011).Decay functions rescaled Euclidean distance variables between zero and one.These functions describe a non-linear relationship be-tween the outcome (i.e., dependent variable) and the predictor variable, where the effect of the predictor variable attenuated to almost nothing after a point specified by the scale.Decay variables were not assessed in our survival analysis.We quantified total linear distances of improved gravel roads, count of energy wells, and the area of total surface disturbance within each spatial scale.
We calculated the percentage surface disturbance by creating a disturbance layer.In developing this layer, we used disturbance data processed by the U.S. Geological Survey that identified well pads and the surrounding area impacted by those pads.Using these data as the baseline for our disturbance layer, we manually digitized additional surface disturbance that was not identified by the U.S. Geological Survey data for accuracy using NAIP imagery in 1-km 2 plots across the study area.Our final disturbance layer consisted of all energy infrastructure including well pads, compressor sites, transfer stations, and haul roads as well as human dwellings.Haul road disturbance was identified by buffering Count of energy wells within analysis region Quadratic transformations assessed.à In addition to its linear form, measured as decay distance in the form e (Àd/a) , where d was the Euclidean distance (m) to an anthropogenic feature, and a was equal to the radius of each scale (decay variables were not used in survival modeling).
§ Time-stamped on the basis of spud date to prevent inclusion of infrastructure into analysis prior to its existence.
v www.esajournals.orgroad line features by 10 m, representing the average road surface disturbance width in the study area.

Experimental design and statistical analysis
Occurrence analysis.-Weemployed a useavailability design to evaluate sage-grouse winter habitat selection (Boyce et al. 2002, Manly et al. 2002, Johnson et al. 2006), where we identified resource use as locations of radio-marked sagegrouse obtained from fixed-wing aircraft flights during 2007-2010 pooled across individuals to represent a population level habitat selection response (a Type 1 Design; Manly et al. 2002, Thomas andTaylor 2006).We evaluated useavailability likelihood with binary logistic regressions using the exponential link function to estimate a RSF (Manly et al. 2002, Johnson et al. 2006, McDonald 2013).The RSF took the following form: where  [Cox 1972]) with the ''counting-process'' method (Allison 2010) to identify relationships between landscape-scale predictor variables and sage-grouse winter survival.The Cox PH model is a robust semi-parametric model that is used to analyze the effect of nominal and continous scale variables on time-to-event data (Cox 1972), such as death (Hosmer and Lemeshow 1999), and is applicable to radio-telemetry data (Winterstein et al. 2001).The counting-process accounted for time-dependence and discontinuous hazard intervals.For the counting-process we distributed our survival data into time intervals for each unique individual.By using this method the baseline hazard was allowed to vary with time (e.g., exposure to variables changing with time; Allison 2010).When using the Cox PH model it is important that survival time relates to the observation interval (Winterstein et al. 2001).For this reason and because our winter survival data were based on radio-telemetry flights spaced approximately 3-4 weeks apart, we estimated weekly survival instead of daily survival.The risk of mortality (hazard ratio [h(tjx t )]) is a function of the non-parametric baseline hazard (h 0 (t)) and the parametric covariates (x's) affecting survival (Hosmer and Lemeshow 1999).The Cox PH model is expressed as: The baseline hazard is unspecified but the effects of the predictor variables are still estimated and interpreted as hazard ratios (exp[b i ]).We left censored individuals entering the study at different times and right censored individuals that did not die during the study.Censoring allowed us to incorporate individuals into the model that were not observed for an entire period or had an unknown event (e.g., a female disapeared from the sample due to a malfuntioning transmitter).The Cox PH model assumes that the hazards are proportional over time (proportional hazard assumption; Hosmer and Lemeshow 1999).Thus, we tested the variables in our top survival models independently and collectively (e.g., top models) for proportionality at a ¼ 0.5 (Le 1997, Hosmer andLemeshow 1999).We assessed if a particular observation was disproportionately influential on a coefficient estimate for each variable by testing for inflated residuals and leverage (dfbetas; Hosmer andLemeshow 1999, Allison 2010).No observations were removed as a result of these diagnostics.We used the covs option in Statistical Analysis Software (SAS) to obtain a robust sandwich estimate for the covariance matrix, which results in robust standard errors for parameter estimates and the aggregate option with the ID statement to provide a summing of the score residuals for each distinct variable pattern (SAS Institute 2011).
We constructed time (t) intervals for each individual sage-grouse by assigning variable information across intervals centered at the flight observation time to the midway point of the next flight observation.Thus time intervals may have included full weeks and partial weeks (e.g., week ¼ 3 to week ¼ 5.5).During each winter, survival time started at t ¼ 0 corresponding to the first flight in November and ended with the last flight in March each year.Therefore, our survival analysis included approximately 17 weeks of exposure time per winter.

Model development
We used variables with the most predictive potential to make population-level inferences about occurrence and survival (Boyce et al. 2002).Thus, initial steps for occurrence and survival modeling focused on pre-model screening to assess correlation and univariate significance to identify variables with a reasonable likelihood of association.Consequently, we removed non-informative variables (85% confidence intervals [CIs] around parameter estimates that included 0; Hosmer andLemeshow 1999, 2000;Arnold 2010).To prevent multicollinearity we omitted variables or variable combinations when r !j0.7j or t (tolerance statistic) j0.40j (Allison 2009, SAS Institute 2011).We also checked for stability and consistency of regression coefficient estimates when variables were moderately correlated (j0.3j r j0.7j).Undetected correlations between variables can cause instability in the signs of coefficients and also result in inflated standard errors (Doherty 2008).If variables were correlated, the variable with the lowest AIC or AIC SUR score (Liang and Zou 2008) was retained.We did not permit correlated variables to compete in the same model at any level of model selection.
For our occurrence modeling, we used Akaike's Information Criterion to assess model support.For all scale-dependent variables, we examined our three spatial scales to determine the scale that was most correlated to occurrence by testing each variable scale individually and comparing AIC scores (Arnold 2010, Carpenter et al. 2010, Doherty et al. 2010).For each variable we retained the scale with the lowest AIC score corresponding to the greatest predictive potential (Burnham and Anderson 2002).The initial variable screening removed unsupported predictor variables, thereby reducing our likelihood of overfitting models in our model selection process (Burnham andAnderson 2002, Arnold 2010).
The model fitting process for our survival model was consistent with that of occurrence modeling with a few variations.We used a derivation of the AIC technique adapted specifically for survival modeling (AIC SUR ) to select the best supported survival models (Liang and Zou 2008).Again, we examined three spatial scales to determine the scale that best explained survival by testing each variable-scale individually and comparing (AIC SUR ) scores (Arnold 2010, Carpenter et al. 2010, Doherty et al. 2010) for scaledependent variables.We retained the variable scale with the lowest (AIC SUR ) score.
For both occurrence and survival, we used a sequential model selection approach (Arnold 2010) by evaluating the relative importance of predictor variables for occurrence and survival within two variable subsets.In the first level of model selection, environmental and anthropogenic model subsets were modeled separately.Within these subsets we explored all variable combinations (Burnham and Anderson 2002).At this stage, we considered models with AIC or AIC SUR scores in the range of two to seven units (Burnham and Anderson 2002) to be competitive with the top model.However, models with AIC scores effectively equivalent (,2 AIC or AIC SUR ) to the null model were not considered competitive (Allison 2010).We assessed variable importance by summing Akaike model weights across models that included the variable of interest (Arnold 2010).We brought forward the variables with the greatest potential as predictors of occurrence or survival within each subset to the final level of model selection.
Top models within each variable subset (e.g., environmental and anthropogenic) that were brought forward were then allowed to compete across subsets to assess whether additional information produced a more parsimonious model (Arnold 2010).We judged improvements in model parsimony or fit by the weight of evidence (w i [Akaike weights]) and difference between AIC or AIC SUR for the top model and AIC or AIC SUR for the ith candidate model (D i; Burnham and Anderson 2002).For example, we explored whether the final model(s) from the environmental subset had the most support by itself, or if a combination of top models from environmental þ anthropogenic subsets produced a model with greater support.When a single top model was not apparent based on AIC or AIC SUR scores ( 7 units considered competitive) we used multi-model inference to calculate final parameter coefficients, 95% confidence intervals, odds ratios, and risk ratios within confidence sets.We determined confidence sets for those models where Akaike weights were within 10% of the top model (Burnham and Anderson 2002).At the final level of model selection we further screened variables with little support having parameter estimates with 95% CIs that overlapped 0 (Hosmer andLemeshow 1999, 2000).
We calculated winter survival estimates with the Kaplan-Meier product-limit estimator (Kaplan and Meier 1958) modified for staggered entry (Pollock et al. 1989).In following with our survival analysis period, Kaplan-Meier adult female winter survival estimates were calculated from November (t ¼ 0 weeks) to March (t ¼ 17 weeks) each year over the 3 winters.

Model validation
We performed 5-fold cross validation to evaluate the predictive performance of our occurrence model (Boyce et al. 2002).For each of the 5 data folds (bins) the withheld set was assessed against the model predictions of that training data set using correlations between bin ranks of the RSF values.A high score corresponded to good predictive performance (Boyce et al. 2002).We used the overall C statistic (C index), designed specifically for survival modeling, to assess the discrimination ability of our final survival models (Pencina and D'Agostino 2004), where C statistic values between 0.7 and 0.8 were considered to have acceptable discrimination, values between 0.8 and 0.9 have excellent discrimination, and values 0.5 indicate the model predicted the outcome no better than chance.For the survival modeling, we also assessed model fit with a robust sandwich estimate (Wald sandwich test) for the covariance matrix (Lin and Wei 1989) and influence statistics for each distinct observation (Langholz andJiao 2007, Dzialak et al. 2011).As an assessment of the survival model predictive performance, we tested for a statistical difference in weekly risk with one-tailed t-tests with unequal variance.Thus, if the survival model was accurately predicting survival risk the death locations should be exposed to significantly greater risk than the live locations.We conducted all statistical analyses with SAS version 9.2 (SAS Institute 2011).We report all Kaplan-Meier survival estimates 6 SE.

Map development
We mapped our final occurrence and survival models with 30-m pixel resolution across our study area.We rescaled the RSF predictions between 0 and 1 to map our model coefficients (Manly et al. 2002, DeCesare et al. 2012): For the occurrence map we distributed predicted probabilities or relative risk into five quantiles on the basis of percentile breaks in predicted probabilities (Sawyer et al. 2006).Similar to the occurrence map, we mapped survival risk on the landscape as the relative risk of weekly survival binned into four quantiles, but we also included a low occurrence bin where minimal use restricted us from predicting risk.Following Aldridge and Boyce (2007) and Kirol (2012), the RSF occurrence and the relative risk predictions were combined to rank habitat quality to form a winter habitat quality map.In forming the habitat quality categories and associated map, the low occurrence bin (see Fig. 2) was left unchanged whereas the top two and bottom two occurrence bins were merged to generate high and moderate occurrence bins, respectively.Similarly, the top two and bottom two relative risk bins were merged to form high and low risk bins, respectively.Occurrence and risk bins combined formed five winter habitat quality bins categorized as: highly effective winter habitat (high occurrence-low risk), effective winter habitat (moderate occurrence-low risk), ineffective winter habitat (moderate occurrence-high risk), highly ineffective winter habitat (high occurrence-high risk), and low occurrence winter habitat (risk was not predicted in the low occurrence bin because of minimal use).

Occurrence
The top environmental model included 6 predictor variables represented by three spatial scales.In the environmental modeling set, no other model was competitive with the top model (DAIC !24.6, w i ¼ 0.99).At the 0.75-km 2 (0.49km radius) scale, TWI had a strong positive correlation with winter occurrence.For every one unit change in TWI the relative probability of sage-grouse occurrence increased by approximately 65%.Standard deviation of Wyoming big sagebrush cover (Wysagesd) was a strong negative predictor of occurrence at the 4.71-km 2 Fig. 2. Predicted probability of sage-grouse winter habitat occurrence in northwest Colorado and south-central Wyoming, USA, winters 2007USA, winters -2008USA, winters , 2008USA, winters -2009USA, winters , and 2009USA, winters -2010.This map displays a resource selection function that was binned into five quantiles of predicted relative probability of occurrence.The white areas within the predictive map are areas with no data values in the Homer et al. (2012) products because they were masked as non-rangeland habitats.
v www.esajournals.org(1.224-km radius) scale with average standard deviation at use locations less than at available locations.A one unit increase in the heterogeneity (i.e., 1 SD) of Wyoming big sagebrush cover resulted in a decrease in relative probability of occurrence by approximately 38.6 times.At the largest scale (18.83 km 2 , 2.448-km radius) total shrub cover and its quadratic form and standard deviation of shrub height were positively correlated with occurrence, whereas standard deviation of shrub cover was negatively correlated with occurrence at the 18.83-km 2 scale.Relative probability of occurrence increased by approximately 1.8 times for every one unit (i.e., 1%) increase in total shrub cover up to 13% cover; for a one unit increase in total shrub cover greater than 13% cover, we predicted a decrease in occurrence by approximately 3.6 times.We predicted that a one unit increase in variability of shrub height would result in 6.3 times higher relative probability of occurrence.Conversely, the relative probability of occurrence decreased by approximately 5.1 times with every one unit increase in shrub cover heterogeneity.
When we combined environmental and anthropogenic models, the inclusion of anthropogenic variables improved model fit over the environmental model by DAIC ¼ 12.068 to 10.079 (Table 2).Anthropogenic variables associated with improved model fit included surface disturbance (Dstbarea) and count of energy wells (wells) at the 0.75-km 2 (0.49-km radius) scale.Model averaging indicated that the 95% confidence interval for the odds ratio estimate for count of energy wells overlapped one; therefore our model with the greatest support was restricted to a single anthropogenic variable (Dstbarea) and included 5 environmental variables (Table 3).The relative probability of occurrence decreased by approximately 3.3% for every 1% increase in surface disturbance within 0.75 km 2 .The interpretations of change in odds ratios (occurrence probabilities) per one unit change in variables presented were calculated  , winters 2007-2008, 2008-2009, and 2009-2010.

Model
Model fit statisticsà  Note: Variable importance, odds ratios and 95% confidence intervals not calculated for intercept parameter.
v www.esajournals.orgas the median change in odds ratios bound by the range of variable values for that variable with the other variables in the model held at their mean value.Cross-validation indicated our final model (top model; Table 2) was a strong, positive predictor of grouse winter habitat selection in our study area (r s ¼ 0.99, P ¼ 0.001, n ¼ 10).Our top model with model averaged coefficients was used to map winter occurrence across the study area as a resource selection function binned into five quantiles based on the relative probability of occurrence (Fig. 2).
The winter survival model had moderate discrimination ability, C statistic of 0.64 and the Wald sandwich test indicated good fit to the data (v 2 ¼ 32.12, df ¼ 2, P , 0.0001).The survival model contained the predictor variables standard deviation in big sagebrush cover and standard deviation in shrub cover.We did not detect an association between anthropogenic variables and female survival risk in the winter.At the 0.75-km 2 (0.49-km radius) scale Bsagesd was positively related to weekly survival, and at the largest scale (18.83 km 2 , 2.448-km radius), Shrubsd was negatively associated with weekly female survival (i.e., riskier habitat; Table 4).Within a 0.75-km 2 area a one unit change in the heterogeneity (i.e., 1 SD) in big sagebrush cover resulted in the relative probability of weekly survival increasing by approximately 1.5 times.A one unit increase in heterogeneity in total shrub cover within an 18.83km 2 area resulted in the relative probability of weekly survival decreasing by approximately 2.3 times.This survival model applied to our study area landscape predicted that death locations (n ¼ 20) were exposed to statistically greater risk (t 20 ¼ 3.44, P ¼ 0.001) compared to live locations (n ¼ 507).The interpretations of change in hazard ratios per one unit change in variable presented were calculated as median change in hazard ratios bound by the range of values for that variable with the other variables in the model held at their mean value.

Winter habitat quality
Our habitat selection model predicted 39.4% of the study area was high occurrence (Fig. 2), and of the high occurrence areas, less than half (43.4%) were low risk (Fig. 3).Highly effective (high occurrence-low risk) winter habitat encompassed only 17.1% of the study area (Fig. 4).Of the remainder of occupied habitat within the study area, the habitat quality map predicted that 21.6% was effective (moderate occurrence-low risk), 19.3% was ineffective (moderate occurrence-high risk), and 22.4% was highly ineffective (high occurrence-high risk) winter habitat.Further, 19.6% of the study area was largely unoccupied (low occurrence).

DISCUSSION
We evaluated the importance of environmental and anthropogenic landscape variables for fe- § Variable importance was calculated by adding Akaike weights for all models a variable occurred in-the closer the value was to 1 the more important the variable was in the set.
} A hazard ratio .1 indicates that per-unit time (e.g., week) the likelihood of death increases with an increase in the predictor variable.
v www.esajournals.orgmale sage-grouse winter occurrence and survival across multiple scales.We further used the most robust predictor variables from these occurrence and survival models to map and categorize winter habitat quality across a large landscape used by wintering sage-grouse.Linking habitat occurrence with demographic parameters is imperative to understanding habitats that promote local population persistence especially within landscapes undergoing human disturbance.Sage-grouse have been shown to alter their winter habitat use relative to the proximity and density of energy field infrastructure (Doherty et al. 2008, Carpenter et al. 2010); however, effects of energy infrastructure on survival have not been studied.Through linking habitat occurrence with female survival we predicted that 17.1% of the landscape in our study area was composed of high occurrence-low risk habitat for wintering sage-grouse.
Our environmental occurrence model was a good predictor of winter habitat selection by female sage-grouse.Wintering sage-grouse selected areas with higher TWI, intermediate total  , winters 2007-2008, 2008-2009, and 2009-2010.Relative risk binned into four quantiles with a low occurrence bin (Fig. 2) where minimal use prevented us from predicting risk.High risk values indicate that female sage-grouse were less likely to survive the winter in these areas.The white areas within the predictive map are areas with no data values in the Homer et al. (2012) products because they were masked as non-rangeland habitats.
v www.esajournals.orgshrub cover, higher variability in shrub height, and less heterogeneity in Wyoming big sagebrush and total shrub cover across three spatial scales.Our results were consistent with Doherty et al. (2008) who found that sage-grouse winter occurrence was more likely in areas of greater sagebrush cover at a scale of 4 km 2 and with others who reported the quadratic form of big sagebrush (Carpenter et al. 2010, Dzialak et al. 2012) or total shrub cover (Dzialak et al. 2012) was a better univariate predictor of winter sagegrouse occurrence than total shrub or big sagebrush cover alone, suggesting that grouse were selecting areas of moderate shrub cover while avoiding the lowest and highest cover available.Selection for intermediate sagebrush cover has also been documented for nesting (Aldridge andBoyce 2007, Doherty et al. 2010) and brood-rearing female sage-grouse (Aldridge andBoyce 2007, Aldridge andBoyce 2008) at both local and landscape scales.Sage-grouse in our study also selected areas with less heterogeneity in Wyoming big sagebrush cover and total shrub cover.Other studies have demonstrated sage-grouse selection for continuous cover of sagebrush in winter (Doherty et al. 2008, Car-Fig. 4. Winter habitat quality map as a function of predicted relative probability of occurrence combined with relative risk of weekly survival for female sage-grouse in northwest Colorado and south-central Wyoming, USA, winters 2007USA, winters -2008USA, winters , 2008USA, winters -2009USA, winters , and 2009USA, winters -2010. .Map displays five categories of winter habitat quality.The white areas within the predictive map are areas with no data values in the Homer et al. (2012) products because they were masked as non-rangeland habitats.
v www.esajournals.orgpenter et al. 2010).Wintering grouse in our study used areas with greater TWI values.TWI is useful in predicting soil moisture, plant productivity, and species density (Zinko et al. 2005, Sorensen et al. 2006).Sagebrush patches in areas with higher TWI are likely more dense and TWI may be positively correlated with individual sagebrush plants that have access to greater resources; thus, providing higher nutrient quality favored by sage-grouse in winter (Remington andBraun 1985, Frye et al. 2013).Consequently, higher TWI may be a proxy for high quality sagebrush that protrudes above snow and is likely of greater importance to grouse winter habitat selection following major snow accumulation events.Finally, our results suggest that selection for greater variability in shrub height may be important for food and thermal cover when faced with varying amounts of snow throughout winter (Beck 1977, Remington and Braun 1985, Hupp and Braun 1989, Homer et al. 1993, Schroeder et al. 1999, Connelly et al. 2000, Frye et al. 2013).
The inclusion of anthropogenic variables in our habitat selection models received overwhelming support during model selection, suggesting that anthropogenic features were negatively influencing habitat selection during winter.We found that sage-grouse in our study were avoiding areas with greater surface disturbance within a 0.49-km radius.Although the same relationship was true for oil and natural gas well density explicitly, it was not as well supported as surface disturbance.Surface disturbance was likely a better predictor because it represented a combination of all energy infrastructures such as compressor stations, roads, and well pads.During winter, Doherty et al. (2008) found that sage-grouse in northeast Wyoming and southcentral Montana were more likely to inhabit sagebrush habitats absent of natural gas wells within a 4-km 2 area.Our results suggested that avoidance of surface disturbance due to fossil fuel development resulted in indirect loss of otherwise suitable and perhaps critical winter habitat.
Our survival analyses results illustrate habitatspecific variations in survival or risk across the landscape at multiple scales.That is, in the winter, habitats being used by sage-grouse fell on a spectrum of high to low survival risk at finer and broader scales.We found that habitats with greater heterogeneity in total shrub cover at a large landscape scale (18.83 km 2 ) were riskier habitats for female sage-grouse in the winter.Thus, large habitat areas with more homogenous shrub cover were safer habitats.However, at a finer scale (0.75 km 2 ) heterogeneity in big sagebrush cover was positively associated with winter survival.This suggests that risk of death was reduced when sage-grouse selected homogenous shrub cover at larger spatial scales interspersed with patches of variable sagebrush cover at smaller scales.
Human disturbance may have important implications for habitat effectiveness, individual fitness, and population productivity for a range of species (see Johnson et al. 2005, Gilbert andChalfoun 2011).Previous research on sagegrouse suggests oil and gas development influences sage-grouse populations through lowered male lek attendance and declining lek persistence (Holloran 2005, Walker et al. 2007, Hess and Beck 2012); lower yearling male recruitment to disturbed leks (Holloran et al. 2010); lower nest initiation rates (Lyon and Anderson 2003); lower annual adult female survival (Holloran 2005, Holloran et al. 2010); and increased chick mortality or brood loss (Aldridge andBoyce 2007, Kirol 2012).We, however, did not detect an association between oil and gas infrastructure variables and winter female survival.Yet, our occurrence model demonstrated that highly disturbed areas were primarily being avoided by sage-grouse.Consequently we expected that any potential survival consequences were unlikely to be realized and not detected or if the females were not using these areas we did not have the information needed to inform our survival model.In a companion study focused on the reproductive period we also did not find a negative relationship between adult female survival and anthropogenic variables at the landscape scale (Kirol 2012) and this finding was corroborated by another recent study also focused on survival of females in summer (Dinkins 2013).Avoidance response to energy infrastructure, manifested as decreased population densities, has been detected in sagebrush obligate songbirds including the Brewer's sparrow (Spizella breweri ), sage sparrow (Amphispiza belli ), and vesper sparrow (Pooecetes gramineus; Gilbert v www.esajournals.organd Chalfoun 2011).Thus, these authors suggested that increased energy development intensity in sagebrush steppe habitats has the potential to further exacerbate regional declines in sagebrush songbirds (Gilbert and Chalfoun 2011).
Current development stipulations intended to minimize impacts on sage-grouse are often implemented at local scales (i.e., surface use restrictions focused on leks) and have been shown to afford limited protection to sage-grouse (Holloran 2005, Walker et al. 2007, Naugle et al. 2011).Sage-grouse are a landscape-scale species (Knick and Connelly 2011) dependent on a variety of habitat characteristics as well as habitat connectivity throughout its life cycle (Aldridge and Boyce 2007, Doherty et al. 2008, Doherty et al. 2010, Aldridge et al. 2012, Fedy et al. 2012, Kirol 2012).Sage-grouse continue to decline throughout its range (Connelly andBraun 1997, Connelly et al. 2004); therefore, methods to effectively reduce impacts of energy development on sage-grouse and conserve critical habitat areas as development proceeds are of critical need.We found that sage-grouse avoided areas of greater surface disturbance that our predicted risk map suggested were low risk areas.Our results suggest that survival of sage-grouse in our study area may be independent of anthropogenic features, but avoidance of energy development may have resulted in loss of otherwise effective habitat.For example, if we make the assumption that the habitat conditions (i.e., environmental predictor variables) would be similar without oil and gas development, our occurrence model predicted an approximate 12.8% increase in high occurrence habitat when we omitted surface disturbance from the model.While our research illustrates the consequences of development on wintering sage-grouse, it also provides a method in which critical wintering habitats can be identified and conserved during development.Thus, our winter habitat quality map provided critical information for conservation planning of the limited amount of effective winter habitat.This habitat represented areas that sage-grouse concentrate in and provide the greatest prospects for survival through the winter.However, due to the relatively short duration of this study, we caution that yearly variability in winter weather patterns (i.e., deep snowpack, limited access to food and cover) may drastically alter sage-grouse occurrence and fitness.For example, deep snows may influence occurrence by altering the availability of shrubs (Dzialak et al. 2013), and female survival may be lower in winters with severe weather events (Moynahan et al. 2006).
During winter, Beck (1977) reported sagegrouse in northern Colorado restricted nearly 80% of habitat use to 7 areas comprising less than 7% of the total study area.Further Hupp and Braun (1989) reported that during a severe winter less than 10% of the sagebrush vegetation in the Gunnison Basin of southwestern Colorado was available to Gunnison sage-grouse (C.minimus).Therefore, areas that provide habitat features needed by sage-grouse in winter, that are conducive to survival, seem to be limited in many landscapes used by sage-grouse.Furthermore, individuals that undergo large movements to wintering areas are likely exposed to greater risk during these movements, which is especially true for juvenile birds (Beck et al. 2006).Because sage-grouse have high fidelity to wintering areas (Connelly et al. 2004), population persistence could be negatively compromised if migrating females arrive at wintering habitats that have been compromised by human activities.As our winter habitat quality map demonstrates, the habitat areas contributing the most to population persistence in winter (highly effective habitat) were often centered within habitats that our survival model identified as riskier (highly ineffective habitat).Consequently, if these highly effective habitats experience energy development, avoidance response could result in females being forced into surrounding habitats.The consequences of avoidance could result in use of riskier habitats creating a ''perceptual trap'' scenario where high-quality wintering habitats are avoided because they become less attractive (Patten and Kelly 2010).Our findings suggest that sage-grouse may be avoiding potentially effective habitats resulting from energy development and identifies undisturbed high quality habitat areas that should receive conservation priorities as development continues.

Fig. 3 .
Fig. 3. Predicted relative risk of weekly sage-grouse survival during winter.Map displays areas of high and low risk in northwest Colorado and south-central Wyoming, USA, winters 2007-2008, 2008-2009, and 2009-2010.Relative risk binned into four quantiles with a low occurrence bin (Fig.2) where minimal use prevented us from predicting risk.High risk values indicate that female sage-grouse were less likely to survive the winter in these areas.The white areas within the predictive map are areas with no data values in theHomer et al. (2012) products because they were masked as non-rangeland habitats.

Table 2 .
Top and competitive (w i !10% of top model w i ) models best explaining sage-grouse winter habitat selection in northwest Colorado and south-central Wyoming, USA Number of parameters (K ), change in Akaike's Information Criterion score from the top model (DAIC), and Akaike weights (w i ).