This paper examines changes in broom snakeweed populations (Gutierrezia sarothrae [Pursh] Britt. & Rusby) from 1979 to 2014 at three prairie grassland sites in New Mexico. Data gathered each fall were used to study broom snakeweed population dynamics and to estimate the probability that the relatively short-lived subshrub will die off or invade blue grama (Bouteloua gracilis [H.B.K. Lag]) rangelands. Annual broom snakeweed standing crop data were used to categorize populations as None (<100 kg ha-1), Light (<300), Moderate (<750), or Heavy (≥750). Ordered logit regression was then used to estimate the frequency of transition between these categories over time depending on environmental and site factors. Significant variables found to influence annual variation in broom snakeweed included the broom snakeweed standing crop and density observed the previous period ( effect for continued broom snakeweed); grass standing crop the previous period (-); rainfall received from April to June ( ); and average temperatures during April ( ) and June (-). The probability of broom snakeweed invading an area that is currently without the plant ranges from about 1% to > 40% depending on environmental conditions and the amount of grass standing crop present. Transition probability estimates were also used in a Monte Carlo simulation model to evaluate the economics of broom snakeweed control. The economics of chemical broom snakeweed control were most strongly related to the rate of snakeweed reinvasion on treated areas and to the probability of natural die-off if infested areas were not sprayed.
Introduction
Some ecologists have considered the presence of moderate to dense stands of broom snakeweed to be an indication of disturbance or overgrazing (Campbell and Bomberger, 1934; Parker, 1939; Jaynes and Harper, 1978; Pieper and McDaniel, 1989). Yet multiyear establishment, survival, and change in plant density has led others to describe broom snakeweed populations as episodic and with changes in plant numbers not necessarily tied to disturbance (Jameson, 1970; Sosebee et al., 1979). Climatic patterns, particularly precipitation amount, are now considered to be the principle factors driving regional broom snakeweed populations in the southwestern United States (McDaniel and Sosebee, 1988; Beck et al., 1996; Beck et al., 1999). Besides weather, other factors shown to influence broom snakeweed establishment and survival include grazing, physical disturbance, fire, insect damage, and other natural events (Pieper and McDaniel, 1989; Davis et al., 2000).
In central New Mexico, the second quarter of the year (April through June) is a particularly important time for broom snakeweed propagation and survival (Pieper and McDaniel, 1989). Broom snakeweed seedlings can emerge any time of the year, but optimal establishment occurs under moist conditions when surface soil temperatures range from 19°C to 25°C (Wood et al., 1997). A 5-yr broom snakeweed population field study on the New Mexico State University Corona Ranch reported highest fecundity rates in April and May and highest mortality rates in June (McDaniel et al., 1997; McDaniel et al., 2000). Pulse establishment of broom snakeweed can be anticipated during the second quarter when optimal environmental conditions occur, especially when soil moisture is sufficient and air temperatures are moderate (Wood et al., 1997; Ralphs and McDaniel, 2011). Broom snakeweed seedlings and mature plants are highly competitive for soil moisture, but both are vulnerable and succumb easily to drought when air temperatures are high in June (Ragsdale, 1969; Osman and Pieper, 1988; Wood et al., 1997). According to a 75-yr long-term data set compiled by Dittberner (1971) on the US Department of Agriculture (USDA) Jornada Experimental Range in southern New Mexico, broom snakeweed seedling survival is < 1% if June rainfall is not at least normal or above. Dittberner (1971) estimated the mean life span of broom snakeweed that survives the first year to be about 4 yr, but some plants may live longer than 15 yr.
Early economic evaluations of broom snakeweed control (Carpenter et al., 1991; Torell et al., 1988) followed traditional procedures and assumed an average annual benefit from broom snakeweed control over a specified treatment life. Net present value (NPV) of snakeweed control was found to be positive with average beef prices if the control treatment lasted with certainty for 5 yr. Uncertainty about treatment life was a noted limitation of both economic studies. Torell et al. (1989) expanded to a stochastic assessment using a first-order stationary Markov chain model (SMM) to incorporate uncertain future outcomes and treatment benefits into the economic assessment. NPV was estimated to be —$0.06/ha (break-even), and uncertainty about future snakeweed infestations was found to be a major factor limiting the economic potential of broom snakeweed control.
In this paper, we examine a 35-yr record of change in broom snakeweed density and standing crop at 3 prairie locations in New Mexico. The analysis is an update and improvement of the SMM procedure used to estimate the transition probabilities in the earlier economic research (Torell et al., 1992). The widely used SMM is based on the rather restrictive assumption that transition probabilities remain constant over time. Probability estimates are obtained by computing the proportion of times that the weed infestation level at time t + 1 was j, given that it was in state i at time t (Hillier and Lieberman, 2010). An emerging modeling alternative used to overcome the time-invariant limitation is to use multinomial logistic regression to estimate probabilities as a function of explanatory variables (ecological examples include Augustin et al., 2001; Boltz and Carter, 2006; Breininger et al., 2010; Bino et al., 2015). We use multinomial regression to estimate how broom snakeweed transition probabilities change in response to weather and vegetation conditions. We then use the probability estimates in a simulation model to determine how broom snakeweed removal by herbicide control influences expected gains in standing grass production over time. The improved estimates of transition probabilities were then used to revisit the economics of broom snakeweed control.
Methods
Study Sites
In May 1979, three study sites were established to examine the effectiveness of various herbicides for broom snakeweed control (McDaniel, 1984). Other sites were also included in the herbicide control study in other areas of the state, but they are not discussed here. The three study sites are located on private ranches near Vaughn, New Mexico in Guadalupe County (34°29′07.48″ N 105°02′16.18″ W); near Yeso, New Mexico in De Baca, County (34°25′51.32″ N 104°41′45.45″ W); and north of Roswell, New Mexico (33°50′45.09″ N 104°54′08.26″ W) in Chaves County. All sites are in the Prairie region of central New Mexico's Highland Major Land Resource Area designation by the Natural Resource Conservation Service (NRCS). A more complete description of each site and the original herbicide study design is found in McDaniel (1984).
Broom snakeweed growing in association with primarily blue grama (Bouteloua gracilis [Kunth in H.B.K] Lag. Ex Griffins) forms the vegetation mosaic on the three study sites. Broom snakeweed is the dominant overstory plant, with occasional scatterings of walking stick cholla (Opuntia imbricata [Harr.] DC.). Common warm season grasses in addition to blue grama include black grama (Bouteloua eriopoda [Torr.] Torr.), sideoats grama (Bouteloua curtipendula [Michx.] Torr.), and sand dropseed (Sporobolus cryptandrus Hitchc.). Annual or perennial broadleaf species were sparse on the sites and most commonly occurred following monsoonal rains. All three sites occur on undulating shallow limestone hills with a restrictive caliche layer at about a 30- to 35-cm depth. The primary climate-soil-plant community classifications of the NRCS (2013) are Shallow at the Vaughn site (Ecological Site ID: R070CY113NM, see https://esis.sc.egov.usda.gov/); Shallow Plains (Cool) at Yeso (Ecological Site ID: R070BY069NM); and Gravelly at North Roswell (Ecological Site ID: R070CY119NM).
Precipitation and temperature data for each site were obtained from nearby online cooperating weather stations managed by the Western Regional Climate Center (WRCC, 2016). Average temperature was recorded as the monthly average of the average daily temperature. Average quarterly rainfall totals were similar across the three study sites, averaging 3.79 ± 2.40 cm (mean ± SD) during quarter 1, 8.21 ± 5.31 cm during quarter 2, 16.04 ± 5.75 cm during quarter 3, and 6.32 ± 4.28 cm during quarter 4. Average annual air temperature ranged from 1.95°C in January to 23.4°C in July. Below-average rainfall totals occurred during 1980, 1983, 1989, 1993, 1995, 2003, 2009, and 2011 — 2012 at all three study sites (Fig. 1).
Vegetation Sampling
Untreated replicated (2) plots (0.1 ha in size) established at each site in the McDaniel (1984) study were monitored to follow the long-term change in the native broom snakeweed population. Chemically treated plots were not used in the assessment. Initially (from 1979 to 1987), every site was visited in the spring (April) and fall (October) and broom snakeweed seedlings and mature plants were counted to determine plant density. Later (1988–2014), sites were mainly visited in the fall and only occasionally in the spring to search for seedlings in those years when early-season rainfall was known to be above normal. In this paper, we include field data taken only in the fall. Broom snakeweed density was determined by counting new (first year seedlings) and mature plants separately within 20 sample frames (30 × 60 cm) placed in each replicated plot (n = 40). When making plant density counts, we considered a plant alive if any green tissue was observed from nodes on basal stems. Broom snakeweed standing crop (kg ha-1) was determined by clipping plants in 10 sample frames per plot (every other frame previously counted) to a 3-cm stubble height. Placement of transect lines (two lines per plot, 38 m in length) was changed each year to prevent repeated harvesting. Grass standing crop (kg ha-1) was determined within the 20 sample frames per plot by ocular estimate (McDaniel, 1984; McDaniel, 1989). Clipped grass material from two random frames per plot was weighed in the field and later oven dried to adjust estimates to a dry weight basis. Grass standing crop at the Yeso site was collected in 20 of the 35 study yr. Missing years were 1982–1990, 1992, 1994–1995, 1998, 2001, 2010, and 2011. Similarly, grass standing crop data were not collected at the other two sites in 2010 and 2011. When grass standing crop was not recorded, a modified version1 of the sigmoid overstory-understory equation estimated by McDaniel et al. (1993) provided the grass standing crop estimate.
Grazing use was not specifically monitored during the study. However, when sites were visited each fall, an ocular estimate of standing crop removed by livestock (%) within the study plots was made. The study plots were located at least 1.5 km from a livestock water source, so typically grazing use was low and most years was rated as none to light (< 10% of herbaceous standing crop removed). As expected, grazing use varied by site and when use was rated above 10% we adjusted grass standing crop upwards to account for grazing biomass removal. Over the 35-yr study this correction was made 6 different yr on the Roswell and Vaughn sites and twice on the Yeso site. Landownership on these sites did not change, and our impression was that the ranchers adjusted herd size and grazing use responsibly according to prevailing weather and standing crop conditions. From 1979 through 1984 cattle and sheep grazed in mixed herds on the three ranches. However, when the wool market declined (roughly 1985), nearly all sheep were removed from the ranches by the late 1980s and replaced with cattle for the remainder of the study. A severe drought in New Mexico from roughly 2009–2012 forced ranchers to remove most livestock, and numbers remained low when our final observations were made in 2014.
Statistical Analysis
Transition Probabilities
We used ordered logit regression to estimate the probability that snakeweed will transition to different levels depending on site and weather conditions. The ordered logit regression model has an ordinal response variable, denoted as y, which is a discrete realization of a continuous random variable, y*. The y* variable has several defined threshold levels (kj), and the value of y depends on whether a particular threshold has been crossed (Borooah, 2002). In this application, we categorize broom snakeweed into four discrete states on the basis of the level of snakeweed standing crop (SWSCt) at time t:2
A visual illustration of these categories at the Vaughn site is given in Figure 2. As shown in McDaniel et al. (1993), when new snakeweed plants invade an area, they immediately suppress grass production but there is a diminishing marginal rate of suppression as the level of snakeweed infestation increases. With average rainfall, McDaniel et al. (1993) estimated an approximate 200 kg ha-1 reduction in grass biomass at the three study sites with movement through each snakeweed category.
The value of y* is defined to be a linear function of K factors such that
There are K betas and M-1 ancillary cut-off parameters to estimate, where M is the number of categories used.
Transition probabilities for each snakeweed state were estimated following procedures detailed by Borooah (2002). The probability that snakeweed standing crop will transition to state K is
Note that the estimated probabilities sum to 1. Further, the estimated transition probabilities are nonstationary and change each year with varying weather and site conditions.
For statistical analysis, we used the ordered logit subroutine in the STATA software package (StataCorp, 2011). With the literature review as the foundation, stepwise selection was used to select the variables to include in the final model with significance levels of 0.15 required to enter and stay in the model. The potential variables considered for possible inclusion are detailed in Table 1 along with averages and standard deviations computed for the variables. Reported statistics were similar by site, and thus they are not reported separately by site. We tested for multicollinearity using the variance inflation factor (VIF) diagnostic option. We explored several alternative models and selected the top candidate model with Akaike's information criterion (AIC).
Simulation Model
State transition probability estimates were incorporated into a Monte Carlo simulation model to predict future changes in broom snakeweed density and standing crop, as well as grass standing crop, under different rainfall and temperature conditions. Weather inputs were assigned in the simulation model using distribution fitting routines within the @RISK 7.5 simulation software (Palisade Corp., 2017). Weather data by site was used to estimate and select a probability distribution for each input based on the lowest Akaike's Information Criterion (AIC).
Key input variables selected for inclusion, the selected probability distributions, and distribution parameters in the stochastic simulation model are shown in Figure 3. Six explanatory variables predict the snakeweed state that will likely prevail at time t (Table 2). These include snakeweed density at t-1 (SWDENt-1), snakeweed standing crop at t-1 (SWSCt-1), grass standing crop at t-1 (GrassSCt-1), rainfall in quarter 2 of yr t (RAINQtr2t), average temperature in April of yr t (tAPRt), and average temperature in June of yr t (tJunt). Four explanatory variables estimate the grass standing crop (GrassSCt) available for grazing at time t (rainfall in quarter 3 at t [RAINQtrt], rainfall in quarter 3 at t-1 [RAINQtrt-1], rainfall in quarter 2 at t [RAINQtrt], and snakeweed standing crop at t [SWSCt]) (see footnote 1). Correlation between inputs was incorporated into the simulation with the following Pearson correlation matrix (estimated from weather data collected at the study sites):
Assigning an initial snakeweed state (SWSTATE) is the first step in the simulation process. Depending on the state chosen, the initial SWSCt-1, GrassSCt-1, and SWDENt-1 are set at mean levels computed for the selected initial state category (see Table 1). Next, using the appropriate distribution (see Fig. 3), the @RISK software randomly assigns values to each weather-related input variable and the probability of transition to each of four possible snakeweed states is estimated, given the simulated weather conditions for a particular year. Once the transition probabilities are estimated for yr 1, an @RISK discrete distribution function selects the snakeweed state that will prevail during yr 1. Depending on the computed SWSTATE for yr 1, SWDEN1 and SWSC1 are set at mean levels on the basis of the state selected and the process continues to compute values for the 15 yr considered in the simulation model.
Table 1
Mean and standard deviation of potential explanatory variables considered in the stepwise regression model
Table 2
Ordered logit regression results
Economics of Snakeweed Control
When broom snakeweed control is undertaken, there is uncertainty about both a natural population decline if the area were left untreated and the amount and length of time that forage and other potential benefits will be gained after control. The ©risk simulation model uses a standard decision tree structure and Monte Carlo simulation to capture these uncertainties in the economic analysis. Expected values are computed using the multinomial logit regression probability estimates to select the level of snakeweed present at each point in time. Expected economic returns with and without snakeweed control were estimated, and expected discounted net return differences, less the cost of treatment, define the net present value (NPV) of the broom snakeweed control alternative.
On the basis of survey results of 62 ranchers and land managers in New Mexico who had aerially sprayed broom snakeweed with the recommended chemical rate (Torell et al., 1989), we assumed that a fall treatment with Picloram would move snakeweed production on the treated area to state N with an 85% probability; state L, 10%; and state M, 3%, and it would remain in state H2% of the time. Snakeweed production proceeds to state k the next year. The “untreated” scenario starts the process in a degraded state j and proceeds to some other state i the next year. Weather variables were defined to be different between years but the same between treated and untreated branches of the decision tree within a year. Thus, differences in transition probability estimates result only because of differences in the starting snakeweed state at time 0. When k > l, grass standing crop differences from control would be economically beneficial and this would be expected to last for several years. When k = l, grass standing crop differences and thus economic benefits are zero. It is also possible to observe k < l, which is the situation when broom snakeweed present on the treated area is higher than that on the untreated area.
The added standing crop from broom snakeweed control was valued at its leasehold value. A 30% rule of thumb is recommended to adjust grazing lease rates to a payment for grass-only in range improvement assessments (Torell et al., 2014). This adjustment recognizes that about 30% of the average forage lease price pays for services provided by the lessor of the lease (e.g., periodic checking of water and livestock, supplemental feeding, maintenance of improvements and facilities). The $ AUM-1 lease price reported by USDA-NASS (2016) for New Mexico was first adjusted for inflation to a 2014 base yr and then adjusted downward by 30%. Using USDA data from 1979 through 2014, the distribution fitting routines within @Risk found a triangular distribution to be the best fit for grazing lease rates (see Fig. 3). The estimated annual $ AUM-1 leasehold price was converted to price per kilogram assuming 363 kg AUM-1. With the 30% discount and the assumed triangular distribution, grass standing crop added from broom snakeweed control was valued at an average of $13 AUM-1 ± $2.47 AUM-1.
An aerial herbicide treatment was considered in the economic analysis. The recommended herbicide treatment for broom snakeweed in the Southwest is a fall application of Picloram at 0.28 kg ae ha-1 or metsulfuron at 0.03 kg ai ha-1 (McDaniel and Duncan, 1987; Ralphs and McDaniel, 2011). The Picloram treatment currently costs about $40.76 ha-1, and this amount was subtracted from the estimated discounted forage value (using a 4% discount rate) to estimate NPV of the investment.
Results
Observed Changes in Broom Snakeweed Populations
Sites included in this study were selected, in part, because broom snakeweed was uniformly abundant across research plots and herbicide control testing was of primary interest (McDaniel, 1984). Rainfall the first yr of the study (1979) was near the long-term regional average for central New Mexico (see Fig. 1). On the untreated plots used in this study, broom snakeweed had healthy foliage and produced a robust flower and seed crop by the end of the 1979 growing season. However, from late autumn 1979 through the 1980 growing season, precipitation across much of central New Mexico was well below average. When sites were visited in 1980, mature broom snakeweed plants were noted as either highly drought stressed or dead. Broom snakeweed standing crop declined from fall 1979 to fall 1980 by 38% and 81% at Roswell and Vaughn, respectively (Fig. 4). All broom snakeweed was dead on the Yeso site when sampled in fall 1980. Precipitation in 1981 was near the long-term average, but favorable early spring and summer rains provided sufficient soil moisture at regular intervals over the growing season. A prolific number of new broom snakeweed plants (seedlings) propagated on all sites in 1981, presumably from the 1979 seed crop. New plants in 1981 that survived the first year comprised the majority of mature broom snakeweed plants that were observed in later years. In general, broom snakeweed seedlings were not common on the study sites except in 1987, 1993, and 2004. Largely, the most significant increase in the broom snakeweed population over the 35-yr study occurred in 1981.
From 1979 to 2014 the yearly difference in broom snakeweed density or standing crop was similar between the Vaughn and Yeso sites (see Fig. 4). Broom snakeweed was present at these locations in all years except in 2002 when no live plants were found at Vaughn and in 2009 when plants were absent at Yeso. At the North Roswell site broom snakeweed density and standing crop was not different from the Vaughn and Yeso sites the first 8 yr of the study (i.e., 1979–1986). However, in 1987–1988 the broom snakeweed population at North Roswell was reduced and then completely eliminated by a substantial population of round-headed root borers (Crossidius puchellus LeConte) (McDaniel et al., 1993). Larvae from these insects destroy the root system and are common in this region (Richman and Huddleston, 1981). In later years, few new broom snakeweed plants were observed on the North Roswell site and these plants were typically short-lived. During 12 of the past 19 yr, no live broom snakeweed was found on the Roswell site.
A significant competitive relationship was observed between the broom snakeweed-overstory and grass-understory, especially in the early years of the study (McDaniel et al., 1993). From 1979 to 1987, broom snakeweed states on the three sites were categorized as either moderate or heavy because broom snakeweed standing crop usually exceeded 600 kg ha-1 (see Fig. 4). Competition from broom snakeweed during this period resulted in highly suppressed grass growth (i.e., grass standing crop usually < 250 kg ha-1). Grass standing crop did not increase dramatically on any study site until broom snakeweed declined by natural causes to the light or none state categories.
Figure 4.
Standing crop of broom snakeweed and herbaceous grasses at each study site (1979–2014). Each stacked bar is the total standing crop measured during the year in kg ha-1. The X is the density of snakeweed measured as the number of snakeweed plants m-2.
Figure 5.
Broom snakeweed transition probabilities estimated from the logit regression model for average conditions. The shaded box is the category of broom snakeweed present at time t-1. Arrows and percentages indicate the probability of transition to that state at time t.
Table 3
Marginal effects of explanatory variable xi on transition probability j evaluated at the mean
Ordered Logit Model
The estimated statistical results for the logit model of best fit, including beta coefficients and significance levels, are presented in Table 2. Of the 15 variables considered for possible inclusion in the stepwise regression, three plant (SWSC, SWDEN, GrassSC) and three weather (tAPR, tJUN, RAINQtr2) variables met the specified rule for inclusion (P < 0.15). All included variables were statistically significant at the P = 0.05 level or higher. The likelihood ratio test for the overall significance of the equation was highly significant (P < 0.0001). The pseudo-R2 was estimated to be 0.43. Multicollinearity was not determined to be a problem with all variance inflation factor (VIF) values estimated to be < 3.0. The beta coefficients of the model are the ordered log-odds (logit) regression coefficients, and while they have no intuitive interpretation, the sign of the estimated parameters can be used to predict population direction. A negative beta coefficient indicates an increased probability of shifting downwards from one state to another when the explanatory variable increases. For example, the negative beta sign for tAPR and RAINQtr2 suggest as these variables increase, there is a greater probability for transitioning from state L or N toward states M or H. Similarly, increased June temperatures decrease the probability of broom snakeweed remaining within heavy or moderate states.
Figure 6.
Probability density function of discounted net present value for the Vaughn-Yeso and Roswell sites.
The DROSW site dummy variable was statistically different (P < 0.0001) from the DYESO and DVAUGHN dummy variables, but these two sites were not different (P > 0.15) (see Table 2). The early and prolonged decline in broom snakeweed populations at Roswell as compared with the Vaughn-Yeso sites is a logical explanation for the observed differences. We ran additional regressions that dropped the data from the Roswell site to see if regression coefficients were significantly changed. Evidence from Wald tests of various linear hypotheses indicated regression coefficients were not different when Roswell data were excluded from the regression model (P > 0.05).
Transition Probabilities
With all explanatory variables set at average levels, the one-step transition probability estimates for Roswell and the Vaughn-Yeso combined data sets are given in Figure 5. Depending on the initial broom snakeweed state specified, the regression model estimates that the snakeweed infestation will stay in the current state with a relatively high probability. For example, when snakeweed started in state H at the Vaughn-Yeso sites, the probability of remaining in this state was estimated to be 96%. If no snakeweed was the initial state, snakeweed infestation levels would remain at that low level with an estimated 30% probability; move to a light infestation with an estimated 37% probability; move to moderate with a 28% probability; and move to heavy with a 5% probability. At the Roswell site, once dense stands of snakeweed were reduced to low levels in 1988, snakeweed returned to a heavy state only one time, in 2005 (see Fig. 4). This is reflected in the probability estimates with < 1% chance that the process will move to a heavy state from any other state (see Fig. 5). The process will remain in a desired none or light state with > 90% probability.
Table 3 gives STATA estimates of marginal effects for all explanatory variables included in the multinomial logit regression model when evaluated at mean levels. Interpretation of the coefficients must be made separately for each variable assigned to each marginal effect category. For example, the coefficient + 0.0687 shown for tAPR under the heavy snakeweed group would be interpreted to mean that when the average temperature during April increases by 1°C, the probability of transitioning to a heavy stand of snakeweed would increase by 6.9%. The probability of transitioning to a light snakeweed stand would decrease by 4.8%. If there were a 100-kg ha-1 increase in grass standing crop (GrassSCt-1), the probability of transitioning toward a heavy category of snakeweed would decrease by about 5% while increasing the probability of having no snakeweed by about 2%. At mean levels, if the site location changed to Roswell, the probability of transitioning to a heavy or moderate snakeweed stand would decrease by an estimated 47% (—0.445 — 0.0270 ≈ —0.47) with an equivalent increase in the probability of moving to light or none (0.312 + .160 ≈ +0.47).
Table 4
Detail of mean net present value calculations.
Economics of Snakeweed Control
The NPV of snakeweed control, calculated over 10,000 iterations, ranged upwards to about +$250 ha-1 when the treatment was simulated to be highly effective and long-lasting. NPV was negative when snakeweed died from natural causes on untreated areas and reinvaded quickly on the treated areas (implying a negative amount of harvestable standing crop added from the treatment) (Fig. 6). Discounted NPV of broom snakeweed control was most strongly correlated with the amount of snakeweed on the treated (Spearman rank correlation, r = 0.64) and untreated areas (r = -0.44). Over 70% of the variation in NPV resulted from variation in snakeweed levels on the treated and untreated areas. Other inputs, including rainfall and assigned forage value, explained about 15% of the variation in NPV.
Differences in the probability of snakeweed transitions explain the PDF distribution differences between the Roswell and Vaughn-Yeso combined sites. Snakeweed production at the Roswell site was estimated to be relatively stable with a high probability of moving to state N with only minimal snakeweed (see Fig. 5). This decreased the variability in NPV for this site relative to the Vaughn-Yeso site (see Fig. 6). Table 4 gives additional detail about how the mean NPV was calculated in the simulation model, and this detail is compared across study sites.
Vaughn-Yeso In-depth Analysis
The average simulated GrassSC at the Vaughn-Yeso combined sites was estimated to be 222, 344, 502, and 654 kg ha-1 (± 111 SD) when broom snakeweed production was in state H, M, L, and N, respectively. GrassSC was not considered to be harvestable if snakeweed was in state H, as described earlier. Thus, average estimated harvestable standing crop (HARSCt) for each snakeweed state was 0, 344, 502, and 654 kg ha-1. Across all iterations, HARSCt was weighted by the probability that each of the alternative states was realized. For example, following snakeweed control, the yr 1 snakeweed state was defined to transition to state H, M, L, and N with probabilities 2%, 3%, 10%, and 85% (based on the assumed success of treatment). With these probabilities, the weighted average HARSC for yr 1 was estimated to be 616 kg ha-1 (see Table 4). By comparison, if left untreated and with the stochastic weather patterns generated in the simulations, 92% of the iterations stayed in undesirable state H, 7% were in state M, and 1% in state L, resulting in a weighted average of only 27 kg ha-1 of additional harvestable standing crop. The productivity shift from the treatment resulted in an average of 589 kg ha of additional harvestable standing crop during yr 1 (see Table 4).
When the Vaughn-Yeso sites have a moderate or heavy stand of broom snakeweed, the area remains essentially unusable for long periods if not treated. Undesirable levels of broom snakeweed were estimated to persist to yr 15 on > 81% (66% in state H, 15% state M) of the 10,000 iterations when untreated (see Table 4). As snakeweed died, harvestable grass standing crop on the untreated area increased to about 160 kg ha–-1 by the end of the 15-yr planning horizon, decreasing to 185 kg ha-1 if the area had been treated. The average NPV for the Vaughn-Yeso site was estimated to be $41.39 ha-1 with a 4% discount rate. The chemical control treatment paid for itself within 3 yr with the assumed prices and costs.
Roswell Site In-depth Analysis
On the basis of population data from the Roswell site, if broom snakeweed control was successful then the treatment was estimated to be long-lived. Only 8% of the iterations remained in an undesirable state 15 yr after control (see Table 4). The corresponding production of harvestable standing crop on the treated area remained nearly constant at about 600 kg ha-1. Added grass standing crop benefits attributable to the chemical control treatment were estimated to be short-lived, however. This is because if the area were left untreated, the snakeweed stand would die from natural causes at a relatively fast pace. After 5 yr, only 22% of the iterations remained in an undesirable state (see Table 4) and within about 10 yr snakeweed production and density on treated and untreated areas would be indistinguishable. Average grass production on the treated and untreated areas would be relatively high and similar. The average NPV for the Roswell site was estimated to be $1.22 ha-1 a breakeven investment. However, because of the sustained natural population decline of broom snakeweed on this site, the discounted NPV was negative 56% of the time (see Fig. 6).
Discussion
Ralphs and McDaniel (2011) note that snakeweed populations often establish in years with above-average spring and winter precipitation, especially if the area was disturbed. Snakeweed is not particularly drought tolerant and tends to die off during dry periods when summer temperatures are hot. Once established, broom snakeweed is competitive and will persist for a number of years. These observations were consistent with our logit regression results, which provided a unique quantification of the probabilities of mortality, stability, and natural snakeweed die-off.
Initially, the level of snakeweed infestation on all three study sites was similar. We believe the transition probability model defined for the Vaughn-Yeso sites represents the normal situation when weather variables drive plant growth on the grasslands of central New Mexico. However, the Roswell model demonstrates the importance of considering an unexpected response, such as that which occurred with the localized infestation of round-headed root borers that eliminated the broom snakeweed population. It was a surprising result to us that the area remained virtually snakeweed free for the next 25 yr. This illustrates the importance of repeated sampling over time and the importance of having multiple sites.
Economic studies have assigned different levels of importance to decreased grazing capacity from snakeweed infestation, animal death losses, and reduced animal production efficiency. The economic valuation of this paper directly considers forage production differences as snakeweed invades and indirectly considers death loss and production efficiency differences. Because no grazing was considered with heavy snakeweed present, a 100% discount of production returns is implied when a heavy stand of snakeweed is present. With a 50% production discount for heavily infested snakeweed areas, the calculated NPV was reduced from $41.39 to $24.70 ha-1 at the Vaughn-Yeso site and from $1.22 to —$4.39 ha-1 at the Roswell site.
The estimated transition probabilities are linked to stochastic weather variation and annual site characteristics. The logit regression model could be further improved if grazing use decisions (i.e., the amount of residual forage that remains after grazing, percent forage use) were incorporated as an explanatory variable. We only superficially considered grazing effects because these data were not available. This is a limitation of our simulation model and rangeland management models in general. Grazing management decisions could potentially be greatly improved if knowledge about grazing impacts on transition probabilities could be incorporated into decision-support models.
The NPV analysis conducted here calculates over a 15-yr treatment life with declining benefits over time (see Fig. 6). Asset replacement theory suggests a more frequent treatment schedule could be more profitable. There is an opportunity cost of letting the range continue to deteriorate to where all benefits from the initial treatment have been exhausted (Torell and McDaniel, 1986). Given the transitory and variable length of time before snakeweed reinvades back to a heavy state, a retreatment rule might be to retreat the area with chemicals whenever it has returned to a heavy snakeweed infestation. Further, as described by McDaniel and Ross (2002), other control methods, such as burning, might be a feasible lower-cost alternative to chemical treatments.
Implications
The long-term monitoring data used in this study provided a unique opportunity to quantify transition probabilities for the infestation and die-off of broom snakeweed on the shortgrass plains of central New Mexico. Further, our study highlights and demonstrates the usefulness of long-term data measured across multiple study sites. Very different conclusions would have been drawn if data had been gathered for only the first 10 yr when snakeweed populations were high versus the more variable snakeweed populations observed over the 35-yr duration of the whole study (see Fig. 4).
Quantitative ecological data are seldom available but considered to be more reliable than professional opinion and observation for predicting ecological dynamics and estimating transition probabilities (Allen-Diaz and Bartolome, 1998; Kachergis et al., 2013). The logit regression analysis conducted here was particularly useful for estimating how selected environmental variables influenced transition probabilities between various states. It can be expected that logit and probit regression models will have an increasing role defining and predicting ecological dynamics.
Acknowledgments
We thank former graduate students Teresa Sedlacek and Jose Mayen Solorzano for their assistance in data collection and initial model development. We also thank numerous former graduate students in the Department of Animal and Range Sciences for their assistance in collecting data over the numerous years of this study. This research was supported by the New Mexico Agricultural Experiment Station and Southwest Regional Climate Hub.
References
Notes
[1]
1 The sigmoid equation of McDaniel et al. (1993) was reestimated after adding quarter 3 rainfall measured the year before as an additional explanatory variable. The nonlinear revised equation was estimated to be 
. The additional lagged rainfall variable was found to be important in subsequent research reported by Torell et al. (2011).
