Density estimation for small mammals from livetrapping grids: rodents in northern Canada
Abstract
Management agencies and quantitative ecologists need robust estimates of population density. The best way of converting population estimates of livetrapped small mammals to population density is not clear. We estimated population density on livetrapping grids with 4 estimators applied to 3 species of boreal forest and 3 species of tundra rodents to test for relative differences in density estimators. We used 2 spatial estimators proposed by Efford (2009) and 2 traditional boundarystrip estimators designed for grid livetrapping. We, analyzed markrecapture data from 104 trapping sessions from the boreal forest at Kluane, Yukon (n = 4,818 individuals), and 56 trapping sessions from tundra areas of Herschel Island and Komakuk Beach in northern Yukon (n = 1,327 individuals). For boreal forest rodents on average both boundarystrip methods produced density estimates larger than Efford's maximumlikelihood (ML) estimator by as much as 50% at all population densities up to 25 animals/ha. For tundra rodents both boundarystrip methods produced density estimates smaller than Efford's ML at low density (<1.5/ha) and larger than Efford's ML density by 36–63% at high density (25/ha). Efford's inverse prediction estimator produced larger density estimates than the ML estimator by 4% for the boreal forest and 32% for the tundra rodents. Relationships were high between all the estimators, such that trends in density could be inferred from all methods. Determining the bias in population density estimators in small mammals will require data from populations spatially closed and completely enumerated. For our small mammals Efford's ML estimator typically provided density estimates smaller than those produced by conventional boundarystrip estimators.
 DENSITY 4
 density estimation
 Dicrostonyx groenlandicus
 Lemmus sibiricus
 Microtus
 Myodes rutilus
 Peromyscus maniculatus
 Yukon
