Sagebrush Steppe Productivity, Environmental Complexity, and Grazing: Insights From Remote Sensing and Mixed-effect Modeling

ABSTRACT Domestic livestock grazing is the primary land use across the planet, but the relationship between grazing and rangeland productivity is difficult to determine because it is influenced by a variety of ecological and management factors. Fine-scale environmental data available through remote sensing are increasingly used to understand land use changes, such as grazing. In this study, we assessed the relationship between a variety of grazing and rangeland productivity metrics while accounting for environmental complexity within the sagebrush steppe ecosystem of Montana. We created mixed-effect generalized linear models using remotely sensed productivity as response variables. Explanatory variables included management and field-based grazing data combined with remotely sensed abiotic and biotic environmental factors. We found point-level field measures of grazing (e.g., cow patties, percentage of dung in Daubenmire plots, and number of plants grazed) showed positive effects, especially on perennial forbs and grasses. Grazing measures at the pasture-level showed a small negative effect on annual forbs and grasses. Grazing metrics tended to have smaller covariate effects on rangeland productivity compared to environmental factors, and interaction effects between grazing and environmental factors were common. This study provides insight into the relationship between grazing and plant productivity in the sagebrush steppe rangeland of Montana and highlights the importance of assessing the effects of grazing using multiple scales while accounting for environmental complexity.


