Introduction

As of the turn of the twenty-first century, sage-grouse occupied roughly half of their former historical range (Schroeder and others, 2004; Miller and others, 2011), and populations have subsequently experienced marked declines in many parts of their current range (Garton and others, 2011; Western Association of Fish and Wildlife Agencies, 2015). A 2021 study led by the U.S. Geological Survey (USGS), in collaboration with the Bureau of Land Management (BLM) and western state wildlife agencies, revealed an approximately 3.1-percent annual average population decline range-wide (Coates and others, 2021). However, variation in trends was described across different spatial and temporal scales.

Decades of literature have attributed sage-grouse population declines to loss and fragmentation of sagebrush communities as well as to a suite of environmental stressors (Connelly and others, 2004; Schroeder and others, 2004; Doherty and others, 2016). Since 1999, sage-grouse have been petitioned for legal protection under the Endangered Species Act (ESA) on nine occasions, and actions to conserve and restore sage-grouse habitats are now central to guiding land-management actions and policies across most of the western United States. Specifically, in recent years, the resource needs of sage-grouse have been used to help guide management actions aimed at improving conditions in sagebrush ecosystems, with resultant practices thought to benefit other sagebrush-dependent species (Rowland and others, 2006; Hanser and Knick, 2011; Dinkins and others, 2021). Sage-grouse are considered an indicator for the function of sagebrush ecosystems and an umbrella for the protection of other sagebrush-obligate or sagebrush-dependent species because of their near complete reliance on sagebrush ecosystems (Rich and Altman, 2001; Rich and others, 2005; Rowland and others, 2006; Hanser and Knick, 2011). However, recent (2018) literature indicates that certain less associated species might not be well covered by the sage-grouse umbrella (Carlisle and others, 2018). Importantly, several federal resource management plan amendments accompanying the U.S. Fish and Wildlife Service (USFWS) “not warranted” 2015 ESA listing determination called for greater integration of sage-grouse management into their land-use planning, specifically identifying how to implement adaptive management. An unprecedented level of conservation effort and planning among federal (Bureau of Land Management, 2015; U.S. Department of Agriculture Forest Service, 2015a, b), state, and private stakeholders was identified as the primary driver for the USFWS decision in the most recent (2015) status assessment (U.S. Fish and Wildlife Service, 2015).

The purpose of this report is to provide updated results on sage-grouse population trends across their geographical range in the western United States. This report reflects previous modeling methodology (Coates and others, 2021, 2022a) and includes additional data to inform population trends through 2023. A detailed description of data collection and compilation, population clustering methods to identify spatial extents, and trend modeling methodologies was provided in the first three objectives described in Coates and others (2021). Additional methodologies developed since Coates and others (2021) are described in Coates and others (2022a). The USGS, in cooperation with the BLM, are providing this scientific information to fulfill a prominent information gap, thus informing the adaptive management of sage-grouse population trends, as well as monitoring and conservation management strategies.

Study Area

Our study extent consisted of the geographic range of sage-grouse in the United States and was described previously in Coates and others (2021). Briefly, this area represents the sagebrush biome occurring across western North America and extending east from the Sierra Nevada and Cascade Mountain ranges to the western regions of the Great Plains of the United States. The vegetation communities vary in precipitation, temperature, soils, topographic position, and elevation (Miller and others, 2011). The most abundant shrub species include sagebrush (Artemisia spp.), with less abundant non-sagebrush species, such as rabbitbrush (Chrysothamnus spp.), horsebrush (Tetradymia spp.), greasewood (Sarcobatus spp.), common snowberry (Symphoricarpos spp.), serviceberry (Amelanchier spp.), fourwing saltbush (Atriplex spp.), and bitterbrush (Purshia spp.). The primary herbaceous species include wheatgrass (Agropyron spp.), fescue (Festuca spp.), bluegrass (Poa spp.), needlegrass (Stipa spp.), and squirreltail (Sitanion spp.), whereas less abundant forb species include phlox (Phlox spp.), milk-vetch (Astragalus spp.), and fleabane (Erigeron spp.).