Small mammals in boreal forest and tundra ecosystems fluctuate dramatically in abundance, and to quantify these fluctuations we need to estimate density. The most common method for small rodents is to use mark–recapture data from livetrap grids. Techniques for the estimation of numbers for closed populations from mark–recapture data have progressed rapidly from the initial simple estimates provided by Lincoln–Petersen estimators to those used in program CAPTURE (Otis et al. 1978) and program MARK (White 2008). However, to get estimates of population density we need to determine the area occupied by the population members. Management agencies and quantitative ecologists need estimates of population density rather than numbers, and this leads to complications in estimation. In particular, the actual area being sampled by a trapping grid extends beyond the trapped area, and we need to estimate the effective area sampled to produce an estimate of density.
Two approaches have been developed to estimate population density on the area being trapped. The classical approach is to add a boundary strip to the edges of the trapping grid (Otis et al 1978). The conventional estimator of the width of the boundary strip is half the mean maximum distance moved (MMDM). The MMDM is identical to the average observed range length (distance between the most extreme captures of an individual). The 2nd approach is to add a boundary strip equal to half the asymptotic traprevealed range length (ARL—Jett and Nichols 1987). When sample sizes are small, both these estimators of the boundary strip are difficult to estimate precisely.
A new method of estimating density from markrecapture of closed populations was developed by Borchers and Efford (2008), Efford (2004), and Efford et al. (2009). Efford's approach was to fit a simple spatial model of animal trapping that estimated the probability that an animal is caught in a trap at a given distance from the center of its home range. This probability is assumed to follow a 2parameter spatial detection function. This approach assumes that animals occupy stationary home ranges whose centers are a realization of a homogeneous random spatial point process with density D. An individualbased model uses the assumption that the probability of detecting an individual is a decreasing function of the distance between its homerange center and a live trap. The simplest detection function used by Efford (2004) has 2 parameters that correspond to a measure of homerange size (sigma [a]), and the probability of capture at the center of the home range (g(0)). These 2 parameters can be estimated from the distances between recaptures of marked animals and the frequency of capture. The approach used in Efford's method assumes that the spatial point process is Poisson and that the shape of the detection function is halfnormal. Two possible methods of estimation of the spatial model have been developed. Efford (2004) developed the 1st approach through inverse prediction (IP) and simulation. A 2nd approach was developed by Borchers and Efford (2008) and uses maximumlikelihood (ML) methods to derive a density estimate. Efford (2009) and Borchers and Efford (2008) used simulation to demonstrate that the ML method gave the most accurate estimate of true density in comparison to other methods using the boundarystrip approach. In program DENSITY 4 Efford (2009) compares the results of his estimators with the enumeration data of Parmenter et al. (2003). The details of the estimation problem are given in Efford (2004) and Borchers and Efford (2008) and implemented in program DENSITY 4.
We compared the boundarystrip and the spatial model approaches to estimate density from grid livetrapping of small mammals in northern Canada. We compared density estimates from 178 trapping sessions of 6 species of small rodents in boreal forest and tundra habitats. We estimated densities of 3 species from the boreal forest—redbacked voles (Myodes [= Clethrionomys]rutilus), deer mice (Peromyscus maniculatus), and field voles (Microtus oeconomus and M. pennsylvanicus)—over 24 years (1986–2009) in the Kluane Lake region of southwestern Yukon, and 3 species from tundra areas of northern Yukon—brown lemmings (Lemmus sibiricus), collared lemmings (Dicrostonyx groenlandicus), and tundra voles (Microtus oeconomus) from Herschel Island and Komakuk Beach—from 2006 to 2009. We asked 3 questions: How closely related are the different density estimators? Is there any evidence of consistent differences in the estimators? Do any consistent differences vary among rodent species?
A basic limitation to developing general answers to these questions exists because they can depend on trap density in relation to small mammal density, trap arrangements, number of trap checks within a single trapping session, and possibly the particular species involved in a study. We were able to develop generality in density estimation with a combination of simulation studies and many specific empirical examples. We note that few data are available on field tests of different methods for estimating density (Efford et al. 2005; Jett and Nichols 1987; Parmenter et al. 2003). All field studies to date can determine only relative differences between estimators because the true density of the population is not known. Simulation studies suggest that the ML spatial estimator is unbiased and that the boundarystrip methods underestimate density (Borchers and Efford 2008; Efford et al. 2009).
Materials and Methods
Study areas.—The major study site was located in boreal forest of southwestern Yukon near Kluane Lake within the Shakwak Trench system (61°01′N, 138°24′W), and lies within the rain shadow of the St. Elias Mountains. Climate and vegetation for this region are described in Krebs and Boonstra (2001) and Turkington et al. (2002).
The 2nd set of study sites was on arctic tundra of northern Yukon, on Herschel Island, which lies 5 km offshore (69°34.2′N, 138°54.1′W), and on the adjacent mainland at Komakuk Beach (69°35.1′N, 140°11.2′W). Mean annual precipitation from the closest weather station at Shingle Point is approximately 161 mm and includes an average annual snowfall of approximately 78 cm (Burn and Zhang 2009). Herschel Island is dominated by 2 vegetation types. Approximately half of the higher ground is covered by a tussock tundra community composed of Eriophorum vaginatum, Salix pulchra, and an assortment of forb, moss, and lichen species. On most of the other half of the higher ground a 2nd vegetation type occurs in which the common plant species are Dtyas integrifolia, Poa arctica, Salix arctica, Lupinus arcticas, other forbs, lichens, and mosses. Along the coast, a few alluvial fans occur, dominated by willows (Salix richardsonii and S. arctica) and sedge and moss species.
At Komakuk Beach each of our trapping grids overlapped 3 habitats: cottongrassmoss tussock, sedgecottongrass meadow, and heath polygons. The tussock habitat is dominated by E. vaginatum and a variety of mosses (notably Aulacomnium sp.). Extensive cover is provided by Rubus chamaemorus, S. árctica, and Vaccinium vitisidaea. Forbs and grasses are uncommon. The sedgecottongrass meadow changes along a gradient of water flow from wet and sloping to wet and flat. The more sloping sites are dominated by a thick growth of Carex aquatilis, with frequent clumps of Salix pulchra. Eriophorum angustifolium, Carex saxatilis, and Carex chordorrhiza are common. On flatter sites Eriophorum scheuchzeri is dominant, with a variety of mosses. The heath polygons combine elements of each of the other 2 habitats. The major contrast between the 2 northern sites is that the vegetation of Herschel Island is largely upland, dry tundra, whereas that of Komakuk Beach is largely coastal plain, wet tundra.
Trapping methods.—The trapping methods are essentially the ones that we have been using for the last 50 years to study small rodents (Krebs 1966). Our approach has always been to have an excess of traps to prevent trap competition and to space the intertrap distance on the trapping grids in such a way as to catch a high fraction of the animals present.
At Kluane Lake redbacked voles, deer mice, and field voles were livetrapped from 1986 to 2009 on 4 unmanipulated sites spaced along a 25km section of the Alaska Highway, just south of Kluane Lake. Trapping grids were separated by a minimum of 2 km in continuous boreal forest, with scattered small patches of grassland. Each grid had 100 stations 15 m apart in a 10 × 10 array with 50 Longworth traps (Little Critter Live Traps, Rogers Manufacturing, West Kelowna, British Columbia, Canada) in alternate rows (i.e., Al, B2, A3, etc.). Actual grid area (1.81 ha) was designed to fulfill the recommendation of BondrupNielsen (1983) that grids be approximately 16 times the average home range of the species being studied. For the 3 rodent taxa at Kluane Lake homerange areas, estimated from body mass using the allometric equation of Harestad and Bunnell (1979), ranged from 0.05 to 0.10 ha. When rodent numbers were high and trap saturation was possible we used 100 traps, 1 at every station. Traps were prebaited with seed oats 1 week before trapping, and, where necessary, traps were placed inside mesh or metal cages to prevent disturbance by squirrels. Traps were covered with a board to protect them from sun and rain. Traps (locked open) were left in place all year. Trapping sessions were conducted over 3 days, typically in May and in August each year (but sometimes in September), and traps were set in the evening of day 1, checked in morning and evening of day 2, and locked open in the last check on the morning of day 3. Small mammals captured were tagged on the right ear with numbered fingerling fish tags (Boonstra and Krebs 2006; Gilbert and Krebs 1981). In this part of Yukon M. pennsylvanicus and M. oeconomus are nearly impossible to tell apart when alive, and we have grouped them together in the analysis. In snaptrap collections in this region M. pennsylvanicus has been the more abundant in a ratio of about 3:1 (Krebs and Wingate 1985).
We trapped brown and collared lemmings and tundra voles on Herschel Island and at Komakuk Beach from May to September from 2007 to 2009 as part of the International Polar Year program of research on climate change. Each livetrapping grid had 256 stations 20 m apart in a 16 × 16 array with 128 Longworth traps in alternate rows (i.e., Al, B2, A3, etc.). Because lemmings have large home ranges, actual grid size (9 ha) was larger than that used in the boreal forest. Brooks (1993) reported home ranges for brown lemmings from 0.4 to 1.3 ha and 0.03 to 6 ha for collared lemmings. Larger home ranges were for males at low population densities. Traps were prebaited with apple for 3–6 days before livetrapping. Traps were left in place and locked open all year. When set, traps were checked every 4–6 h, and in most cases a trapping session involved 5 or 6 checks over 2–3 days. Traps occasionally were closed in inclement weather.
Density estimators.—We investigated 4 density estimators: program CAPTURE population estimates adjusted with a boundary strip of ARL/2; program CAPTURE population estimates adjusted with a boundary strip of MMDM and MMDM/2; Efford's IP spatial model (Efford 2004); and Efford's ML spatial model (Borchers and Efford 2008). We investigated both the IP and the ML estimators because we have found some cases in which one or the other of these 2 estimators could not be solved from the available data. Efford et al. (2005) have shown that the IP estimator is affected by the choice of the program CAPTURE model used for estimating population size. We calculated markrecapture population estimates for all small mammals from the CAPTURE module in program DENSITY 4 using 3 heterogeneity (M_{h}) estimators: the jackknife estimator, Chao's estimator (Chao 1987; Chao et al. 1992), and the betabinomial estimator of Dorazio and Royle (2003). We also calculated population estimates from the heterogeneitytime (M_{th}) estimator of Lee and Chao (1994). The betabinomial estimator did not perform well with our data, providing confidence limits that were biologically impossible, so we did not use it. The remaining 3 estimators were highly correlated with one another, and we chose the jackknife estimator for all our data (r = 0.969, P < 0.01 between the jackknife and Chao estimates, and r = 0.970, P < 0.01 between the jackknife and the Lee and Chao estimates).
We do not know the true population density for any of the species studied in any of the trapping sessions, and our purpose here is to compare alternative methods to see how much they differ, rather than to conclude which method provides accurate density estimates. We used Efford's ML estimator as the standard of comparison because in simulations Efford (2009) showed that it was an unbiased estimate of true density. We do not assume here that the ML estimator gives the true density, which was unavailable in our study.
Relative differences among estimators were estimated as an average deviation of the fitted regressions by comparing the regression estimate with the value expected on the line of equal estimates (x = y). If the relative difference for any 2 estimators is 0, the predicted regression estimate will fall on the 45° line as shown in the figures. Relative differences among estimators are thus a statistical average, and any individual estimate of density could have more or less relative difference. Five of the 6 linear regressions between the 3 estimators and the Efford ML estimator were judged to be nonlinear from an analysis of regression residuals. A squareroot or log transformation rectified this statistical problem, and all further analyses were carried out with the appropriate transformation on both of the regression variables.
All density calculations were carried out in DENSITY 4.4 (http://www.otago.ac.nz/density). The detection function was assumed to be half normal. The buffer width (the spatial model analogue of the boundary strip) for boreal forest and tundra rodents was set at 100 m, approximately 3–4 × a. It was rare to get individual movements above these distances, and the estimates of density were robust to increasing the buffer width beyond these limits (Efford et al. 2009). The spacing for the integration mesh for the ML estimator was set to 64 × 64 points on the trapping grids. We used the jackknife estimator in DENSITY 4 to estimate population size. We used full likelihood to fit all models, and each trapping session was treated as an independent sample for estimation. The withinsession models of SECR (Spatially Explicit Capture Recapture) were the “dot” models, as defined in Efford (2009). In general, we took the default values for all the computations in DENSITY 4, except for the buffer width. All statistical analyses were done in NCSS (Number Crunching Statistical System, Kay, Utah; www.ncss.com). All correlations reported here are Pearson product moment correlations (r). All regressions computed are orthogonal regressions (Ricker 1984). All livetrapping of rodents was carried out in accord with the animal care guidelines of the American Society of Mammalogists (Gannon et al. 2007), and all our protocols were approved by the University of British Columbia Animal Care Committee.
Results
Descriptive statistics for the 6 sets of rodent species at the boreal and tundra sites were derived (Table 1). Not all the estimators could be obtained for every trapping session. In particular low densities are a problem with all markrecapture estimators. For our data none of Efford's estimators could be solved when only 2 animals were caught in a trapping session and only 23% of cases where 3 or 4 animals were caught (n = 26 cases). All of these small sample size cases were excluded from the data analysis. In virtually all cases where the minimum number alive was >4 we could solve all the estimators, provided some individuals moved between trap locations. This result reinforces the wellknown problem of trying to estimate numbers when capture rates are low.
Efford's spatial models assume no competition for traps, and except for 2 trapping sessions we never approached trap saturation (Table 2). When rodent numbers were high we increased the numbers of traps available to avoid trap saturation.
For each regression between estimstorss and the Efford ML estimator we asked whether the different species fitted the same regression line. We tested this hypothesis with a general linear model, with species as a covariate, and in none of the 6 regressions was species significant. Sample sizes for each species are given in Table 1, and for all of these regressions P 0.50 for a test of the species covariate. We have no evidence that the relationships described here are speciesspecific in these 2 ecosystems. Consequently for each ecosystem all species were combined in the analysis that follows.
For boreal forest rodents the boundarystrip ARL/2 estimator was closely related to Efford's ML estimator, but yielded on average a predicted higher density (Fig. 1a).This ranged from 31.5% greater at 1/ha (ML estimator) to 27.4%, 26.5%, and 26% at 5/ha, 10/ha, and 25/ha, respectively. For tundra rodents, a similar trend existed (Fig. 1b); although densities from the ARL/2 estimator initially were smaller than Efford's ML estimator (5% less at 1/ha), they were 14%, 23%, and 36% at 5/ha, 10/ha, and 25/ha, respectively. For both groups boundary strips estimated by ARL/2 were too small, and hence, density estimates too large (relative to the ML estimator), except for tundra rodents at low densities.
The 2nd boundarystrip estimator MMDM/2 also was closely related to Efford's ML estimator for boreal rodents (Fig. 2a). It produced on average density estimates that were larger than ML estimates, so that at an ML density of 1/ha, MMDM/2 estimated density that was 49.0% higher and at 25/ha 50.2% higher, a nearly constant average difference of 50% relative to ML estimates. For tundra rodents a similar trend existed. The MMDM/2 boundarystrip estimator (Fig. 2b) produced density estimates that were initially smaller than ML estimates below an ML density of 0.6/ha and then larger than ML density estimates above 0.6/ha. At 1/ha MMDM/2 estimates were 7% higher, at 5/ha 32% higher, at 10/ha 44% higher, and at 25/ha 63% higher than those produced by the ML estimator. Thus, for both boreal and tundra rodent species boundary strips estimated by MMDM/2 were too small relative to the ML estimates, and the difference in density estimation from the ML estimator for both boreal and tundra rodents was larger than that produced by the ARL/2 estimator. Thus the ARL/2 estimates were closer to the ML estimates than were the MMDM/2 estimates.
Because they are both based on the spatial prediction approach, Efford's IP model for density estimation would be expected to be the estimator most closely related to Efford's ML. The IP estimator for boreal rodents produced density estimates below those of the ML estimator below an ML density of 10.9/ha (Fig. 3a), so that at a density of 1/ha the IP estimate was 23.3% lower than the ML estimate, and at 10/ha 0.5% lower. Above the ML density of 10.9/ha the IP estimator estimated densities slightly higher than the ML estimator, so that at 15/ha the IP estimates were 1.6% higher, at 25/ha 3.7% higher, and at 40/ha 5.2% higher than ML estimates. Consequently, except at very low densities, the IP estimator produced only slightly different density estimates relative to the ML estimator, so IP provided a better estimate relative to ML than the boundarystrip estimators. For tundra rodents Efford's IP model also was closely related to Efford's ML (Fig. 3b). The IP estimator produced smaller estimates below a density of 0.9/ha, so that at a density of 0.2/ha the average IP density estimate was 11% smaller than the ML estimate. Above the density of 0.9/ha IP estimates were higher than ML estimates by 16% at 5/ha, by 23% at 10/ha, and by 32% at 25/ha. Consequently, for tundra rodents Efford's IP estimator produced estimates that differed from ML estimates by about the same amount as the ARL/2 boundarystrip estimator. We had few data points above 15 animals/ha, and more data from high densities would be useful for defining this tundra regression more precisely.
We investigated whether the MMDM might produce density estimates closer to ML estimates in comparison with MMDM/2 because it makes the boundary strip wider for grid trapping. For our data MMDM boundary strips underestimated ML density for both boreal forest and tundra rodents and did not duplicate ML estimates.
We investigated whether any of the capture parameters listed in Table 1 were correlated with population density. There was a weak correlation between density estimated by Efford's ML and a: r = −0.51, n = 96, P < 0.01 for Kluane rodents, and r = −0.38, n = 51, P < 0.01 for tundra rodents. This result is consistent with the general finding that rodents move less as density rises, although we recognize that a is a composite estimate that is partly related to movements but can be affected by the bait used and any aspects of trapping that affect capture rates. We had expected the probability of capture to decline with density, but we found only a slight tendency for this to occur in the Kluane rodents (r = −0.30, n = 96, P < 0.01) and none for the tundra rodents (r = −0.07, n = 51, P>0.10). This result is consistent with our attempts to always have an excess of traps relative to rodent numbers.
Discussion
An estimate of small mammal density that is both precise and unbiased has been the desired goal of many ecologists. Our results from these 6 species of small rodents from northern Canada are specific to our particular trapping methods and grid designs, and we do not know how general our conclusions are. A fundamental limitation in this study is that we do not know the true population density for any of our grids, and we only can compare among 4 estimators to evaluate relative differences. For the assessment of estimator accuracy we require data from all of these estimators for populations of known size, a luxury currently unavailable. Nor were we able to carry out trapping webs to test estimators (Jett and Nichols 1987; Parmenter et al. 2003). Trapping web designs, although conceptually elegant, are laborintensive, and our studies have focused on livetrapping procedures that can be carried on week after week, yearround to assess temporal dynamics in longterm studies.
One of the advantages of using the Efford ML estimator of population density is that trapping need not be done in a grid. We have used 2 transect lines of live traps for lemmings on Herschel Island and Komakuk Beach and found that Efford's ML density estimates obtained from these transect lines corresponded closely with grid estimates that used many more live traps, as long as a sufficient number of movements of individuals existed up and down or between the transect lines (D. Reid, pers. obs.). Unfortunately we do not have enough data (n = 9 sessions) for a formal analysis, but these preliminary results are sufficient to encourage a more detailed exploration of transect lines as possible sources of reliable density data for small rodents.
Obbard et al. (2010) used Efford's spatial SECR models to estimate density of black bears (Ursus americanus) in Ontario and reached conclusions similar to ours. They sampled 11 areas with densities varying about 3 to 4fold (coefficient of variation [CV] of density of 21–32%) and found that SECR methods produced density estimates that were 20–200% lower than those obtained by boundarystrip methods. They did not have the advantage we had of a series of populations existing at widely different densities (CV of density of 135–187%) to determine the change in potential bias in MMDM/2 and MMDM measures of boundarystrip size over a range of animal densities. They highlighted the management problems associated with estimates of large carnivore density that rely on traditional boundarystrip measures of density.
We have 3 important conclusions that can be hypotheses for future research. First, all of the species studied both in the boreal forest and in the tundra ecosystems fit the same regressions. We had thought that different estimators might be applicable to different species of rodents, but this was not the case. We need more information to determine if this conclusion applies to other rodent populations in other ecosystems. The rodent species we have reported on here are similar in body size in habitats of relatively low productivity so that a single grid size is sufficient to encompass many of their movements. Very large trapping grids provide one method of reducing edge effects (Gurnell and Gipps 1989; RajskaJurgie 2001), but they also require significant personeffort that is difficult to maintain in a longterm project. Similar problems occur with trapping webs and nested grids.
Second, the boundarystrip methods showed density estimates that are larger on average than those obtained by Efford's ML estimator (except at low densities). Our results are consistent with the suggestion of Efford (2009) that both the boundarystrip measures ARL/2 and MMDM/2 would produce higher estimates of population density than the ML estimator because the boundarystrip measures are expected to be underestimates of the effective trapping area.
Third, the 4 potential estimators of population density are very highly related in our data. Thus population trends will be similar no matter which estimator is used, as long as it was not changed during the course of a longterm study. The critical point, however, is that estimates of rates of population increase and quantitative estimates of standing crop and potential offtake by predators may be strongly biased unless one has an unbiased density estimator.
We do not know whether Efford's ML estimator is at present an unbiased density estimator for small mammal populations, and we do not assume this here. Efford's ML spatial estimator assumes a homogeneous habitat and circular home ranges, and we do not know how much these simple assumptions might affect the accuracy of the resulting estimates of density. If further work suggests that these spatial estimators are unbiased, we would be able to avoid the problems of estimating edge effects, highlighted by Gurnell and Gipps (1989) who suggested that we might require speciesspecific, seasonspecific, and habitatspecific correction factors to estimate the effective trapping area of a livetrapping grid. We still would be left with the problem of density estimation when rodent numbers are very low (minimum number alive < 4), and in these cases we might have to rely on counts of minimum numbers known to be alive, adjusted by an average ARL/2 boundary strip, or if possible to pool data across sampling units. Our main conclusion is that for our small mammal species Efford's ML estimator typically provided density estimates smaller than those produced by conventional boundarystrip estimators. Further work is required to determine which of the available density estimators for markrecapture trapping of small mammals are unbiased.
Acknowledgments
Without the efforts of M. Efford in developing the new spatial density estimators and coding them in the DENSITY 4 program, none of this work would have been possible. We thank him for his work. We also thank D. Moore, M. Leung, G.O. Cimon, A. Blachford, A. Fehr, M. Müller, D. Fehr, E. McLeod, S. McLeod, S. Ale, and B. Halliday for their assistance in field work. Comments by M. Efford, J. Hone, and 4 referees helped us to sharpen the presentation. Research funding was provided by the Natural Science and Engineering Research Council of Canada (RB and CJK), the International Polar Year Program (DR and CJK), and the EJLB Foundation (RB and CJK). The facilities of the Kluane Lake Research Station of the Arctic Institute of North America and Herschel IslandQikiqtaruk Territorial Park were essential to this research. We thank A. Williams at Kluane and R. Gordon and the park wardens on Herschel Island for their assistance with logistics.
Footnotes

Associate Editor was Dirk H. Van Vuren.
 © 2011 American Society of Mammalogists