Introduction
Domestic livestock grazing (hereafter, "grazing") is the most prevalent land use globally: 22% to 26% of all ice-free land surfaces are rangelands used for grazing ( Phelps and Kaplan, 2017 ).Rangelands are vast, inherently spatially and temporally variable areas ( Bastin et al., 2012 ) that provide vital ecosystem services, including ✩ Funding has been provided by Federal Aid in Wildlife Restoration Grants W-158-R and W-181-R to Montana Fish Wildlife and Parks, US Bureau of Land Management Agreement #L15AC0 0 097 and #L20AC0 0508 , US Geological Survey Montana Cooperative Wildlife Research Unit at the University of Montana , Safari Club International Foundation , USDA Natural Resources Conservation Service (NRCS) Conservation Innovation Grant Agreement #69-3A75-10-151 , NRCS and US Fish and Wildlife Service Interagency Agreement # 60181BJ653 facilitated by the Intermountain West Joint Venture and Pheasants Forever, USFWS Cooperative Agreement F14AC01224 , NRCS Conservation Effects Assessment Project , MFWP Upland Gamebird Enhancement Program , and Big Sky Upland Bird Association .Teague and Barnes, 2017 ;Hewins et al., 2018 ).A consensus is that long-term, high-intensity grazing (i.e., overgrazing) degrades rangeland quality over time, and more moderate grazing levels can maintain or improve rangelands ( Milchunas and Lauenroth, 1993 ;Fuhlendorf and Engle, 2001 ;Bates and Davies, 2014 ;Davies et al., 2014Davies et al., , 2018 ) ).
Many environmental and management factors affect the outcome of grazing and operate within the context of the evolutionary history of grazing across a rangeland ( Heady, 1974 ;Mack and Thompson, 1982 ;Milchunas and Lauenroth, 1993 ).These factors include moisture availability ( Milchunas et al., 1988 ), vegetation composition and functional traits ( Adler et al. 20 01 ,20 05 ), pasture size ( Barnes et al., 2008 ;Allred et al., 2011 ), livestock characteristics (e.g., breed, age, sex, herding behavior;McNaughton, 1984 ), and livestock breeding status ( Connelly et al., 2004 ).The influence of these factors and their interactions depends on sitespecific processes and spatiotemporal scale ( Milchunas and Lauenroth, 1993 ;Adler et al., 2001 ).As a result, every livestock operation is unique and operates within specific constraints to achieve production goals and ecological objectives ( Budd and Thorpe, 2009 ).While established principles and metrics often guide grazing management decisions (e.g., Butler et al., 2003 ;O'Brien et al., 2003 ), measuring the subsequent outcomes at scales relevant to ranchers and landscape processes can be challenging.
Remote sensing is a valuable tool for assessing landscape changes, including those induced by land use practices ( Kerr and Ostrovsky, 2003 ;Lawley et al., 2016 ).Since rangelands are large, heterogeneous areas, assessing how grazing influences rangelands requires broad-scale studies that account for site-specific processes ( Milchunas and Lauenroth, 1993 ;Bastin et al., 2012 ;Lawley et al., 2016 ).Collecting observations of all grazing activities, rangeland response, and abiotic complexity across large scales requires extensive field data and is logistically challenging.However, remote sensing tools can facilitate broad-scale grazing studies that would otherwise be unfeasible by offering an increasingly wide range of data across large spatial areas.Recent advances provide more fine spatial resolutions for identifying site-specific differences in vegetation response ( Seto et al., 2004 ;Rafique et al., 2016 ).Specifically, remote sensing data and techniques can capture vegetation changes, including metrics of productivity, across rangelands at scales that translate to pasture-or landscape-level management ( Jansen et al., 2021 ).
Here, we integrate field data with remote sensing variables to determine rangeland productivity response to grazing by isolating its influence from the effects of environmental factors.For our model system, we chose rangelands within the sagebrush steppe, an ecosystem of conservation concern where grazing is the predominant land use ( Knick et al., 2003 ).Various productivity metrics were used to measure rangeland response.Productivity is an ecologically important component of rangeland systems ( Linderman et al., 2005 ;Pettorelli et al., 2005 ;Watson et al., 2014 ;Jung et al., 2019 ) that can be manipulated by grazing ( Adler et al., 2001 ;Ferraro and Oesterheld, 2002 ;Davies et al., 2016 ).Previous empirical evidence shows grazing is capable of increasing productivity ( Dyer et al., 1986 ;Augustine et al., 2003 ;Davis et al., 2015 ).We explore if grazing will stimulate plant growth in the sagebrush steppe of Montana at the spatial and temporal scales of this study.We base our investigation on the grazing optimization hypothesis that states grazing activity increases vegetation productivity until it reaches a maximum at moderate grazing levels, then is reduced at the highest levels of grazing use ( McNaughton, 1976 ;Hilbert et al., 1981 ;Dyer et al., 1986 ).As grazing intensity and duration increase, we predict a corresponding increase in remotely sensed rangeland productivity up to a point where grazing pressure will overwhelm the ability of the vegetation to regrow, and the relationship will become negative.We collected grazing data in the field and pasture management data from ranchers and managers.Remotely sensed data was used for environmental variables and rangeland response.We used Bayesian generalized linear mixedeffect models (GLMMs) to assess the relationships with productivity as the response and all other factors as the predictors.The results of this study will improve the ecological understanding of the relationship between grazing and sagebrush steppe rangelands.

Study area
Our study area was located in the Northern Great Plains near Roundup, Montana (46.4488 • N, 108.5438 • W, Fig. 1 ).It was approximately 2,160 km 2 of rolling hills dominated by sagebrush shrubs ( Artemisia tridentata ssp.Wyomingensis ) and perennial grasses such as Agropyron smithii , Stipa comata , Poa sandbergi , and Bouteloua gracilis interspersed with conifer and riparian areas.The study area elevation ranged between 938 and 1425 m ( USGS, 2022 ).The climate is characterized as cold and semi-arid with an average monthly temperature and precipitation varying from -3.8 • C to 21.8 • C and 0.94 cm to 7.3 cm, respectively, in Roundup from 1991 to 2020 ( NOAA, 2022 ).

Data collection
We grouped the factors we considered most influential to rangeland productivity in the sagebrush steppe into three categories: grazing, absorbed photosynthetically active radiation (APAR), and light use efficiency (LUE; Fig. 2 ).Grazing was represented in models using field and management data measuring the main components of grazing regimes (i.e., intensity, duration, and timing).APAR represents the amount of available light energy absorbed by a plant, and LUE is a measure of how efficiently plants convert absorbed light into biomass.Both APAR and LUE are environmental variables expressed by remotely sensed metrics.APAR was incorporated using estimates of shortwave radiation and fraction absorbed photosynthetically active radiation (FPAR).LUE was accounted for using measures of temperature, moisture, and plant composition.A full description of the metrics used to represent each factor can be found in Table S1.

Grazing data
Grazing data were collected from private landowners from 2011 to 2020.This study is part of a broader program focused on assessing the efficacy of grazing regimes incentivized by the Sage Grouse Initiative as part of ongoing greater sage-grouse ( Centrocercus urophasianus ) conservation effort s.The grazing data were collected at two spatial scales: pasture-and point-level.Management data were collected at the pasture-level, and field data was collected at the point-level.
Pasture-level grazing data were collected annually and included 1) changes in individual pasture boundaries over the course of the study, 2) the total number of days livestock were grazed in the calendar year, and 3) grazing intensity measured via stocking rate using the number of animal units per month (AUMs; Figs.1-2 ).AUMs were calculated by multiplying the total number of livestock by the total number of months grazed (i.e., number of days divided by the average number of days in a month) and then dividing by the pasture area in square meters.We accounted for differences in livestock type by multiplying the total number of livestock by the appropriate animal use equivalent (e.g., 1 0 0 0 lb cow and her calf is equal to 1.0 animal use equivalent).Our data were compiled spatially by georeferencing 97 pasture boundary maps  and digitized pasture boundaries, which resulted in ∼523 digitized pastures per year (i.e., 5,234 total pastures over the study, mean pasture size = 2.9 km 2 , median pasture size = 2.1 km 2 , minimum = 0.0061 km 2 , maximum = 33.8km 2 , standard deviation = 3.3 km 2 ) and covered ∼1 515 km 2 of the study area.We then linked the digitized pasture polygons to their corresponding grazing information by year and grazing occurrence.
Point-level grazing field data were collected from April to August each year using 30 m circular plots located randomly or at sage-grouse nesting sites ( Fig. 1 ; Smith et al., 2018 ).Three different metrics of grazing intensity were collected at this spatial scale: average percentage dung, patty count, and utilization frequency.The average percentage of dung in a plot was determined by estimating the percent cover of dung in Daubenmire frames.We took the mean percent cover that cow patties occupied within twelve 20 × 60 cm quadrats placed three meters apart in each cardinal direction from the center of the plot ( Daubenmire, 1956 ;Smith et al., 2018 ).The total number of cow patties in the plot was counted and recorded along with the age of each patty (i.e., from previous years or current year).Utilization frequency was measured as the number of randomly inspected plants out of 100 that had visible signs of grazing within the 30 m plot.

