Although statistical population reconstruction (SPR) provides a flexible framework for estimating demographics of harvested populations using age-at-harvest data, that information alone is insufficient. Auxiliary data are needed to ensure all model parameters are estimable and to improve the precision and accuracy of the estimates. We examined the influence of two types of auxiliary information, independent estimates of annual abundance and annual harvest mortality from radio-telemetry studies, on the stability and precision of abundance estimates from SPR. Further, we evaluated whether the timing of auxiliary studies in the reconstruction affected the precision of abundance estimates. Monte Carlo studies simulated auxiliary data with precision levels defined by the coefficients of variation (CV) of 0.05, 0.125, 0.25 and 0.50 corresponding to the three levels of precision suggested by Robson & Regier (1964) for accurate research, accurate management and rough management and a minimal information scenario. For comparable levels of precision, radio-telemetry studies used to estimate harvest mortality stabilized the reconstructed population trends better than independent abundance surveys. However, independent abundance surveys were superior at improving the precision of reconstructed abundance estimates. We found that the timing of auxiliary studies did not influence the stability of SPR estimates, which has important implications for managers designing studies to collect auxiliary data. Our research highlights that different types and quality of auxiliary studies affects the precision and stability of SPR models differently.

In recent years, statistical population reconstruction (SPR) has emerged as a feasible method for estimating the demographics of harvested wildlife over large geographic areas using age-at-harvest data which are commonly collected by wildlife agencies (Gove et al. 2002, Skalski et al. 2007, 2011, Broms et al. 2010). Its origin can be found in the vast history of quantitative stock assessment in fisheries (Quinn & Deriso 1999). In wildlife science, where methods such as the Downing (1980) method and sex-age-kill (Millspaugh et al. 2009) are still the norm, there are advantages in using SPR which offers flexibility and robustness not available in traditional techniques. Further, SPR allows for the simultaneous estimation of survival, abundance, recruitment and harvest mortality whereas traditional reconstruction techniques estimate only total abundance.

Age-at-harvest data provide the primary source of information in SPR models. However, age-at-harvest data alone are insufficient to reconstruct population demographics using SPR. In addition to age-at-harvest data, one or more sets of auxiliary data are needed to estimate one or more of the parameters from the age-at-harvest likelihood, either survival rates, harvest rates, recruitment or abundance. The joint likelihood structure of SPR models is flexible enough to incorporate almost any form of auxiliary data. In the past, it has been common to include catch-effort data as auxiliary information for SPR models. Skalski et al. (2007) calibrated a black-tailed deer *Odocoileus hemionus* SPR model using catch-effort data. That was possible because hunter effort was deliberately manipulated to produce a strong contrast in harvest rates with alternative levels of effort. Similarly, Skalski et al. (2011) used a five-fold change in trapping effort over time to construct a catch-effort relationship and reconstruct the abundance of American martens *Martes americana* in the upper peninsula of Michigan, USA. In other populations where hunter or trap effort may be relatively constant over time, catch-effort data will likely be an insufficient form of auxiliary information. Thus, although catch-effort data has been a staple in SPR analyses other types of auxiliary data may be more useful.

Several other types of auxiliary data have been used in previous SPR analyses including radio-telemetry information on harvest rates and independent estimates of abundance. Broms et al. (2010) used radio-telemetry to help reconstruct greater sage-grouse *Centrocercus urophasianus* abundance in Oregon, USA. Radio-telemetry data were used to estimate vulnerability coefficients associated with harvest mortality of greater sage-grouse in Oregon. Fieberg et al. (2010) also used radio-telemetry to estimate harvest rates and reconstruct the abundance of a Minnesota black bear *Ursus americanus* population. Other researchers have used independent estimates of total abundance as auxiliary data. For example, Gast et al. (submitted) used independent mark-recapture estimates of total abundance to help calibrate an SPR model of elk *Cervus elaphus* in the upper peninsula of Michigan. Alternatively, Fieberg et al. (2010) chose not to use independent DNA mark-recapture estimates of abundance when reconstructing the Michigan black bear population. Instead, they chose to use that information as an independent source of confirmatory information. In contrast, Skalski et al. (2007) found that an independent browse damage index of deer abundance had little or no benefit in reconstructing a black-tailed deer herd. Such index data can help characterize the trend of a population but not its absolute abundance. These examples illustrate the flexibility of SPR, but also raise questions about the potential utility of various types and quality of auxiliary data. Uncertainty remains as to the effect of auxiliary data on model stability and how much auxiliary data are needed. Also, it is unclear how the precision of the auxiliary data (i.e. coefficients of variation (CV)) changes the precision of the reconstructed abundance estimates. Additionally, it is unclear if the timing of auxiliary studies changes their effectiveness.

