Trend of the Yellowstone Grizzly Bear Population

Yellowstone’s grizzlies ( Ursus arctos ) have been studied for more than 40years. Radiotelemetry has been used to obtain estimates of the rate of increase of the population, with results reported by Schwartz et al. (2006). Counts of females with cubs-of-the-year “unduplicated” also provide an index of abundance and are the primary subject of this report. An exponential model was ﬁtted to n = 24 such counts, using nonlinear leastsquares. Estimates of the rate of increase, r , were about 0.053. 95% conﬁdence intervals, were obtained by several di ﬀ erent methods, and all had lower limits substantially above zero, indicating that the population has been increasing steadily, in contrast to the results of Schwartz et al. (2006), which could not exclude a decreasing population. The grizzly data have been repeatedly mis-used in current literature for reasons explained here.


Introduction
Yellowstone's grizzlies have been studied extensively for over 40 years.A recent monograph [1] provided a detailed analysis of much of the accumulated data.Two methods have been used to estimate trends of the population.One, the "unduplicated" counts of females with cubs-of-the-year, was developed by Dr. R. R. Knight and has been maintained since 1973.The second approach uses radiotelemetry data and was pioneered by Craighead [2].The "unduplicated" counts along with some earlier observations were converted to an index described as a "minimum population" estimate [3].That index used a three-year moving average to smooth the annual numbers of adult female bears and a calculated proportion of adult females in the population of 0.274.A long series of "study committees" [4] continued use of the index because a "minimum population" estimate appeared to be something readily understood by administrators and the public.
A weakness in the "minimum population" index is the presence of a strong serial correlation, induced by the fact that each annual estimate contains the counts from 2 other years.Difficulties caused by serial correlation have been assessed by Watt [5,6], Chapman [7], and Freckleton et al. [8].Eberhardt [9] indicated the sizable degree of correlation thus induced when the underlying data are constructed only from random numbers and showed that a correlation of −0.707 could be induced in analysis of two independent sets of random numbers.Royama [10] reported the same result.The "minimum population" index is not used here, and the Durbin-Watson test [11] does not show evidence of serial correlation in the data used here.
An alternative approach using radiotelemetry was applied by Eberhardt, et al. [12] using an approximation to Lotka's equation proposed by Eberhardt [13].That analysis yielded an estimate of λ of 1.046 with bootstrap 95% confidence limits of 1.00 to 1.09.That study depended on bears caught and marked in a set of "backcountry" traps, located away from centers of human activity.The major, and essentially the only, cause of adult grizzly mortality comes from conflicts with humans.When such conflicts arise, efforts are made to capture and move bears.Much experience has shown that bears involved in such a situation have a low survival rate.Consequently, bears first captured in "conflict" situations were not used in the analysis, but bears firstcaptured in the so-called "research" trapping continued to be used in the analysis after they were involved in conflicts with humans.Subsequently, Pease and Mattson [14] proposed that there are "wary" and "unwary" bears so that the "unwary" bears first-captured in "conflict" situations should be included in the telemetry sample.This, of course, reduced the lower confidence limit on population trend below 1.0.The issues thus involved have been extensively explored in the monograph by Schwartz et al. [1] and need not be further considered here.However, probability distributions generated by stochastic simulations of Yellowstone grizzly numbers [15] have lower limits below λ = 1.00, thus not excluding the prospect that the population has not increased.Again, these results have been amply examined in the monograph of Schwartz et al. [1].The purpose of the present study is to utilize the original "unduplicated" counts as an independent measure of population trend and to correct some erroneous uses of the data.

Study Area
The study area constitutes Yellowstone and Grand Teton National Parks, 6 adjacent National Forests plus some state and private lands, and totals about 34,500 km 2 .It is known as "The Greater Yellowstone Ecosystem" (GYE).The GYE constitutes the Yellowstone Plateau, and surrounding mountain ranges above 1500 m.Long cold winters and short summers characterize the climate.Low elevations are covered by grasslands or shrub steppes.Douglas fir (Pseudotsuga menziesii) and lodgepole pine (Pinus contorta) are the dominant tree species.A detailed description of the study area is available in Schwartz et al. [1].

