ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 7 THE IMPORTANCE OF INDIVIDUAL VARIATION IN DEFINING HABITAT SELECTION BY MOOSE IN NORTHERN BRITISH COLUMBIA Michael P. Gillingham and Katherine L. Parker Natural Resources and Environmental Studies Institute, University of Northern British Columbia, 3333 University Way, Prince George, British Columbia, Canada V2N 4Z9, email: michael@unbc.ca ABSTRACT: Understanding resource use and selection has been central to many studies of ungu- late ecology. Global positioning satellite (GPS) collars, remote sensing, and geographic information systems (GIS) now make it easier to examine variation in use and selection by individuals. Resource selection functions, however, are commonly developed for global (all animals pooled) models and important information on individual variability may be lost. Using data from 14 female moose (Alces alces) collared in the Muskwa-Kechika Management Area of northern British Columbia, we examined differences among global and individual resource selection models for 5 seasons (winter, late winter, calving, summer, and fall). The global models indicated that moose selected for mid-elevations, and for deciduous burns and Carex sedge areas in all seasons. Resource selection models for individuals, however, indicated that no individuals selected the same attributes as the global models. We also ex- amined selection ratios among seasons with individual moose as replicates, and within individuals with bootstrapping techniques. We discuss the importance of considering individual variation in defining resource selection and habitat use by moose and contrast the results of selection ratios and resource selection models. We also use these data to illustrate some of the pitfalls that can be encountered using the 2 methodologies. ALCES VOL. 44: 7-20 (2008) Key words: Alces alces, habitat selection, home range, individual variation, resource selection, selec- tion ratio To better inform management strategies, wildlife research has long focused on under- standing use of habitats and, when combined with the availability of resources, what animals select and avoid on the landscape. Early studies using radio telemetry examined population- level habitat selection (e.g., Neu et al. 1974, reviewed by Thomas and Taylor 1990) usually concentrating specifically on habitat types (reviewed by Alldredge and Griswold 2006). The availability of global positioning satellite (GPS) collars and advances in remote sensing and geographic information systems (GIS) now enable researchers to more easily exam- ine variation in selection among individuals (Thomas and Taylor 1990, 2006). Numerous recent studies on selection by ungulates including moose (e.g., Osko et al. 2004, Dussault et al. 2005, Poole and Stuart- Smith 2005, 2006, Poole et al. 2007) have used resource selection functions (RSF, sensu Boyce and McDonald 1999, Boyce et al. 2002, Manly et al. 2002), although other multivariate approaches have also been employed (Nikula et al. 2004). The RSF models provide a broad- scale perspective of general selection patterns on the landscape (Boyce and McDonald 1999, Manly et al. 2002). They also accommodate any type of habitat variables (categorical and continuous) and easily incorporate spatial data acquired from GIS or remote sensing (Boyce and McDonald 1999). Our best understanding of variation in resource selection by ungulates comes from study designs in which use and availability of resources are measured for individual animals VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 8 (Design III in Thomas and Taylor 1990, 2006). GPS locations can provide relatively accurate estimates of use by ungulates (but see D’Eon et al. 2002, Frair et al. 2004). Multiple as- sumptions, however, are inherent in estimating resource availability for individual animals. Various studies have assumed that random points drawn from the landscape or on an overlay of all home ranges (e.g., Poole and Stuart-Smith 2006), from individual home ranges (e.g., Nielson et al. 2002, Gillies et al. 2006), from a buffer of potential movement radius (e.g., Arthur et al. 1996), or from a ra- dius of movement between consecutive fixes (i.e., matched-case; e.g., Johnson et al. 2002) constitute a sample of individual resource availability. When availability is estimated for all animals as a group (Design II in Thomas and Taylor 2006), or when individuals are combined (pooled) in the analysis of data, an estimate of variation in individual selection potentially is lost. Despite wide use of RSFs, there continues to be debate about their appropriateness (i.e., Keating and Cherry 2004, but see Johnson et al. 2006). Thomas and Taylor (2006) identi- fied several problems with researchers not meeting the statistical assumptions of selection modeling, but the main concerns include pool- ing of unequal sample sizes when individual animals are combined into global models, and a variety of issues around the estimation of unused points. Some of the statistical concerns regarding selection models are difficult for researchers to accommodate. For example, the inherent attributes of telemetry (mortal- ity, premature collar failure, fix-acquisition bias) are problematic because researchers should exclude animals or fixes in analyses (i.e., throw away data) in order to balance a sample design before building global models. Until new statistical techniques emerge in the literature, however, RSFs will continue to be used by wildlife biologists. Understanding selection of habitat attributes, in addition to measures of habitat use only, allows for a better understanding of the relative values of specific habitats in different landscapes. While studying resource selection by moose in northern British Columbia, we ob- served that many individual collared animals completely avoided specific habitats. This avoidance may have occurred because of low animal densities or juxtaposition of particu- lar habitats within seasonal ranges, but also because of avoidance of specific attributes such as low forage availability or predation risk (Gillingham and Parker 2008). The avoided habitats differed among individual moose such that when taken as a group, all available habitats were used by some collared animals in all seasons. In statistical model- ing of resource selection, habitats that are completely (or nearly always) avoided must be dropped from individual resource selec- tion models because of issues of complete, or near-complete separation (i.e., no or very low use of some levels of categorical variables). Simpler analytical techniques such as selec- tion ratios (Manly et al. 2002) are not subject to the same constraints, but they cannot deal with continuous variables such as elevation and distance to specific features. To highlight the potential importance of individual variation in moose behaviour and selection, we constructed both individual and global RSF models and contrasted the results. As a baseline for habitat selection, we also examined individual and pooled selection ratios for collared moose in different seasons. We suggest that some of the biological infor- mation that may be lost in pooling animals warrants equal consideration with some of the statistical arguments that apply to the study of resource selection. STUDY AREA The study area was located between 57°11’ and 57°15’ N, and 121°51’ and 124°31’ W, south of the Prophet River and including the Besa River, within the Muskwa Ranges and Rocky Mountain Foothills. It covered a total ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 9 area of ~740,887 ha within the Muskwa-Ke- chika Management Area (MKMA) in northern British Columbia. The Besa-Prophet study area is characterized by numerous east-west drainages and south-facing slopes. The under- lying sedimentary rock formations are folded and faulted, and as is common along the east- ern slopes of the Rockies, potentially contain significant oil and gas reserves. At this time there is relatively little access into the Besa- Prophet region other than several permanent outfitter camps and 1 government-designated, all-terrain vehicle trail. The majority of human activity occurs during the summer and fall with the start of hunting seasons; some snowmobile activity occurs during winter. There are primarily 3 biogeoclimatic zones in the Besa-Prophet study area: boreal white and black spruce (Picea glauca and P. mariana) at lower elevations, spruce-willow- birch (Salix spp., Betula glandulosa) at mid- elevations (~1300-1600 m), and alpine tundra above ~1600 m (Meidinger and Pojar 1991). Valleys at ~800-1300 m are lined with white spruce, some lodgepole pine (Pinus contorta) and trembling aspen (Populus tremuloides) on dry sites, and black spruce, willow-birch communities on poorly drained sites. There also are slopes that have been burned by the British Columbia Ministry of Environment and local guide outfitters to enhance ungulate populations, primarily Stone’s sheep (Ovis dalli stonei). The spruce-willow-birch zone of the subalpine area is characterized by an abundance of willow and scrub birch, as well as some balsam fir (Abies lasiocarpa) and white spruce often in krummholz form, and various grasses, sedges, and fescues (Festuca spp.). Alpine areas consist of permanent snow- fields, rock, mat vegetation, and grasslands (Demarchi 1996). METHODS Fifteen adult female moose were fitted with GPS collars (GTX, Advanced Telemetry Systems, Isanti, MN) in March 2003. Collars were programmed to record locations every 6 hours for a 1-year sampling period. We defined 5 seasons that were distinguished by biological and ecological characteristics for our analy- ses of habitat selection by moose: winter (1 November–28 February) that corresponded to the formation of sex-specific groups follow- ing rut; late winter (1 March–15 May) when movement rates were lowest (Gillingham and Parker, unpublished); calving (16 May–15 June) when parturient females became solitary and the onset of plant greening occurred; sum- mer (16 June–15 August) from plant green-up through peak vegetation biomass to the start of plant senescence; and fall (16 August–31 October) when senescence of vegetation oc- curred, males and females formed mixed sex groups, and females came into estrus. RSF Model Inputs Vegetation Classification – The vegeta- tion classification system for the Besa-Prophet study area was developed using remote-sens- ing imagery and 227 field training sites (Lay 2005). Fifteen general vegetation associations were classified with a 2001 Landsat Enhanced Thematic Mapper (TM) image with 25-m resolution. We amalgamated several of these associations into 10 habitat classes to ensure that we had sufficient samples sizes for our analyses and an overall classification accuracy of >80%. Classes were lumped according to similarity of vegetation and elevation, and associations relevant to moose (Table 1). The 2 burn classes may also include some other small disturbed areas such as avalanche chutes, which could not be distinguished separately with remote-sensing imagery. We used the normalized difference veg- etation index (NDVI) derived from Landsat TM and Enhanced Thematic Mapper (ETM) to describe seasonal changes in vegetation (model described in Gustine et al. 2006a). The TM (4 June and 22 July 2001) and ETM (15 August 2001) images were used as a measure of vegetation biomass (June, July, and August) VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 10 and vegetation quality (change in NDVI from June to July and from July to August). We assumed that the vegetation classes and rela- tive differences in biomass and quality among classes were comparable among years in our study area. Other GIS Inputs – We obtained eleva- tion, slope, and aspect layers from a digital elevation model (DEM) in the 1:20,000 Brit- ish Columbia Terrain and Resource Inventory Management program (British Columbia Ministry of Crown Lands 1990). To minimize issues of perfect separation between used and available points, we modeled aspect as 2 continuous variables (i.e., northness and eastness; Roberts 1986); we did not assign an aspect to pixels with a slope ≤1°. Northness (the cosine of aspect) values of 1.00 and -1.00 suggest selection for north and south aspects, respectively, whereas values near 0.00 suggest selection for east and west aspects. Eastness (the sine of aspect) values show selection for east (i.e., 1.00) and west (i.e., -1.00) aspects; values of 0.00 show selection for northern/ southern exposures. We also defined potential risk of predation to moose using resource selection functions with logistic regression models by season from GPS-collared wolves (Canis lupus) and grizzly bears (Ursus arctos) in the Besa- Prophet area (details of predator models are in Gustine et al. 2006a, b, Walker et al. 2007). Grizzly bears and wolves are assumed to be the most significant, large mammalian predators in the Muskwa-Kechika Management Area (Bergerud and Elliott 1998). The predation- risk models included slope, aspect, elevation, habitat class, fragmentation (an index of vegetation diversity), and distance to linear features (e.g., seismic lines). We generated a risk surface to define which areas have the highest selection values for grizzly bears or wolves in each season by applying the coef- ficients from models to each 25 x 25-m pixel in the Besa-Prophet, based on its topographic and vegetation features. We scaled values from 0 to 1 to standardize selection surfaces, and then assumed that the risk of predation to moose from grizzly bears and wolves was directly related to selection values from the RSFs of those species. Habitat Class Description Non-vegetated Rock and rock habitats; permanent snowfields or glaciers and water bodies. Elymus burn Recently burned and open disturbed sites dominated by Elymus innovatus. Deciduous burn Older burned and disturbed areas containing Populus tremuloides and Populus balsamifera shrubs (<2 m) and trees (≥2 m). Subalpine Deciduous shrubs ≥1600 m in elevation; and spruce-shrub transition zone at middle to upper elevations with white and hybrid spruce (Picea glauca and P. glauca x engelmanni), and dominated by birch and willow. Stunted spruce Low productivity sites often on north-facing slopes with Picea glauca of limited tree height and percent cover. Pine-spruce White and hybrid spruce-dominated communities; and lodgepole pine- dominated communities. Riparian Low-elevation, wet areas with black (Picea mariana) and hybrid spruce; often with standing water in spring and summer; exposed gravel bars adjacent to rivers and creeks. Alpine Dry alpine tundra habitat ≥1600 m characterized by Dryas spp.; and wet alpine tundra habitat ≥1600 m dominated by Cassiope spp. and sedge (Carex spp.) meadows. Low shrub Deciduous shrubs <1600 m dominated by birch and willow. Carex Wetland meadows dominated by sedges (Carex spp.) at elevations <1600 m, with intermittent Salix shrubs. Table 1. Description of the 10 habitat classes used to describe habitat selection by moose in the Besa- Prophet area of northern British Columbia. ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 11 Determining Use and Availability GPS locations from telemetered moose were screened for fix quality (points with Positional Dilution of Precision >25 were dropped) and for improbable fixes (Spatial Viewer, unpublished Visual Basic program; M. P. Gillingham). To determine availability of resources for individual animals, we used all movement rates from consecutive 6-hour GPS fixes for each animal in a season, and de- termined the 95th percentile distance traveled during 6 hours. Our reasoning was that 95% of the time an animal typically moves within this movement potential (Arthur et al. 1996). The movement potential, therefore, generally represents how far an animal could have moved and the movements shorter than the potential represent choices that the animal made. The remaining 5% of movements included longer distances traveled during a 6-hour period within a season and were likely evoked by less common conditions (e.g., migratory movements and transitional movements prior to calving; Gustine et al. 2006a, 2006b). A circular buffer with the corresponding 95th percentile radius was then placed around each GPS location (used point) and we randomly selected 5 points from that area to represent availability. The circular buffer was defined by distance only and did not exclude physical constraints or barriers to movements such as cliffs. Nonetheless, we believe this is a better representation of what was available in the vicinity of the moose versus selecting points from a very large minimum convex polygon (MCP) of home-range size or a kernel based on density of use. To examine possible issues of lack of independence among animals (e.g., 2 or more collared animals spending large amount of time in close proximity), we calculated the minimum distance between every collared animal whenever a location was obtained. Finally, we checked to ensure that no 2 points were used twice and that there was no overlap between used and available points (Manly et al. 2002). We then used a raster GIS (Image- works XPACE; PCI Geomatics Corp. 2001) to query attributes in all GIS layers for used and random points. We did not consider the used and available points to be matched (i.e., case- controlled) because buffers for random points were selected based on seasonal movement potential and not the distance moved from the last fix (i.e., we used a Design III rather than a Design IV; Thomas and Taylor 2006). Resource Selection Modeling We developed 11 a priori, ecologically plausible models (Table 2) to describe resource selection in an information theoretic frame- work using Akaike’s Information Criterion (AIC; Burnham and Anderson 2002), and evaluated the relative importance of each of the variables in the models using selection coefficients (β) from logistic regression. We used tolerance scores (threshold of <0.20) to assess all model variables for collinearity, which can inflate selection coefficients and error terms (Menard 2002). The same suite of models was used for individual moose and all moose (pooled data), but not all models were run in all seasons (i.e., no risk of predation by bears during hibernation; some variables were dropped because of collinearity in models; Table 2). We ranked the suite of models using AIC values corrected for small sample size (AICc; Burnham and Anderson 2002) and then vali- dated all top models for individuals and pooled animals using k-fold, cross-validation (Boyce et al. 2002) and an averaged Spearman’s rank correlation coefficient. Within each model set (i.e., by individual and season), we calculated Akaike weights (wi), which are an estimate of the relative weight of evidence that the top model is the best within a model set. In cases for which there was not a single model with wi ≥0.95, we considered competing models until the sum of wi was ≥0.95 (Burnham and Anderson 2002). For each model set, the se- lection coefficients (β) in competing models VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 12 were averaged based on their relative wi. All statistical analyses were run in STATA (version 9.2; StataCorp 2006) and we used the add-in Desmat (Hendrickx 1999) for deviation coding of categorical variables. Logistic regression models do not provide reliable estimates if there is either complete or near-complete separation (few cases of presence or absence) in levels of categorical variables (Menard 2002). In our study, this occurred whenever individual moose com- pletely avoided an available habitat or used it very infrequently in a season. To avoid issues of separation, for each individual we dropped both used and available points in habitats for which either use or available points were rare (i.e., <5 points). Therefore, strong avoidance of an available habitat by individual moose is not reflected in the final individual resource selection functions if that avoidance was complete or near-complete. Estimates of variation around selection for both individual and global models were obtained directly from fitting logistic regressions. Selection Ratios Because we dropped several habitat classes for each moose in our RSF modeling, we also calculated selection ratios (Manly et al. 2002) so that rarely used habitats were not ignored in our analyses. We took the ratio of used (GPS locations) to available (random) points; available habitat types were divided by 5 (because we chose 5 random points per fix) before calculating each ratio. Each individual moose was treated as a replicate and selection ratios by habitat class were averaged across individual animals within each season. To compare selection ratios to the β coef- ficients from the individual RSF models, we estimated variation in selection ratios for each individual by bootstrapping the used and avail- able points. For each animal and season, we randomly selected 100 replicates; in each rep- Model Calving Summer Fall Winter Late Winter Elevation1+Aspect+Habitat2 Yes Yes Yes Yes Yes Elevation+Slope+Aspect+Habitat Yes Yes Yes Yes Yes Wolf 3+Habitat Yes Yes Yes Yes Yes Elevation+Slope+Aspect+Wolf+Bear3+Biomass+Habitat Yes Yes Yes Elevation+Slope+Aspect+Wolf+Bear+Quality+Habitat Yes Yes Aspect+Wolf+Bear+Biomass+Habitat Yes Yes Yes Aspect+Wolf+Bear+Quality+Habitat Yes Yes Elevation+Slope+Aspect+Wolf+Biomass+Habitat Yes Yes Yes Elevation+Slope+Aspect+Wolf+Quality+Habitat Yes Yes Aspect+Wolf+Biomass+Habitat Yes Yes Yes Aspect+Wolf+Quality+Habitat Yes Yes Table 2. Candidate models for all animals pooled (global) and individuals, developed a priori to de- scribe resource selection by moose by season in the Besa-Prophet area of northern British Columbia. Vegetation biomass for calving, summer, and fall were based on NDVI values from June, July, and August, respectively. Vegetation quality, assessed by the change in NDVI between summer months, was used only in calving and summer models. No risk of predation by grizzly bears was included during hibernation (winter and late winter seasons). 1 Elevation was modeled as a quadratic with both a linear and squared term. 2 Habitat classes are described in Table 1. 3 Wolf and Bear represent risk of predation by wolves and grizzly bears, respectively; see text for details. ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 13 licate 80% of the available used points (along with their corresponding random available points) were sampled. These data were then used to estimate a within-animal and season variance for each selection ratio. RESULTS We obtained 14,534 GPS locations from 14 of the collared moose. The fix rates (i.e., the number of GPS fixes recorded as a percentage of the number of attempted GPS fixes) averaged 76.7 ± 0.03% (X ± SE, range among individuals = 56-90%). Animals were assumed to be independent. Average distances between individuals at any one point in time ranged from 8.6 km in late winter to 12.2 km in summer, with a maximum value of 54.5 km between individuals during fall. Closest loca- tions among collared animals also occurred in fall, but fewer than 1.3% of locations were within 250 m and <1.9% were within 500 m of another collared individual. There was relatively good agreement in the signs of significant selection coefficients for continuous variables when individual and global models were compared (Table 3). In no seasons, however, did the sign of the coef- ficients that were significant in the individual models correspond completely with significant attributes in the global model. For example, moose appeared to always select for mid- elevations (positive linear term and negative quadratic term) based on global models, but this was reversed for at least 1 individual animal during both the calving and winter seasons. Given that there were 14 individual models for each season (except in winter when n = 12), many of the significant global coef- ficients corresponded to similar selection in less than half of the individual models (Table 3). Individual variation in selection associated with the continuous variables also indicated that there were several parameters that were significant in some individual models that were not supported by the global models. For example, slope was important in 9 of the 14 individual summer models, but there was no significant selection for slope in the global summer model. Pooling of animals to build the global seasonal models resulted in RSF models in which all habitat classes were included (Fig. Parameter Calving Summer Fall Winter1 Late Winter + – All + – All + – All + – All + – All Elev (km) 7 1 7.44 10 0 30.15 7 0 24.07 6 2 4.12 7 0 8.75 Elev (km2) 1 7 -2.39 0 10 -10.7 0 7 -8.23 1 7 -1.52 0 7 -2.85 Slope 1 6 -0.03 1 8 0 8 -0.03 2 4 2 4 -0.01 Northness 0 1 0 1 1 0 0 1 Eastness 1 0 0 1 0 1 Wolf Risk 1 2 2 0 1.06 4 0 2 0 1 2 Bear Risk 1 1 -1.18 2 2 Biomass 0 1 1 3 -1.7 1 5 -0.76 Quality 2 1 1.99 3 0 0.45 Table 3: Comparison of significant selection coefficients by season for continuous variables from in- dividual and global (all animals) resource selection models for 14 female moose in the Besa-Prophet area of northern British Columbia. For each season, the number under the + indicates the number of individual final models that showed significant selection for that parameter; the number under the – indicates the number of individuals that significantly avoided that attribute. The significant β coefficients in each seasonal global model for all animals are shown under ‘All’. 1 Only 12 animals were used in the winter models because of collar failure. VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 14 1A). Even though all habitat classes were available to all moose, individual animals frequently avoided habitat classes such that those classes could not be included in indi- vidual models. In winter, 2 moose used all 10 habitat classes, but in all other seasons each of the collared moose completely avoided at least 1 habitat class. As many as 7 (calving and fall), 6 (late winter and summer), and 4 (winter) habitats were dropped from all analy- ses because of near-complete separation for individual animals. The global models were able to incorporate all habitat classes without issues of complete separation, but there were no ‘average’ individuals in our sample that exhibited this habitat selection. Because many of the individuals completely avoided the same habitat classes in a given season, inference in the global seasonal models for some habitat classes was actually based on as few as 1 - 3 individuals. The global RSF models suggested selec- tion for, and avoidance of, many more habitat classes than did the selection ratios (Fig. 1A and 1B), in part because of the small sample sizes for the selection ratios (number of individu- als). Most habitat classes that were identified as important by the selection ratios, however, were also important in the RSF models, with the exception of Subalpine and Riparian habitats during calving, and Stunted spruce, Pine-spruce, and Riparian habitats during fall (Fig. 1). There were no instances in which a habitat that was significantly selected in the RSF models or selection ratios was signifi- cantly avoided in the other. Although the bootstrapped estimates of individual selection ratios were only an approximation of selection, they were not nearly as affected by near-complete avoid- ance of specific habitat classes in particular seasons as RSF models were (i.e., selection ratios could be computed as long as there was at least some use of a habitat type; Table 4). Because of the large number of habitat classes that were dropped from the RSF modeling due to separation issues, the majority of individual RSF models for moose did not show strong selection or avoidance of any habitats (n = 14 models in all seasons except winter for which n = 12; Table 4: RSF β+ and β–). In contrast, the bootstrapped estimates of selection ratios provided much stronger evidence of selection and avoidance of most habitat classes in most seasons (although these estimates were not be- ing corrected for the continuous variables that also were incorporated in the RSF models). For example, selection ratios indicated that twice as many individual moose selected for Subalpine in summer and against Alpine in fall as when determined by RSF models. In the 0 1 2 3 4 5 -3 -2 -1 0 1 2 Calving Summer Fall Winter Late Winter ***** * * * * * * ** ** * * * *** * * * ** * * * * * ** * ** * * * *** ** * * * * * * * * * A B No n-v ege tate d Ely mu s B urn De cid uou s B urn Su bal pin e Stu nte d S pru ce Pin e-S pru ce Rip aria n Alp ine Low Sh rub Ca rex No n-v ege tate d Ely mu s B urn De cid uou s B urn Su bal pin e Stu nte d S pru ce Pin e-S pru ce Rip aria n Alp ine Low Sh rub Ca rex R S F S el ec tio n  C oe ffi ci en ts S el ec tio n R at io Fig. 1. Comparison of seasonal selection of habitat classes by female moose based on selection coef- ficients (β) of global resource selection models (A) and selection ratios (B) in the Besa-Prophet study area in northern British Columbia. Er- ror bars represent 1 SE; * indicate significant β coefficients (A; significantly different from 0) or selection ratios (B; significantly different from 1). ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 15 extreme, 10 of 14 moose selected for Carex associations during calving when analyzed by selection ratios, yet no significant selection by any individual occurred in the RSF models. DISCUSSION Resource selection functions are a power- ful tool for incorporating both continuous and categorical variables in studies of resource selection and they are widely used for many wildlife species. Our results, however, sug- gest that care should be taken in interpret- ing global (cross-animal) RSF models even when all of the statistical assumptions of the analyses are met. In particular, although most models are developed to make population- level inferences, the variation in individual selection may be important to researchers as well as to outcomes of global models. The work of Thomas and Taylor (1990) focused considerable attention on examining indi- vidual variation in selection (Design III) and led to the development of new methods of modeling resource selection (e.g., Manly et al. 2002), but many studies inherently avoid examining individual variation by building global models. Recently, new analytical approaches allow for individual effects. Thomas et al. (2006), for example, presented a Bayesian random-effects model to assess resource selection. In addition, approaches for explic- itly accounting for individual animals have been developed in discreet-choice models (e.g., Buskirk and Millspaugh 2006). These models provide simultaneous estimation of both individual- and population-level selec- tion. Individual effects can also be included in more traditional RSF approaches (Gillies et al. 2006). There is a difference, however, between accounting for individual effects in models and examining important individual Habitat Calving Summer Fall Winter1 Late Winter RSF SR RSF SR RSF SR RSF SR RSF SR β+ β– R+ R– β+ β– R+ R– β+ β– R+ R– β+ β– R+ R– β+ β– R+ R– Non- vegetated 0 1 0 3 0 4 0 6 1 3 1 6 0 4 0 8 0 2 0 7 Elymus burn 1 1 7 2 2 1 4 7 4 0 3 7 1 1 3 8 2 0 4 5 Deciduous burn 3 1 11 2 5 0 10 2 8 0 8 3 6 1 7 4 2 0 8 5 Subalpine 3 3 5 1 10 3 11 0 12 0 9 0 10 2 2 0 5 7 Stunted spruce 1 0 7 4 5 0 4 9 0 1 1 10 1 0 4 6 6 1 7 2 Pine-spruce 2 2 5 9 5 0 9 5 1 3 2 10 2 6 3 8 2 1 5 8 Riparian 1 1 3 9 3 1 3 7 0 2 3 7 1 0 3 3 1 0 3 9 Alpine 1 0 0 7 0 5 0 10 0 4 0 10 0 1 1 5 Low shrub 1 0 7 5 3 4 10 3 9 0 12 0 6 0 9 3 2 0 6 5 Carex 10 3 1 1 6 2 6 1 1 0 3 3 2 4 Table 4. Comparison of significant categorical habitat classes using selection coefficients (β) from individual resource selection (RSF) models and individual bootstrapped selection ratios (SR; see text) for 14 female moose in the Besa-Prophet area of northern British Columbia. For each season, the number under the β+ indicates the number of individual final RSF models that showed selection for that parameter; the number under the β– indicates the number of individuals that avoided that attribute in RSF models. Similarly, R+ and R– correspond to the number of individuals that showed significant selection ratios for and against the habitat class, respectively. 1 Only 12 animals were used in the winter models because of collar failure. VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 16 differences in selection. We argue that the challenge (both statistical and methodologi- cal) is to understand what individual animals select and avoid. Although techniques are available for analyzing presence-only data (see Pearce and Boyce 2006), we cannot examine the choices animals make without comparing used and available points. From a management perspective, there is a demand for global RSF models. In prin- ciple, global models derived from appropri- ately pooling data (see review of Thomas and Taylor 2006) from individuals should provide the average response of the population if the sample of individuals is representative of the selection strategies within the population. In addition, many studies of resource selection that use RSF models result in resource selec- tion probability functions (RSPF; Manly et al. 2002) that are mapped as surfaces in a GIS. Despite some problems interpreting these surfaces (e.g., Keating and Cherry 2004), they can spatially and concisely depict the results of global models, but unfortunately without individual variation. In our study, there was generally good agreement between global and individual models, although there were some important differences. In order to meet statistical as- sumptions, data were dropped for categorical (habitat class) variables, and in some cases, these necessary statistical procedures may have affected the biological interpretation of selection by moose on the landscape. Only the Pine-spruce habitat was used by all individu- als in all seasons. When several individuals almost completely avoided a particular habitat class in a given season, inference in the global seasonal models was based on a few individu- als and the resultant selection coefficients were sometimes misleading. For example, during late winter, use of Carex habitats by 11 of 14 female moose was so limited that Carex could not be included in individual models, even though it was almost completely avoided by those individuals. When all animal locations and habitat availability were incorporated, however, the global model for this season indi- cated that moose selected for Carex because of the behaviour of a few individuals. The same situation occurred in fall when only 4 moose used Carex habitats enough to be included in individual models, but again the global model indicated selection for Carex. When comparing all animal-use locations to all available locations without regard to sample sizes per individual, the global model appears to be biologically misleading. In all, we recorded 5 instances for which at least half the moose avoided a habitat class so exten- sively that it had to be dropped from individual models, but the global model indicated selec- tion for that class. In addition, there were 5 instances in which habitats were dropped for more than half of the moose, but for which selection coefficients in the global model did not indicate avoidance. There were also 6 season-habitat models for which more than half of the animals completely avoided the habitat, and the ‘correct’ conclusion was drawn from the RSF models because the remaining animals showed significant avoidance of those habitats even though that inference was being drawn from as few as 1 individual. It is important to note that selection by individuals reflects the choices that they have to draw from, which may include different seasonal range configu- rations separating foraging from resting areas, for example, or individual demands related to physiological condition. These differences among animals are real and, therefore, it is important to understand the variability within the population. Selection ratios can be calculated from a single use point in a habitat class, unlike RSF selection coefficients, and seemingly may pro- vide a better measure of habitat selection. They are also much simpler to calculate. Selection ratios and similar selection indices, however, can only accommodate use and avoidance of categorical variables (Alldredge and Gris- wold 2006), and the influence of continuous ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 17 variables is not incorporated. Theoretically, continuous variables such as elevation or slope could be partitioned into categories and then combinations of multiple categorical variables using selection ratios or similar indices could be assessed, but the interpretation of those results would be difficult. Therefore, in most cases researchers are probably dependent on analytical techniques that result in the exclusion of little-used habitat classes when examining the influence of both continuous and categorical variables on resource selection. In those instances when continuous variables such as slope, elevation, and aspect do not vary much on the landscape, selection ratios could be effective in quantifying selection of rarely used habitat classes. We are not advocating the substitution of selection ratios for RSF modeling. We have used these ratios to demonstrate what infor- mation was lost by RSF models alone in this study. If managers and biologists are interpret- ing RSF coefficients as the average response of animals in a population, they should also examine use and availability to ensure that the responses of all individuals are reflected in the coefficients and their measures of variation. Categorical variables (e.g., habitat class) that are never or rarely used by individuals should be reported as measures of avoidance if those resources are available. The purpose of our analyses was to show the extent and influence of individual variation in defining habitat selection by moose. These findings are from relatively few individuals (n = 14) from a relatively short period of time (1 year). They provide initial insights into habitat selection by moose in northern British Columbia, but more importantly, they appear to show a large range in variation among individuals. It is possible that with a GPS sampling rate of >4 fixes per day and over multiple years, fewer habitat classes would have been dropped in seasonal analyses, although the landscape of the Besa-Prophet area is spatially heterogeneous and animals could easily have used and had access to all habitat classes within our 1 year of measure- ments. Longer-term data sets would certainly lend themselves to more robust RSF analyses even though year effects might be introduced. Within-individual variation probably also oc- curred in response to reproductive status (e.g., with or without a calf) and age, both of which can influence habitat selection in relation to nutritional demands and predation risk (e.g., Bowyer et al. 1999). Unfortunately, we did not have exact age and reproductive information for the GPS-collared animals in this study. Nonetheless, if the individual variation that we observed within these few animals was representative of many more animals at the population level, knowledge of the different selection strategies should be important to wildlife biologists and managers. Our findings indicate that these strategies may be masked using global selection models. In this study, the acquisition rates of the collars were low and may reflect features of the terrain (D’Eon et al. 2002, Frair et al. 2004), vegetation (Rempel et al. 1995, D’Eon et al. 2002), or leaf cover (Dussault et al. 1999, D’Eon 2003). The habitat class that would likely have the poorest signal attenuation would be the Pine-spruce habitat, but it was the only habitat in which use was recorded for all animals. We do not know which habitat classes may have been under-represented in our samples. By using stationary collars we potentially could have developed corrections for each habitat class (Frair et al. 2004), but even the behaviour of individual animals has been shown to influence fix rates (D’Eon 2003, Graves and Waller 2006). Nonetheless, be- cause we used the same data, our comparisons of individual and global models and between RSF models and selection ratios per se would not be affected by fix acquisition biases. Unless habitats are categorized into very broad types and all habitats are used exten- sively by all individuals, issues of complete, or near-complete separation will continue VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 18 to pose statistical problems when analyzing resource selection data. When using RSF models, researchers should be explicit about which habitat classes (or levels of any categori- cal variables) are dropped during analyses so that this information (i.e., total avoidance of specific attributes) is made known in addi- tion to model results. If we had not set out to examine individual variation in selection by moose, we would not have realized that inference about habitat selection was based on very few individuals for several habitat classes in different seasons. Therefore, we believe that caution must be taken when pool- ing individuals (in addition to stated statistical limitations) not only because of the potential loss of important individual variation, but also because all animals probably do not exhibit the average responses predicted by global models. Researchers with access to long-term data sets with numerous individuals and a high frequency of sampling in heteroge- neous environments should attempt to define whether incorporating more individuals in a global model can ever encompass the range of individual strategies for a given population, or whether knowledge of the different selection strategies within a population is more impor- tant to effective management of habitats. ACKNOWLEDGEMENTS Support for this project was provided by the Muskwa-Kechika Trust Fund, the British Columbia Ministry of Environment and the University of Northern British Columbia. We extend particular thanks to G. Williams who helped familiarize us with the Besa-Prophet landscape. R. B. Woods of the British Co- lumbia Ministry of Environment captured and collared all of the animals monitored in this study. We are also grateful to those who assisted in data visualization and analyses (J. B. Ayotte, S. G. Emmons, D. D. Gustine, R. J. Lay, B. Milakovic, A. B. D. Walker), and to external reviewers who provided very con- structive comments on the manuscript. REFERENCES Alldredge, J. r., and J. griswold. 2006. Design and analysis of resource selection studies for categorical resource variables. Journal of Wildlife Management 70:337- 346. Arthur, s. M., B. J. MAnly, l. l. McdonAld, and g. w. gArner. 1996. Assessing habi- tat selection when availability changes. Ecology 77:215-227. Bergerud, A. t., and J. P. elliott. 1998. Wolf predation in a multiple-ungulate system in northern British Columbia. Canadian Journal of Zoology. 76:1551-1569. Bowyer, r. t., V. VAn BAllenBerghe, J. g. Kie, and J. A. K. MAier. 1999. Birth- site selection by Alaskan moose: ma- ternal strategies for coping with a risky environment. Journal of Mammalogy 80:1070-1083. Boyce, M. s., and l. l. McdonAld. 1999. Relating populations to habitats using resource selection functions. Trends in Ecology and Evolution 14:268-272. _____, P. r. Vernier, s. e. nielsen, and F. K. schMiegelow. 2002. Evaluating resource selection functions. Ecological Modelling 157:281-300. British coluMBiA Ministry oF crown lAnds. 1990. Terrain and resource inventory management in British Columbia speci- fications and guidelines for geomatics: digital baseline mapping at 1:20,000. Ministry of Crown Lands for the Govern- ment of British Columbia, Victoria, British Columbia, Canada. BurnhAM, K. P., and d. r. Anderson. 2002. Model Selection and Multi-model Infer- ence: A Practical Information–theoretic Approach. Second edition. Springer- Verlag, New York, New York, U.S.A. BusKirK, s. w., and J. J. MillsPAugh. 2006. Metrics for studies of resource selec- tion. Journal of Wildlife Management 70:358-366. d’eon, r. g. 2003. Effects of a stationary ALCES VOL. 44, 2008 GILLINGHAM AND PARKER - VARIATION IN HABITAT SELECTION 19 GPS fix-rate bias on habitat selection analysis. Journal of Wildlife Management 67:858-863. _____, r. serrouyA, g. sMith, and c. o. KochAnny. 2002. GPS radiotelemetry error and bias in mountainous terrain. Wildlife Society Bulletin 30:430-439. deMArchi, d. A. 1996. Introduction to the ecoregions of British Columbia. British Columbia Wildlife Branch, Ministry of Environment, Lands and Parks, Victoria, British Columbia, Canada. dussAult, c. r., r. courtois, J. ouellet, and J. huot. 1999. Evaluation of GPS telemetry collar performance for habitat studies in the boreal forest. Wildlife So- ciety Bulletin 27:965-972. _____, J. ouellet, r. courtois, J. huot, l. Breton, and h. Jolicoeur. 2005. Linking moose habitat selection to limiting factors. Ecography 28:619-628. FrAir, J. l., s. e. nielsen, e. h. Merrill, s. r. lele, M s. Boyce, r. h. M. Munro, g. B. stenhouse, and h. l. Beyer. 2004. Removing GPS collar bias in habitat selec- tion studies. Journal of Applied Ecology 41:201-212. gillies, c. s., M. heBBlewhite, s. e. nielsen, M. A. KrAwchuK, c. l. Aldridge, J. l. FrAir, d. J. sAher, c. e. steVens, and c. l. Jerde. 2006. Application of random effects to the study of resource selection by animals. Journal of Animal Ecology 75:887-898. gillinghAM, M. P., and K. L. PArKer. 2008. Differential habitat selection by moose and elk in the Besa-Prophet area of northern British Columbia. Alces 44: 41-63. grAVes t. A., and J. s. wAller. 2006. Under- standing the causes of missed global posi- tioning system telemetry fixes. Journal of Wildlife Management 70:844-851. gustine, d. d., K. l. PArKer, r. J. lAy, M. P. gillinghAM, and D. C. heArd. 2006a. Calf survival of woodland caribou in a multi-predator ecosystem. Wildlife Monograph 165. _____, _____, _____, _____, and _____. 2006b. Interpreting resource selection at different scales for woodland caribou in winter. Journal of Wildlife Management 70:1601-1614. hendricKx, J. 1999. Using categorical variables in STATA. STATA Technical Bulletin 52:2-8. Johnson, c. J., K. l. PArKer, d. c. heArd, and M. P. gillinghAM. 2002. A multiscale behavioral approach to understanding the movements of woodland caribou. Eco- logical Applications 12:1840-1860. _____, s. e. nielsen, e. h. Merrill, t. l. McdonAld, and M. s. Boyce. 2006. Resource selection functions based on use-availability data: theoretical motiva- tion and evaluation methods. Journal of Wildlife Management 70:347-357. KeAting, K. A., and s. cherry. 2004. Use and interpretation of logistic regression in habitat selection studies. Journal of Wildlife Management 68:774-789. lAy, r. J. 2005. Use of Landsat TM and ETM+ to describe intra-season change in vegetation with consideration for wildlife management. M.Sc. Thesis, University of Northern British Columbia, Prince George, British Columbia, Canada. MAnly, B. F., l. l. McdonAld, and d. l. thoMAs. 2002. Resource Selection by Animals: Statistical Design and Analy- sis for Field Studies. Second edition. Chapman-Hall, London, U.K. Meidinger, d., and J. PoJAr. 1991. Ecosystems of British Columbia. British Columbia Ministry of Forests, Victoria, British Columbia, Canada. MenArd, s. 2002. Applied Logistic Regres- sion Analysis. Second edition. Sage, Thousand Oaks, California, U.S.A. neu, c. w., c. r. Byers, and J. M. PeeK. 1974. A technique for analysis of utilization- availability data. Journal of Wildlife Management 38:541-545. VARIATION IN HABITAT SELECTION – GILLINGHAM AND PARKER ALCES VOL. 44, 2008 20 nielsen, s. e., M. s. Boyce, g. B. stenhouse, and r. h. M. Munro. 2002. Modeling grizzly bear habitats in the Yellowhead ecosystem of Alberta: taking autocorrela- tion seriously. Ursus 13:45-56. niKulA, A., s. heiKKinen, and e. helle. 2004. Habitat selection of adult moose Alces al- ces at two spatial scales in central Finland. Wildlife Biology 10:121-135. osKo, t. J., M. n. hiltz, r. J. hudson, and s. M. wAsel. 2004. Moose habitat preferences in response to changing avail- ability. Journal of Wildlife Management 68:576-584. Pci geoMAtics corPorAtion. 2001. PCI Works Version 7.0. Richmond Hill, On- tario, Canada. PeArce, J. l., and M. s. Boyce. 2006. Mod- elling distribution and abundance with presence-only data. Journal of Applied Ecology 43:405-412. Poole, K. g., r. serrouyA, and K. stuArt- sMith. 2007. Moose calving strategies in interior montane ecosystems. Journal of Mammalogy 88:139-150. _____, and K. stuArt-sMith. 2005. Fine- scale winter habitat selection by moose in interior montane forests. Alces 41:1-8. _____, and _____. 2006. Winter habitat selec- tion by female moose in western interior montane forests. Canadian Journal of Zoology 84:1823-1832. reMPel, r. s., A. r. rodgers, and K. F. ABrAhAM. 1995. Performance of a GPS animal location system under boreal forest canopy. Journal of Wildlife Management 59:543-551. roBerts, d. w. 1986. Ordination on the basis of fuzzy set theory. Vegetatio 66:123-131. stAtAcorP. 2007. STATA Version 9.1. Col- lege Station, Texas, U.S.A. thoMAs, d. l., d. Johnson, and B. griF- Fith. 2006. A Bayesian random effects discrete-choice model for resource se- lection: population-level selection infer- ence. Journal of Wildlife Management 70:404-412. _____, and e. J. tAylor. 1990. Study designs and tests for comparing resource use and availability. Journal of Wildlife Manage- ment 54:322-330. _____, and _____. 2006. Study designs and tests for comparing resource use and avail- ability II. Journal of Wildlife Management 70:324-336. wAlKer, A. B. d, K. l. PArKer, M. P. gill- inghAM, d. d. gustine, and r. J. lAy. 2007. Habitat selection and movements of Stone’s sheep in relation to vegeta- tion, topography and risk of predation. Ecoscience 14:55-70.