The purpose of our paper is to provide game managers with guidance on how best to incorporate auxiliary studies in SPR. We evaluate whether precision of the auxiliary study is the sole consideration, or whether the types of parameters being estimated are also important. To this end, we compare the performance of SPR where independent estimates of abundance or harvest mortality are available. We also consider the timing of the auxiliary studies in relation to the duration of an SPR and the relative benefit of more than one auxiliary study in the precision of SPR.

## Methods

### Overview of SPR

SPR is based on age-at-harvest data collected over time by game management agencies. The observed counts h_{ij} (i = l,...,Y; j = 0,..., A) are modeled as a function of the initial abundance of a cohort and the subsequent natural survival and harvest over time and perhaps probabilities of reporting and age determination. Skalski et al. (2007) modeled the diagonals of this age-at-harvest matrix as independent multinomial distributions (L_{ij}) where the joint likelihood can be written as

This likelihood is often accompanied by a likelihood model describing the annual probabilities of harvested animals being report and/or aged, i.e.

where L_{i}are binomial sampling models for the fractions of animals harvested in year i being reported and/or aged. This likelihood can be omitted if there is 100% reporting and aging of all harvested animals.

Together, these two likelihoods are incapable of estimating the demographic parameters of interest, annual abundance, recruitment, natural survival and harvest mortality. Gove et al. (2002) proved that at least one demographic parameter must be estimated independent of the age-at-harvest data for SPR to be possible. Sometimes there may be, say, k independent auxiliary studies contributing to the reconstruction; hence, the joint likelihood model may be written as

Skalski et al. (2007) suggested using a catch-effort likelihood as an auxiliary where the annual harvest numbers are modeled as binomial, random variables as a function of the unknown total abundance in year i (i.e. N_{i}) and hunter effort (f_{i}) where

### Auxuliary likelihoods

We considered two alternative forms of auxiliary likelihoods in our evaluation. In one case, an unbiased annual abundance estimate _{i} was assumed to be available with estimated standard error _{i} The auxiliary likelihood then assumed the estimate was asymptotically normally distributed with the likelihood

We took this approach to allow the auxiliary likelihood to be independent of the form of the abundance survey and simply reflect survey precision (i.e. CV = σ_{i}/N_{i}). The second auxiliary approach was based on a hypothetical radio-telemetry study to estimate harvest mortality during the hunting season. Here we used a binomial of the form

_{i}= number of animals tagged and at risk of harvest, d

_{i}= number of tagged animals harvested, and where harvest mortality was reparameterized in terms of a vulnerability coefficient (c) and annual year-specific hunter effort (equation 1). Precision in the case of the radio-telemetry study was expressed in terms of

### Monte Carlo simulations

A Monte Carlo simulation study was used to determine the precision of population reconstruction estimates based on the amount, type and timing of auxiliary studies. A stochastic Leslie matrix model was used to generate age-at-harvest data for populations with different levels of natural survival rates and harvest rates. Recruitment levels were adjusted to produce populations with stationary abundance of approximately 6,000 animals in expectation. Recruitment was generated using a Poisson process, and natural survival and harvest were modeled as binomial processes.

In each simulation, 20 years of data were generated to establish demographic trends with years 21-44 used in the population reconstruction analysis. The full age-class data were generated and used in standard population reconstruction models. The same data were also reanalyzed after pooling the adult age-at-harvest data (i.e. 2.5+ year olds) using the pooled adult reconstruction of Skalski et al. (2012).

Demographic scenarios were performed to represent a range of scenarios expected for harvested large mammal populations. Natural survival probabilities were simulated at 0.75 or 0.90 and harvest rates at 0.10 or 0.25. To minimize the number of scenarios investigated, survival and harvest rates were assumed constant across all age classes. Auxiliary data were simulated to estimate either the annual abundance (_{i}) or a harvest probability (_{i}) with CV equal to 0.05, 0.125, 0.25 or 0.50. The CVs of 0.05, 0.125 and 0.25 correspond to precision levels described by Robson & Regier (1964) as appropriate for accurate research, accurate management and rough management. A fourth CV of 0.50 was simulated to represent a minimum information scenario. At this level of precision, a parameter is estimated within ± 100% of the true value 95% of the time. The effect of timing of the auxiliary data was tested for each parameter combination by staging the auxiliary study at either the beginning, middle, near end (i.e. year 23 of the reconstruction) or end (i.e. year 24) of the reconstruction.

Average measurement error of the reconstructed abundance estimates was estimated from the variance component expression

where =/(10,000 – 1) is the empirical variance among the abundance estimates in the jth year of reconstruction and =/(10,000 − 1) is the empirical variance among the true abundance values in the jth year of the reconstruction. A total of 10,000 simulations per scenario were used to obtain precise estimates of and. This approach provides a model-independent estimate of measurement error. Across the 24 years of reconstruction, reported precision was calculated in terms of median CV of measurement error.### Black-tailed deer sensitivity analysis