Remotely sensed data
We extracted variables representing the rangeland productivity response and environmental factors from the Google Earth Engine catalog as raster data ( Fig. 2 , Table S1).Several measures of productivity in the sagebrush steppe were extracted to determine how grazing affects different aspects and measurements of productivity.Variables with a temporal scale finer than one year (e.g., daily) were summarized using the mean and standard deviation for the growing season (i.e., the last killing freeze in the spring to the first killing freeze in the fall [-4.44 °C]).We used gross primary production (GPP), annual net primary production (NPP; Robinson et al., 2018 ), and normalized difference vegetation index (NDVI; USGS, 2020 ) as measures of productivity.We also considered biomass as a metric of productivity represented by the annual biomass accumulation of annual forbs and grasses, perennial forbs and grasses, or both ( Allred et al., 2021 ).

Data cleaning
Data manipulation and analyses were completed using program R ( R Core Team, 2021 ).First, we stacked all pasture-level grazing and raster data (i.e., environmental factors and productivity).We then extracted the value of those variables at the center point of the point-level grazing data using the Raster package ( Hijmans et al., 2020 ).The resulting dataset included pasture-level grazing management, point-level field grazing metrics, and raster data for each location where point-level grazing data were collected.We visually assessed and calculated the coefficient of variation to determine dispersion of each variable.We removed extreme outliers judged to have resulted from erroneous data entry or remote sensing model predictions.We removed variables with a coefficient of variation below 10% or above 100% from consideration for models except for the grazing metrics, which we retained irrelevant of dispersion.Finally, we centered and scaled all remaining variables to allow for a direct comparison of their effects.