Data Compilation and Inputs

All digitized field observations of sage-grouse lek counts since 1953 were compiled from all state wildlife agencies that monitored sage-grouse populations and entered into a single unified database. We worked with each state agency to ensure the fullest understanding of the data to maximize the number of appropriate records kept in the database. We addressed spatial errors, and we reviewed all data products with wester state wildlife agency biologists. Data compilation rules and the detailed methodology used in this analysis was published in Coates and others (2021) and O’Donnell and others (2021). The additional years of data (up to 2023) used to update this trend analysis followed the exact quality assessment and quality control (QA/QC) measures described in the previous publications. Furthermore, state-specific summaries of number and percentage (relative to maximum value in any given year) of leks and observations retained followed the sequential application of Rules 1–6 described in appendix 2 of Coates and others (2021). Using all the rules for selecting data appropriate for population modeling, we retained 108,759 observations across 5,413 leks. The sage-grouse lek data used in this update either are not publicly available or have limited availability owing to unique restrictions held by each state (data are managed by 11 western states and are not public due to the sensitivity of the species and state regulations, policies, or laws). Contact individual western state wildlife agencies (see the “Acknowledgments” section) for more information.

Changes in population abundance are affected by environmental factors that operate on multiple spatial and temporal scales, which follow ecological rather than geopolitical boundaries. Hence, we examined population trends across biologically relevant and hierarchically nested units to improve the detection of factors driving change across multiple spatial scales. We grouped sage-grouse lek sites into hierarchically nested clusters or populations using least-cost minimum spanning trees, a clustering algorithm, and a suite of relevant biotic and abiotic spatial products described in O’Donnell and others (2022). Briefly, we selected two cluster levels to represent a fine (neighborhood cluster, NC) and a broad spatial scale (climate cluster, CC) in trend analyses. Model output included estimated abundance ($\widehat{N}$) and intrinsic rates of population change ($\widehat{r}$) at the lek level.

Range-wide Sage-Grouse Population Model

Detailed formulation of the sage-grouse population model was described in Coates and others (2021). Coates and others (2022a) describe analytical updates that (1) identify population nadirs (lowest points within cycles) at spatial scales of the lek (breeding ground) and neighborhood cluster (group of leks) and (2) truncate prior distributions on rate of change in apparent abundance values to more realistic boundaries for leks with missing data. Briefly, we used a state-space model (SSM) that relied on Markov chain Monte Carlo sampling to derive posterior probability distributions (PD) of $\widehat{N}$ and $\widehat{r}$ using lek count data across the sage-grouse geographical range. This approach allowed inferences for each lek, as well as higher-order, nested spatial extents, such as NC and CC, during each year of the time series. Advantages and assumptions inherent to the SSM approach were described in Coates and others (2021, 2022a).

Range-Wide Population Trends