Material and Methods
The model used here is the exponential r is the annual rate of increase and t is time in years, while a denotes an initial value.This equation was fit using nonlinear least-squares [16], as implemented in the R-language [17], and ordinary least-squares regression after a logtransformation of the trend index.Utility of the logtransformation was demonstrated by Eberhardt [13].The data were also examined using a smoothing technique, "lowess" [18,19].Confidence intervals on r were obtained by bootstrapping [19,20] for the nonlinear model and from the usual linear regression model for the log-transformed data.Several different confidence intervals were calculated in the R-language program boot.ci, as expounded by [20, chapter 5] A generalized linear model (GLM) [21] with Poisson errors was also fit to the data using program glm in the R-language.All bootstraps used 5,000 calculations.The Durbin-Watson test [11] was used to test for serial correlation.This test depends on the fact that the squared difference between successive deviations can be used to approximate the variance of the deviations if the pattern of deviations is random.The test is where the observations are deviations from a fitted model.Values of d lie between 0 and 4, with zero correlation at d = 2. Bias can be estimated in bootstrapping, using the following: where θ * is the mean of the bootstraps and tF * is the original estimate [19].Data used here on the "unduplicated" females with cubs in the GYE are those given in

Results
The unduplicated data set was initially assessed by the lowess program, which is essentially a data smoothing program that calculates values from a weighted linear regression of a subset of points adjacent to the current location.It was used to produce Figure 1.
The J-shaped nature of the plot is in consequence that measures to protect and enhance grizzly populations in Yellowstone did not take effect until about 1983 [1,3].All subsequent analyses given here use the data from 1983 onwards.Fits of the data from 1983 onwards using nonlinear least-squares with (1) appear in Figure 2.
The estimate of rate of increase (r) for the exponential fit (1) was 0.052 with a standard error of 0.0071.The Durbin-Watson test when applied to the exponential fit yielded a D-statistic of 2.32, a nonsignificant result, indicating no evidence of serial correlation.Bias calculated from (3) was negligible in all instances.Results of the various analyses appear in Table 1.

Discussion
The fact that the difference between r from the exponential fit and 0 is over 7 times the standard error provides strong evidence that it is significantly different from zero.This result is clearly supported by the several confidence interval calculations (Table 1), so it seems clear that the "unduplicated"  counts show a significant upwards trend.Two different models are considered, one the linear model with additive errors resulting from the logarithmic transformation of the data and the second a multiplicative error structure for the nonlinear least-squares fit to the exponential model.Confidence intervals from both models demonstrate that the rate of increase, r, is substantially greater than zero.We use r, the "intrinsic rate of increase" here, rather than λ = exp(r)as calculated from reproductive and survival rates, but, for small r, the difference is trivial.
The data indicate that the GYE bear population does not yet show signs of approaching an asymptote.This is a bit surprising, inasmuch as Schwartz et al. [1].detected signs of decreasing survival in the population surrounding the primary inhabited area.One possible explanation is that the expansion zone, being crudely an annulus of substantial area, can still provide a substantial increment to the population  in spite of having lower overall survival rates.Clearly, the population will ultimately level off, and its numbers may then be fit with a logistic curve.We have favored a "generalized logistic" model in which the rate of increase remains high until the asymptote (K) is approached [22], but believe the present data set is too variable for that model to be useful until the population begins to approach an asymptotic value.
Harris et al. [15].proposed use of the Chao2 estimator [23] to "estimate the total number of females with cubs present from the estimated number observed" [15, page 17].Using their estimates [24,Table 4] we fitted (1) by nonlinear least-squares giving the results shown in Figure 3.This fit had a residual standard error of 8.36 as contrasted to a value of 6.11 from the data of Figure 2, suggesting that the Chao2 estimator gives substantially more variable results than using the unadjusted estimates of females with cubs of the year.

International Journal of Ecology
Conducting live-trapping, marking and radiotelemetry operations on grizzlies are very expensive operations.Hence, the main lesson from the present study is that index studies should also be attempted as a check on results.Transforming observed index data with the Chao2 estimator produced slightly more variable results, so the method may not be useful unless some improvements are possible.The radiotelemetry results provide essential data for management so should also be conducted whenever possible.
The time span of the reported data (approaching 40 years) has led to its use in a long list of applications [25][26][27][28][29][30][31][32].All of these ignore the fact that the initial data in the series were constructed from different data sources not compatible with the 1973 development of the "unduplicated" index.As noted above, the population was evidently decreasing until about 1983.This lead to various forecasts of ultimate extinction [27,28] a somewhat decreased risk [31, Table 1], a lower limit on lambda of 0.99 [26], and a very small probability of extinction [29].Staudenmayer and Buonaccorsi [31] indicate that some of their estimates indicate a rate of increase less than zero.They also used the data to estimate extinction rates.Brook and Bradshaw [25] used the data to study densitydependence and did not consider extinction of the population.

Figure 3 :
Figure 3: Plot of Chao2 data with fitted exponential function.

Table 1 :
Estimates of rate of change (r), bias (3), and 95% confidence intervals for various methods of estimation using the exponential model.