Analysis
We assessed relationships between each rangeland response variables and grazing and environmental metrics using Bayesian GLMMs.Bayesian GLMMs offer the advantage of including additional model complexity, such as adding components with different distributions or categorical variables.These models also allow estimates of covariate effects that can be easily interpreted.To create the final models for each combination of rangeland response and grazing metric, we 1) assembled the variables considered for each pairing, 2) filtered the variables based on univariate effect and collinearity, and 3) performed backward stepwise selection on the full model (i.e., all variables were incorporated in the initial model run).All models used only the subset of the data that had observations for all variables considered in each model (i.e., list-wise deletion).
Our first step for every rangeland response variable was to analyze univariate models for each environmental explanatory variable.We only considered those variables for future models if their univariate coefficient value was significant (i.e., not overlapping zero) for that rangeland response.Furthermore, we examined the variables with a significant univariate effect for collinearity using Pearson's R ( Cohen, 1960 ), correlation matrices, and variance inflation factors.We removed variables until none were highly cor-related (i.e., variance inflation factors < 5, Everitt and Skrondal, 2010 ;|r| > 0.7, Dormann et al., 2013 ).The retained variables were selected from among correlated variables based on their value to rangeland productivity as determined by the literature and the univariate coefficient values.
For the next step, we created a GLMM for each rangeland response (Table S1: "Response Variables") and individual grazing metric (Table S1: "Grazing Variables") pairing, including 1) a quadratic term for that grazing metric to represent the grazing optimization hypothesis, 2) non-correlated environmental factors with a significant univariate effect (Table S1: "Environmental Variables"), 3) interaction terms between the grazing, temperature, and moisture variables, and 4) a random effect for year to capture any other annual variability.The final model was chosen through backward stepwise selection.We performed mixed effect generalized linear regression using Jags ( Depaoli et al., 2016 ) through the jag-sUI package ( Kellner 2021 ).We used vague normal priors for the intercept and covariate effects and vague uniform priors for standard deviation of rangeland response and year effect.All models used 10 0 0 0 iterations, 5 0 0 0 burn-in, no thinning, and three chains.We performed backward stepwise selection by running a full model that included all possible covariates.We removed variables that were insignificant (i.e., covariate effect overlapping zero) one at a time, beginning with the variables with the smallest coefficient effect in the univariate models, until our model fit metrics indicated the best fit.We also used DIC, WAIC, and leave-one-out AIC from the loo package ( Vehtari et al., 2017 ) + year e f f ect( year i ) (1) We visualized our results by creating coefficient density plots using the Bayesplot package ( Gabry et al., 2019 ) and partial effects plots using the ggplot2 package ( Wickham, 2016 ).Coefficient density plots allowed straightforward comparison of coefficient effects by displaying the posterior distribution of the coefficients for each explanatory variable in the final models.Partial effect plots show the isolated effects of the grazing metrics by predicting the rangeland response values for the full range of the grazing variable and assigning the other variable values in the model to be fixed at their median values.

Results
We collected grazing data from 5 234 pastures (510-529 per year) and 3,588 survey points (202-566 per year) over the 9-year study period.The data subsets developed through list-wise deletion for the final models had 1 483 or 1 778 observations depending on whether the grazing variable in that model was a pasturescale or point-scale grazing metric, respectively.
The coefficient values of grazing produced by the final GLMMs for each unique pairing of productivity and grazing metric varied; patterns were dependent on whether the grazing metric was a pasture-or point-scale grazing metric.Overall, grazing explanatory variables did not have a higher coefficient value than the environmental variables in any final GLMM (Tables S2, S3, Fig. S1).We found the quadratic effect of grazing on rangeland response statistically insignificant, with an effect overlapping zero, for all rangeland response metrics.S2 and S3 for all statistically significant coefficient values in the final models and Table S1 for covariate abbreviations.S2 and S3 for all statistically significant coefficient values in the final models and Table S1 for covariate abbreviations.