Our model fit the observed data well (Bayesian p-value=0.50). Median $\widehat{N}$ from the SSM revealed seven distinct range-wide population nadirs across the 64 years of data, which were 1966, 1975, 1986, 1996, 2002, 2013, and 2021. The range-wide, across-years, mean male count was 17.4 per lek (95-percent confidence interval=17.3–17.6). The number of years for complete oscillation periods (nadir-to-nadir) was relatively consistent across periods (average=9.2; 95-percent credible interval [CRI]=4.9–15.0; table 1). Model estimates revealed evidence of range-wide decline, on average, from every historic abundance nadir to 2021 (fig. 1; table 2); for example, the average annual $\widehat{\lambda }$ for Period 5 (19 years, two oscillations), Period 3 (35 years, four oscillations), and Period 1 (55 years, six oscillations) was 0.971 (median; 95-percent CRI=0.969–0.973), 0.970 (median; 95-percent CRI=0.968–0.972), and 0.972 (median; 95-percent CRI=0.969–0.974), respectively. These trends imply declines of 41.1, 64.5, and 78.4 percent, relative to population sizes observed 19, 35, and 55 years earlier, respectively. In an earlier version of the model (Coates and others, 2021), we specified the final nadir for all populations using the final year of the dataset (2019), which was the lowest point of $\widehat{N}$ for most populations at that time. A more recent (2023) report provided an update of the original analysis and included years 2020–22 (Coates and others, 2023). In that report, researchers concluded that not all populations had reached a nadir in 2019. The additional year of count data (2023) provided in this report supports previous conclusions of 2021 being the most recent nadir for populations at the range-wide scale. However, some populations (CC-B; Washington area) continue to decline, making 2023 a new tentative nadir at finer scales. Despite these updates, trends estimated at CC (fig. 2) and NC (fig. 3) scales are generally consistent with the previous report (Coates and others, 2023). We estimated median $\widehat{\lambda }$ to be less than 1.0 for 87.1, 93.0, and 97.9 percent of NCs during Period 5, Period 3, and Period 1, respectively, throughout the sage-grouse range.

Climate Cluster Population Trends

Climate cluster A (CC-A; Bi-State area) consisted of 11 NCs that encompassed 726,907 hectares (ha). Two NCs did not have sufficient lek data to estimate trends. Climate cluster A consisted of 95 leks, representing approximately 1 percent of the range-wide database. After QA/QC, 58 leks met criteria for use in the SSM (table 2), totaling 1,720 field observations. Mean male count was 20.1 (95-percent confidence interval=18.8–21.5). For CC-A, we estimated six population abundance nadirs that dated back to 1960 and included nadirs of 1969, 1978, 1983, 1995, 2002, 2008, and 2019 (table 1). We estimated $\widehat{\lambda }$ for Period 5 (2002–19, two oscillations, 17 years), Period 3 (1983–2019, four oscillations, 36 years), and Period 1 (1969–2019, six oscillations, 50 years) as 0.976 (95-percent CRI=0.964–0.986), 0.990 (95-percent CRI=0.977–0.999), and 0.979 (95-percent CRI=0.969–0.987), respectively (table 2). During the past 17, 36, and 50 years, sage-grouse populations have experienced declines in abundance equal to 32.2, 29.7, and 64.7 percent, respectively (fig. 4A). We estimated median $\widehat{\lambda }$ to be less than 1.0 for 88.9, 100.0, and 100.0 percent of NCs during Period 5, Period 3, and Period 1, respectively.

Climate cluster B (CC-B; Washington area) consisted of four NCs that encompassed 1,139,954 ha. One NC did not have sufficient lek data to estimate trends. Climate cluster B consisted of 54 leks, representing approximately 0.6 percent of the lek database. After QA/QC, 43 leks met criteria for use in the state-space trend model (table 2), totaling 1,205 field observations. Mean male count was 14.0 (95-percent confidence interval=13.2–14.8). We estimated six population abundance nadirs that dated back to 1960 and included 1964, 1976, 1986, 1994, 2008, 2017, and 2023 (table 1). We estimated $\widehat{\lambda }$ for Period 5 (2008–23, two oscillations, 16 years), Period 3 (1986–2023, four oscillations, 37 years), and Period 1 (1964–2023, six oscillations, 59 years) as 0.957 (95-percent CRI=0.945–0.969), 0.976 (95-percent CRI=0.966–0.983), and 0.972 (95-percent CRI=0.963–0.981), respectively (table 2). During the past 16, 37, and 59 years, sage-grouse populations have experienced declines in abundance equal to 48.3, 58.3, and 80.7 percent, respectively (fig. 4B). We estimated median $\widehat{\lambda }$ to be less than 1.0 for all NCs during Period 5, Period 3, and Period 1, respectively.