The previous Monte Carlo simulation studies looked at the relationship between the precision of abundance estimates from population reconstruction and the use of auxiliary data. This section examines the effect of auxiliary data on the stability of reconstructed population trends for one particular realized data set. The black-tailed deer reconstruction of Skalski et al. (2007) was selected for illustration because no auxiliary likelihood was incorporated in the original population reconstruction. Over a 24-year period, the abundance ranged between 1,500 and 3,500 does. Only catch-effort data were used to calibrate the model. The example is therefore convenient for illustrating the relative merits of population reconstruction without and with auxiliary data of varying degrees.

Skalski et al. (2012) recommended using point-deletion techniques to determine the stability of population reconstruction to varying amounts of historical information. For a reconstruction to be reliable, the estimated abundance trends should be relatively insensitive to the amount of historical data used in the demographic analysis. They recommended determining how stable the reconstruction abundance estimates were when 0, 1, 2, … years of the historical data were sequentially eliminated for the analysis. Following the advice of Skalski et al. (2012), simulated survey data to estimate abundance (_{i}) and harvest probability (_{i}) were added to the original population reconstruction with Cvs = 0.05, 0.125, 0.25 or 0.50. One such survey was assumed to have occurred either at the middle (i.e. 1991) or at the end (2002) of the 24-year population reconstruction (1979-2002). Stability was measured by the relative absolute deviation (RAD) in abundance defined as

_{i}= abundance estimate in year i from original population reconstruction using all years of data, N

_{ik}= abundance estimate in year i from a reconstruction with k historical years of age-at-harvest data deleted and y = number of years in the truncated reconstruction. The number of years deleted ranged from k = 0, 2, 4, 6, 8, 10 and 12 of the original 24 years of reconstruction.

## Results

### Monte Carlo simulations

As the precision of the auxiliary studies increased, precision of the reconstructed abundance estimates increased roughly proportional. With minimal information from an auxiliary study with a CV of 50%, the median CVs of the abundance estimates were intolerably large, usually > 100% (Table 1). As the CV of an auxiliary study used in estimating abundance went from 0.25 to 0.05, the median CV of the reconstructed abundance estimates was reduced by more than half and ranged from 0.354-0.798 to 0.082-0.354, respectively. The use of auxiliary data in conjunction with pooled age-class data (i.e. age classes 0.5, 1.5 and 2.5+) had the same pattern of improvement in precision as occurred for full age-class reconstruction (see Table 1). The only difference was a slight additional reduction in the anticipated CVs.

Auxiliary abundance studies had a greater influence on the precision of reconstructed abundance estimates than auxiliary harvest probability studies for equal precision (see Table 1). For example, when the auxiliary harvest study had a CV of 12.5% (i.e. with ± 25% of the true value 95% of the time), the resulting median CVs for the reconstructed abundance estimates ranged from 19.8 to 58.3% (see Table 1). For a similar level of precision in an auxiliary harvest mortality study, the population reconstruction estimates had CVs in the range 57.4-92.8%. In general, reducing the CV of an auxiliary study produced a commensurate reduction in the CVs of the reconstructed abundance estimates.

Timing of the auxiliary studies was generally not important. The same improvement in the precision of the population reconstruction estimates occurred regardless of whether the auxiliary study was conducted at the beginning, middle or near the end of the investigation. The only exception was when the auxiliary study was performed in the last (i.e. current) year of reconstruction. Under this circumstance, precision of the population reconstruction will not be measurably improved until a year thereafter.

Incorporating multiple auxiliary studies, 1/3 and 2/3 of the way through the reconstruction period, increased precision of the resulting abundance estimates (see Table 1). The effect is most substantial with auxiliary studies of minimal precision (see Table 1). However, a single auxiliary study with a CV of 0.125 results in better precision than two auxiliary studies with a CV of 0.250 each (see Table 1).

### Black-tailed deer example

Augmenting the original black-tailed deer data with an auxiliary abundance survey or a telemetry study, to estimate the vulnerability coefficient, greatly improved the stability of the population reconstruction when the amount of historical data was reduced. Without any auxiliary data, the black-tailed deer reconstruction was very sensitive when four or more years of harvest data were omitted from the analysis (Table 2 and Fig. 1). With six of 24 years of historical data removed, the RAD exceeded 150%. With even more years of data deleted, the reconstruction virtually disintegrated (see Table 2). However, the presence of rough auxiliary studies with a CV = 25% (i.e. ± 50% of the true value 95% of the time) resulted in reasonable stability of the reconstructed population trends. When the auxiliary data estimated the harvest probability (;_{i}) with a CV ≤ 0.25, the RAD ≤ 7.48% with as many as 12 years of historical data eliminated (see Table 2, Fig. 2). When the auxiliary data provided an abundance estimate within a CV = 0.25, the RAD ≤ 22.25% with as many as 12 years of data deleted (see Table 2, Fig. 3). By the time the auxiliary studies had a precise level suitable for accurate management purposes (CV = 0.125; Robson & Regier 1964), the RADs ≤ 6.32% with 12 years of data deleted (i.e. 50%) when auxiliary telemetry data were available.