Pasture-scale metrics
The grazing coefficient values tended to be slightly negative for all grazing metrics measured at the pasture-scale: days grazed and AUMs ( Fig. 5 ).AUMs had the strongest negative coefficient value on the productivity of annual forbs and grasses of all productivity and grazing metric pairings (-0.077 [95% CRI -0.11 --0.041]) for the interaction between AUMs and minimum temperature standard deviation), but did not have a statistically significant coefficient value for any other rangeland response variables ( Figs. 3 and  4 , Tables S2 and S3, Figs.S1 and S2).There was no significant effect of the number of days grazed on any rangeland response variable besides NDVI and NDVI SD.The number of days grazed had a significant negative interaction effect with temperature and moisture levels on only NDVI and NDVI standard deviation, which was the strongest covariate effect found for those productivity response metrics ( Fig. 5 , Tables S2 and S3, Fig. S1).The number of days grazed also had the highest positive effect of the pasture-level vari-ables on NDVI with an estimated coefficient value of 0.024 (95% CRI 0.0048-0.044).

Point-scale metrics
Overall, the grazing metrics measured at the survey points (i.e., percentage dung cover, total cow patties, new cow patties, old cow patties, and utilization frequency) had a positive coefficient effect on rangeland productivity ranging from -0.055 (95% CRI -0.076 to -0.035, interaction effect of new cow patties and moisture [VPD] on perennial plant biomass) to 0.13 (95% CRI 0.1 -0.16, effect of total cow patties on perennial plant net primary production, Fig. 5 , Tables S2 and S3).Percentage dung cover was a statistically significant positive predictor for several rangeland response metrics, including perennial plant net primary production and biomass, overall biomass, and the standard deviation of GPP ( Fig. 5 ).Percentage dung cover had a slightly negative coefficient value for GPP, and otherwise had a statistically insignificant effect ( Fig. 5 , Tables S2  Figure 5. Covariate effects produced from the final GLMMs for each grazing covariate and rangeland response biomass variable.Yellow circles indicate negative effects of the covariate on the response variable, and teal circles indicate positive effects of the covariate on the response variable.Gray circles indicate that the 95% credible intervals overlapped zero and the effect of the covariate was neutral.The size of the circles represents the effect size.Each covariate is grouped by grazing predictor.Interactions between the grazing predictor and the covariate are indicated by an 'x' preceding the covariate.

and S3
).The total number of cow patties was a significant predictor more frequently than all other grazing metrics, with a relatively strong positive coefficient value for almost all rangeland productivity metrics ( Fig. 5 , Tables S2 and S3).The exceptions were GPP standard deviation and annual grass NPP and biomass, with a negligible or slightly negative coefficient for its interaction effects ( Fig. 5 , Tables S2 and S3).The other grazing metrics for cow patties (i.e., old and new cow patty count) showed similar, generally positive coefficient value patterns ( Fig. 5 , Tables S2 and S3).However, old cow patties did not interact with environmental effects ( Fig. 5 , Tables S2 and S3).The utilization frequency was less often statistically significant and tended to be positive when significant ( Fig. 5 , Tables S2 and S3).