Climate cluster C (CC-C; Jackson Hole, Wyoming area) consisted of two NCs that encompassed 66,733 ha. Climate cluster C consisted of 17 leks, representing approximately 0.2 percent of the lek database. After QA/QC, 14 leks met criteria for use in the SSM (table 2), totaling 340 field observations. Mean male count was 13.7 (95-percent confidence interval=12.1–15.4). For CC-C, we estimated population abundance nadirs during 1963, 1969, 1984, 1999, 2003, 2011, and 2019 (table 1). We estimated $\widehat{\lambda }$ for Period 5 (2003–19, two oscillations, 16 years), Period 3 (1984–2019, four oscillations, 35 years), and Period 1 (1963–2019, six oscillations, 56 years) as 0.964 (95-percent CRI=0.944–0.981), 0.969 (95-percent CRI=0.942–0.993), and 0.966 (95-percent CRI=0.950–0.985), respectively (table 2). During the past 16, 35, and 56 years, sage-grouse populations have experienced declines in abundance equal to 42.3, 65.7, and 85.1 percent, respectively (fig. 4C). We estimated median $\widehat{\lambda }$ to be less than 1.0 for all NCs across this temporal scale.

Climate cluster D (CC-D; Eastern area) consisted of 169 NCs that encompassed 25,920,530 ha. There were 24 NCs that did not have sufficient lek data for trend estimates. Climate cluster D consisted of 3,069 leks, representing approximately 35.3 percent of the lek database. After QA/QC, 1,910 leks met criteria for use in the SSM (table 2) and totaled 37,491 field observations. Mean male count was 16.2 (95-percent confidence interval=16.0–16.4). For CC-D, we estimated six population abundance nadirs that dated back to 1960 and included 1966, 1981, 1986, 1997, 2004, 2013, and 2019 (table 1). We estimated $\widehat{\lambda }$ for Period 5 (2004–19, two oscillations, 13 years), Period 3 (1986–2019, four oscillations, 31 years), and Period 1 (1966–2019, six oscillations, 51 years) as 0.961 (95-percent CRI=0.957–0.964), 0.969 (95-percent CRI=0.964–0.972), and 0.967 (95-percent CRI=0.962–0.971), respectively (table 2). During the past 13, 31, and 51 years, sage-grouse populations have experienced declines in abundance equal to 38.0, 61.1, and 81.3 percent, respectively (fig. 4D). We estimated median $\widehat{\lambda }$ to be less than 1.0 for 91.4, 96.4, and 99.3 percent of NCs during Period 5, Period 3, and Period 1, respectively.

Climate cluster E (CC-E; Great Basin area) consisted of 241 NCs that encompassed 34,627,182 ha. There were 17 NCs that lacked sufficient data to estimate trends. Climate cluster E consisted of 4,189 leks, representing approximately 48.2 percent of the lek database. After QA/QC, 2,377 leks met criteria for use in the SSM (table 2) and totaled 45,317 field observations. Mean male count was 16.0 (95-percent confidence interval=15.8–16.2). For CC-E, we estimated six population abundance nadirs that dated back to 1960 and included 1967, 1975, 1985, 1996, 2002, 2013, and 2021 (table 1). We estimated $\widehat{\lambda }$ for Period 5 (2002–21, two oscillations, 19 years), Period 3 (1985–2021, four oscillations, 36 years), and Period 1 (1967–2021, six oscillations, 54 years) as 0.967 (95-percent CRI=0.964–0.969), 0.974 (95-percent CRI=0.970–0.977), and 0.972 (95-percent CRI=0.969–0.974), respectively (table 2). During the past 19, 36, and 54 years, sage-grouse populations have experienced declines in abundance equal to 45.3, 60.2, and 77.8 percent, respectively (fig. 4E). We estimated median $\widehat{\lambda }$ to be less than 1.0 for 83.3, 90.5, and 96.8 percent of NCs during Period 5, Period 3, and Period 1, respectively.