Smaller RADs were obtained using radio-telemetry auxiliary data to estimate harvest probabilities rather than auxiliary abundance surveys for equal levels of precision (see Table 2, Table 3). Stability of the population reconstruction was not affected by whether the auxiliary study was conducted in the middle or near the end of the time series (see Table 3).

## Discussion

Our simulation and sensitivity analyses illustrate several important trade-offs in the quality and type of auxiliary studies used in SPR modeling. First, the value of even rather imprecise auxiliary data (e.g. CVs = 25%) on the precision and stability of population reconstruction were evident. This finding demonstrates the general utility of auxiliary information in SPR and can be used by managers as a guide on the required quality of future auxiliary studies. Second, different types of auxiliary data have different benefits to SPR estimates and the choice of auxiliary data components ultimately depends on the goal of the resource manager. If the primary concern is precision of the abundance estimates, auxiliary abundance studies are more beneficial than auxiliary radio-telemetry studies for comparable levels of precision. However, managers might be more interested in the inter-annual stability of abundance estimates when designing harvest regulations. If so, auxiliary radio-telemetry studies are more beneficial than abundance auxiliary studies for comparable levels of precision.

Another important question in SPR modeling relates to the timing of auxiliary data collection, and our findings suggested that improvements in precision are expected when auxiliary data are collected at any point during the reconstruction except in the last year. This finding is important for wildlife managers who might have historical age-at-harvest data and are considering conducting a contemporary auxiliary study. These results suggest the continued value of collecting auxiliary data after the collection of age-at-harvest data has begun. In other words, age-at-harvest and auxiliary data do not need to be collected simultaneously from the start of the study. Further, because auxiliary data can be collected at nearly any point during the reconstruction, a manager can be less concerned if data collection during an auxiliary study proves unsuccessful. For example, if radio-collars malfunction, a manager can attempt a telemetry study at a later date or alter plans and collect another type of auxiliary data without losing the ability to use SPR.

Our analyses also indicate that multiple auxiliary studies can further enhance the precision of population reconstruction estimates. In an adaptive management framework where decisions are updated as more and more information becomes available through time, this finding is relevant. SPR is not only flexible enough to handle multiple auxiliary data types, but we can expect precision of the resulting estimates to improve. However, if a manager is considering whether to complete one or two auxiliary studies that estimate abundance, it is important to consider that one study with high precision will improve precision of SPR estimates compared with two auxiliary studies with low precision. Although precision of demographic estimates is generally not considered in traditional models of population reconstruction (Millspaugh et al. 2009), we encourage managers to carefully consider the quality of auxiliary data and how it ultimately affects the results of SPR.

Our assessment, however, is void of any cost-precision comparison for radio-telemetry vs auxiliary abundance surveys. The feasibility for each type of survey will vary by species, geographic factors and labor cost. Given that precision of population reconstruction estimates improve roughly proportional to the improvement in precision (i.e. reduction in CV) of auxiliary studies, it should be fairly straightforward to perform a cost-benefit analysis. Field investigators should consider their end goals and perform a cost-benefit analysis when planning auxiliary studies to augment SPR models. Such an approach is likely to produce the most useful population estimates and gain support from administrators and stakeholders holding managers accountable for costs and reliability of modeling results.

## References

## Table 1.

Median coefficient of variation (CV) of measurement error in simulated population reconstruction models when including either abundance or harvest probability auxiliary studies at varying levels of precision in the middle of the reconstruction. Further, we considered one (single abundance auxiliary) and two abundance auxiliary studies (double abundance auxiliary). Populations were simulated at high and low levels of harvest and natural survival probabilities.

## Table 2.

Relative absolute deviation in annual abundance estimates from point-deletion sensitivity analyses performed on a statistical population reconstruction of female black-tailed deer (Skalski et al. 2007). Models had either no auxiliary data or auxiliary data that estimated abundance or the vulnerability coefficient. Auxiliary studies had a coefficient of variation (CV) of 0.05, 0.125, 0.25 or 0.50 and were simulated in the final year of study.

## Table 3.

Relative absolute deviation of abundance estimates from a point-deletion sensitivity analysis of female black-tailed deer comparing auxiliary studies simulated at the end (i.e. 2002) or the center of the reconstruction (i.e. 1991).

**The utility of auxiliary data in statistical population reconstruction**," Wildlife Biology 19(2), 147-155, (1 June 2013). https://doi.org/10.2981/12-076