Discussion
Our prediction that rangeland growth is positively correlated with grazing activity was partially supported.After accounting for environmental factors, we saw a relatively small positive influence of point-level grazing metrics (i.e., cow patty counts, utilization frequency, and percentage dung in the Daubenmire sampling plots) on rangeland productivity at lower grazing intensities.However, the slight negative effect of grazing for pasture-level metrics and the lack of evidence of a quadratic term does not support that vegetation productivity would increase to a certain threshold of grazing intensity, beyond which it would decline (i.e., the grazing optimization hypothesis).This may be because sampling was opportunistic and landowners in this area already practice moderate grazing and adaptive management based on pasture conditions ( Budd and Thorpe, 2009 ;Davies et al., 2014 ).Thus, existing grazing practices were potentially not intense enough in our study area to reach the inflection point in the grazing optimization curve.More variability in grazing regimes, especially a complete lack of grazing or higher use of the available forage, may show a more pronounced or quadratic short-term effect of grazing.
Our small grazing effect sizes and conflicting results between the pasture-and point-level grazing metrics also support the idea that rangelands respond to grazing differently at varying temporal and spatial scales ( Adler et al., 2001 ).Grazing metrics in our study were limited to annual temporal resolutions, which may have reduced our ability to detect grazing effects.Measurements of grazing intensity should be taken during time periods when grazing is influencing productivity to make more accurate assessments of the relationship between grazing and productivity (i.e., match remote sensing temporal resolutions with grazing occurrence; Fig. 7 ).For example, the field grazing data were collected April through August, but if cattle were rotated into the pasture after the field data were collected, they could be directly reducing the remotely sensed productivity we detected across the whole growing season.The temporal mismatch in grazing field data collection would not reflect this grazing pressure or the resulting reduced productivity.Furthermore, this study focused on the short-term effects of grazing, while some effects on rangelands may not be detectable for decades or over a century ( Fuhlendorf and Smeins, 1997 ).Similarly, the differences in results between the point-and pasture-level grazing metrics is likely due to some combination of three spatial scale considerations: 1) the point-level analysis of this study did not capture the effects of grazing at the pasture-level, 2) the pasture-level grazing influences rangeland productivity differently because grazing does not usually occur uniformly across a pasture ( Chapman et al., 2007 ;Raynor et al., 2021 ) or 3) metrics between point-and pasture-level are measured differently and are capturing a different effect of grazing.
Although the effects of grazing metrics varied based on spatial scale, every final GLMM showed environmental factors were more influential than grazing based on their covariate effects in our study area.Specifically, measures of moisture (e.g., VPD), temperature (e.g., minimum temperature), plant composition (e.g., percentage of annual forbs and grasses in the pixel), and FPAR often had a more significant effect on rangeland productivity than grazing.This supports previous work that suggests the influence of environmental factors to be more influential to rangelands than grazing ( Milchunas and Lauenroth, 1993 ) and for organisms in the sagebrush steppe (e.g., Smith et al., 2018 ;Reintsma et al., 2022 ).The covariate effects of the interactions between temperature, moisture, and grazing were also often significant (e.g., increased positive effects of cow patties as minimum temperature variability increases, decreased negative effects of AUMs for higher minimum  S1-S3) which indicates the effect of grazing on productivity changes in different environmental conditions.These results also align with other studies that demonstrate similar interactions between environmental factors and grazing across the globe (e.g., Yan et al., 2013 ).Environmental factors, especially lack of precipitation and high temperatures (i.e., drought conditions), may also influence grazing regimes outside of a strict experimental study due to landowners changing grazing regimes adaptively to the conditions of the grazed pastures.Therefore, it is necessary to account for environmental factors to isolate the effect of grazing on rangeland productivity ( Cibils and Coughenour, 2001 ).
We suggest further investigation into the complex relationship between grazing and rangelands in different systems and scales while accounting for environmental factors using remote sensing.This study is one of the few examples using remotely sensed estimates of productivity to understand the effects of grazing while accounting for both environmental factors and plant functional groups ( Díaz et al., 2007 ).It was feasible to conduct this study because access to remote sensing databases with appropriately fine resolutions improved.Most metrics necessary to examine the effects of grazing on productivity are now available in remote sensing databases, including other response metrics or variables more relevant to other ecosystems.However, fine-scale information that could be important to those relationships, such as vegetation species and palatability ( Adler et al., 2005 ), is not yet available.Furthermore, while remote sensing affords access to broad temporal and spatial measures of vegetation indices, our results reiterate the importance of spatiotemporal alignment between data sources.Spatial and temporal scale mismatches between remotely sensed vegetation data and livestock grazing may dilute rangeland responses.Future studies may benefit from using remotely sensed and field measures of vegetation metrics collected before, during, and after active grazing to further understand rangeland response to grazing ( Fig. 7 ).