Climate cluster F (CC-F; western Wyoming area) consisted of 56 NCs that encompassed 8,899,755 ha. Three NCs lacked sufficient data to estimate trends. Climate cluster F consisted of 1,275 leks, representing approximately 14.7 percent of the lek database. After QA/QC, 1,011 leks met criteria for inclusion in the SSM (table 2) and totaled 22,686 field observations. Mean male count was 22.5 (95-percent confidence interval=22.1–22.8). For CC-F, we estimated six population abundance nadirs that dated back to 1960 and included 1965, 1975, 1987, 1996, 2002, 2013, and 2021 (table 1). We estimated $\widehat{\lambda }$ for Period 5 (2002–21, two oscillations, 19 years), Period 3 (1987–2021, four oscillations, 34 years), and Period 1 (1965–2021, six oscillations, 56 years) as 0.976 (95-percent CRI=0.972–0.981), 0.971 (95-percent CRI=0.967–0.975), and 0.976 (95-percent CRI=0.972–0.980), respectively (table 2). During the past 19, 34, and 56 years, sage-grouse populations have experienced declines in abundance equal to 35.4, 62.1, and 73.7 percent, respectively (fig. 4F). We estimated median $\widehat{\lambda }$ to be less than 1.0 for 90.6, 92.5, and 98.1 percent of NCs during Period 5, Period 3, and Period 1, respectively.

Neighborhood cluster trends for different temporal scales that are not listed in this report can be viewed in Coates and others (2022b). Coates and others (2022b) is the accompanying data release to this report.

Probability of Future Extirpation

Expanding on the population trend analysis, we projected $\widehat{N}$ for each lek and NC across three temporal scales that reflected two, four, and six oscillations using the same model and dataset. We used the mean duration of one complete oscillation (9.2 years) based on estimated $\widehat{N}$ from SSM results. We then calculated the proportion of the posterior probability distribution of $\widehat{N}$ that was less than two males (minimum number to represent a lek) for the last prediction year of each temporal scale. Although this value is not true extirpation (zero birds), we refer to it as extirpation to align with state definitions of lek inactivity. Extirpation of leks within an NC was thought to reflect a loss of a meta-population resulting from reduced demographic rates.

The mean extirpation probability calculated across all leks was 0.47 (SD=0.29), 0.56 (SD=0.23), and 0.62 (SD=0.19) for two (fig. 5A), four (fig. 5B), and six (fig. 5C) oscillations, respectively. However, distributions of predicted extirpation rates for leks appeared bimodal with peak probabilities near 0.12 and 0.81 (two oscillations, fig. 5A), 0.35 and 0.82 (four oscillations, fig. 5B), and 0.48 and 0.84 (six oscillations, fig. 5C). Approximately 47, 58, and 70 percent of leks possessed an extirpation probability greater than 0.5 for two, four, and six oscillations, respectively. The spatial distribution of lek extirpation probabilities demonstrated higher concentrations of higher probabilities at the edges of the species range (fig. 6A–C). The mean extirpation probability calculated across all NCs was 0.15 (SD=0.25), 0.22 (SD=0.27), and 0.26 (SD=0.29) for two (fig. 5D), four (fig. 5E), and six periods (fig. 5F) of future oscillations, respectively. Approximately 13, 19, and 24 percent of NCs possessed an extirpation probability greater than 0.5 for two, four, and six oscillations, respectively. The highest probabilities of NC scale extirpation often occurred at the periphery of the species range (fig. 6D–F). Conversely, interior NCs consisted primarily of low extirpation probabilities ranging from 0 to 0.25.

Neighborhood cluster extirpation probabilities for different temporal scales that are not listed in this report can be viewed in Coates and others (2022b). Coates and others (2022b) is the accompanying data release to this report.

Watches and Warnings from a Targeted Annual Warning System

The Targeted Annual Warning System (TAWS) is a hierarchical monitoring strategy that contrasts estimates of $\widehat{r}$ across nested spatial scales on an annual basis. Comparisons can act as a powerful analytical tool to help target when and where to adjust monitoring efforts, or where to carry out management actions. Methodology of TAWS was described in detail in Coates and others (2021) and Prochazka and others (2023). TAWS produces signals referred to as watches and warnings, which signify progressively greater degrees of evidence for aberrant decline. Evidence of aberrant decline is assessed using independent sets of standardized thresholds that seek to stabilize the range-wide population (range-wide thresholds) versus individual climate clusters (climate cluster thresholds). The original published report describing TAWS results by Coates and others (2021) described each set of thresholds to provide managers with multiple strategic options. Here, we chose to report TAWS results using CC thresholds based on feedback of implementation from state and federal agency personnel. The primary reason for preferential use of CC thresholds is the superior performance in targeting peripheral and sparsely distributed populations. Like the trend analysis, thousands of historic lek surveys underwent QA/QC, as described in the results of Objective 1 for the TAWS analysis (Coates and others, 2021).

During 1990–2023, we estimated 0.725 and 0.578 of 4,808 sage-grouse leks experienced at least one watch or warning (table 3), respectively, across the species’ range. We calculated a mean annual proportion of leks that experienced a first watch or warning (no watches or warnings during any preceding year) to be 0.022 (105 leks) and 0.018 (84 leks), respectively. The mean annual proportion of leks that experienced a repeat watch or warning (at least one watch or warning occurred during a preceding year) was 0.060 (291 leks) and 0.060 (291), respectively. The CC with greatest proportion of activated watches at leks across the 34 years was CC-A (Bi-State area), where 0.941 proportion of leks activated one or more times (table 3). Conversely, CC-E (Great Basin area) consisted of the greatest number of watches compared to other clusters (number of first watches=1,461; number of repeat watches=3,999; table 3). The CC with the least proportion of watches was CC-C (Jackson Hole, Wyoming, area), at 0.643 (number of first and repeat watches=9 and 43, respectively; table 3). The CC with the greatest proportion (0.726) of warnings was CC-F (western Wyoming area), whereas CC-E had the greatest number (first=1,164 and repeat=4,080). The second highest proportion of watches (0.901) was CC-F, where we estimated 840 (repeat=3,027). The second highest proportion of warnings (0.706) was CC-A, where we estimated 36 (repeat=192).

During 1990–2023, we estimated 0.606 and 0.448 of 424 NCs experienced at least one watch or warning (table 3), respectively, across the species’ range. An average of 0.018 (repeat=0.051) and 0.014 (repeat=0.041) proportion of clusters activated per year, which was approximately 7.8 (repeat=21.8) and 5.8 (repeat=17.5) clusters. We reported CC-A (Bi-State area) and CC-B (Washington area) had the greatest proportion (1.000) of watches, whereas CC-E (Great Basin area) had the greatest number (first=176 and repeat=539) of watches across the 34-year timeframe. For warnings, CC-A had the greatest proportion (0.889; table 3), and CC-E had the greatest number (first=122 and repeat=400).

During 2023, we estimated 0.028 and 0.025 proportions of leks experienced a first watch and warning, respectively, range-wide (table 4), which resulted in 133 (repeat=464; fig. 7) and 122 (repeat=576; fig. 7) lek activations. During 2023, the greatest proportion of first watches (0.057) and warnings (0.033) were within CC-D (Eastern area), which were 94 (repeat=172) watches and 55 (repeat=159) warnings (table 4). During 2023, we estimated 0.019 and 0.031 proportions of neighborhoods experienced a first watch and warning, respectively, range-wide (table 4), which resulted in 8 (repeat=34; fig. 8) and 13 (repeat=47; fig. 8) neighborhood activations. Climate cluster A experienced the greatest proportion of NC warnings (0.111) during 2023 followed by CC-D (0.036) and CC-E (0.027). Warnings and watches for each neighborhood cluster are available in Coates and others (2022b).