Implications
This study highlights the complex interplay between grazing, environmental factors, and vegetation productivity in a sagebrush steppe rangeland system.The varying effects of different grazing metrics in our results imply that no single metric can fully represent the grazing-productivity relationship.Our results also suggest that, in our study area, environmental factors like moisture and temperature affect productivity more than grazing.Overall, this work demonstrates the utility of remote sensing for examining grazing and productivity dynamics over broad spatiotemporal scales.However, it also underscores the need for grazing data at appropriate resolutions, as mismatches may mask grazing effects.Thoughtful integration of field measurements and remotely sensed data will continue to provide landscape-scale insights into sustainable rangeland management under variable climatic conditions.
In the future, rangelands, on average, will likely experience increased temperatures, changes in water availability, and higher frequencies of severe weather events such as drought ( Bernstein et al., 2008 ).These factors will intensify stress on native vegetation and render the ecosystem more vulnerable to additional threats such as invasive annual grasses, increased fire frequencies, and conifer encroachment ( Neilson et al., 2005 ;Ziska et al., 2005 ;Miller et al., 2011 ).Thus, maintaining the integrity of rangeland ecosystems may become increasingly difficult, especially given their complex interactions with grazing ( Li et al., 2023 ).To effectively mitigate the impacts of these climate-related threats, grazing management must adopt an adaptive management approach tailored to the specific ecological context of each rangeland system ( Burkhardt and Sanders, 2012 ;Clark et al., 2018 ;

Figure 1 .
Figure 1.Map showing the pasture-and point-level data collection in our study area, near Roundup, Montana.Pastures are the black-outlined polygons and the point-level surveys are shown as yellow points.

Figure 2 .
Figure 2. Graphic showing the grazing management and environmental factors we considered influential to rangeland productivity.Each element is colored according to its data extraction method as indicated in the legend: yellow indicates data sourced from remote sensing products, green indicates data collected from the field, and gray is a category associated with the other components.Components considered to influence the light use efficiency of the plant (LUE) are temperature, moisture, and vegetation composition.Factors affecting absorbed photosynthetically active radiation (APAR) are shortwave radiation (SWrad) and fraction absorbed photosynthetically active radiation (FPAR).The effects of grazing can be measured using many factors, categorized into timing, duration, or intensity.Grazing intensity in this study was measured through the animal units per month (AUMs) divided by the number of meters squared in the pasture, dung patty counts, or the number of grazed plants.
to determine the goodness of fit for each model iteration.Each model has the basic structure shown below where mu i is the average rangeland response value for each individual survey point I , α is an intercept value, ˆ β 1 is the estimated coefficient for the grazing metric at each survey point, ˆ β 2... n are the estimated coefficients for environmental covariates n, and environment represents the remotely sensed values of the environmental variables in the final model at each survey point.

Figure 3 .
Figure 3. Examples of the centered and scaled coefficient density plots found showing the posterior distributions for the coefficient values of the highest grazing predictor coefficient value on the left (perennial grasses and forbs net primary production [PFGNPP] ∼ cow patties [CP]) and the lowest on the right (annual forbs and grasses above ground biomass [AFGAGB] ∼ animal units per month [AUMs]) in the final GLMMs for each rangeland productivity response and grazing metric.See Figures 5 , 6 , S1, and S2 or TablesS2 and S3for all statistically significant coefficient values in the final models and TableS1for covariate abbreviations.

Figure 4 .
Figure 4. Examples of partial effects of grazing predictor centered and scaled coefficient values on rangeland response variables based on final GLMMs for each rangeland productivity response and grazing metric.The most positive grazing coefficient is shown on the left (perennial grasses and forbs net primary production [PFGNPP] ∼ cow patties [CP]) and the lowest on the right (annual forbs and grasses above ground biomass [AFGAGB] ∼ animal units per month [AUMs]).The predicted rangeland response is shown by the black line, the 50% credible intervals shown by the dark gray shading, and the 95% credible intervals shown in the light gray shading.The data values of the grazing predictor and associated rangeland response values are shown using transparent points.See Figures 5 , 6 , S1, and S2 or TablesS2 and S3for all statistically significant coefficient values in the final models and TableS1for covariate abbreviations.

Figure 6 .
Figure 6.Covariate effects produced from the final GLMMs for each environmental covariate and rangeland response biomass variable for a given grazing predictor.Yellow circles indicate negative effects of the covariate on the response variable, and teal circles indicate positive effects of the covariate on the response variable.Gray circles indicate that the 95% credible intervals overlapped zero and the effect of the covariate was neutral.The size of the circles represents the effect size.

Figure 7 .
Figure 7. Diagram showing the hypothetical effects of grazing on rangeland productivity at different times of the growing season.Remote sensing facilitates frequent productivity and environmental factor sampling across broad spatial scales.Pairing remotely sensed data and seasonal grazing intensity metrics with measures of vegetation productivity before, during, and after grazing would better elucidate the response of vegetation to livestock grazing.