UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Habitat suitability modeling from empirical data : application to mule deer in the interior of British.. Simons, V. Brock 2005

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
[if-you-see-this-DO-NOT-CLICK]
ubc_2005-0638.pdf [ 10.98MB ]
Metadata
JSON: 1.0075008.json
JSON-LD: 1.0075008+ld.json
RDF/XML (Pretty): 1.0075008.xml
RDF/JSON: 1.0075008+rdf.json
Turtle: 1.0075008+rdf-turtle.txt
N-Triples: 1.0075008+rdf-ntriples.txt
Original Record: 1.0075008 +original-record.json
Full Text
1.0075008.txt
Citation
1.0075008.ris

Full Text

H A B I T A T SUITABILITY M O D E L I N G F R O M E M P I R I C A L D A T A : APPLICATION TO M U L E D E E R IN THE INTERIOR OF BRITISH C O L U M B I A  by  V . B R O C K SIMONS  B.Sc, University of British Columbia, 2000  A THESIS SUBMITTED IN P A R T I A L F U L F I L M E N T OF THE REQUIREMENTS FOR THE D E G R E E OF  M A S T E R OF SCIENCE  In  THE F A C U L T Y OF G R A D U A T E STUDIES (Forestry)  THE UNIVERSITY OF BRITISH C O L U M B I A April 2005 © Brock Simons, 2005  Abstract Habitat suitability modeling has both strengths and weaknesses as a land management tool. Its utility is highly dependent on the ecological interactions and spatial and temporal scales that are pertinent to land management concerns and the species of interest. To maximize the usefulness of mathematical habitat suitability models, it is important that they are constructed using all the reliable a priori information available, and selected using a method that consistently selects models of an appropriate level of complexity. Application is exemplified here to produce winter and summer habitat suitability models for mule deer (Odocoileus  hemionus) in young,  intensively managed lodgepole pine (Pinus contorta) stands in the interior of British Columbia. The building of all models likely to have good explanatory power was informed by a comprehensive literature review of mule deer habitat requirements. After models were built, multivariate correlations between predictor variables and the dependent variable of standardized pellet-group densities were analyzed to ensure that no strong and sensible relationships suggested by the data were left out of the model set. Akaike's Information Criterion (AIC) was used for model selection, as it is currently the best readily available model selection criterion when 'truth' is of near-infinite complexity. To improve robustness of inference and prediction error estimates, final models are produced as AIC weighted averages of the models most strongly supported by the data. Although models should ideally be validated using independent data, error was estimated here based on the same data used for model fitting.  ii  Table of Contents ABSTRACT  II  TABLE OF CONTENTS  Ill  LIST OF T A B L E S  V  LIST OF FIGURES  VII  ACKNOWLEDGMENTS  :  IX  CHAPTER 1 - INRODUCTION 1.1  THE M A N A G E M E N T OF FORESTS  1.2  M O D E L I N G F O R E C O S Y S T E MM A N A G E M E N T  1 1 :  ..2  1.2.1. Habitat Suitability Models 1.3  3  M A N A G I N G FOR M U L E DEER  5  1.3.1 Forest Management Effects on Mule Deer 1.3.2 Mule Deer Effects on Forest Management 1.4  7 8 8  THESIS O V E R V I E W  CHAPTER 2 - THE HABITAT REQUIREMENTS DEER: A REVIEW  OF ROCKY MOUNTAIN  MULE 10  2.1  INTRODUCTION  10  2.2  METHODS USED FOR DETERMINING HABITAT PREFERENCE  11  2.3  HABITAT REQUIREMENTS: GENERAL OVERVIEW  12  :  2.3.1 Cover 2.3.2 Forage 2.4  WINTER RANGE  15  2.4.1 Winter Cover 2.4.2 Winter Forage 2.5  17 18  S U M M E R RANGE  20  2.5.1 Summer Cover 2.5.2 Summer Forage 2.6  20 20  S U M M A R YA N D CONCLUSION  22  C H A P T E R 3 - HABITAT SUITABILITY M O D E L I N G F R O M D A T A USING MULTIMODEL SELECTION 3.1  12 14  INTRODUCTION  25 25  3.1.1 Model Construction 27 3.1.2 Coefficient Estimation 28 3.1.3 Model Selection 29 3.1.4 Multimodel Selection and Inference 31 3.1.5 Validation 33 3.1.6 Habitat Suitability Modeling Using Multimodel Selection: Application to Mule Deer in the Interior of British Columbia 34 3.2  METHODS  35  3.2.1 The Data 3.2.2 Multimodel Construction and Inference 3.3  35 42  RESULTS..!  49  3.3.1 Multimodel Selection and Inference 3.3.2. Data Mining with Canonical Correlation Analysis  49 55  3.4  DISCUSSION  60  3.5  CONCLUSION  68  CHAPTER 4 - GENERAL CONCLUSION 4.1  70  FORESTS, COMPLEXITY, A N DM O D E L I N G  70  iii  4.2  THESIS OVERVIEW  71  71  4.2.1 Literature Review of Mule Deer Habitat Requirements 4.2.2 Modeling and Model Results 4.3  71  APPLICATIONS A N D FUTURE W O R K  73  73 74 75  4.3.1 Application Without Simulation Over Time 4.3.2 Stand-Level Simulation Models 4.3.3 Landscape-Level Simulation Models 4.4  CONCLUSION  76  REFERENCES CITED  77  A P P E N D I X 1: E X T E N D E D M O D E L A N D D A T A A N A L Y S I S R E S U L T S  96  A P P E N D I X 2: H A B I T A T S U I T A B I L I T Y C O M P O N E N T D E S I G N F O R F O R E C A S T .  100  A2.1  GENERAL STRUCTURE  A2.2  WILDLIFE HABITAT SUITABILITY INTERFACE  100  A2.3  WILDLIFE VARIABLE SELECTION  101  A2.4.  100  BUILD HABITAT SUITABILITY EQUATIONS  104  A P P E N D I X 3: U N I V A R I A T E R E S P O N S E S O F M O D E L V A R I A B L E S STANDARDIZED  TO  P E L L E T - G R O U P DENSITIES  A3.1.  G R A P H S O F UNIVARIATE RESPONSES  IN S U M M E R  A3.2.  GRAPHS O F UNIVARIATE RESPONSES  IN WINTER  A3.3.  SIMPLE CORRELATIONS  107 107 109 112  A P P E N D I X 4: G R A P H S O F P R E S S R E S I D U A L S V S . S T A N D A R D I Z E D GROUP DENSITIES  PELLET114  iv  List of Tables Summary of literature for important winter habitat variables and their probable order of importance  23  Summary of literature for important summer habitat variables and their probable order of importance  24  Details of replicate sites from a study on the response of wildlife to various thinning and fertilizer regimes in young lodgepole pine stands, conducted by T.P. Sullivan (unpublished data)  36  Means and standard deviations of pellet-group counts per site per season, separated by fertilizer treatment. Values are calculated from averages of 5 annually repeated observations  41  Table 3.3  Summary of AICc scoring results for summer data  50  Table 3.4  Summary of AICc scoring results for winter data  51  i aDie J . S  § f q d errors (SSE) and Root mean squared error (JMSE) of the fit between observed pellet-group densities and the values predicted by models derived as weighted averages of the full and reduced (for management applications) model sets. Also included are estimates of  Table 2.1  Table 2.2  Table 3.1  Table 3.2  Table 3.6  Table 3.7  Table A l . l  Table A1.2  Table A1.3  u  m 0  S  u a r e  prediction error using PRESS and root mean PRESS (IMPRESS)  53  Summary of canonical correlation analysis results for fertilized and unfertilized treatments combined  56  Summary of canonical correlation analysis results for fertilized and unfertilized analyzed separately  58  Full A I C output for summer data, including AIC weights and evidence ratios calculated for the full model set  96  Full A I C output for winter data, including AIC weights and evidence ratios calculated for the full model set  97  Final models. Listed values are weighted averages of parameter estimates for multivariate linear regression models of summer and winter data; for the full models created using all a priori selected variables (models cut off at a Aj value of 10), as well for as a reduced (management) model list composed of models with a A value less than .10 that do not include shrub or herb species richness  98  C  C  ;  v  T a b l e A1.4  T a b l e A3.1  T a b l e A3.2  Summary of canonical correlation analysis results for fertilized and unfertilized plots analyzed separately, by site  99  Simple Pearson's correlation coefficients between available variables for summer data  112  Simple Pearson's correlation coefficients between available variables for winter data  113  vi  List of Figures Map of British Columbia showing the location of the three study sites at Summerland, Kelowna, and Gavin Lake (Cariboo)  35  Figure 3.2  Aerial view of the Kelowna site  38  Figure 3.3  Segmenting of a sigmoidal response curve, a) A n approximately Type I exponential shape, b) Approximately linear segment, c) A n approximately Type III exponential shape  43  Figure 3.4  Exponential increase of AIC evidence ratios with increasing A; values....  49  Figure 3.5  Summation of AIC weights over all the models containing each variable, representing a subjectively influenced estimate of the relative importance of variables  54  Figure 3.1  Figure A2.1  Wildlife interface main page  Figure A2.2  The Wildlife Variable Selection component  Figure A2.3  The Build Habitat Suitability Equations component  104  Figure A2.4  Example of a traditional HSI component model  106  Figure A3.1  Standardized pellet-group density vs. shrub volume index for fertilized and unfertilized plots in summer  107  Standardized pellet-group density vs. herb volume index for fertilized and unfertilized plots in summer  107  Standardized pellet-group density vs. shrub species richness for fertilized and unfertilized plots in summer  108  Standardized pellet-group density vs. herb species richness for fertilized and unfertilized plots in summer  108  Standardized pellet-group density vs. stand density for fertilized and unfertilized plots in summer  109  Standardized pellet-group density vs. shrub volume index for fertilized and unfertilized plots in winter  109  Standardized pellet-group density vs. herb volume index for fertilized and unfertilized plots in winter  110  Figure A3.2  Figure A3.3  Figure A3.4  Figure A3.5  Figure A3.6  Figure A3.7  101  vii  102  Figure A3.8  Standardized pellet-group density vs. shrub species richness for fertilized and unfertilized plots in winter  110  Standardized pellet-group density vs. herb species richness for fertilized and unfertilized plots in winter  Ill  Standardized pellet-group density vs. stand density for fertilized and unfertilized plots in winter  Ill  Figure A4.1.  Residuals vs. standardized observations for the full summer model  114  Figure A4.2  Residuals vs. standardized observations for the summer management  Figure A3.9  Figure A3.10  (reduced) model  114  Figure A4.3  Residuals vs. standardized observations for the full winter model  115  Figure A4.4  Residuals vs. standardized observations for the winter management (reduced) model  115  viii  Acknowledgments  In addition to abundant intellectual and moral support, funding for this thesis was generously provided by Dr. Hamish Kimmins. I owe special and heartfelt thanks to Dr. Kimmins for giving me an opportunity to prove myself in a Masters program, when no one else would. I would like to thank my supervisory committee of Dr. Hamish Kimmins, Dr. Brad Seely, Dr. Thomas Sullivan, and Dr. Peter Arcese for the useful advice and comments that they provided throughout the research and writing process. I would like to especially thank Dr. Thomas Sullivan for generously providing the data used for model development. I would also like to thank Dr. Valerie Lemay, Dr. Christina Staudhammer, Rasmus Astrup, and Dr. David Huggard for their efforts and stoic patience in helping me come to terms with the difficult concepts integral to this thesis. Last but certainly not least I would like to thank my family for their unyielding support.  ix  Chapter 1 Introduction  1.1 The Management of Forests Forestry has been defined as the "science, art, and business of managing forests for human benefit" (Seymour and Hunter 1999). Recreation, food, water, aesthetic beauty, timber, medicine, biodiversity, and carbon dioxide regulation are but a fraction of the products and services produced by forests that are of value to humans. However, manipulating a forest to produce desired products is an incredibly difficult task, for two main reasons. First, the absolute and relative values of forest products are highly variable, changing with time, economic conditions, public demand, legislation, and ecological conditions (Thomas et al. 1979). Second, forest management requires planning over long periods of time for large, heterogeneous areas made up of complex ecosystems of which relatively little is known, and which are constantly changing (Christensen et al. 1996, Bunnell et al. 1999). These difficulties apply to management of the forest in general, and also equally well to management of biodiversity, and wildlife in particular. At one time, timber production was broadly viewed as the overwhelming priority of forest managers, while other values were viewed as constraints (Kessler et al. 1992, Seymour and Hunter 1999). In some circles, this paradigm likely still holds. Unfortunately, the most efficient way to produce timber is through intensive silviculture, which involves simplifying naturally diverse plant communities, and harvesting forest stands before they become senescent (Kuusipalo and Kangas 1994, Seymour and Hunter 1999). When biodiversity is seen as a constraint, timber production and biodiversity conservation appear to be in conflict. More and more, however, the importance of biological complexity and connectedness is becoming understood and appreciated, both for the intrinsic value of life, as well as for the provision of ecosystem services critical to human welfare (Christensen et al. 1996). As this understanding and appreciation has accrued in public opinion, the dominant paradigm of forest management has evolved into ecological forestry and ecosystem management (Kessler et al. 1992, Christensen et al. 1996, Seymour and Hunter 1999). Under the new paradigm of ecosystem management, the conservation of biodiversity is not a constraint to be minimized, but rather an integral consideration. 1  Managing an ecosystem is, needless to say, extremely challenging. Ecosystems are highly complicated and complex systems, exhibiting non-linear behaviour, self-organization, and emergent properties (Boyce 1992, Holling 1992, Levin 1998, Levin et al. 1998, Reynolds 2002, Wu and Marceau 2002). Every novel management activity, or a common management action applied in a novel environment, will produce results with a high degree of uncertainty (Walters and Holling 1990). As the process of science is currently the best tool available for learning about complex systems, it is important that the scientific method of testing hypotheses formed with the best knowledge available, monitoring results, quantifying error, and then re-evaluating hypotheses is integrated into the forest planning process (i.e. adaptive management; Walters 1986, Walters and Holling 1990). Management goals and strategies must then be responsive to the improvements in knowledge that is acquired (Christensen et al. 1996). In this way, forest ecosystem management can adapt effectively to the changes in values demanded by society, and minimize the risk of doing persistent environmental damage. Unfortunately, it appears that old paradigms die hard, and many land management decisions for biodiversity conservation continue to be made based on experience, and follow traditional practices (Pullin et al. 2004).  1.2 Modeling for Ecosystem Management Ecology has been defined as "the scientific study of the interactions that determine the distribution and abundance of organisms" (Krebs 1994). Those interactions are, for all practical purposes, of infinite dimension (Burnham and Anderson 2002). Understanding ecosystems therefore requires the aggregation and simplification of available knowledge, retaining what is essential and disregarding that which is not essential at the particular scale of interest (Levin 1992). The concept of scale is imperative to this discussion, because relationships at large-scale cannot be described simply by aggregating information about small scale relationships, and vice versa (Levin 1992, Hobbs 2003). Ecologists need tools to help them simplify and understand the natural world, and to help them make predictions on the effects of changes to it. Modeling is such a tool, and its advantages are so convincing that it is considered by some to be necessary to the process of producing feasible resource management hypotheses and highlighting the uncertainties within them (Walters 1986). The advantages of modeling are numerous and compelling. For one, modeling allows the investigation of ecological systems in ways that experimentation cannot because of practical, political, or financial limitations (Jackson et al. 2000). For another, modeling aids in achieving an understanding of complex processes, as the procedure requires that just enough detail is 2  included to approximate observed patterns (Levin 1992). The very act of scientifically rigorous modeling forces researchers to be explicit about assumptions and to consider processes and interactions that had previously been overlooked (Hilborn and Mangel 1997, Jackson et al. 2000). Also, models provide a structure for incorporating what is known with best guesses of what is not, at the same time highlighting gaps in knowledge (Hilborn and Mangel 1997, Jackson et al. 2000). In addition, models often give insight into systems, showing that it is combinations of variables acting together that control system behaviour rather than variables acting individually, thus cementing the importance of recognizing complexity for both model builders and users alike (Hilborn and Mangel 1997). Finally, through the use of computers, models give structure to the potentially convoluted calculations necessary to provide precise, quantified predictions and estimates of error.  1.2.1. Habitat Suitability Models One of the simplest and perhaps the most frequently used form of ecological model is the habitat suitability model, which is based on the concepts of habitat and carrying capacity (Schamberger and O'Neil 1986). Habitat has been defined as "the range of environments or communities over which a species occurs", and can often be effectively approximated by observing the response of species to a range of environmental variables (Whittaker et al. 1973). Carrying capacity refers to the maximum density of animals that a habitat can support (Krebs 1994, Morris and Davidson 2000). It is assumed that the measure of habitat suitability is proportional to carrying capacity. Habitat suitability models have been used extensively to predict the range of habitat variability that will sustain a particular species, and through that predict the potential impact of habitat alteration (Turner et al. 1995, Kliskey et al. 1999, Marzluff et al. 2002). These models assume that habitat is an important factor in deciding the presence and relative abundance of the species in question (Farmer et al. 1982). There is a strong theoretical foundation for this approach, as habitats are variable in their contributions to fitness, behaviour has a heritable component, and natural selection should therefore promote behaviour that steers animals into higher quality habitat (Krebs, 1994). It seems intuitive that higher quality habitat should receive greater use (Schamberger and O'Neil 1986). However, it is an often unrealistic simplification to assume that population densities are overwhelmingly controlled by habitat quality, and that animals are distributed according to an ideal free distribution (Fretwell and Lucas 1970).  3  Habitat suitability modeling has its limitations, as the presence and abundance of species are not always tightly coupled with particular habitat features on a particular site (Van Home 1983, Maurer 1986, Rotenberry 1986, Hobbs and Hanley 1990). Presence and abundance is a function of many factors that may operate independently of habitat quality, such as history, weather, disease, parasites, predators, and human harvest (Lancia et al. 1982, Schamberger and O'Neil 1986, Levin 1998). Also, relationships between population density and habitat variables are likely to vary with habitat type (e.g. Stauffer and Best 1986), scale of measurement (e.g. Hamel et al. 1986), and populations density itself (Hobbs and Hanley 1990). To make matters even more complicated, animals have occasionally been documented showing actual preference for low quality habitat over higher quality habitat, in a phenomenon known as an 'ecological trap' (Dwernchuk and Boag 1972, Donovan and Thompson III 2001, Battin 2004). The risk and influence of confounding factors can be minimized, though not eliminated, by basing specieshabitat relationships on the mechanisms linking demographic performance to habitat features, rather than on simple correlations (Hobbs and Hanley 1990). Habitat quality is a function of the population density, birth rate, death rate, and social interactions particular to an area (Van Home 1983). Data on population presence and abundance, or indices of density such as pellet-group counts, do not necessarily equate to measures of habitat suitability, but they are often all that is available (Hansen et al. 1993). The many other factors affecting the distribution and abundance of species, and habitat quality itself, are usually excluded from habitat suitability models because they are either difficult to measure, manage, or predict (Farmer et al. 1982, Schamberger and O'Neil 1986). Of course, removal of mechanisms critical to the distribution and abundance of species reduces the explanatory and predictive power of habitat suitability models, especially when applied to very different conditions or for long periods into the future (Schamberger and O'Neil 1986, Conroy et al. 1995). This simplification is necessary, however, to make model construction and application a practical possibility. Even i f the majority of critical factors were included in habitat suitability models, predictive power would still be limited. For one thing, less obvious factors may exert great influence, such as the occurrence of an infrequent event with strong influence, or numerous weak interacting factors combining to have large influences over system behaviour. Also, change in any ecosystem, like all complex adaptive systems, occurs mainly at small scales through random events such as mutation and environmental fluctuation (Levin 1998). As the system changes and evolves, the local rules of interactions also change (Levin 1998). This inherent unpredictability 4  of ecosystems makes the number of potential outcomes of any management action extremely high, and works to hamstring even the very best efforts at prediction. Formation of useful habitat suitability models is extremely challenging even with excellent data. However, i f the data are of poor quality, i f spatial scale of measurement was inappropriate, or i f the data were not collected with habitat quality characterization in mind, any attempt at robust habitat suitability model creation may be utterly futile (Guisan and Zimmerman 2000). Clearly, there are numerous risks and challenges with the use of habitat suitability models for forest management. It should be emphasized at this point that the previous discussion of the conditional weaknesses of habitat suitability modeling is not meant to discourage its use, but instead to carefully outline the considerations necessary to evaluate when they will be useful, and how to maximize their usefulness. "Models, of course, are never true, but fortunately it is only necessary that they be useful. For this it is usually needful only that they not be grossly wrong" (Box 1979). With the exception of exceedingly simple or impossibly complex situations, it seems logical that the process and products of mathematical habitat suitability modeling will be generally superior to an educated guess. In forest management, incorrect assumptions on the outcome of proposed management actions can have significant social, economic and ecological impacts. It is imperative that forest managers have the best information available for decision making, and habitat suitability modeling can be a powerful aide.  1.3 Managing for Mule Deer Mule deer (Odocoileus hemionus), an important game species, have a range that extends 2500 km in length, from central Arizona and New Mexico to northern Alberta and British Columbia, and 1500 km in width, from the east side of the Rocky Mountains west to the Coast range (Wallmo 1981). Mule deer reach the northern limits of their continuous distribution in the Cariboo Forest Region in the central interior of British Columbia (Armledder et al. 1994, 1998). The population trends of mule deer, when generalized over their entire range, show that populations have been in a state of decline since an overabundance peak that occurred through the 1940's to 1960's (Wallmo 1981, Alberta Forestry, Lands and Wildlife, Fish and Wildlife Division 1989). Mule deer populations in B.C. also appear to have followed this general trend (Edwards 1956). Periods of decline appear to have no consistent relationship with over-harvest, but there has been a general relationship with poor fawn survival (Wallmo 1981). Further  5  exploration of the possible reasons for this decline, as well as finding strategies to mitigate it, requires a basic understanding of mule deer biology and ecology (Loveless 1964). Throughout their distribution, mule deer are extremely variable in whether they migrate (Nicholson et al. 1997). In the interior of British Columbia, they do migrate, moving from winter to summer ranges to pursue nutritious forage and to escape the heat, and from summer to winter range to escape deep snow (Willms et al. 1976). Winter is generally accepted to be a harsh time for ungulates in temperate forests. In winter, forage is low in digestible energy and digestible components, such as protein, starches, sugars, and hemicelluloses (Wallmo et al. 1977, Hanley and McKendrick 1985). This is due largely to leaf fall of deciduous plants (Short et al. 1966), maturity of herbs (Willms et al. 1976), and snow covering low lying forage (Robinette et al. 1952). Not only is the intake of energy low in winter, but energy expenses are high, as snow makes  locomotion more difficult  (Parker et al. 1984), and low temperatures make  thermoregulation costly (Hobbs 1989). In late winter, available forage may be insufficient to meet maintenance energy requirements (Wallmo et al. 1977), at which point deer will draw upon stores of fat and protein (Torbit et al. 1985, Anderson et al. 1990). If energy stores are sufficiently depleted, fawns the following spring will be born with low vitality and a low rate of survival (Wallmo 1981). In more severe situations, extensive adult starvation may occur (Wallmo 1981). It seems clear that deer health and survival over winter is affected by snow on the ground and the quality and abundance of available forage. However, it has been suggested that winter ranges rarely meet maintenance requirements, and that management emphasis should therefore be spread more equitably through the year (Edge et al. 1990). Over-winter health and survival may have more to do with the body condition of deer entering winter than the condition of winter range (Hobbs 1989). Be that as it may, much of the focus on mule deer management has focussed on winter range requirements, due largely to the conspicuous relationship between mule deer population declines and winter severity (Wallmo 1981). As forest management activities have the potential to have considerable impact on deer habitat, mule deer winter range requirements have become necessary considerations for forest harvesting operations throughout the Southern and Central Interior of British Columbia (see the Ministry of Sustainable Resource Management's Cariboo-Chilcotin Land Use Plan 1995 and Okanagan-Shuswap Land and Resource Management Plan 2001).  6  1.3.1 Forest Management Effects on Mule Deer Forest management alters a forest by changing the amount and distribution of particular . age classes of forest stands, as well as the species composition of trees and understory plants within those stands. Economic incentives in sustained yield timber management are optimized by felling trees when they reach their maximum mean annual increment in monetary value (Kuusipalo and Kanga 1994). This provides an impetus that has tended to result in substantially reduced amounts of old-growth forest' and old forest attributes (i.e. large course woody debris, patchy canopy structure) on the landscape, and increased amounts of young forest and early serai associates (Kuusipalo and Kanga 1994). There is a strong possibility that the presence of oldgrowth forest improves mule deer winter habitat, as old-growth stands, particularly old-growth Douglas-fir (Pseudotsuga menziesii) stands, intercept large amounts of snow, allow forage production in canopy gaps, and provide additional forage through winter branch breakage (Armleder et al. 1994). Forest management also can affect stand density through thinning, as there is a well established trend that understory production, and therefore mule deer forage, generally increases with decreasing canopy cover (Jameson 1967). This trend is certainly not absolute, however, as in some instances stand thinning activities have not resulted in a change to understory volumes (Lindgren and Sullivan 2001, Sullivan et al. 2002). This is likely to be because thinning can improve overall stand vigour, releasing stagnant stands and causing the canopy to close quickly (Homyack et al. 2004). In addition, forest management can also affect the speed that new stands are re-established after harvest or natural disturbance by planting seedlings and removing competing vegetation. In truth, forest management can affect the forest, and mule deer habitat, in near countless ways. There is great variation in stand and landscape harvest patterns, regeneration techniques, and application of thinning, pruning, fertilization, slashburning, herbicide, and scarification treatments. Forest management, by its very nature, alters stand overstory characteristics, as well as understory vegetation productivity, abundance, and species composition, all of which are components of mule deer habitat.  The term 'old-growth', as used here, refers to a structural description of an over-mature forest stand. Overstory trees will be senescent, typically exhibiting large crowns, as well as signs of decay, disease, and deformation. The stand exhibits canopy gaps, as well as standing and fallen dead trees that are relatively large for the species present (Spies and Turner 1999). 1  7  1.3.2 Mule Deer Effects on Forest Management While the impact of forest management on mule deer is generally of primary concern to forest managers, the impact of mule deer on forests can be significant. A l l large mammalian herbivores alter vegetation and, potentially, vegetation communities (Huntly 1991, Mclnnes et al. 1992). Mule deer affect vegetation directly by consuming leaves, stems, fruits and flowers of palatable vegetation (Cote et al. 2004). Changing plant biomass and community composition may have a cascading affect on other species that are affected by their presence, abundance, and availability (Huntly 1991, Cote et al. 2004). Some mule deer impacts on forests and forest management are most apparent when deer become overabundant. For example, mule deer browsing on tree seedlings may over time affect overstory tree species composition and reduce timber yield, and mitigation measures (e.g. protective cones) may be expensive (Cote et al. 2004). Although outside the realm of traditional forest management, mule deer may also damage gardens, nearby agricultural crops, and, perhaps most importantly, impact with cars and thus endanger human life (Cote et al. 2004). Finally, deer may transmit disease to livestock, other wildlife, and humans (Cote et al. 2004). When deer are not overabundant, their impacts on forest management can still be considerable. For example, considerations for preserving mule deer winter range in the interior of British Columbia have effectively precluded harvest of numerous large, valuable stands of old-growth Douglas-fir. If attempts to integrate mule deer habitat needs into harvesting plans are not successful and population declines continue, large areas may be set aside as reserves, further limiting economic opportunities (Armledder and Dawson 1992).  1.4 Thesis Overview Habitat suitability models generally fall into one of two categories; qualitative and quantitative synthesis of expert opinion, knowledge and understanding, or statistical models derived from empirical data, and incorporating varying degrees of expert influence (Pearce et al. 2001). Quantitative models are generally preferred to qualitative ones, because weaknesses of quantitative models are more easily analyzed and evaluated (Maurer 1986). This thesis is concerned with the production of statistical mule deer habitat suitability models from empirical data, while at the same time exploring in detail the methods necessary to do so. The responses of mule deer pellet-group densities to experimental manipulations of thinning and fertilizer treatments in young lodgepole pine (Pinus contorta) stands in the interior of British Columbia provide the empirical data used. Chapter 2 (The Habitat Requirements of 8  Rocky Mountain Mule Deer: A Review) is a comprehensive review of mule deer habitat requirements. This review was necessary to identify which predictive variables to use in habitat suitability equations, as well as their likely ranking of importance. In Chapter 3 (Habitat Suitability Modeling from Data Using Multimodel Selection and Inference), habitat suitability models are developed for mule deer using an effective and increasingly popular approach. This approach, most commonly known as the information-theoretic approach and championed by Burnham and Anderson (2002), is based on Chamberlin's (1890) theory of multiple working hypotheses and Akaike's information criterion (AIC; Akaike 1973). AIC, in turn, is based on the theory of maximum likelihood and the Kullback-Leibler measure of information (Kullback and Leibler 1951). The intent is to apply this approach in such a way as to produce statistical habitat suitability models for mule deer in the interior of British Columbia that are as robust as the data will allow.  9  Chapter 2 The Habitat Requirements of Rocky Mountain Mule Deer: A Review  2.1 Introduction Habitat is defined as the range of physical and chemical environmental gradients that a given species occupies, and is described in nature by demographic responses to those gradients (Whittaker et al. 1973). As habitats are variable in space and time, natural selection will favour individuals that utilize the portions of their habitat that are most suitable for the successful production of progeny (Krebs 1994). This selection pressure has resulted in many species developing habitat preference behaviours to actively seek out habitat that maximizes net benefits to fitness (Parker and Robbins 1984). Behaviours are never a perfect reflection of natural selection pressures, however, and therefore preference does not necessarily equate with requirement (Krebs 1994). Preference may be shown for a habitat feature that does not contribute positively to fitness, or no preference may be apparent for a feature that is overabundant (Peek et al. 1982). As it is, behaviour is an imperfect indicator of preference, but even identifying specific behaviour can be challenging. Population density, commonly used as an index of preference, may become uncoupled from habitat quality, causing correlations between density and habitat to be potentially misleading (Van Home 1983). This may occur, for example, when territoriality or predator avoidance behaviours create high population densities in sub-optimal habitat. In addition, apparent preference is likely to vary with scale and time of year, making the meaningful measurement and interpretation of preference behaviour even more difficult (Krebs 1994, Apps etal. 2001). Although identifying and interpreting animal behaviour for determination of habitat preference and requirements is challenging, it is a worthwhile endeavour. Land management activities that operate in ignorance or indifference to the needs of wildlife seem more likely to cause significant habitat damage. Such an outcome is considered to be unacceptable, for example, in modem, first-world forestry policy and practice. The following literature review is a detailed look at what is known regarding the habitat preferences and requirements of mule deer, with an emphasis on the Rocky Mountain sub10  species and associated ecotypes that are found in the interior of British Columbia (Wallmo 1981). Methods of determining habitat preference are also reviewed, as familiarity with them is necessary for interpreting literature results, and fully understanding their inferential limitations. The purpose of this review is to have the information on hand necessary to make the most informed decisions possible on which predictor variables to utilize, from a selection of many available, for habitat suitability modeling of mule deer in the interior of British Columbia.  2.2 Methods Used for Determining Habitat Preference Numerous methods have been used to measure mule deer habitat preference. The most common method involves using pellet-group densities as an index of time spent in a habitat and of relative population densities (e.g. Robinette et al. 1952, Julander 1955, White 1960, Loveless 1964, Gilbert and Wallmo 1970, Lyon and Jensen 1980, Tomm et al. 1981, Wambolt and McNeal 1987, Altendorf et al. 2001). Habitat preferences have also been inferred by directly observing mule deer (e.g. Martinka 1968, Constan 1972), or by snow track surveys (e.g. D'Eon 2001). Data collected by tracking radiotelemetry-collared deer have been used to find correlations with habitat features (e.g. Carson and Peek 1987, Kie et al. 1991, Kie 1996, Coe et al. 2001), occasionally using tip-switches to identify habitat used for foraging (e.g. Kie et al. 1991, Kie 1996). Recently, location tracking using GPS collars has also been utilized (D'Eon and Serrouya 2005). Preferences for specific forages have also been established in numerous studies, using several methods. The most common methods employed have been the observation of bites taken by tame deer (e.g. Wallmo et al. 1972, Urness et al. 1975, Carpenter et al. 1979, Deschamp et al. 1979, Hobbs et al. 1983), and microhistological analysis of fecal material to determine forage composition (e.g. Hansen and Reid 1975, Uresk and Uresk 1982, Campbell and Johnson 1983, Hanley and McKendrick 1985, Waterhouse et al. 1994, Irby et al. 2002). Analysis of rumen samples (e.g. Martinka 1968, Willms et al. 1976) and observations of foraging signs (e.g. Constan 1972, Keay and Peek 1980) have also been employed. There are many potentially serious sources of error with all habitat and forage preference assessment methods. For example, increased defecation rate with increased forage digestibility and succulence makes comparisons of mule deer preference between habitats using pellet-group densities difficult (Longhurst 1954, Rogers et al. 1958, Wallmo 1981). Also, results from analysis of digested plant material have an unavoidable bias towards plants with low digestibility or distinctive morphological features (Kufeld et al. 1973). Plants with low digestibility are 11  represented in feces disproportionately to their contribution to total diet because more of their mass remains intact and they are easier to identify relative to more thoroughly digested plants. Plants with distinctive morphological features are also relatively easy to identify in fecal samples. Even taking sources of error into account, however, a sufficient number of studies have been performed to clearly establish general trends.  2.3 Habitat Requirements: General Overview Within the geographic range of mule deer, habitat requirements can be separated primarily into the interacting categories of cover and forage. Cover influences understory production and community structure, but also has important effects on mule deer habitat that are independent of its direct influence on forage production. For example, overstory trees act to intercept snow in winter, thereby reducing energy required for locomotion, and increasing access to forage (Parker et al. 1984, Bunnell et al. 1985). Other categories of habitat may be considered, but are unlikely to be of enough importance to be worthy of focussing on separately. For instance, coarse woody debris over 0.5m deep has been shown to greatly impede mobility (Lyon and Jensen 1980). Access to water is also important (e.g. Nicholson et al. 1997, Stewart et al. 2002), although it is unlikely to be directly limiting.  2.3.1 Cover Forest cover has both positive and negative repercussions for mule deer. The benefits of cover include security from predators, reduced energy costs (Peek et al. 1982, Armleder and Dawson 1992), and increased access to winter forage (Peek et al. 1982). The negative repercussions of cover.include reduced predator detectability, escape obstruction (Mysterud and Ostbye 1999), and reduced forage production (Peek et al. 1982). Cover is composed of the interacting effects of topography and vegetation (Loveless 1964), and can be separated into the interacting functional categories of snow interception cover, security cover, and thermal cover. Snow interception cover, which is of importance only in winter, will be reserved for discussion in the 'Winter Range' section. Security cover, on the other hand, is of importance to deer year-round. Security cover helps deer escape predation and reduce energy expenditures caused by fleeing from human disturbances (Armleder and Dawson 1992). Security cover has been defined as any combination of vegetation and topography capable of hiding 90% of a standing adult deer from view of a human approximately 60 m away (Thomas et al. 1979). Deer 12  have shown strong tendencies to remain close to cover when it is available, with avoidance behaviour beginning at between 40 and 80 m into the open, and use declining sharply 100 m from cover (Tomm et al. 1981). A relationship between stress due to harassment and use of security cover has an intuitive appeal, and is somewhat supported by the literature. For example, use of cover has been shown to remain unchanged in relation to intensity of exposure to human activities (Tomm et al. 1981), but has been shown to increase with focussed harassment of does by an all-terrain vehicle (Yarmoloy et al. 1988). Female mule deer show a particular preference for cover when fawning, probably to reduce the likelihood of insect attack, predation, and extreme weather effects on young (Peek et al. 1982). Thermal cover is sought out by deer to avoid stress due to temperature extremes in a process known as behavioural thermoregulation (Sargeant et al. 1994). Thermoregulation in general is the ability of an animal to maintain body temperature within acceptable limits despite large variations in ambient conditions (Bunnell et al. 1986, Parker and Gillingham 1990). The temperature experienced by an animal, or operative temperature, is affected by air temperature, conduction with the ground (Bunnell et al. 1986, Parker and Gillingham 1990, Sargeant et al. 1994), the cooling effects of wind, and the heating effects of incoming or reflected solar radiation (Parker and Gillingham 1990). Heat, cold and the interacting effect of humidity all elicit reactions in deer (Loveless 1964). Behavioural thermoregulation can involve any behaviour that affects an animal's operative temperature. Body posture, for example, can affect thermoregulation by altering exposure of surface area. Standing or lying stretched out will increase exposure to air flow and help with cooling, while lying and folding legs underneath the body will minimize air flow and retain heat (Parker and Gillingham 1990, Sargeant et al. 1994). Deer may also attempt to simply move to areas of more moderate temperatures. Thermal cover, although perhaps not necessary for maintaining high populations of mule deer (Peek et al. 1982), is often strongly preferred by them for relief from extreme high and low temperatures (Loveless 1964, Sargeant et al. 1994). Deer will seek shade and northerly aspects to cool themselves, while seeking exposure and southerly aspects for warmth (Loveless 1964, Sargeant et al. 1994). When wind causes operative temperatures to drop below lower critical limits, however, deer will leave exposed areas to seek shelter (Loveless 1964). Physiological thermoregulation is also important, and functions by expending metabolic energy when temperatures exceed critical limits (Parker and Gillingham 1990). In mule deer, this is exhibited primarily through panting to reduce heat and shivering to create heat, with cooling 13  by perspiration playing only a minor role (Parker and Robbins 1984, Bunnell et al. 1986). Critical limits will change with season, as Parker and Robbins (1984) found that shivering occurred below -20°C in winter and below 5°C in summer, while panting occurred when temperatures exceeded 2.5°C in winter and 23.5°C in summer. Mule deer produce a winter coat with excellent insulative properties (Parker and Gillingham 1990), and this is at least part of the reason for this seasonal difference.  2.3.2 Forage Presence and availability of forage, although affected by cover, is likely of dominant importance for mule deer habitat. Food is required for energy, and the quality and availability of food affect both birth and death rates (Nicholson et al. 1997). Mule deer have been shown to be highly selective feeders (Wallmo 1981, Collins and Urness 1983), generally preferring forages that have a high nutrient content and high digestibility, when available (Swift 1948, Willms et al. 1976). This selectivity is clearly to the species' advantage, as a restricted diet is essential for efficiency, and efficiency is necessary for competitive success (MacArthur 1955). While diets vary with geographic location, year and season (Morris and Schwartz 1957, Uresk and Uresk 1982, Campbell and Johnson 1983), mule deer generally feed primarily on browse (leaves and twigs of woody plants) and forbs with some component of grasses, sedges, rushes, and mushrooms (Kufeld et al. 1973, Wickstrom et al. 1984, Shackleton 1999). With regard to the twigs of browse, only the terminal current annual growth is considered forage (Hanley and McKendrick 1985). Nuts, fruits, and berries are also likely to be taken when available (Wallmo 1981). Mule deer select forages not only for particular quality and species, but also for variety. When available, a large number of species are generally consumed, particularly.with regard to forbs (Kufeld et al. 1973). This may be in part due to the fact that herbivores are forced to eat a variety of plant material to avoid toxification from excessive exposure to any one plant secondary defence compound (Freeland and Janzen 1974). A varied diet also seems more likely to provide the required amounts of trace minerals and nutrients. The quality of plants as forage varies between seasons and between species. Between seasons, the changes in forage quality follow the general growth patterns of dormancy in winter, re-growth in spring, succulence in summer and tissue hardening in autumn (Short et al. 1966). Succulence refers to the cell content to cell wall ratio, as well as general nutrient content (Short et al. 1966). From summer to winter, the digestible components (Wallmo et al. 1977, Hanley and 14  McKendrick 1985) and protein content (Willms et al. 1976) of forage generally decreases considerably. The difference is most acute in deciduous plants due to the annual shedding of leaves (Short et al. 1966), as well as in herbs and grasses due to the decrease in digestibility that comes with the thickening of the cell wall with age (Wallmo 1981, Spalinger et al. 1986). Also of potential importance is the increasing fat concentration in deciduous plants during winter due to the storage of lipids in branches during dormancy (Short et al. 1966). The changes in digestibility, in particular, can have a dramatic effect on mule deer, as decreased digestibility may increase forage handling time to the point that digestibility limits energy intake below maintenance requirements (Ammann et al. 1973, Wallmo et al. 1977, Torbit et al. 1985). Deer succumbing to starvation with full rumens have been observed frequently (Wallmo 1981). As digestibility increases beyond the point that there is a neutral gain in energy at maximum forage intake, less forage is required to meet energy needs and consumption rate is decreased (Ammann et al. 1973). Finally, as this is potentially of importance in forest management, it should be noted that fertilizer application may increase the digestibility and nutrient concentration of plant parts (Puoli et al. 1991, Johnson et al. 2001). For example, regarding Douglas fir foliage, fertilizer has been shown to increases the growth rate, crude protein level and preference by mule deer (Oh et al. 1970).  2.4 Winter Range It has been well established that mule deer generally exhibit a preference between habitat types within seasons (e.g. Loveless 1964, Wallmo et al. 1972, Deschamp et al. 1979, Nicholson et al. 1997, Altendorf et al. 2001), between seasons (e.g. Loveless 1964, Schoen and Kirchoff 1985, Nicholson et al. 1997, Plante et al. 2004), and between years (e.g. D'Eon 2001). Mule deer display a preference for winter range location that is largely determined by the interacting effects of snow depth, aspect, and topography (Gilbert at al. 1970). Relative to summer ranges, steeper slopes, more southerly exposures, and lower elevations are generally preferred in winter (e.g. Martinka 1968, Wallmo et al. 1977, Telfer 1978, Thomas et al. 1979, Wambolt and McNeal 1987, Nicholson et al. 1997, Shackleton 1999, D'Eon 2001, D'Eon and Serrouya 2005), because of the lesser snow depths in those areas (Loveless 1964, Thomas et al. 1979, Bunnell et al. 1985, D'Eon 2004). However, there is some discussion over whether steep slopes are in fact preferred (Armleder et al. 1994). Snow depth has been implicated in driving deer migration (Gilbert at al. 1970, Bartmann 1984, Schoen and Kirchoff 1985), limiting the amount of winter range available 15  (Wallmo et al. 1977, Schoen and Kirchoff 1985), and driving population dynamics (e.g. Edwards 1956, Wallmo 1981, Bartmann 1984, Parker et al. 1984, Picton 1984, Hanley and McKendrick 1985, Ballard etal. 2001). Snow has a dominant effect on mule deer habitat and habitat selection because snow covers low lying forage (Robinette et al. 1952, Loveless 1964, Martinka 1968, Willms et al. 1976, Hanley and McKendrick 1985) and inhibits locomotion (Robinette et al. 1952, Parker et al. 1984, Hanley and McKendrick 1985, Schoen and Kirchoff 1985, Bunnell et al. 1990, Shackleton 1999). Snow will tend to cover forbs and grasses first, increasing the shrub component of the diet in winter (Martinka 1968, Willms et al. 1976, Carpenter et al. 1979, Hanley and McKendrick 1985). This loss of herbaceous forage results in lower diet quality, most notably through lowered digestible energy (Hanley and McKendrick 1985). If enough snow accumulates, considerable proportions of shrub communities may become covered (Loveless 1964). As winter severity increases, preferred areas become those that receive the most sunlight, are driest in summer, and therefore support less vegetation than alternative areas (Wallmo 1981). These emergency winter refuges may be so small and so intensely used that vegetation has little opportunity to recover, and would be inadequate to prevent starvation over time regardless of population size (e.g. Gilbert etal. 1970). Snow not only affects energy intake, but also energy expenditures by making locomotion more difficult and costly (Parker et al. 1984, Hanley and McKendrick 1985). This is likely to be of greatest importance in situations where frequent encounters with people cause deer to attempt to flee (Parker et al. 1984), and where snow impedes escape from predators (Shackleton 1999). However, the increased energy cost of foraging may also take a heavy toll. As the depth that a deer sinks into snow increases, energy costs of locomotion increase exponentially (Parker et al. 1984). Sinking depth is a function of not only snow depth, but also snow density and surface hardness or supportability (Bunnell et al. 1990). Deer will tend to avoid areas where they sink approximately halfway to their chest (Shackleton 1999). Despite the limitations on the relevance of actual snow depth measures, deer use of habitats has often been related to snow depth. In general, snow depths in interior habitats greater than about 50 cm seem to make an area inaccessible (Loveless 1964, Gilbert at al. 1970), and snow depths of about 25 cm to 30 cm are actively avoided (Loveless 1964, Martinka 1968, Telfer 1978). Although unlikely to be as important as snow, low temperatures also affect energy expenditures through thermoregulation costs (Hobbs 1989), and influence habitat selection towards areas that provide exposure to sunlight or shelter from the wind (Loveless 1964, 16  Sargeant et al. 1994).  With increasing winter severity, energy costs of thermoregulation  increase, energy intake decreases (Hobbs 1989), and mortality rate increases (Robinette et al. 1952, Bartmann 1984, Ballard et al. 2001). The moderating effects of cover on winter conditions experienced by mule deer can therefore be of tremendous importance.  2.4.1 Winter Cover Cover in winter is important to mule deer, to the degree that edges are generally avoided (Altendorf et al. 2001). Forests act as barriers to falling snow, altering snow depth and density (Kittredge 1953, Mysterud and Ostbye 1999). Snow interception generally increases with stand density (Kittredge 1953, Kirchoff and Schoen 1987), while forage production decreases (Mysterud and Ostbye 1999). Although total quantity of forage may be greater in the open, available forage will likely be greater under a forest canopy during severe winter conditions. The most important factors affecting the snow interception capabilities of a forest canopy are mean crown completeness and canopy structure (Bunnell et al. 1985, Nyberg et al. 1986, Kirchoff and Schoen 1987). Snow interception efficiency also varies greatly with magnitude of snowfall (McNay 1985), as well as wind and slope to a lesser degree (Bunnell et al. 1985). Mean crown completeness, which measures the average proportion of sky blocked out by tree crowns in a stand, is considered the best measure available of how snow responds to a canopy (Bunnell et al. 1985). In general, lower mean crown completeness results in less snow intercepted (Kittredge 1953, Bunnell et al. 1985). Crown closure, defined as the proportion of the ground surface in a stand encompassed by vertical projections of the outer edges of tree crowns, is considered an inferior correlate of snow interception (Bunnell et al. 1985). Nevertheless, mule deer appear to prefer stands with high crown closure during periods of deep snow (Armleder et al. 1994, D'Eon and Serrouya 2005). Canopy structure is a function of tree species composition, age, density, and site productivity. Snow interception efficiency increases with roughness of canopy, length and width of crown, inter-whorl distances, as well as flatness and width of foliage (Bunnell et al. 1985). Interception efficiency decreases with branch slope and branch flexibility (Bunnell et al. 1985). Evidence suggests that Douglas-fir is the most effective at intercepting snow of the naturally occurring tree species in British Columbia (Nyberg et al. 1986). Indeed, mule deer appear to show greater preference for Douglas-fir dominated stands as snow depths increase (Armleder et al. 1994). However, there is some debate over whether tree species is of dominant importance in affecting mule deer winter range (Nyberg et al. 1986). Requirements for snow interception cover 17  will of course vary with area and year, as low snowfall would make considerable canopy closure unnecessary, while very deep snow may make areas completely inaccessible (Nyberg et al. 1986). Although the energy cost of thermoregulation is likely to be only a small portion of total energy costs, it may still have a substantial impact on energy reserves (Hobbs 1989). Thermal cover in winter provides relief from extreme high and low temperatures by reducing incoming and outgoing radiation and wind speed (Mysterud and Ostbye 1999). Any vegetation cover impedes air flow, reduces wind velocity, and increases air turbulence, thereby reducing windchill in winter, and this effect increases with vegetation height (Bunnell et al. 1986). Any forest canopy that provides good snow interception will also provide adequate thermal cover, although adequate thermal cover can also be supplied by tall shrubs with negligible snow interception properties (Bunnell et al. 1986). Optimal shelter from severe winter conditions would have a patchy canopy structure to allow ease of locomotion under dense canopy and increased forage production within small canopy gaps. As this structural heterogeneity is generally absent in natural young stands, old-growth stands have traditionally been considered optimal mule deer winter range (McNay 1985). Indeed, mule deer appear to prefer old-growth stands in winter in the central interior of B.C., with preference increasing during periods of deep snow (Armleder et al. 1994).  2.4.2 Winter Forage In winter, particularly in northern and high elevation habitat types, shrubs are generally the most important forage for mule deer, with the remainder of the diet made up mostly of forbs, coniferous trees, and some grasses (Julander 1955, Martinka 1968, Constan 1972, Kufeld et al. 1973, Willms et al. 1976, Keay and Peek 1980, Campbell and Johnson 1983). Lichen and mushrooms may be present in small quantities (Kufeld et al. 1973, Willms et al. 1976, Campbell and Johnson 1983). Lichen supply in winter is likely to be positively correlated with stand age (Waterhouse et al. 1991). Consumption of forbs and grasses tends to decrease sharply, and consumption of tall shrubs, trees, and lichens tends to increase sharply, as snow accumulates and low-lying plants are covered (Willms et al. 1976, Hanley and McKendrick 1985, Nyberg et al. 1986). Richness of plant species in the diet may be of importance in winter, as Waterhouse et al. (1994) found that in the central interior of B.C.; although only a few species made up the bulk of the diet, 87 species were consumed in total.  18  Very few published studies presented actual forage preferences, as most simply presented diet composition. Of those that did demonstrate preference, fewer still presented preferences for forages relevant to the interior of B.C. Willms et al. (1976), working in the Kamloops Forest District, found that sweet clover (Melilotus alba) and other members of the legume family were preferred forbs. For shrubs and trees, false box (Pachistima myrsinites) and Oregon grape (Mahonia aquifolium) were preferred over Douglas-fir (Pseudotsuga menziesii) foliage. This is likely due in part to the fact that essential oils from Douglas-fir have an inhibitory effect on rumen microbes, reducing digestibility (Oh et al. 1970). Working in the Selkirk Mountains of B.C., D'Eon (2001) found that Saskatoon berry (Amalenchier alnifolia), redstem ceanothus (Ceanothus sanguineus), Douglas maple (Acer glabrum), nootka rose (Rosa nutkana), ocean spray (Holodiscus discolor), wild rose (Rosa gymnocarpa), and ninebark (Phyoscarpus malvaceus) were preferred relative to availability. Western hemlock (Tsuga heterophylla) and beaked hazelnut (Corylus cornuta) were avoided, while white birch (Betula papyrifera) and willow (Salix spp.) were browsed in proportion to their availability (D'Eon 2001). Because of the belief that Douglas-fir blow-down is of great importance as forage on critical winter range, at a time when most other forage has either been buried or consumed, a number of studies have examined preference for Douglas-fir foliage. In the Kamloops Forest District, for example, Douglas fir foliage was found to make up 63% of the diet during the month of December (Willms et al. 1976). Foliage from Douglas-fir trees grown in the open has been shown to be preferred over that grown in gullies, unshaded needles have been shown to be preferred over shaded needles, and dry foliage seems to be preferred over moist foliage (Tucker et al. 1976). Foliage from the tops of trees is preferred over low-grown foliage, and older trees are preferred over younger trees (Tucker et al. 1976, Dawson et al. 1990). Seedlings of Douglasfir are generally avoided in uneven aged stands but new growth is likely to be browsed lightly, with a preference for taller seedlings and seedlings grown from cuttings of mature trees (Silen et al. 1986, Dawson et al. 1990). When assessing or predicting the value of winter habitat for mule deer, it is important that the researcher keep in mind that animal population densities may not necessarily be correlated with forage abundance. For example, although mule deer distributions have, in certain situations, been shown to be positively correlated with abundance of shrubs (Loveless 1964, Anderson et al. 1972, Peek et al. 2002), distributions have also been shown to be independent of shrub abundance (Telfer 1978). In one case, deer population density was found to be negatively correlated with shrub production, and deer were concentrating on unproductive areas (Anderson 19  et al. 1972). This is likely due to the fact that in severe winter conditions, snow dictates deer distributions (Martinka 1968, Gilbert at al. 1970, Wallmo et al. 1977, Thomas et al. 1979, Schoen and Kirchoff 1985, D'Eon 2001). As winter severity increases, preferred areas are those that receive the most sunlight, are driest in summer, and therefore have the disadvantage of supporting less vegetation than alternative areas (Wallmo 1981).  2.5 Summer Range Summer range, although generally not implicated directly as a major factor in deer population regulation, plays an important role in physically preparing deer for winter hardships (Julander et al. 1961, Wallmo 1981, Torbit et al. 1985, Hobbs 1989, Edge et al. 1990, Peek et al. 2002). Relative to winter, in summer mule deer prefer more northerly exposures and higher elevations (Gilbert at al. 1970, Thomas et al. 1979, Kie et al. 1991, Nicholson et al. 1997). Deer that move to higher elevations in the summer benefit from a more digestible diet than deer remaining at lower elevations (Hanley and McKendrick 1985). Mule deer also show a preference for edge in summer months (Tomm et al. 1981), likely due to the triple benefits of access to forage, thermal cover, and security cover.  2.5.1 Summer Cover Cover for mule deer in summer generally appears to be a simpler topic than it is in winter. As discussed earlier, deer require security and thermal cover in summer. In addition to normal security cover requirements, deer require cover in summer to camouflage fawns from predators (Shackleton 1999). However, one study that examined correlations with fawn mortality found no relationship between hiding cover and fawn mortality, but did find a strong relationship between forb production and decreased mortality (Hamlin at al. 1984). The authors suggested that this could be due to either increased nutritional intake, or increased alternate prey (e.g. microtines) for coyotes. The requirements for thermal cover in summer were also discussed earlier, and involve providing shade (Sargeant et al. 1994) and cooling airflow (Parker and Gillingham 1990). While shade is desirable, a forested area with too much canopy closure reduces carrying capacity due to reduced forage production (Bennett et al. 1940).  2.5.2 Summer Forage In summer, forbs are utilized much more, and shrubs much less, than in winter (Julander 1955, Martinka 1968, Kufeld et al. 1973). This is a generalization of course, and there are 20  exceptions (e.g. Campbell and Johnson 1983, Mt. Chopaka, Washington). Forbs are preferred when available (Julander 1955, Deschamp et al. 1979), and are clearly an important summer forage (Pederson and Harper 1978). Trees are utilized considerably less in summer than in any other season (Kufeld et al. 1973). Grass is generally consumed least in summer (Kufeld et al. 1973), although it may be an important source of nutrition when it is succulent and actively growing, usually in spring or early summer (Willms et al. 1976, Collins and Urness 1983). Lichen and mushrooms may be present in the diet in small quantities (Kufeld et al. 1973, Willms et al. 1976, Deschamp et al. 1979, Campbell and Johnson 1983). Actual amounts of forage types in mule deer diets vary considerably among studies, likely due largely to differences in availability (Julander 1955). Some studies found that browse dominated summer diets (e.g. Smith 1953, Utah; Wallmo et al. 1972, Colorado; Hansen and Reid 1975, southern Colorado; Austin and Urness 1985, western Utah), while others found that forbs made up the bulk of summer diets (e.g. Martinka 1968, Utah; Willms et al. 1980, southern interior of B.C.). In a review of 17 papers examining summer diets, Kufeld et al. (1973) found that shrubs and forbs were consumed in roughly equal proportions on average, but the proportion of herbs in the diet appears to increase somewhat when regions more arid than those typical to the southern interior of British Columbia are excluded from consideration. There is a clear trend of decreasing consumption of forbs and increasing consumption of shrubs from early to late summer, regardless of proportional contributions (Smith 1953, Willms et al. 1980, Austin and Urness 1985), probably due to the curing of forbs under hot and dry conditions (Coe et al. 2001). Richness of plant species in the diet may be of importance in summer. For example, Wallmo et al. (1972) in Colorado found that while essentially one forage species made up the bulk of the diet per habitat type (Vaccinium spp. in forested areas, Salix spp. on roads), the remainder of the diet was made up of 85 species, 71 of which contributed less than 1% each. Although none of the studies investigated this point in particular, the possibility certainly exists that deer do not require species rich diets, and the large number of species making small volume contributions to the diet may simply be a result of forage sampling behaviour. With regards to general summer forage, deer appear to primarily select foods with high digestibility and low lignin content (Vangilder et al. 1982). This is followed in importance by forage that is abundant, high in soluble material and high in calcium content. Forage that is high in gross energy and lignin content appears to be the third most important component of selection (Vangilder etal. 1982).  21  As with winter forage, few published studies have presented preferences for specific summer forages. Mushrooms appear to be a preferred source of food when available (Deschamp et al. 1979). Julander (1955), working in Utah, observed numerous forage preferences likely to be relevant to the interior of B.C, and Julander et al. (1961) produced a very similar list of preferred species proposed from numerous unpublished studies in Utah. Of browse species, aspen (Populus tremuloides), antelope-bush (Purshia tridentatd), Saskatoon berry, choke cherry (Prunus virginiand), elderberry (Sambucus spp.), falsebox, ninebark, Oregon grape, and wild rose (Rosa spp.) are considered highly preferred. O f forbs, white geranium (Geranium richardsonii), asters (Aster spp.), and lupines (Lupinus spp.) were highly preferred. Ligusticum spp., larkspur (Delphinium spp.), Valerian spp., columbine (Aquilegia sp.) and paint brush (Castilleja sp.) were also identified as highly preferred, although species other than those occurring in the interior of B.C. were specifically referred to. In the spruce-fir zone of western Montana, White (1960) also found preferences for plant species likely to be relevant. O f browse species, black huckleberry (Vaccinium membranaceum), and western mountain ash (Sorbus scopulina) were preferred. Of forb species, golden-aster (Chrysopsis villosa) was preferred to the point of being 'chewed to the ground' wherever it was found (White 1960).  2.6 Summary and Conclusion Clearly, an abundance of studies have established the essential habitat requirements of mule deer. The approximate rankings in importance of habitat variables per season are summarized in Tables 2.1 and 2.2. However, a detailed understanding of the complexity of mule deer ecology remains elusive. The abundance, availability and quality of one habitat requirement are likely to affect preferences and requirements for other habitat components. Although the presence of such interactions has frequently been acknowledged, their strength and importance have proven difficult to establish (Edwards 1956, Loveless 1964, Peek et al. 1982, Parker and Robbins 1984).  22  Table 2.1. Summary of literature for important winter habitat variables and their probable order of importance. Arrows represent either positive (f) or negative (I) support for the importance or relative ranking of a given variable. Two opposing arrows side by side (i.e. f i ) represents conflicting support, likely due to confounding factors. Arrows are followed by a brief statement describing the variable's functional role, i f not obvious or self explanatory, except in regards to forage variables, in which case the percent contribution to the diet (if described), and study location is stated. The variable 'operative temperature' refers to the temperature experienced by an animal, and includes both ambient temperature and wind-chill. When considering forage volume variables, it is considered given that the volume of preferred forages will always be a preferred variable over gross volume. Studies thought to be irrelevant to the southern interior of British Columbia are excluded from mention. Irrelevant studies in the literature review by Kufeld et al. (1973) do not contribute to values presented here (values from Colorado, California, and Arizona excluded; values from Montana, Idaho, N . E . Washington, and Wyoming included). Relevant studies already included in the review by Kufeld et al. (1973) are not mentioned separately. Rank  Variable  Studies Providing Supporting Evidence for Rank  l  snow depth  Robinette et al. 1952; Edwards 1956; Loveless 1964; Martinka 1968; Gilbert at al. 1970; Willms et al. 1976; Wallmo et al. 1977; Telfer 1978; Carpenter et al. 1979; Schoen and Thomas et al. 1979; Bartmann 1984; Parker et al. 1984; Picton 1984; Bunnell et al. 1985; Hanley and McKendrick 1985; Kirchoff 1985; Bunnell et al. 1990; Shackleton 1999; Ballard et al. 2001; all | , forage and locomotion  2  mean crown completeness  Bunnell et al. 1985 (f, snow cover); Kirchoff and Schoen 1987 (f, snow cover)  3  crown closure .  4  stand age  5  stand density  6  overstory tree species  7  sinking depth in snow  Bunnell et al. 1985 (f, snow cover); Armleder et al. 1994 (f, snow cover); D'Eon 2004 (|, snow cover at low elevations) Tucker et al. 1976 (f, forage); Bunnell et al. 1985 (|, snow cover); McNay 1985 (t, snow cover); Dawson et al. 1990 (f, forage); Armleder et al. 1994 (f, snow cover) Kittredge 1953 (f, snow cover); Kirchoff and Schoen 1987 (|, snow cover); Mysterud and Ostbye 1999 (f, snow cover) Willms et al. 1976 (|, forage, 43% Douglas-fir in winter, S. interior B.C.); Bunnell et al. 1985 (f, snow cover); Nyberg et al. 1986 (J., snow cover); Armledder and Dawson 1992 (|, snow cover); Armleder et al. 1994 (f, snow cover) Bunnell et al. 1990 (|); Shackleton 1999 (|) Kufeld et al. 1973 (f, avg. 69%, stdev. 18%, shrubs and trees combined, extensive review); Willms et al. 1976 (f, 40%, S. interior B.C.); Telfer 1978 (|, SW, central AB); Keay and Peek 1980 (|, 85%,W. Montana and N . Idaho); Campbell and Johnson 1983 (t, 47%, N . central Washington); Peek et al. 2002 (1, S. central Oregon) Kufeld et al. 1973 (f i , avg. 30%, stdev. 18%, extensive review); Willms et al. 1976 (t, 16%, S. Interior B.C.); Keay and Peek 1980 (f, 12%, W. Montana, N . Idaho), Campbell and Johnson 1983 (t, 32%, N . central Washington) Loveless 1964 (|); Hobbs 1989 (f); Sargeant et al. 1994 (f); Mysterud and Ostbye 1999(f)  8  shrub volume  9  herb volume  10  operative temperature  11  vegetation height  Bunnell et al. 1986 (t, wind reduction)  12  distance from edge  Altendorf et al. 2001 (edges avoided)  13  forage species richness  Kufeld et al. 1973 (f, herbs more than shrubs); Waterhouse et al. 1994 (f)  23  Table 2.2. Summary of literature for important summer habitat variables and their probable order of importance. Details are as listed for Table 2.1. Rank  Variable  Studies Providing Supporting Evidence for Rank  1  herb volume  Kufeld et al. 1973 (ti, avg. 56%, stdev. 24% extensive review; Willms et al. 1980 (t, 52%, S. interior, B.C.); Collins and Urness 1983 (t, > 50%, N . Utah); Hamlin at al. 1984 (|, N . Montana)  2  shrub volume  Kufeld et al. 1973 ( | | , avg. 42%, stdev. 25%, shrubs and trees combined, extensive review); Willms et al. 1980 (|, 45%, S. interior, B.C.); Collins and Urness 1983 (t, < 50%, N . Utah)  3  distance to edge  Tomm etal. 1981 (f)  4  canopy closure  Bennett et al. 1940 (f, shade reduces forage production), Sargeant et al. 1994 (|, shade for cooling)  5  forage species richness  Wallmo et al. 1972 (|, Colorado); Kufeld et al. 1973 (|, herbs more than shrubs, extensive review); Pederson and Harper 1978 (T, herbs more than shrubs, SE Utah)  24  Chapter 3 Habitat Suitability Modeling from Data Using Multimodel Selection  3.1 Introduction Reality is essentially infinite in complexity, particularly in relation to the limitations of the human mind (Bronowski 1978). Ecosystems display all the typical traits of complex systems, including interactions across scales, nonlinearity, self organization, and emergent properties (Boyce 1992, Holling 1992, Levin 1998, Levin et al. 1998, Wu and Marceau 2002). Any ecosystem model, be it qualitative or quantitative, represents an attempt to create a simplified approximation of reality (Burnham and Anderson 2002, Kristov 2004). However, while reducing a process to its most dominant factors is necessary in order to improve understanding and make prediction possible, the complexity of ecosystems means that they are not completely reducible to their components, and generalities are difficult to obtain (Boyce 1992, Wu and Marceau 2002). Ultimately, a model must incorporate that which is known with the best guesses possible regarding that which is not (Hilborn and Mangel 1997). A l l this simplification and guesswork, while necessary, inevitably leads to error (Farmer et al. 1982, Burnham and Anderson 2002). Great rigor is necessary in order to minimize this error. There appear to be two essential concepts that must be carefully considered before any formal model building should begin. First, modeling is necessarily and unavoidably subjective, incorporating the prior knowledge, bias, beliefs and objectives of the modeller (Farmer et al. 1982, Marcot et al. 1983, Hodges 1987, Sokal and Rohlf 1995, Guisan and Zimmerman 2000, Jackson et al. 2000, Burnham and Anderson 2002). While some may think this to be a flaw, models that are constructed with an understanding of the system of interest are actually more likely to produce useful insights and reliable predictions (Conroy 1993, MacNally 2000, Burnham and Anderson 2002). As stated clearly by Krebs (1994), "Ecological controversies are biological in nature, and will be resolved by biologists rather than mathematicians." The second essential concept in modeling is that the optimal level of simplification for approximation is a delicate balancing act. The frequently invoked Ockham's razor embodies the principle of parsimony, and has been interpreted to mean, "Among competing hypotheses, favour the simplest one" (Jefferys and Berger 1992). To be more precise, Ockham's razor is 25  intended to shave away only excessive complexity, as some detail will of course be required for any description or prediction. The drive for simplicity is largely due to the fact that, when formulated properly, a simpler hypothesis can be more of a generalization about a process or phenomenon, and therefore more likely than not to be correct in its description. A simple hypothesis has the potential to encompass greater variation than a more complex and specific hypothesis within the same conceptual hierarchy. For example, while increasing parameters in a statistical model will improve the fit to data (decreasing error due to bias), it will also increase the variances of parameter estimates and prediction errors from one data set to another (Gorman and Toman 1966, Brieman 1995, Neter et al. 1996, Buckland et al. 1997, Hilborn and Mangel 1997, MacNally 2000, Burnham and Anderson 2002, Cavanaugh 2004, Johnson and Omland 2004). The difficulties with overly complex ('overfit') statistical models go beyond the problem of increasing parameter estimate variance. With enough parameters an excellent fit can be attained even i f selected variables have little or no correlation with the dependent variable (Brieman 1995, MacNally 2000, Burnham and Anderson 2002). The inclusion of such variables would of course reduce the reliability of a model's predictions (Neter et al. 1996). Also, the more complicated a model the more likely it is to deviate from assumptions (Chatfield 1995), the more difficult it is to understand, and the less likely it is to provide useful insight into the system of interest (Neter et al. 1996, Hilborn and Mangel 1997). Finally, model complexity has a practical cost, as data for each variable must be collected and processed. This increases the financial cost for model testing and application, and decreases the chance that the model will be used properly or at all (Chatfield 1995, MacNally 2000). While excessively complex models may be problematic, excessively simple ones may be even more so (Burnham and Anderson 2002). As models become too simplistic, the risk increases that important explanatory variables or processes will be left out (Neter et al. 1996, Hilborn and Mangel 1997, Burnham and Anderson 2002). Also, the fit of a simple model to the data used for parameter estimation will often be poor (i.e. a large bias, or 'underfit'). Optimal model complexity is a compromise between parsimony and goodness of fit (Sokal and Rohlf 1995), and increases in complexity are desirable only when they result in dramatic improvements to predictive accuracy (Jefferys and Berger 1992). The proper construction of statistical models, in general, is a process that involves four main steps: careful a priori model construction, parameter estimation, model selection, and robust error rate estimation (Burnham and Anderson 2002). While the procedure 26  for  mathematical modeling discussed in this chapter is applicable for virtually any subject, the emphasis here is on the production of habitat suitability models, which relate the distribution and abundance of organisms to measurable habitat attributes for which they display a preference (Farmer etal. 1982).  3.1.1 Model Construction Model construction consists of the selection of relevant variables, variable interactions, and an appropriate functional form (Sokal and Rohlf 1995, Neter et al. 1996). Before beginning model construction, modeling objectives and assumptions must be clearly stated to maximize the model's usefulness and minimize misinterpretation (Farmer et al. 1982). The selection of variables is generally a major difficulty in model building, as correlations between predictor variables often make the identification of causal variables through data analysis extremely difficult (MacNally 2000). Multicollinearity is said to occur when variables are present that are strongly correlated with causal variables, but have no actual relationship with the dependent variable (Hocking 1976, Sokal and Rohlf 1995, MacNally 2000, Tabachnick and Fidell 2001). The inclusion of redundant variables serves to increase error terms and produce less reliable predictive equations, especially when predicting outside the range of sample data (Hocking 1976, Tabachnick and Fidell 2001). Therefore, it is important that variable selection not be based strictly on data analysis (Flack and Chang 1987). Choices of variables and functional forms can only be properly performed by someone with considerable knowledge of the population from which data are drawn (Gorman and Toman 1966). Ideally, habitat variables should have mechanistic linkages to a population's performance, as simple correlative relationships are likely to produce inferences that are either unreliable or reliable only for a very limited set of conditions (Hobbs and Hanley 1990, Guisan and Zimmerman 2000). Modeling causal factors is not always practical however, and predictor variables may be limited to indirect variables that are easy to measure or manipulate (Marcot et al. 1983). Habitat suitability modeling has the added limitation that habitat itself may not actually be limiting population distribution or abundance (Farmer et al. 1982). It seems to be relatively undisputed in the statistics literature that careful, well informed thought is an integral part of model construction, serving to reduce the incidence of spurious (highly variable from one data set to the next, biologically irrelevant) relationships (Anderson et al. 2001, Pearce at al. 2001, Burnham and Anderson 2002) and models that are more complex than the data can reasonably support (Stauffer and Best 1986, Buckland et al. 1997). These a 21  priori considerations lead to a more powerful strength of inference (MacNally 2000, Burnham and Anderson 2002). As stated by Chamberlin (1890), "Laudable as the effort at explanation is in itself, it is to be condemned when it runs before a serious inquiry into the phenomenon". What does appear to be in dispute is at what point in the modeling process a subjective filter should become explicitly involved, and to what degree. For example, Chatfield (1995) and Burnham and Anderson (2002) state that models should be constructed entirely a priori using all available information, such as personal experience, expert knowledge, and accepted theory. In contrast, Cox (1958) and Hilborn and Mangel (1997) believe that it is important to first become intimately knowledgeable of the information contained in a data set before beginning model construction. However, the point is made well by Burnham and Anderson (2002) that variables, relationships, and model form suggested by data can be considered after a priori considerations have been exhaustively explored, but a priori based modeling cannot be performed after the data set has been explored. A middle ground appears to be advocated by Sokal and Rohlf (1995) and Neter et al. (1995), who suggested that prior knowledge be used as the foundation for model construction, but not as the exclusive influence. This more moderate approach is attractive when a priori knowledge is insufficient, or when there is the possible influence of a novel but important factor. Certainly the risk of modelling spurious effects will increase (Burnham and Anderson 2002). However, it needs to be kept in mind that i f preconceived theories are held too tightly, facts will inevitably be misinterpreted and potentially valuable insights will be missed (Chamberlin 1890).  3.1.2 Coefficient Estimation Once variables and form are selected, model parameters must be estimated. Occasionally, a priori information may suggest specific coefficient values or value limits, and that information may be used to produce a more reliable model (Burnham and Anderson 2002). More commonly, however, parameters are based on empirical data, in which case the quality of the data will limit the quality of the model (Turner et al. 1995). Before parameters can be estimated from data, the data must be checked for outliers. Also, multicollinearity and serious departures from the assumed distribution of residuals must be assessed within a given model (Gorman and Toman 1966, Burnham and Anderson 2002). Collinearity diagnostics exist to help screen for multicollinearity (Neter et al. 1996, Tabachnick and Fidell 2001). Careful examination of residual plots can reveal deviations from assumptions of equal variance and the shape of the relationship between the dependent and independent 28  variables (Gorman and Toman 1966). Normal probability plots can help assess i f residuals follow a normal distribution. A global model, including all variables included in less complicated models, has been suggested as the model used for such diagnostic assessments (Burnham and Anderson 2002). Transformation of variables maybe necessary to satisfy assumptions (Hodges 1987), and prior knowledge and experience should be used to suggest which transforms are appropriate (Box and Cox 1982, Neter et al. 1996). When there is no a priori reason for selecting a particular transformation, a method such as the Box-Cox procedure might be used to suggest one (Box and Cox 1964, Sokal and Rohlf 1995). With this approach, however, there is the risk of finding spurious transformations that may be sub-optimal when applied to a different sample. With regard to outliers, any that are the clear result of error should be removed regardless of the circumstances. Deletion of valid outliers, however, will lead to underestimation of predictive uncertainty i f uncertainty is estimated with the same dataset used for parameter estimation (Chatfield 1995).  3.1.3 Model Selection Model selection may be thought of as the process of making inferences from the data on what relationships they are capable of supporting (Burnham and Anderson 2002). A priori information should be incorporated as much as possible, preferably within a Bayesian framework (Akaike 1981). However, a point will often come where knowledge is too vague to assign prior probabilities to models or parameters with confidence (Akaike 1981). When this occurs, statistical model selection techniques become appealing. Some common single model selection methods and criteria include the likelihood ratio method, stepwise selection, adjusted R , and Mallow's C . Likelihood ratio tests operate on the 2  p  principle of parsimony, and select a more complex model only when the maximum likelihood of that model is significantly greater than that of the simpler model (Johnson and Omland 2004). Stepwise selection, which includes forward and backward methods by single variables or subsets, operates by making decisions to add or remove variables based on whether their presence or absence elicits a statistically significant change in the coefficient of multiple determination (Sokal and Rohlf 1995, MacNally 2000). The adjusted coefficient of multiple determination includes a penalty to the calculation of R for the number of parameters in a model 2  (Neter et al. 1996). Mallow's C statistic is used to graphically compare the bias and error of p  fitted models (Gorman and Toman 1966). 29  These model selection procedures are simple to operate, particularly with modern statistical analysis programs, but they are plagued with difficulties. First, likelihood ratio tests and stepwise selection are based on statistical significance, which frequently has little or no relation to biological significance (Hilborn and Mangel 1997). This is because p-values are arbitrary, fluctuating with sample size and relying on assumptions that are only approximately met (Cox 1958, Preece 1990, Johnson 1999, Anderson et al. 2000). Second, multicollinearity and limitations of the sample may make reliable predictor variables impossible to recognize in the data by any automated method (Flack and Chang 1987, Neter et al. 1996). Third, stepwise selection methods can be highly variable in the subset of variables selected, both between forward and backward methods and between samples (Maurer 1986, James and McCulloch 1990). Stepwise procedures were not designed for, nor are they effective at ranking variables (Hocking 1976, James and McCulloch 1990, Sokal and Rohlf 1995). Results of statistical methods for variable selection should be used as tentative suggestions at best, and selected variables should not be used unless confirmed by theory (Flack and Chang 1987). Model selection by adjusted R  is a well established technique and a distinct  2  improvement over simply modelling to maximize fit. However, as it is based on R , its 2  applications are limited. Comparisons of R values between models become problematic when 2  there is no intercept, when comparing between different transformations, and for weighted least squares (Anderson-Sprecher 1994). Also, although this practice would be ill advised, R values 2  appear to be biased upwards when models are constructed using stepwise regression (Rencher and Pun 1980). The C statistic has been found to have a one-to-one correspondence to the p  adjusted R (Kennard 1971). 2  It must be noted at this point that any data driven method that results in the selection ofa single model leads to a pervasive underestimation of error (Draper 1995). When the model is not known and is selected based on the data itself, which model is 'best' will likely vary from sample to sample (Burnham and Anderson 2002). However, calculations of parameter estimation error are made assuming that the model is known. When this is not the case, error arising from model selection uncertainty is being left out of calculations and estimates of error will be overly conservative (Chatfield 1995, Draper 1995, Buckland et al. 1997). Model selection uncertainty has its greatest impact when sample size is small, when there is a lack of true replication, and when there are many possible models (Burnham and Anderson 2002). A l l this concern with underestimation of error does not apply; however, when estimates of error are based on a data set  30  completely independent from that used in parameter estimation and model selection (Chatfield 1995).  3.1.4 Multimodel Selection and Inference In order to help resolve the difficulty of model selection error, Chamberlin's (1890) theory of multiple working hypotheses is incorporated into model selection (Anderson et al. 2000). Science operates by challenging more than one description of nature with data, and it makes sense to stay true to this philosophy when modeling (Hilborn and Mangel 1997). The best model is very unlikely to be 'truth', and therefore it is illogical to treat it as such, especially when alternative models may contain important information (Stoica et al. 2004). Including that information allows for more robust inference as well as an estimation of model selection uncertainty, which is often ignored in traditional model selection procedures (Anderson et al. 2000, Burnham and Anderson 2002). In addition, the act of producing multiple working hypotheses tends to stimulate non-linear thinking, which helps in viewing the problem from multiple perspectives (Chamberlin 1890). Once multiple models are constructed, they must be ranked according to how well they fit the data while taking model size and complexity into account (Cavanaugh 1999, Cavanaugh 2004). In this way the principle of parsimony is incorporated into model selection (Schwarz 1978, Atilgan 1996). Differences between rank scores can then be transformed into weights that provide measures of relative support for each candidate model (Johnson and Omland 2004). These weights can then be used to calculate weighted averages of parameter estimates and predictive output, which will be relatively robust in terms of model selection bias and model selection uncertainty (Johnson 1999, Anderson et al. 2000, Johnson and Omland 2004). A number of criteria have been developed for ranking models, including Mallow's C  p  (Gorman and Toman 1966, Mallows 1973), Schwarz Criterion (SC; Schwarz 1978), K I C (Cavanaugh 1999), and Akaike's Information Criterion (AIC; Akaike 1973). These criteria act as objective and consistent methods for ranking and weighting candidate models (Johnson 1999, Anderson et al. 2000). A l l but Mallow's C operate by providing approximately unbiased p  estimates of Kullback-Liebler (Kullback and Leibler 1951) information loss, which classifies them as information-theoretic approaches (Akaike 1973, Hurvich and Tsai 1989, Cavanaugh 1999, Anderson et al. 2000). The performance of a criterion appears to be largely determined by how well its penalty term approximates the appropriate bias adjustment (Kim and Cavanaugh 2004). Sample size is 31  incorporated into the bias adjustment such that the ranking of complex models improves with sample size (Buckland et al. 1997). It was discovered quite some time ago that AIC tends to overfit in situations where sample size is small, or when the number of parameters is a large fraction of the sample size (Sugiura 1978, Hurvich and Tsai 1989). To resolve this difficulty, a version of AIC with a corrected penalty term was developed (AICc; Sugiura 1978, Hurvich and Tsai 1989, Bedrick and Tsai 1994). There appear to be important differences in the performance of different criteria. Fortunately, thorough investigations to clarify those differences have been undertaken. In circumstances that involve a true model of finite complexity where the dimension of the true model is represented in the candidate set, SC is clearly preferred (Hurvich and Tsai 1989, Hurvich and Tsai 1990, Cavanaugh 1999). AIC in this situation is not consistent in model ranking and tends to overfit (Hurvich and Tsai 1989). In contrast to SC, AIC operates under the assumption that the model representing 'truth' is infinitely complex (Buckland et al. 1997). Also, the true model does not need to be present in the set of candidate models for optimum results (Anderson et al. 2000). When the 'true' model is highly complex, SC tends to select underfit models (Cavanaugh 1999) and AICc outperforms all competing criteria (Hurvich and Tsai 1989). While KIC has been shown to compare favourably with AIC, it does not outperform it, and is not optimal for as broad a class of models as AIC (Cavanaugh 1999, 2004, K i m and Cavanaugh 2004). AIC appears to be consistent in outperforming Mallow's C (Hurvich and Tsai 1989, Fujikoshi and Satoh 1997). p  More alterations of A I C have recently been suggested, including a 'modified' AIC (MAIC; Fujikoshi and Satoh 1997), and an 'improved' AIC ( A i d ; Hurvich et al. 1990). M A I C appears to be optimal only in situations where the dimension of the true model is the same as the largest dimension represented in the candidate model set (Cavanaugh 2004). M A I C is therefore not likely to be relevant to application in ecology. In contrast, A I Q appears to be a genuine improvement over AICc. A I Q provides only slightly better model selections when the number of model parameters is no more than half the sample size, but performs much better when some models have higher dimension (Hurvich et al. 1990). This is because AICi was designed under the assumption that the candidate models are either correct in dimension or overfitted. Therefore, it seems likely that in ecological modeling applications, where the ratio of model parameters to sample size should remain low and true model complexity is very high, the slight increase in performance of AICi over AICc is not worth the considerable increase in calculation difficulty of AICi. 32  Although a substantial improvement over single model approaches, a multimodel ranking approach does have weaknesses that must be recognized. First, model ranking is still based on data, and therefore the reliability of results is sensitive to the limitations and characteristics of those data (Chatfield 1995). Second, the reliability of results is sensitive to the construction and selection of candidate models to be considered. If important variables or better models are left out of the candidate set, they cannot be identified by selection criteria (Chatfield 1995, Burnham and Anderson 2002). Similarly, i f spurious models are included that closely fit quirks of the data, selection criteria may rank them higher than is appropriate (Johnson and Omland 2004). Third, if none of the proposed models is a useful description of reality, then resultant inferences will be weak regardless of robustness to model selection uncertainty (Buckland et al. 1997).  3.1.5 Validation  Although a priori model construction and multimodel consideration minimizes the underestimation of error that results from fitting and evaluating models, error in predicting the behaviour of any open system is inevitable (Oreskes et al. 1994). Models can never be completely true, but fortunately they only ever need to be useful (Box 1979). To evaluate usefulness, habitat suitability model output should be tested by comparison with observations in a process known as validation (Hansen et al. 1993, Rykiel 1996, Hilborn and Mangel 1997, Huggard 1999, Burnham and Anderson 2002). Validation is necessary to evaluate model reliability and establish credibility (Marcot et al. 1983, Conroy et al. 1995, Rykiel 1996, Vanclay and Skovgaard 1996). Validity itself is a relative measure that must take into account model objectives, ecological context of planned application, as well as specific precision and accuracy requirements (Farmer et al. 1982, Schamberger and O'Neil 1986, Rykiel 1996, Vanclay and Skovgaard 1996, Guisan and Zimmerman 2000). The desired end result is a sufficient and quantified understanding of the model's attributes such that the limits of reliability can be anticipated within the planned range of application (Farmer et al. 1982). It is generally accepted that accurate estimates of prediction errors require an analysis of how well model output fits data that are independent of those used in the construction of the model or the estimation of model parameters (Hocking 1976, Farmer et al. 1982, Lancia et al. 1982, Conroy 1993, Chatfield 1995, Conroy et al. 1995, Vanclay and Skovgaard 1996, Fielding and Bell 1997, Huggard 1999, Yang et al. 2004). Occasionally, truly independent data may not be available. Unfortunately, error estimates made with data that are not independent are likely to 33  be biased and overly optimistic (Chatfield 1995, MacNally 2000). When lacking independent data, some researchers suggest splitting the dataset into a portion for model fitting and the remainder for error estimation - a technique known as cross-validation (e.g. Anderson-Sprecher 1994, Draper 1995). However, it seems likely that cross-validation is a poor substitute for independent data (Chatfield 1995). More recently, researchers have suggested re-sampling the dataset with replacement (bootstrapping) to form numerous new datasets that are treated as i f they were independent (Draper 1995, Buckland et al. 1997, Hilborn and Mangel 1997, Guisan and Zimmerman 2000). Bootstrapping appears to have the effect of 'fleshing out' the sample distribution. The dependence of the 'new' data to the original data, however, means that derived error estimates will be unreliable i f the original data is not a reasonable representation of the population of interest (Vanclay and Skovgaard 1996). It must be kept in mind that agreement between model output and independent observations does not in itself confirm model assumptions (Oreskes et al. 1994). A faulty model may appear to be correct due to sheer coincidence, as many models may arrive at the same prediction under very different assumptions (Conroy et al. 1995). Also, models that adequately fit historical data may not necessarily be reliable in predicting the future, as model errors will propagate and modeled relationships may become invalid under changing conditions (Conroy 1993, Oreskes etal. 1994).  3.1.6 Habitat Suitability Modeling Using Multimodel Selection: Application to M u l e Deer in the Interior of British C o l u m b i a  Mule deer are fairly amenable to habitat suitability modeling, as they show strong habitat preferences (e.g. Loveless 1964, Wallmo et al. 1972, Deschamp et al. 1979, Nicholson et al. 1997, Altendorf et al. 2001), and population dynamics are often driven by the interaction between severe winter conditions and habitat (e.g. Wallmo et al. 1977, Schoen and Kirchoff 1985, Edwards 1956, Wallmo 1981, Bartmann 1984, Parker et al. 1984, Picton 1984, Hanley and McKendrick 1985, Ballard et al. 2001). As an important game species, mule deer have requirements for winter range that are necessary considerations in forest harvesting operations throughout the Southern and Central Interior of British Columbia (e.g. Ministry of Sustainable Resource Management's Cariboo-Chilcotin Land Use Plan 1995 and Okanagan-Shuswap Land and Resource Management Plan 2001). Therefore, modeling habitat suitability of mule deer has both a sound ecological basis, as well as a management incentive to find more efficient methods of satisfying winter range requirements. The goal of this project was to produce separate mule 34  deer habitat suitability models for the summer and winter seasons, while at the same time exploring the methods necessary for doing so. To be considered acceptable, final models must make good use of available data, fit the data well, and must be functional for management applications.  3.2 Methods 3.2.1 The Data Empirical data used for fitting and validating habitat suitability models were provided by T.P. Sullivan (Department of Forest Sciences, Faculty of Forestry, University of British Columbia), and were obtained from a split-split plot experiment (Hicks and Turner, 1999) he conducted at three replicate sites (see Figure 3.1) in the interior of British Columbia (this and all details on the experimental design from T.P. Sullivan, personal communication).  Gavin Lake  Summerland  Kelowna  Figure 3.1. Map of British Columbia showing the location of the three study sites at Summerland, Kelowna, and Gavin Lake (Cariboo). Study sites (see Table 3.1) consisted of forest stands with an overstory dominated by young lodgepole pine. The location, proximity, and sizes of study sites were selected based on a combination of minimizing spatial autocorrelation (Hurlbert 1984, Legendre 1993), and practical considerations regarding logistics and access for operational-scale stand manipulations. Harvesting of lodgepole pine on the Summerland site began in 1978 in response to mountain pine beetle infestation (Dendroctonus ponderosae). Harvesting was done mostly by 35  clearcutting, but with some dispersed and aggregated reserves of Douglas-fir seed-trees. The number of residual Douglas-fir trees ranged from zero to two trees per ha, depending on original stand condition and post-harvest windthrow losses. Lodgepole pine now present arose from natural regeneration to become the dominant tree species. In addition to Douglas-fir, interior spruce, subalpine fir, Ponderosa pine (Pinus ponderosa), willow (Salix sp.), Sitka alder (Alnus sinuatd), and trembling aspen also made up minor components of the forested stands at the Summerland site.  Details of replicate sites from a study on the response of wildlife to various thinning and fertilizer regimes in young lodgepole pine stands, conducted by T.P. Sullivan (unpublished data). Montane Spruce B E C zone, dry mild subzone. Sub-Boreal Spruce B E C zone, dry mild subzone. Table 3.1.  ^^3^®.... (SRrarftir  'Smsmm^mJj  ** M S  d  m  *  M S  d  m  *  SBSdm  1450-1520m  1220-1240m  850-870m  gently rolling  gently rolling  gently rolling  cool continental climate; cold winters and short, warm summers  cool continental climate; cold winters and short, warm summers  continental; severe snowy winters, relatively warm, moist, short summers  . 0.5-4.7°C  0.5-4.7°C  1.7-5.0°C  380-900mm (ms)  380-900mm (ms)  440-900 mm  9980 to 11150 stems per ha  8686 stems per ha  3333 stems per ha  5 +/- 0.2 cm to 9.5 +/- 0.2 cm  . 6.2 +/- 0.2 cm to 10.8 +/- 0.2 cm  8.5 +/- 0.2 cm to 12.7 +/- 0.3 cm  4.1+/- 0.1m to 5.1 +/- 0.1m  5.3+/- 0.1m to 6.4+/- 0.1 m  6.1+/-0.1m to 7.6 +/- 0.2 m  17-19 yrs.  17-18 yrs.  18 yrs.  4.4 to 11.3 ha  9.5 to 12.6 ha  1.5 to 4.5 ha  36  Harvesting of the Kelowna site by clearcutting occurred from 1979 to 1980, and a stand dominated by lodgepole pine was naturally regenerated. Minor tree species components were the same as the Summerland site, with the addition of western larch (Larix occidentalis). Unlike the other two sites, the Kelowna site was virtually split by a riparian area with a steep ravine that varied in width from 75 to 300 m. The Cariboo site was dominated by lodgepole pine that had regenerated naturally after wildfire, with minor components of interior spruce, subalpine fir, and some Douglas-fir. The area was clearcut in 1976, and the current lodgepole pine dominated stand arose from a combination of natural regeneration, and some planting that took place in 1983. In the understory of the young forest stands at each site, the main herb species included yarrow (AchUlea millefolium), rosy pussytoes (Antennaria microphylla), racemose pussytoes (A. racemosa), field pussytoes (A. neglecta), heart-leaved arnica {Arnica cordifolia), fireweed (Epilobium angustifqlium), bunchberry (Cornus canadensis), wild, strawberry  {Fragaria  virginiana), assorted grasses, white-flowered hawkweed (Hieracium albiflorum), Arctic lupine (Lupinus arcticus), and common dandelion (Taraxacum officinale). At the Cariboo site, fringed aster (Aster ciliolatus), northern bedstraw (Galium boreale), creamy peavine (Lathyrus ochroleucus), red-clover (Trifolium pratense), white clover (T. repens), and American vetch (Vicia americana) were also prominent (Sullivan 2005). The main shrub species consisted of Sitka alder (Alnus sinuata), twinflower (Linnaea borealis), black twinberry (Lonicera involucrata), red twinberry (L. utahensis), falsebox (Pachistima myrsinites), black gooseberry (Ribes lacustre), red raspberry (Rubus idaeus), kinnikinnick (Arctostaphylos uva-ursi), willow (Salix spp.), birch-leaved spirea (Spiraea betulifolia), and various Vaccinium species. At the Cariboo site, Saskatoon berry (Amelanchier alnifolia),  tall Oregon-grape (Mahonia aquifolium), prickly rose (Rosa acicularis), and  soopolallie (Shepherdia canadensis) were also notable (Sullivan 2005). Each site was divided into eight treatment stands and one control stand, with treatments consisting of one fertilized and one unfertilized stand at each of four thinning densities; 250 stems per ha, 500 stems per ha, 1000 stems per ha, and 2000 stems per ha (see Figure 3.2 for an aerial view of the Kelowna site for an example). Thinning of all stands was done in the fall of 1993, and fertilizer was applied in 1994, 1997, 1998, 2000, and 2002. In addition to fertilization and thinning, 250, 500 and 1000 stems per ha treatment stands were pruned to a height of 3 m at each site in September to October of 1998.  37  Stand density measurements taken in 1998, 5 years after thinning, showed that stand densities of crop trees were substantially different from stated treatment densities. The densities measured in 1998 were thus used in all analyses. Although stand densities would certainly have changed over the 5 years of the study (1998 to 2003), no better measurement was available. Each replicate site had one unthinned unfertilized stand per replicate site. However, in the unthinned stands, only total stems per ha of trees were measured, compared to the separate measurement of stems per ha over 3 m in height that was performed in all thinned plots. It seems unlikely that trees less than 3 m in height in unthinned stands would contribute much to deer habitat, as the low crown volumes of small trees would provide little shading, environmental moderation and snow interception at stem densities low enough to allow deer movement. I felt that including data from the unthinned unfertilized stands would confound more than clarify, and they were excluded from consideration.  TFL 49 Riverside Forest Product? Ltd.  Figure 3.2. Aerial view of the Kelowna site.  Trees were sampled by nested circular plots, accommodating 10 sample trees per plot. Understory vegetation was sampled each summer from 1999 through 2003 using three 25 m transects with five plots each located systematically within each stand following the methods of 38  Stickney (1980). Each plot consisted of three nested subplots; a 5 m x 5 m plot for sampling trees, a 3 m x 3 m plot for sampling shrubs, and a 1 m x 1 m plot for herbs. Ground cover was estimated visually for each plant layer by treating the two horizontal dimensions of plant crowns as ellipse axes and solving for area. Each plant layer (trees, plants, herbs) was also divided into height classes of 0 - 0.25 m, 0.25 - 0.5 m, 0.5 - 1.0 m, 1.0 - 2.0 m, 2.0 - 3.0 m, and 3.0 - 5.0 m, and multiplied by crown area estimates to calculate a crown volume index which was then averaged over each stand (m / O.Olha; Stickney 1980). Species richness was calculated for each vegetation layer as the number of species sampled per stand, adjusted for sample size using the rearefaction method (Krebs 1999). In addition, species diversity was also calculated for each vegetation layer, along with total structural diversity of all layers combined. Although species richness is considered a measure of species diversity (Krebs 1994), the term diversity as used here refers exclusively to heterogeneity and proportional distribution. The Shannon-Wiener index (Shannon 1948) represents the evenness of the distribution between categories (species for species diversity and height class for structural diversity), and was used to calculate measures of diversity. Structural diversity was also calculated using the Shannon-Wiener index, with plant species represented by crown volume index of vegetation per height class, and coniferous trees represented by the stem density of trees per height class. Within each stand, pellet groups were counted for mule deer (Odocoileus hemionus), moose (Alces alces) and snowshoe hare (Lepus americanus) in the first two weeks of May to represent winter habitat use, and the first two weeks of October to represent summer habitat use. Pellet-group densities are a potentially problematic index of habitat use, however the method has a long history of use, and its strengths and weaknesses are well known. A thorough discussion of those strengths and weaknesses can be found in the Discussion section (Section 3.4). Overall, the method should provide a reasonable index of deer habitat use given the experimental design. Pellet groups were first counted for the winter of 1998/1999 after first clearing all plots in the fall of 1998, and were repeated twice a year for five years. The sampling design utilized permanent 5.0 m circular sample plots located systematically in five-plot arrays placed at the same permanent plots used to sample crop trees in each stand. In total, there were 55-145 plots per stand in Summerland, 60-140 plots per stand at Kelowna, and 35-100 plots per stand at Gavin Lake. This range of plots per stand was due to uncontrolled variation in fertilizer boundaries. Effort was made to count all pellet groups within each plot, and vegetation was moved i f it concealed the ground. Pellets that had become a part of the duff or litter layer were 39  excluded. After counting, pellets were again cleared off plots. Clearing plots takes time but provides complete assurance that the same pellet group will not be counted twice (British Columbia Ministry of Environment, Land and Parks 1998, Neff 1968). Pellet group counts were averaged over permanent sample plots to obtain counts per ha per stand per replicate per year. In total, there were 120 observations (4 treatments x 3 replicate sites x 5 years x 2 seasons) from unfertilized stands and 120 observations from fertilized stands. Observations were split by season due to the differing habitat requirements of deer during summer and winter months. However, as observations were repeated in the exact same location every year for five years, all observations were not independent. Repeated observations are, in fact, spatially autocorrelated (Hurlbert 1984, Legendre 1993). Autocorrelation causes model parameters to be estimated with variances that are not minimized or unbiased (Dutilleul and Legendre 1992). To correct this problem, observations repeated over time were simply averaged to obtain one observation per treatment stand per replicate. This resulted in a final count of 24 independent observations per season. Although this dramatically reduces degrees of freedom, it will prevent error estimates used in AIC calculations and confidence intervals from being irrationally conservative, which they would be if each observation was treated as independent. As a further consideration, pellet group densities themselves could not be used as the dependent variable because of the large expected differences in mule deer population densities between sites. For example, mule deer are known to reach the northern limits of their continuous distribution in the Cariboo Forest Region (Armleder et al. 1994, 1998). Because of this, mule deer will generally be at lower population densities at the Cariboo site relative to the other two sites, and will therefore deposit lower pellet group densities there. Means of pellet-group counts per replicate site support this assumption (see Table 3.2). Using unaltered pellet-group counts would have been appropriate only i f modeled habitat variables are capable of determining regional population densities. This would not be correct, however, as mule deer distribution and abundance is likely to be constrained more strongly by different factors at the regional scale, such as climate. Therefore, all pellet-group counts were standardized by region (each pellet count observation was divided by the highest pellet group density per region) combining fertilized and unfertilized treatments. For separating by region, Summerland and Kelowna were grouped together, and Gavin Lake was treated individually. In this way, pellet group counts were transformed into measures of preference that can be compared directly between areas of different population densities. Fertilized and unfertilized treatments were combined during standardization because fertilized treatments are expected to be of higher habitat suitability, and when fertilized 40  sites are being considered it is important to keep all observations relative to the pertinent maximum. Unfortunately, although standardizing resolves much of the variation due to unequal population densities, it creates the problem of seemingly setting an optimum habitat suitability level at each site. However, no combination of variables captured in the observations at any site is necessarily optimal habitat. Instead, standardizing is setting a particular combination of variable levels as 'best observed'.  Table 3.2. Means and standard deviations of pellet-group counts per site per season, separated by fertilizer treatment. Values are calculated from averages of 5 annually repeated observations. Itasca Summerland  Summer  Kelowna  Gavin Lake  Summerland  Winter  Kelowna  Gavin Lake  unfertilized  201.49  189.25  fertilized  494.18  563.88  unfertilized  176.67  117.82  fertilized  488.02  165.51  unfertilized  21.64  16.21  fertilized  39.39  23.66  unfertilized  268.39  106.78  fertilized  553.75  186.99  unfertilized  277.01  93.65  fertilized  674.63  75.09  unfertilized  106.73  54.67  fertilized  76.68  40.98  If fertilization was not being considered as a management option, it would be best to renormalize the data to the maximum of unfertilized pellet-group counts per site, and exclude fertilized observations. Unfortunately, this would also have the effect of cutting the degrees of freedom for regression in half, and the number of parameters considered in the more complex models would have to be reduced (i.e. only shrub volume, herb volume, and lodgepole pine stem density would likely be considered).  41  3.2.2 Multimodel Construction and Inference For candidate model set construction, the functions that define the relationship between independent variables and the response variable should be assigned as much as possible prior to data analysis (Johnson and Omland 2004). Available variables were: crown volume, species richness, and species diversity measurements for each vegetation layer (herb, shrub, tree), as well as total species richness, total species diversity, total structural diversity and stand densities after thinning (T.P. Sullivan, personal communication). Based on available information in the literature (as summarized in Chapter 2), variables selected as likely to be important for modeling were: herb volume, shrub volume, herb species richness, shrub species richness, stand density, and interactions between stand density and each of herb and shrub volume. For the models produced for management applications, herb and shrub species richness were excluded, because, unlike understory volumes (e.g. Jameson 1967, Peek at al. 2001) and stand density, they would likely be too difficult to estimate and too expensive to measure (Palmer 1990). No available information suggested that 'evenness'  of species abundances (species diversity) was of  importance to mule deer, nor were any of the other calculated variables that were excluded from the a priori model set. Tree crown volume was excluded because although it is likely a superior indicator of canopy closure relative to stand density, the difficulty and expense involved in acquiring the measure prohibits its application in models intended for management use. Total structural diversity would likely have had an impact on thermal properties of winter cover, security cover, as well as represent robustness of the understory relative to the canopy; however, it seems that the connection would not be sufficiently important to include in models. Once variables were selected, they were then evaluated for multicollinearity in the global model. This refers to the situation where predictor variables are so highly correlated that they become redundant and inflate error terms (Tabachnick and Fidell 2001). Variables were examined using collinearity diagnostics in SAS (SAS Systems 2001). If the reciprocal of the variance inflation factor (the tolerance) for any variable is too low (say less than 0.01 or 0.001), the variable in question is considered to be highly related to the others (Neter et al. 1996, Tabachnick and Fidell 2001). In this instance, no variables appeared to be redundant. Functional, shape then had to be assigned, which requires the modeller to describe the likely shape(s) taken by the response variable (in this case, pellet-group, counts as an index of habitat suitability) as it varies with the independent variables. Although linear forms are often assumed, ecological processes are inherently nonlinear (Boyce 1992). Over the entire range of habitat conditions, it seems most likely that habitat suitability would take a sigmoidal form, as 42  suggested for numerical response of predators to changing prey densities by Hassell and May (1974). At very low habitat quality, deer would likely not differentiate between slight improvements that keep the habitat within an "unprofitable" range, and habitat suitability would increase very slowly. Above the "profitable" threshold in habitat suitability, increases in habitat quality would likely translate directly into increases in preference, and habitat suitability would increase rapidly. At high habitat quality, the habitat suitability curve should flatten out, as deer would eventually reach their maximum point of habitat exploitation due to limitations imposed by handling time, space, disease, predation, or other factors.  Figure 3.3. Segmenting of a sigmoidal response curve, a) A n approximately Type I exponential shape, b) Approximately linear segment, c) A n approximately Type III exponential shape. Although habitat suitability is likely to take on a sigmoidal shape over the entire range of habitat variables, the form may not be sigmoidal over the range captured within the observations (see Figure 3.3). If the observations are centered near the "unprofitable" area, the shape may appear to be a Type I exponential (see Figure 3.3a). If the middle of the total range is observed, the relationships would appear to be approximately linear (see Figure 3.3b). And i f the range captured is mostly in the high habitat quality area, the shape may appear to be a Type III exponential curve (see Figure 3.3c; British Columbia Ministry of Forests 1994). It is possible that the sample sites in summer would be well represented by a Type III exponential curve due to the low canopy cover of young thinned stands (i.e. high forage production). Winter sites may intermediate habitat until the snow begins falling, at which point it would become poor (Type I exponential). While a non-linear form is the most intuitive, limited preliminary data analysis did 43  not suggest strong non-linear relationships. In addition, modeling of non-linear forms adds another level of difficulty to parameter estimation, as well as a potential increase in the number of parameters (increased model complexity). Therefore, only linear models were explored. Once variables and model form(s) have been selected, model construction can proceed. Burnham and Anderson (2002) suggested that the number of candidate models be kept below 20 in number, and that one model in the set should encompass all predictor variables and interactions included in less complicated models (a global model). This makes parameter estimation and model selection easier, but is not absolutely necessary. For the example presented here, 19 linear models were produced by listing all plausible combinations of the dependent variables selected for modeling, with the global model having 10 parameters. The full candidate set of models is listed in Table A l (summer data) and A2 (winter) of Appendix 1. A reduced set of models was produced for management application by removing models that included herb and shrub species richness, which are very expensive to measure, and very difficult to estimate (Palmer 1990). Error structure was then defined for parameter estimation (Johnson and Omland 2004). For count data, the error distribution would generally be expected to take on a Poisson distribution (Sokal and Rohlf 1995). However, as the values for each treatment are in fact averages over numerous sample plots and repeated observations, then according to the Central Limit Theorem errors should approach a normal distribution, assuming reasonable model fit (Sokal and Rohlf 1995). As only linear forms with normally distributed residuals were being explored, least squares linear regression was applied for parameter estimation using Proc G L M in SAS (SAS Systems 2001). Fertilizer was used as a classification variable, which allows the use of fertilized and unfertilized observations together to provide degrees of freedom for estimation of parameter coefficients, while separate intercepts are calculated. Although it seems logical that fertilized sites may exhibit steeper slopes of increasing habitat suitability with increasing habitat quality, to incorporate that into models would involve adding interaction terms between fertilizer and each of the other independent variables. It seems unlikely that the available data would have a high enough resolution to support a model with such a proliferation of parameters. Parameters were then estimated, once using all three sites for model fitting, and once using only the Summerland and Cariboo sites. In the latter case, the Kelowna data was reserved to be used as an independent dataset for error estimation. After parameters were estimated, residuals from the global model were examined for serious outliers or deviations from 44  assumptions regarding functional form or error distribution (Burnham and Anderson 2002). When the assumptions of the global model appear to be satisfied, more parsimonious models will also fit the data (Burnham and Anderson 2002). To test for deviations from normality, the Shapiro-Wilk W test was referred to due to its greater power relative to the KolmogorovSmirnov test (Dutilleul and Legendre 1992). However, as any statistical hypothesis test is sensitive to sample size, reliable determination of true departures from normality is very difficult (Cox 1958, Johnson 1999, Neter et al. 1996). Therefore, normal probability plots were used as the primary assessment of normality. Plots of residuals over predicted values were used as the primary evaluation of fit and error variance homogeneity. In residual plots, a model with errors of equal variance, and a modelled relationship shape (e.g. linear) that fits the data well will display residuals that are evenly spread around zero (Sokal and Rohlf 1995, Neter et al. 1996). Residuals will show a highly uneven spread around zero i f an incorrect functional form was used, or if error variances are simply heterogeneous (Neter et al. 1996). After models were constructed, parameters were estimated, and agreement of models to structural assumptions was evaluated; models were then ranked. As the model representing 'truth' is assumed to be essentially of infinite complexity, model selection criteria from the AIC family were selected due to their competitive advantage over all other readily available methods in the circumstances encountered here (Hurvich and Tsai 1989). The number of parameters per model was calculated by adding together the number of variables (including the classification variable fertilizer) plus one parameter each for the intercept and the residual variance (Burnham and Anderson 2002). Because the ratio of sample size to number of parameters is considerably less than 40 (the ratio varied from 2.4 to 6.0), AICc was used to correct for small sample bias in AIC (Burnham and Anderson 2002). Equations for AIC and AICc are shown below, and follow Hurvich and Tsai (1989). Here, n refers to the number of independent observations, __ £,-  2  refers to the sum of squared errors of the regression, and K refers to the number of parameters in a given model:  Ts  2  AIC = n l o g O ^ - * - +1) + 2(K +1) n  45  A  I  C  c  = „ , o g ( ^  «  C  +  1)  +  t  ^ +2  = AIC H-  2  (  *  +  1  )  (  *  +  2  )  n-K-2  n The more the data support a model, the lower the AICc score is for that model. Note that the first term of AIC will decrease with sample size and fit. The second term acts as a penalty for complexity by increasing as the number of parameters increases, thus encouraging the selection of parsimonious models (Burnham and Anderson 2002). Once AICc values had been calculated, further analysis followed Burnham and Anderson (2002). Absolute values of AICc are fairly unimportant, as it is the AICc values relative to the best ranked candidate model in the set that is of interest. Aj values were calculated as: A= t  AIC -min  AIC  t  Akaike weights must then be calculated for each model as:  exp(-^A ) ;  Wv =  R  1  __exp(--A ) r  r=\  1  The Akaike weights represent the amount of relative support for each model considered, as provided by the data (Burnham and Anderson 2002). That relative support was then also expressed, perhaps more clearly, using evidence ratios, wherein the weight of the strongest model is divided by the weight of each successive model in turn. The relative importance of predictor variables was estimated by summing Akaike weights across all models that contained that variable. Models with Aj greater than 10 were then eliminated from further consideration, as they explain relatively little of the variation in the data (Burnham and Anderson 2002). Akaike weights were then recalculated for the reduced list of models, and used to produce a weighted average of the parameter coefficients of those models. This produces a multiple regression model that is a weighted average of the short-listed model set (Burnham and Anderson 2002). Parameter estimates are averaged by:  46  i=\ where 6 is the coefficient being averaged across models. Estimates produced by the weighted averaged model are identical to the weighted averages of estimates produced by each model separately. Once model construction, parameter estimation, and model scoring with AIC had been completed, the standard error of parameter estimates was estimated following Burnham and Anderson (2002). First, the standard error of the estimate for each parameter (represented by  )  was taken from SAS Proc G L M output (SAS Systems 2001). The variance of the model averaged parameter estimate was then calculated by:  where 6 is the model averaged parameter estimate, and 0 is the parameter estimate of each i  weighted model. The variance of the model averaged estimate is equal to the standard error of the model averaged estimate squared. The error in the fit between the final weighted average models and the observed data was calculated as the root mean squared error:  where SSE is the sum of squared errors, Y represents the observed values, Y represents the t  t  predicted values, and m is the number of independent variables in the model. However, as the same observations are being used for model fitting and error estimation, the root mean squared error is expected to produce an overly conservative estimate of prediction error. A n alternative approach involving calculation of the root mean PRESS statistic was also performed. The PRESS (prediction sum of squares) criterion estimates prediction error by calculating the sum of squared differences between each observation and the estimation of that observation with, a model whose parameters were estimated in the absence of that same observation (Neter et al. 47  1996). As each final model here is in fact a weighted average of numerous models, this involved, for each observation; removing one observation, estimating parameters for each candidate model, re-scoring and re-weighting each model with AICc, producing new weighted averaged models, and finally calculating new estimates. The root mean PRESS criterion was then calculated by:  root mean PRESS  The relative importance of predictor variables was first estimated by summing Akaike weights across all models that contained that variable (Burnham and Anderson 2002). After a priori modeling and parameter estimation was completed, the underlying relationships between standardized deer pellet-group counts and the habitat variables, as described by the data, were examined more objectively using the canonical correlation analysis procedure in SAS (Proc C A N C O R R ; SAS Systems 2001). A strong a priori filter is necessary in model building to avoid the inclusion of spurious relationships, but a dataset's potential for enhancing learning about a system is limited i f the data are not eventually examined for potentially valuable insights. One must be cautious, however, not to be tempted to pursue a highly iterative approach to model building (Burnham and Anderson 2002). Canonical correlation analysis creates canonical variables for each side of an equation that are linear dimensions of the correlations between dependent and independent variables (Ter Braak 1986, Tabachnick and Fidell 2001). This method has a distinct advantage over comparisons of univariate correlations, where relationships can become more easily distorted by covariates. Canonical correlation analysis also has a descriptive advantage over standard multiple linear regression, where coefficients are more difficult to interpret because they do not arise from standardized independent variables (Ter Braak 1986), and the importance of variables must instead be judged by the black and white determinant of statistical significance. For examining predictive relationships, canonical correlation analysis also has an advantage over principal component analysis due to its emphasis on correlations rather than variance. Multivariate normality is not required when using this method descriptively, although a linear relationship in the data is. Correlations between independent variables and their associated canonical variable are used as the units of measurement because the information they provide is insensitive to  48  correlations among the independent variables, unlike standardized canonical coefficients (Ter Braak 1986).  3.3 Results 3.3.1 Multimodel Selection and Inference For the global model that was fit using data from all three sites, analysis of normal probability plots and the Shapiro-Wilk W test suggested that the distribution of residuals did not deviate substantially (or statistically significantly) from normal for either summer or winter data. The model that was fit using only Summerland and Cariboo data combined, however, revealed a substantial lack of fit when compared to observations at the Kelowna site. Investigation of models fit only to Summerland and Cariboo data was therefore discontinued; the significance of this will be discussed later. AICc scores and weights, along with model variables and numbers of parameters, are listed in full in Table A l and Table A2 of the Appendix 1, for summer and winter data respectively. Models with A; values greater than 10 were eliminated from further consideration, as they explain relatively little of the variation in the data (Burnham and Anderson 2002). For example (see Table A2), a model with a Aj value of 10.25 has approximately 170 times less support from the data than the best model. Although 10 is an arbitrary cut-off point, it is a reasonable one. As shown in Figure 3.4, evidence ratios increase exponentially with A;, and begin to increase rapidly around Aj=10.  Figure 3.4. Exponential increase of AIC evidence ratios with increasing A; values. 49  Models with Aj values less than 10 for the full model set are listed in Table 3.3 and Table 3.4, for summer and winter data respectively. AIC weights have been re-calculated for the shortlisted model set, while evidence ratios remain as calculated for the full model set to provide contextual reference. There is a clear trend of decreasing AIC support for models of increasing complexity, with the simplest (2 variable) models receiving the best (lowest) scores for both winter and summer data. The trend is certainly not absolute, however, as there were a number of instances where more complex models possessed lower A; values than simpler models (see Table A l . l and A1.2). It also seems apparent from Tables 3 and 4 that models containing only the variables herb volume, shrub volume, and density were dominant in terms of support by the data relative to the other models considered.  Table 3.3. Summary of AICc scoring results for summer data. AIC weights are calculated considering only this short-listed model set. K represents the number of parameters per model, and A is the difference between the AICc values of the model in question and the best model (lowest AIC value). Evidence ratios express the relative support for a model by the data, and are calculated by dividing the weight of the strongest model by the weight of each successive model. ;  jm  K 17 fertilizer, shrub volume  4  5.38  0.00  0.37281  1.00  19 fertilizer, density  4  5.52  0.13  0.34880  1.07  18 fertilizer, herb volume  4  7.38  2.00  0.13736  2.71  16 fertilizer, shrub volume, herb volume fertilizer, shrub volume, herb volume, 15 density fertilizer, shrub volume, density, 9 density x shrub volume fertilizer, shrub volume, herb volume, 14 shrub species richness fertilizer, shrub volume, herb volume, 13 herb species richness fertilizer, herb volume, density, 10 density x herb volume  5  8.90  3.51  0.06441  5.79  6  10.86  5.47  0.02417  15.43  6  10.94  5.55  0.02319  16:07  6  12.47  7.08  0.01080  34.52  6  12.62  7.24  0.01000  37.30  6  12.95  7.57  0.00847  44.03  50  Table 3.4. Summary of AICc scoring results, for winter data. Details are the same as for Table 3.3.  M€  e  E  19 fertilizer, density  4  7.27  0.00  0.30538  1.00  17 fertilizer, shrub volume  4  7.31  0.05  0.29822  1.02  18 fertilizer, herb volume  4  7.48  0.22  0.27398  1.11  16 fertilizer, shrub volume, herb volume  5  10.31 3.05  0.06661  4.58  6  13.11 5.84  0.01646  18.55  6  13.59 6.32  0.01294  23.60  6  14.29 7.02  0.00911  33.53  6  14.37 7.11  0.00874  34.95  6  14.41 7.15  0.00857  35.65  fertilizer, shrub volume, herb volume, herb species richness fertilizer, shrub volume, herb volume, 15 density fertilizer, shrub volume, herb volume, 14 shrub species richness fertilizer, herb volume, density, 10 density x herb volume fertilizer, shrub volume, density, 9 density x shrub volume  13  The weighted averages of estimated parameter coefficients for models with A; values less than 10 are was used to produce final models, and are listed in Table A1.3 of Appendix 1. Also listed in Table A1.3 are calculated standard errors for each model averaged coefficient. The Akaike weights used for averaging were those recalculated to consider only the models in that reduced set.  These coefficients form multiple regression models that take into account the  numerous models most strongly supported by the data in each season, and the relative variation in parsimonious predictive power captured by each. However, as one of the objectives for modelling mule deer habitat suitability is to provide a useful land management tool, the practical utility of the final model must be taken into account. Given the weak relative support by the data and prior theory of shrub and herb species richness as predictor variables, as well as the expense of their measurement in the field, a reduced set of models was produced considering only those models with A; values less than 10 that did not include those two variables. Akaike weights were again re-calculated to consider only this further reduced set of models, and weighted averages were calculated for estimated parameters (see Table A l .3). Results of calculations describing the estimated prediction error and model fitting error in weighted averaged models from the full and reduced (excluding of species richness measures) management model set are listed in Table 3.5. Given the methods used for data handling, 51  parameter estimation, and model averaging, validation options were greatly reduced. After averaging across repeated measures, the number of independent observations remaining was too few to put aside a portion for model validation (i.e. cross-validation). Therefore, the sum of squared errors (SSE)  and root mean sum of squared errors  ) were used to quantify the  error between observed and predicted values of standardized pellet-group densities. The P R E S S (prediction sum of squares) criterion and root mean P R E S S (IMPRESS  ) were used to estimate  prediction error. The root mean sum of squared errors can be thought of as an estimator of the standard deviation of residuals (Neter et al. 1996). As shown in Table 3.5, this mean error estimate represents approximately 25% of the range of 'habitat suitability' for both models in the summer, and only slightly more in the winter. The closer P R E S S is to SSE,  the better  is as an  indicator of actual prediction error (Neter et al. 1996). As the P R E S S values are similar to SSE in summer,  will be a reasonably good estimate of prediction error for that season. In winter,  calculated P R E S S values are substantially different from SSE  and therefore -JMSE  will  generally be a poor indicator of prediction error for that season. The root mean P R E S S can be thought of as an estimate of the standard deviation of predictions from observed values. As the same observations were not used in fitting and validation in the production of the root mean P R E S S statistic, it should produce a less conservative estimate of error, and does (see Table 3.5). The root mean P R E S S indicates a prediction error in summer in excess of 30% of the range of 'habitat suitability' observations, and 40% in winter. However, as the remaining observations used in model fitting were from the same data set, root mean P R E S S should itself be a conservative estimate of the errors that would likely occur in the presence of a truly independent validation data set. For an indication of how model prediction errors (derived from P R E S S residuals) vary with standardized pellet-group densities, see Figures A4.1 to A4.4 of Appendix 4. In summer, errors show a negative residual bias (overestimation, as residuals arise from observed minus predicted values) at low deer use and a strong positive residual bias (underestimation) at high deer use. Also, residuals from summer models reveal a general trend of increasing degree of error as observed standardized pellet-group densities increase. In winter, errors are considerable both at low pellet-group densities, where the model tends to overestimate, as well as at high pellet-group densities, where the model tends to underestimate deer preference. The greatest estimation accuracy in winter is apparent at moderate levels of standardized pellet-group densities. 52  Table 3.5. Sum of squared errors (SSE) and Root mean squared error (4MSE ) of the fit between observed pellet-group densities and the values predicted by models derived as weighted averages of the full and reduced (for management applications) model sets. Also included are estimates of prediction error using PRESS and root mean P R E S S (IMPRESS ). ^nnrmoorsrp  mM  M D  mM  SSE  1.03  1.03  1.27  1.29  JMSE  0.26  0.25  0.29  0.28  PRESS  1.48  1.49  2.93  2.93  IMPRESS  0.31  0.30  0.44  0.42  As mentioned previously, the relative importance of predictor variables in the a priori model set of mule deer pellet-group counts can be estimated by summing Akaike weights across all models that contain that variable (Burnham and Anderson 2002). It is important to keep in mind that, when examining A I C indices of the relative strength of explanatory variables (see Figure 3.5), much of the variation in apparent variable strength arises not from the data, but from the importance placed on variables during a priori model construction. The summed relative strengths of variables over an a priori constructed model set therefore include a pervasive subjective component. For example, no models were put forward with herb or shrub species richness as the sole independent variables because of the weakness of support for those variables in the literature. However, exploratory analysis performed after a priori model set fitting and ranking revealed that models with only fertilizer (the classification variable) and herb or shrub species richness were among the most strongly supported models, particularly in winter. Had all variables been entered in every possible combination without any subjective weighting, the support for each variable by the data would be much closer than shown in Figure 3.5. This subjective component is prudent to apply in predictive modelling, because it reduces the probability of capturing spurious relationships, and increases the chance of producing reliable models (MacNally 2000, Burnham and Anderson 2002).  53  0.6  * 2 oj TS  o  0.5  Vt Vl  o 0.4 TJ a) 0.3 £ 01  sumrnerData  Winter Data  Vt  0.2  £  o.i  <  + + shrub  herb  volume  volume  density  ++  —+  + -  shrub  herb  densityx  densityx  species  species  shrub  herb  richness  richness  volume  volume  Figure 3.5. Summation of AIC weights over all the models containing each variable, representing a subjectively influenced estimate of the relative importance of variables. Positive and negative signs represent the calculated direction of influence on mule deer pellet-group densities. 'Density' refers to stand density in stems per ha. The value of A I C indices of the relative strength of explanatory variables becomes apparent when comparing relative strength between variables of equal a priori assigned importance (i.e. shrub volume, herb volume, and stand density), and referring to model coefficients (i.e. Table A1.3) for discerning direction of the apparent relationship (represented by positive and negative signs in Figure 3.5). Also, potential flaws in the a priori model set may be identified when variables that were appointed low a priori importance still come out as very strong relative to others. In the example presented here, no variables that were assigned as weak a priory overwhelmed the results, but comparisons of variables assigned equal importance are informative. For summer data, shrub volume was the most strongly supported variable, density was the second most important, and herb volume was third. In summer, shrub volume was 1.2 times more strongly supported by the data than density, and about twice as strongly supported as herb volume. In winter, shrub volume, herb volume, and density are similar in their estimated strength of support by the data, with density the weakest of the three. The remaining variables appear be extremely weak in comparison (see Figure 3.5). Differences between the weaker variables are small, and certain to be meaningless after sampling error and bias are considered. 54  3.3.2. Data M i n i n g with Canonical Correlation Analysis Any analysis of the data without incorporating prior biological knowledge runs the risk of providing misleading results, because the symptomatic 'noisiness' of ecological data is likely to obscure any true underlying relationships. Nonetheless, canonical correlation analysis was performed to analyze the data for any potentially important relationships that may have been missed, and to further investigate the data set's usefulness for providing insight into mule deer habitat preferences. For the sake of interest, graphs showing univariate responses of selected model variables are displayed in Appendix 3. These graphs occasionally show an approximately linear response or little response at all. Also in Appendix 3 are univariate Pearson's correlation coefficients between all available variables, which occasionally show near independence between predictor variables and standardized pellet-group counts. In general, however, these simple graphs and correlations simply reveal substantial noise and variance. However, it is risky to read too much into a univariate analysis of data that is expected to be affected by multivariate interactions. Univariate perspectives  conceal the effects  of co-varying factors, therefore  distorting  relationships. Although simple to understand, in complex situations such a reductionist approach is likely to mislead more than clarify. For this reason, only multivariate analyses are considered in depth here. Overall results of canonical correlation analysis are fairly consistent with the expectation of substantial noise, as many variables appear to have moderate multivariate correlations with pellet-group densities, and few displaying truly weak or strong relationships (see Table 3.6). For summer data, the canonical variable explained 58% of the standardized variance of deer pelletgroup counts. Of all the variables included in the analysis, total structural diversity contributed the most, with a fairly strong correlation coefficient of 0.72. The possible effects of total structural diversity on mule deer habitat seem to be convoluted and difficult to interpret. The most likely explanation for this strong correlation may be due to relatively strong relationships between total structural diversity and variables that are more directly causal to mule deer habitat preference. For example, there are fairly strong negative correlations between total structural diversity, density and tree volume, and a moderate positive correlation with shrub volume. A stand with high structural diversity appears to be fairly open, with associated good understory development; two stand traits known to be of value to mule deer in summer.  55  Table 3.6. Summary of canonical correlation analysis results for fertilized and unfertilized treatments combined. Correlations listed are between available independent variables and their canonical variable. The canonical variable here is a linear component that maximizes the correlation between the independent variables and deer pellet-group counts standardized by region.  Herb Volume  0.52  0.56  Shrub Volume  0.69  -0.07  Herb Species Richness  0.22  -0.01  Shrub Species Richness  0.13  -0.25  Stand Density  -0.45  0.29  Tree Volume  -0.24  0.31  Total Structural Diversity  0.72  -0.25  Herb Species Diversity  -0.47  -0.08  Shrub Species Diversity  -0.40  -0.08  Tree Species Richness  0.06  -0.19  Tree Species Diversity  0.54  -0.20  Total Species Richness  0.21  -0.12  Density x Shrub Volume  0.25  0.14  Density x Herb Volume  0.05  0.51  Standardized Variance of Pellet-Group Counts Explained by Canonical Variable of Independent Variables  0.58  0.53  Shrub volume was also correlated fairly strongly to the canonical variable, with herb volume showing a more moderate strength (see Table 3.6). The fact that these variables are strong relative to the others is well in line with theory. The importance of herb volume and shrub volume relative to each other has often been shown to vary greatly from location to location (Kufeld et al. 1973). For summer data, the relationships of the remaining variables to pelletgroup counts appear to be generally quite weak. For winter data (see Table 3.6), there were no strong correlations between any independent variables and the canonical variable, which explains only 53% of the standardized 56  variance of pellet-group counts. Somewhat notable within the winter data (see Table 3.6) is the low to moderate positive correlation of stand density and tree volume to the canonical variable. This is in contrast to the negative correlation of those variables in the summer data. In addition, total structural diversity shows a low to moderate negative correlation with the canonical variable, in contrast to the much stronger and positive correlation in the summer data. Clearly, in winter mule deer pellet-count densities are greater in stands that are more closed than in summer, although apparent strength of preference in winter is not as strong. The canonical correlation results for winter, which show shrub volume to be almost independent of standardized deer pellet-groups, are in sharp contrast to model fitting and ranking results, which show shrub volume to be slightly stronger in importance than herb volume (see Table 3.4 and Figure 3.5). To investigate this discrepancy, canonical correlation analysis was repeated, this time analyzing fertilized and unfertilized plots separately. As there are fewer independent observations available when the data set is split in half, further analysis cannot support as many independent variables. Therefore, only the variables included in the a priori model set were included in this additional analysis. Results are summarized in Table 3.7, where the reason for the discrepancy mentioned above becomes immediately apparent. When fertilized and unfertilized observations are combined, shrub volume is nearly independent of pellet-group counts, but when dealt with separately, shrub volume has somewhat more influence than herb volume on fertilized sites. A number of other interesting insights can be gained by separating observations by fertilization treatment (see Table 3.7). First, it is worth noting that canonical variables generally explain more of the standardized variance in pellet-group counts when fertilized and unfertilized plots are examined separately, with the exception of fertilized plots in winter. In that particular case, it appears that deer preference is being driven by factors that were not considered here. The amount of variance explained by a canonical variable is a necessary value for translating correlation strengths between predictor variables and the canonical variable into multivariate correlations between predictor variables and standardized pellet-group counts. It is immediately apparent (see Table 3.7) that both fertilization and season have a substantial impact on the relationships between pellet-group densities and the identified habitat variables. For example, in summer and winter data, shrub volume is nearly independent of deerpellet counts on unfertilized sites, yet more strongly and positively correlated on fertilized sites in summer, and more strongly and negatively correlated on fertilized sites in winter. Also in  57  winter, herb volume is nearly independent on unfertilized sites, but more strongly correlated to pellet-group counts on fertilized sites.  Table 3.7. Summary of canonical correlation analysis results for fertilized and unfertilized analyzed separately. Correlations listed are between available independent variables and their canonical variable  Herb Volume  0.42  0.20  0.09  0.40  Shrub Volume  0.06  0.69  0.09  -0.62  Herb Species Richness  0.54  0.30  0.55  -0.85  Shrub Species Richness  0.46  0.67  0.16  -0.49  Stand Density  -0.49  -0.47  0.63  0.23  Density X Shrub Volume  -0.28  0.22  0.72  -0.26  Density X Herb Volume  -0.41  -0.22  0.62  0.50  Standardized Variance of PelletGroup Counts Explained by Canonical Variable of Independent Variables  0.77  0.61  0.65  0.36  Finer details of the results summarized in Table 3.7 are interesting, but may simply typify noisy ecological data, with an abundance of weak and moderate relationships. For example, herb volume appears to be slightly more strongly correlated to deer pellet-group counts on unfertilized sites than on fertilized sites in the winter. This appears to be counterintuitive and may simply be a 'quirk' in the data, but it may be because forage quality is likely to influence habitat selection more strongly than abundance when supplies are not limiting. If this is the case, preference shown for fertilized habitats may be much greater than differences in volume measurements may suggest (Fryxell and Sinclair 1988). Data on herb nutritional content on each plot would be required to test this hypothesis. In both summer and winter (see Table 3.7), herb and shrub species richness appear to have moderate correlations to pellet-group counts. In summer, species richness measures are moderately important for all but herb species richness on fertilized sites. In winter on unfertilized sites, herb species richness has a moderate positive strength of correlation to pellet-group counts, 58  while shrub species richness shows little relationship. On fertilized sites, however, both species richness counts have a negative multivariate relationship to pellet-group counts, with herb species richness showing a relatively strong negative relationship. Also of interest are differences in density and density by volume interaction strengths. In summer, multivariate correlations between density and standardized pellet-group counts are moderate and negative, regardless of fertilizer treatment. However, in winter the multivariate correlation is positive; relatively strong on unfertilized plots, but weak on fertilized plots. When fertilizer treatment differences are ignored (see Table 3.6), the density by shrub volume interaction is very weak in winter, yet it appears to be relatively strong on unfertilized sites (see Table 3.7). Also, the density by herb volume interaction is almost independent of the canonical variable in summer when fertilizer treatment differences were ignored. However, fertilized and unfertilized sites both show negative relationships to pellet-group counts when examined separately. For an even deeper investigation, canonical correlation analysis was repeated a third time, this time separating observations by site as well as fertilizer treatment. Results suggest substantial differences between sites (see Appendix 1, Table A4). First of all, it important to note that the amount of variance in standardized pellet-group counts explained by each canonical variable varies widely by site, season, and fertilizer treatment. The predictor variables explain variance in standardized pellet-group counts at the Gavin Lake site the most effectively, with a fairly strong relationship on average. The variances explained in Summerland and Kelowna pellet-group counts are generally similar to each other, and are weak to moderate in strength. Finer details of differences by site are pervasive and readily apparent. For example, in summer, herb volume is moderately negatively correlated with standardized pellet-group counts on the Gavin Lake site on both unfertilized and fertilized plots. However, in Kelowna, herb volume is weakly negatively correlated on unfertilized plots, and strongly positively correlated on fertilized plots. In Summerland, herb volume is weakly correlated with pellet-group counts on unfertilized and fertilized sites. For another example, shrub volume in winter on fertilized sites is moderately negatively correlated with its canonical variable in Summerland, moderately positively correlated in Kelowna, and very strongly positively correlated in Gavin Lake. What makes this even more significant is the fact that the canonical variable for the Gavin Lake data on fertilized plots in winter actually explains a large proportion of the variance in pellet-group counts for that site, unlike the canonical variables for the Summerland and Kelowna sites. This suggests that unlike 59  the other two sites, shrub volume is actually quite important to mule deer in winter at Gavin Lake. However, when sites are analyzed together (see Table 3.7), shrub volume is strongly negatively correlated with a canonical variable that explains very little of the variance in pelletgroup counts. Although only a couple of examples are presented here, differences in relationships by sites between predictor variables and standardized pellet-group counts are pervasive, and occasionally fairly strong. These site differences suggest that analyzing sites together may be distorting or obscuring regionally relevant mule deer relationships with the habitat variables, resulting in weaker predictive models. However, given the paucity of independent observations, lumping sites together was necessary for model development.  3.4 Discussion With the completion of careful and well informed a priori model construction and model averaging using AIC-based weighting, habitat suitability models have been assembled and evaluated. However, habitat suitability modeling, even in the best of circumstances, is an exercise in rough approximations. Population density in a habitat, or an index of population density, is often the best available measure of habitat quality. However, habitat quality is a function of more than just population density, but also of the birth rate, death rate, and social interactions particular to an area (Van Home 1983). Unfortunately, these extra factors are much more difficult and expensive to measure. To further complicate matters, presence and abundance may be driven by additional factors that operate independently of habitat quality, such as history, weather, disease, parasites, predators, and human harvest (Lancia et al. 1982, Schamberger and O'Neil 1986). Compounding this complexity is the fact that ecosystems are constantly changing, due to both autogenic succession (Kimmins 2004) as well as random, unpredictable events (Levin 1998). Still, scientifically rigorous habitat suitability modeling is likely to be a dramatic improvement over educated guesses, and land management decisions must be made. While habitat suitability models are rough approximations in the best of circumstances, i f the data are of poor quality, i f the spatial scale of measurement was inappropriate, or i f the data were not collected with habitat quality characterization in mind, any attempt at robust model creation may be utterly futile (Guisan and Zimmerman 2000). In the case of the example presented here, data used were taken from an experiment testing select mammal responses to various thinning and fertilization regimes in young lodgepole pine stands. Data quality (i.e. measurement error and consistency) and spatial scale of measurement are unlikely to be serious impediments, although certainly the data were collected for a different purpose than that which it 60  is used for here. The greatest problem in this case is likely to revolve around the concept that wildlife-habitat relationships may become confounded when sampling is not complete along habitat resource gradients (Best and Stauffer 1986). In the best of situations, models based on empirical observations are likely to only represent the realized ecological niche, and not the fundamental niche, which limits application over varying environmental conditions (Guisan and Zimmerman 2000). Certainly any predictions made outside the range of conditions represented in the sample will have a highly uncertain range and distribution of error. There are many other inherent difficulties with the data as well, and these will require some extensive discussion to explore them fully. Data for use in habitat suitability modeling obviously should include measurements of animal response over the natural range of habitat conditions that are of interest. However, the range of habitat conditions captured in the dataset used here is of such small extent that the habitat suitability models produced are of similarly low applicability. First of all, canopy structure is a function of tree species composition, age, stand density, and site productivity (Bunnell et al. 1985), but only young lodgepole pine stands on sites of similar productivity were sampled. Also, as stands were pruned and thinned, the canopy structures of the lodgepole pine stands measured are unlikely to be represented anywhere but in intensively managed plantations. Canopy structure is an important habitat feature because it influences ground shading, which affects understory communities and therefore herbivore forage year-round, as well as snow interception, a potentially critical factor for mule deer winter habitat. In addition to canopy structure, snow interception varies with magnitude of snowfall, and, to a lesser degree, wind and slope which would all, in turn, vary with elevation, topography, and climate (Bunnell et al. 1985, McNay 1985). Understory productivity, species composition, and overall quality of forage are also likely to vary with elevation, B E C zone and B E C sub-zone. As all measured stands were composed of young lodgepole pine at mid-elevation in only two B E C sub-zones (and only one subzone each), only a very small range of mule deer habitat is captured in the data. As such, habitat suitability models based on stand density would likely not be very reliable when applied in other conditions. To make matters more challenging, the stand density measurements, which are used as an index of canopy structure, are not entirely accurate. Measurement of stand densities, during the same period that pellet-group densities were being measured, occurred only once, in 1998. As stands will naturally decrease in density as they age (Oliver and Larson 1990), stand densities each year of the study will certainly be different. This discrepancy will have affected the 61  regression equations that were derived here. The resulting error, however, will not be reflected in the error estimated using the present dataset, as all density measures are similarly affected. This issue could have been resolved by using only observations from 1998, but this would have resulted in losing representation of the seemingly random variation in deer pellet-group densities between plots over time. Under normal circumstances, stand density may be a reasonable indicator of canopy closure, as snow interception generally increases with stand density (Kittredge 1953, Kirchoff and Schoen 1987), and forage production decreases (Mysterud and Ostbye 1999). In the example presented here, however, where stand density measurements are inaccurate, and canopies have been altered by pruning, density may be a somewhat weaker index of canopy structure. Still, the importance of canopy structure suggests that it should be included in models in some way, and density has the advantage of being a stand characteristic for which data is easy to obtain. In addition to the problems with the measures of stand density, the nature of the other independent variables available for modeling is potentially problematic. For example, shrub and herb crown volume measurements naturally lump together palatable and unpalatable plant species. While there may be a strong correlation between deer habitat preference and palatable vegetation, the relationship with plant volume may be weak if palatable vegetation is not a large component. Of course, measuring total volume without separating by species makes field measurements for validation in a new area much more economical, which would increase the appeal of models for management application. There are several other characteristics of the data set used here that limit model predictive power and applicability. One factor is the reliance on pellet-group densities, which are recognized as a potentially problematic index of relative population densities. Pellet-group densities have received frequent use for assessing the relative preference of mule deer for different habitats (e.g. Robinette et al. 1952, Julander 1955, Loveless 1964, Gilbert and Wallmo 1970, Lyon and Jensen 1980, Tomm et al. 1981, Wambolt and McNeal 1987, Altendorf et al. 2001). The general theory behind the technique appears to be sound, as mule deer tend to defecate fairly continuously whenever grazing or travelling (Collins 1981). However, upon closer examination, the shortcomings of the technique are many, and have been well recognized for some time. For example, defecation rates have been shown to increase with increasing forage intake and nutrient concentration (Smith 1964), as well as forage digestibility and succulence (Arnold and Reynolds 1943, Longhurst 1954, Rogers et al. 1958, Wallmo 1981). Comparing results of pellet-group surveys between habitats may therefore be of questionable validity, as 62  intake, nutrient content, succulence and digestibility of forage will vary. This makes direct comparison between sites something that must be done cautiously. Also, as fertilizer application may increase the digestibility and nutrient concentration of plant parts (Puoli et al. 1991, Johnson et al. 2001), preferences of deer for fertilized sites, as indexed by pellet-group counts, may be exaggerated due to the resultant laxative effect. This is unlikely to be a complicating factor in the data set used here, however, as fertilized plots are not so large or so segregated that the deer defecating upon them are likely to be different individuals than those defecating on the unfertilized plots. Strong variations in individual preference and aversion for fertilized forage would cause enough of a segregation to affect results, and this is possible, although fairly unlikely. Defecation rates have been shown to decrease with the age of deer (Smith 1964), and so differing population structures between areas may cause incorrect inferences regarding habitat preference. This also may affect the utility of combining sites for regression, although the degree of difference between populations' age structures is likely to be too small to be of concern. The final concern seems to be comparing pellet-group densities as an index of deer use directly across seasons, because defecation rates are generally lower in winter than in summer (Rogers 1987). This makes sense, as forage intake by deer is known to decrease in winter (Wallmo 1981) and available forage is less digestible (Hanley and McKendrick 1985). This concern is accommodated by the separation of the data set by season. There have been few studies that have actually tested the validity of relating habitat preference to pellet-group densities, and those that exist have produced varied results. For example, Collins (1981) found that although mule deer spent three times as much time in an oldgrowth forest habitat as in a wet meadow habitat, pellet-group densities were three times as high in the wet meadow. This appears to support the hypothesis that defecation rates increase with forage succulence. However, Gilbert and Wallmo (1970) found a reasonable correlation between pellet-group counts and distribution inferred from actual mule deer observations. On the other hand, Anderson et al. (1970) found mean shrub utilization and mean fecal groups per plot to be only weakly related. The authors suggested that this may have been due to measuring preference for winter range during a mild winter, when preference is expected to be less. After considering the evidence put forth in the literature, it seems logical that using pellet-group densities as an index of habitat preference is fairly reliable as long as comparisons are made within a season, and different habitats do not vary too much in nutrient concentration, digestibility, or succulence of forage. 63  Regarding the winter habitat suitability models derived here, there is an additional factor limiting model application, namely, normal winter snow depths. On the experimental sites, these are approximately 0.5 to 1 m at the Summerland and Kelowna sites and 1.5 m at the Gavin Lake site (T.P. Sullivan, personal communication). However, the available literature suggests that snow depths in interior habitats greater than about 0.5 m generally make an area inaccessible to mule deer (Loveless 1964, Gilbert at al. 1970), and snow depths of about 0.25 m to 0.30 m are actively avoided (Loveless 1964, Martinka. 1968, Telfer 1978). Based on this information, as well as casual observations made during random mid-winter visits to some sites (T.P. Sullivan, personal communication), it seems unlikely that deer are using the experimental sites during the regularly occurring periods of deep snow. Instead, it is likely that pellet-groups counted at the end of winter were deposited before and after periods of deep snow. Although the winter habitat dealt with here is not without value, it would not be considered worthy of special management attention. This is in contrast to winter range that provides shelter from severe winter conditions, especially deep snow, which is very important to forest managers in the interior of British Columbia (e.g. Ministry of Sustainable Resource Management's Cariboo-Chilcotin Land Use Plan 1995 and Okanagan-Shuswap Land and Resource Management Plan 2001). The prevailing theory is that such habitat is very likely of critical importance to mule deer populations for access to forage at a time when accessible forage is extremely scarce (Robinette et al. 1952, Loveless 1964, Martinka 1968, Willms et al. 1976, Hanley and McKendrick 1985), and for minimizing energetic costs of locomotion at a time when nutritional intake is low (Robinette et al. 1952, Parker et al. 1984, Hanley and McKendrick 1985, Schoen and Kirchoff 1985, Bunnell et al. 1990, Shackleton 1999). These factors combine to suggest that the winter habitat suitability models produced using the data at hand may be of limited management value. Another factor that may reduce the predictive power of the models produced here is the standardization of observed mule deer pellet-group densities, which runs the risk of altering relationships between pellet-group densities and the independent variables. The choice to standardize by region (i.e. Summerland and Kelowna standardized together, Gavin Lake standardized separately) seems the best option, given the use of general linear models for parameter estimation, and the similarities in expected population densities and observed mean pellet-group densities between the Summerland and Kelowna sites relative to the Gavin Lake site (see Table 3.2). Nonetheless, it is likely that standardizing will have biased relationships in some" way.  64  The issue of differences between sites appears to go even deeper than differences in population densities. Using all three sites in the same regression equations may greatly decrease the effectiveness of the resultant models due to differences in correlative relationships between sites. This is most apparent when comparing canonical correlation results from the analysis of each site separately (see Table A4) to the analysis done with sites lumped together (see Table 3.7). While a large component of the apparent difference between sites is bound to be due to random variation, it is reasonable to assume that a substantial component is also due to legitimate differences. For example, the dominant positive relationship of shrub volume with pellet-group densities in winter on fertilized plots in the Cariboo is in sharp contrast to the weak influence of shrub volume at the other sites. When the sites are analyzed together, the result is a strong negative relationship between shrub volume and a canonical variable that explains very little of the variation in pellet-group densities. It appears that grouping sites together may very well be distorting and obscuring regionally-relevant relationships. This is really an issue, once again, of utilizing data that were never intended for habitat suitability modeling, as well as a potentially sub-optimal method of parameter estimation (which will be discussed shortly). The difference between sites acted to foil efforts to validate models properly. As mentioned previously, habitat suitability model output should be validated to evaluate usefulness (Hansen et al. 1993), ideally with independent data (Hocking 1976). Attempts to utilize independent data for validation involved using two sites for fitting (i.e. Summerland and Gavin Lake) and one site for validation (i.e. Kelowna). However, residuals from the fitted model displayed substantial non-normality and lack of fit, which would have made any efforts at predicting error estimates based on the normal distribution completely pointless. The same result occurred regardless of the combination of fitting and validation sites. It is impossible to say exactly why this occurred, but the possibilities are numerous. For example, it may be that there were too few observations in the fitting data sets. A l l sites used together for fitting provided more independent observations and resulted in distributions of residuals that were approximately normal. However, while three sites used for fitting obviously allowed more numerous independent observations, the differences should not have been as dramatic as they were. Instead, it seems possible that the normality of residuals obtained when using all sites together for parameter estimation is simply coincidence. Unfortunately, a good independent data set would be required to test this hypothesis. In addition, while data from one additional site, i f available, could be used for validation, resultant error estimates may be more or less extreme than the overall population of such sites, depending on random chance. Only with numerous 65  validation attempts or one very extensive data set could a realistic estimate be obtained. Therefore, the validation performed here, which was done without independent data, is likely as good an estimate as can realistically be made, given the circumstances. Ideally, with more numerous independent observations, region-specific regression equations could be derived, and problems with aggregating sites together,  either for  standardizing or to establish independent observations for regression, would be greatly reduced. Through simulation and experimental results, it has been suggested that sample size is the most important factor in producing accurate wildlife-habitat models (Pearce and Ferrier 2000, Hirzel and Guisan 2002). What sample size is appropriate for a given study, however, depends on the degree of wildlife response to changes in habitat variables, the amount of variation in response, and the risks of inaccuracy (Hicks and Turner 1999). However, the option of gathering a larger sample is not available, and in truth, describing the problem as simply a lack of independent observations is somewhat misleading. Given the circumstances that are present, it is quite possible that parameters should have been estimated using a different regression technique that would resolve some of the difficulties imposed by aggregating repeated observations and standardizing pellet-group counts by region. General linear modeling, which was used for parameter estimation in this thesis, is not capable of representing the random interacting effects of site or time without a proliferation of parameters (Laird and Ware 1982, Schabenberger and Pierce 2001). Therefore, pellet-group densities were standardized by region (aggregating sites) and all observations were averaged over time. Unfortunately, the aggregation of sites is likely to distort relationships, while averaging over time means that all the variation between repeated observations is discarded prior to analysis, and sample size is reduced considerably. If included, that variation could alter the relationships as described in the final model, and increase estimates of error variance (in the absence of independent validation data) to more realistic values (Raudenbush and Bryk 2002). There is occasionally a concern with repeated measures experiments that the effects of initially applied treatments will wear off after time, altering relationships (Schabenberger and Pierce 2001). This concern is not addressed by standard repeated-measures models, which assume that each effect is constant over a particular time period (Laird and Ware 1982). In the data used here, fertilizer was applied repeatedly (in 1994, 1997, 1998, 2000, and 2002), and therefore it seems unlikely that the effects of fertilization would diminish. Nonetheless, the effects of stand thinning and site may change over time, along with other relationships as the stand matures over the five years of repeated observations. 66  As an alternative to general linear models, mixed models could have been used on pelletgroup counts that were not averaged over time nor standardized by site. Mixed modeling is a flexible approach that incorporates at least two random variables (one of which is model error) with standard fixed effects (Schabenberger and Pierce 2001). Determining exactly what constitutes a random variable appears to be a topic of some confusion (see Schabenberger and Pierce 2001 for a thorough discussion). As a general rule it seems effective to state that effects are considered random if they come from a probability distribution, and they should otherwise be considered fixed (Robinson 1991). Schabenberger and Pierce (2001) state further that it is likely reasonable to specify that all experimental effects other than treatments effects should be considered random as long as the outcome is stochastic in nature. In the case of the data used here, it seems likely then that all explanatory variables would be fixed effects, and time and site would be random effects. There are clearly numerous ways that the habitat suitability models derived here could have been improved. In addition to requiring more independent observations, models would also have been improved i f more of the natural range of habitat variation was sampled equitably in each region (Guisan and Zimmermann 2000). With more natural variation sampled, maximum pellet-group density levels could be set with more confidence, and models could be applied to a larger range of conditions. Numerous sampling methodologies have been thoroughly explored, and their relative efficiencies mathematically established (see Lohr 1999 for a comprehensive discussion). The most basic approach is simple random sampling, which is self explanatory, and often inefficient unless a population is thoroughly mixed. Slightly more efficient, as it is not constrained by the difficulty of achieving true randomness, is systematic sampling. Systematic sampling involves taking samples at equal intervals, and produces results similar to simple random sampling as long as population units are not in a periodic order of variation (Lohr 1999). Stratified sampling is somewhat more complex, and involves purposefully dividing sampling effort among different strata so that each stratum has a higher probability of being evenly represented by the samples. Effort can also be allocated differently per strata to further increase efficiency in certain circumstances. In stratified random sampling with proportional allocation, sampling effort is allocated to be proportional to the size of each stratum. This increases the likelihood of representative sampling as long as the within-stratum variation is similar across strata. If the variation differs considerably, optimal allocation will be preferred, which involves sampling a greater percentage of more variable strata (Lohr 1999). A variation on stratified sampling with optimal allocation is known as gradsect sampling, and operates by 67  selecting transects that fall along the steepest environmental gradients present in the study area (Austin and Heyligers 1989). Gradsect sampling has been suggested to be the optimal methodology for most habitat suitability modeling requirements due to its efficiency (Wessels et al. 1998).  3.5 Conclusion Despite difficulties due to incompatibilities between the available data and the use to which they were applied here, models produced may yet be of use for management application. Application however, must be by necessity over a very limited set of conditions. This is mostly due to the fact that the actual sampled range of natural mule deer habitat conditions was very narrow, and therefore any application of the models produced here would be advisable only in young, thinned, pruned stands of lodgepole pine in the Montane Spruce and Sub-Boreal Spruce Biogeoclimatic Ecosystem Classes (dry mild subzones). Calculated relationships are only correlative, and stands with different overstory and understory species composition, age structure, site productivity, and management history will likely exhibit different relationships between mule deer pellet-group densities and the predictor variables that were utilized here. That being said, within the range of habitat that was sampled, the estimated standard deviation of prediction error is not too severe, especially in summer (see Table 3.5). Models tend to overestimate habitat preference at sites of low general preference in both seasons, but much more so in winter. Models tend to underestimate habitat preference at sites of high general preference, and more so in summer. However, prediction error overall appears to be more of a concern in winter. This is particularly represented by the much larger error estimates that resulted from calculation of the PRESS statistic through cross-validation relative to model fitting errors in winter. Validation with a representative and independent data set is necessary for a more reliable estimate of prediction error variance in both seasons. Problems with compatibility of the available data for their application here are pervasive. This raises an interesting issue that is sure to increase in prevalence as predictive ecosystem modeling becomes more commonplace; increased multiple uses of data sets for purposes that were not considered during data collection. To maximize scarce research resources, the design of applied ecology experiments will have to consider the breadth of demands that will come upon resultant data. That is not a criticism of the experimental design that produced the data set used in this thesis, however. In order to investigate the original research question posed (e.g. can old growth attributes be produced through intensive management of young lodgepole pine stands?), 68  the operative scale was necessary, and further replication was simply not a logistical possibility. If anything, the addition of mule deer pellet-group censuring to the existing sampling scheme was a noble attempt to maximize the utility of the resultant data in anticipation of demand. In the end, the compatibility problems between the available data and their use in this thesis are important to note, but they in no way invalidate the efforts applied here to produce predictive models with them. Indeed, the predictive models produced in Chapter 3 still may be useful for informing management of young lodgepole pine stands for deer habitat. Models may be further improved by utilizing a mixed models approach for parameter estimation. Mixed modeling has the potential to incorporate all the variation that was lost by the aggregation of observations over space and time that was done to facilitate general linear modeling in this project. Still, producing mule deer habitat suitability models that make maximal use of the data and fit the data well was not the only objective of this project. The other main objective was to gain familiarity with one of the better methods available for producing useful mathematical habitat suitability models from data, and produce a thorough, easy to follow example of it. The general procedure of a priori model set construction and model averaging with A I C is sound, and has been extensively explored and exemplified herein. The subjective filter of an a priori information emphasis in variable selection and model construction aids in coping with the symptomatic noise of ecological data in a logical manner. Using A I C (or AICc for small samples) allows an approximately unbiased estimation of the Kullback-Liebler information loss between models, thereby providing a logical foundation for model weighting and averaging that consistently and evenly incorporates the principle of parsimony. When used rigorously, the model selection and weighted averaging methods demonstrated here will aid greatly in producing robust habitat suitability models from empirical data.  69  Chapter 4 General Conclusion  4.1 Forests, Complexity, and Modeling The management of forest ecosystems is a daunting task. Forest managers must not only manage for the numerous, highly variable demands that society places upon forests (Thomas 1979), but they must also plan over long periods of time and over large, complex, heterogeneous, dynamic, and relatively enigmatic ecosystems (Christensen et al. 1996, Bunnell et al. 1999). Due to the complex nature of ecosystems, any novel management activity, or a common management action applied in a novel environment, will produce results with a high degree of uncertainty (Walters and Holling 1990). Indeed, it has been suggested that true generality may not even exist in ecosystems, and that novelty may be the norm (Fielding and Haworth 1995). Ecologists need tools to help them effectively simplify and gain understanding of natural systems, and to help them make predictions on the effects of changes to it. As the process of science is currently the best method available for learning about complex systems, it is important that the scientific method be integrated into the forest planning process (i.e. adaptive management; Walters 1986, Walters and Holling 1990). A n important tool within science is modeling, which is considered by some to be so advantageous as to be necessary for producing feasible resource management hypotheses and highlighting the uncertainties within them (Walters 1986). B y incorporating scientifically rigorous modeling, it is likely that forest ecosystem management can adapt as quickly as possible, therefore lessening the risk of causing undesirable ecological changes through clumsy manipulations. When attempting to maintain viable wildlife populations, effective management requires that species' habitat needs are well understood and predictable (Clark et al. 1993). Also, the effective management of wildlife populations is certainly aided by, i f not outright requiring specific, accurate, predictive models of wildlife-habitat relationships (O'Neil and Carey 1976). Although fraught with limitations (i.e. see Lancia et al. 1982, Van Home 1983, Schamberger and O'Neil 1986), formal habitat suitability models act as an aide to mental models, helping to facilitate critical analysis, the exchange of ideas, and comparison of the consequences of multiple possible management scenarios (Boyce and McNab 1994). Also, the associated numeric outputs  70  of simulation models are helpful because managers are often interested in meeting specific targets, and may need to defend their decisions in court (Turner et al. 1995).  4.2 Thesis Overview 4.2.1 Literature Review of M u l e Deer Habitat Requirements  In Chapter 2 (The Habitat Requirements of Rocky Mountain Mule Deer: A Review), a thorough literature review of mule deer habitat requirements was conducted in order to inform the selection of variables for habitat suitability modeling as well as is realistically possible. Winter habitat has often been a focus for management attention, largely due to the conspicuous relationship between mule deer population declines and winter severity (Wallmo 1981). The variables suggested by the literature for describing winter habitat, in probable order of importance are; snow depth due to its effect of covering forage and increasing the energy cost of locomotion; mean crown completeness, crown closure, stand age, stand density, and overstory tree species for their effects on snow interception and ground shading; sinking depth in snow for its effect on increasing energy costs of locomotion; shrub volume followed by herb volume for their importance in providing maintenance energy; operative temperature for its effect on energy costs of thermoregulation; vegetation height for wind reduction and lowering energy costs of thermoregulation; distance from edge for proximity to security cover and access to productive forage; and forage species richness for access to all minerals and nutrients required in the diet. The variables suggested by the literature for describing summer habitat, in probable order of importance, are: herb volume followed by shrub volume for sustenance; distance to edge for proximity to security cover and access to productive forage; canopy closure, for shade from the heat of summer as a benefit, and affecting forage production; and forage species richness for access to all required minerals and nutrients in the diet. Limited by the availability of variables, those selected for habitat suitability model construction were overstory stand density, shrub volume, herb volume, shrub species richness and herb species richness  4.2.2 Modeling and M o d e l Results  In Chapter 3 (Producing the Most Robust Habitat Suitability Models Possible from Data Using Multimodel Selection and Inference), the philosophy and methodology of one approach towards rigorous mathematical habitat suitability modeling was briefly reviewed, and an example was presented. It appears that one of the best methods currently available for producing 71  habitat suitability models from empirical data involves a process of building multiple models based on information known a priori to data analysis, scoring and averaging those models using Aikaike's information criterion (AIC), and then estimating error in the best way possible given available data and resources. Methods were outlined in detail. In the end, two sets of habitat suitability models were produced: one model for each for summer and winter seasons involving all selected variables, and one model for each of summer and winter seasons excluding forage species richness measures. The latter models were built with management application in mind, as species richness measures would be difficult to obtain over. large areas, and are difficult to estimate from generally available forest cover data. Resultant models will hopefully be found to be of practical use for management application, but will only be appropriate in a fairly limited set of conditions. Due to the relatively small extent of mule deer habitat captured in the data, it is likely that predictive error will increase greatly beyond that estimated i f models were to be applied outside of young, thinned, pruned, lodgepole pine stands in the Montane Spruce and Sub-Boreal Spruce Biogeoclimatic Ecosystem Classes (dry mild subzones) in the interior of British Columbia. As modeled relationships were primarily correlative and not causative, it is to be expected that they will vary with environmental conditions from one region to another. Error sum of squares estimates made with cross-validation suggest that error variance is not too severe, particularly in summer. This is despite data analysis results that suggest substantial differences in habitat preference relationships between replicate sites. However, cross-validation error variance estimates are still likely to be overly conservative, as fitting and validation observations arise from areas likely to be more similar to each other that they would be to a randomly sampled observation from the larger population of relevant sites. Reliable estimates of prediction error therefore require validation with a completely independent data set. Unfortunately, insufficient data were available for site specific prediction models or an independent test of model predictive accuracy because repeated observations were averaged over time. A mixed models approach to parameter estimation may result in enough observations for site specific models, representation of potentially meaningful changes over time and space, as well as sufficient independent observations for proper evaluation of the models' predictive accuracy. Should the models produced here not find application in a management context, the exercise was still of tremendous value. Modeling is a necessary tool for ecosystem management, and habitat suitability modeling, despite its shortcomings, holds great value for providing insight 72  into the potential repercussions of land management activities. The methods of model selection outlined in Chapter 3 are broadly applicable, and are among the best currently available for producing habitat suitability models from empirical data. As such, an explicit and concise review of those methods, with step-step-by step instructions, and all reinforced by the detailed example provided, should be of great value as a learning tool and guide for application in other modeling endeavours.  4.3 Applications and Future Work 4.3.1 Application Without Simulation Over Time As they currently are, parameterized models derived in Chapter 3 (see Appendix 1, Table A3 for parameter values of multiple linear regression models) should be applied to young, intensively managed lodgepole pine stands in the Montane Spruce and Sub-Boreal Spruce Biogeoclimatic Ecosystem Classes (dry mild subzones) of the interior of British Columbia. Optimum habitat suitability ratings (i.e. habitat suitability scores near 1) do not necessarily represent optimum habitat, but rather the habitat structure that is estimated to be most preferred in young, intensively managed lodgepole pine stands in the sampled ecotypes. Application of these models is likely to provide useful comparison or contrast with expert opinion, as well as the quantified estimates of mule deer response that expert opinion is less well suited for. Wherever possible, it would be of great advantage to validate the models parameterized in Chapter 3 with truly independent data. Data from a wide variety of locations with similar stand and biogeoclimatic characteristics would provide a fairly robust evaluation of the predictive accuracy of these models in the circumstances that they are best suited for. Diverse data from areas with different stand and environmental conditions would be of great use for demonstrating how predictive accuracy varies when models are applied outside the limited range of the calibration data. A diversity of locations for validation data collection is purposely stressed above, as a narrow set of independent data is more likely to be skewed by atypical relationships between mule deer and their habitat. If enough additional data are collected or available, it is possible that they may be pooled with the data used in Chapter 3 and models can be re-parameterized. Ideally, region specific models could be produced, which I believe have the potential to substantially improve predictive accuracy. Without collecting new data, as previously mentioned, a mixed modeling approach has the potential to increase the number of observations available for parameter estimation, while  73  eliminating complications that have arisen due to averaging over repeated observations and standardizing by region. A static estimate of a forest stand's habitat suitability for certain species is bound to be useful in certain circumstances. More often, however, estimates of habitat suitability are likely to be needed in the context of forest planning; over long periods of time and large areas of land. Because of the dynamic nature of forests, effective estimates of how habitat suitability in an area will change over time requires linking habitat suitability models to models of forest growth and succession (Brand et al. 1986, Holt et al. 1995). To be compatible with forest planning, forest simulation models must be able to represent the impacts of alternative management activities, and must be at a spatial scale and resolution of detail necessary to represent critical relationships (Smith 1986).  4.3.2 Stand-Level Simulation Models Just exactly what level of detail is required in stand-level vegetation simulation models depends on the factors being considered for manipulation by forest management, the habitat attributes of importance to wildlife, and the manner in which the two interact. Models will likely have to be spatially explicit for accurate predictions i f the habitat at the scale of interest is not homogenous, i f organisms are not continuously well mixed, i f the organisms are not highly mobile, or i f organisms interact with each other and their environment over long distances (Law et al. 2000). For large mammals especially, it seems likely that aspatial approximations at a stand-level scale should be fairly effective, except perhaps under non-uniform cutting patterns (i.e. aggregated retention). This would of course depend on the stand structures represented in the calibration data set. FORECAST is an aspatial, stand-level simulation model that operates using a hybrid of bioassay and process-based approaches to predicting forest stand change and development (Kimmins et al. 1999). It is a powerful and flexible decision support tool, allowing representation of numerous adjustable forest management activities (e.g. thinning, harvesting, slashburning, fertilization, etc.) and numerous adjustable levels of complexity (e.g. soil, bryophyte, plant, and tree layers, single or multiple species). In addition, FORECAST incorporates economics. As ecological and economic considerations are highly interconnected in forest management (Lippke and Bishop 1999), models that can calculate both ecological and economic factors are certainly useful (Turner et al. 1995, Hansen et al. 1995), and occasionally required (Liu 1993). Also, as proper sustainable ecosystem management must be ecologically viable, economically feasible, 74  and socially responsible (Salwasser 1992), F O R E C A S T is additionally functional by also estimating the impacts of management scenarios on employment. To further improve the functionality of FORECAST, a habitat suitability component has been designed and directly integrated (see Appendix 2 for description). This component allows for the flexible definition of habitat variables and habitat suitability equations. Equations can be designed in a traditional Habitat Suitability Index fashion (HSI; U.S. Fish and Wildlife Service 1981), or as. empirically derived linear or non-linear equations. Although economic considerations are of pivotal importance to ecosystem management, economic repercussions of habitat management options are rarely incorporated into habitat suitability models (Marzluff et al. 2002). Through FORECAST, however, both can be examined concurrently.  4.3.3 Landscape-Level Simulation Models  For mule deer, which are highly mobile animals, planning at the landscape level probably does not require explicit spatial considerations. Exceptions are likely to occur during heavy snowfall period in winter, when the arrangement of stands providing a good balance of snow interception and forage opportunities may be of critical importance. Available data for the habitat suitability equations developed in Chapter 3 did not allow for spatial analysis, and resultant equations operate under the assumption that the importance of spatial arrangement is either immaterial or implicit. For many other species, however, aspatial approximations at the landscape level are less likely to be appropriate. The distribution of habitat has the potential to have powerful biological repercussions (Fahrig and Merriam 1985). Forest fragmentation, which produces a series of residual habitat patches surrounded by a matrix of altered land, has been recognized as a potentially significant force of change to wildlife habitat for some time now (Saunders et al. 1991). Fragmentation causes changes to microclimate (i.e. increased edge effects) and the isolation of habitat patches, both of which are modified by the size, shape, and position of the fragments (Saunders et al. 1991). In general, whether or not fragmentation is of importance in affecting habitat depends largely on a particularly species' ability to disperse, and the scale at which it perceives the landscape (Fahrig and Paloheimo 1988, D'Eon et al. 2002). At the very least, a landscape perspective is helpful for beginning to anticipate how activities in one area may influence the properties of nearby areas (Turner et al. 1995). Often, geographic information system (GIS) layers are used to project habitat suitability models over a landscape (e.g. Dale et al. 1998, Woolf et al. 2002, Store and Jokimaki 2003). As most GIS data 75  is very course in resolution, GIS-based habitat suitability analysis may be most effective for species with general habitat requirements (Clark et al. 1993). On its own however, GIS is not capable of dynamic modeling over time, but instead produces a static map. For projection over time, GIS layers must be combined into a landscape simulation program, such as Possible Forest Futures (PFF) or SIMFOR (Daust and Sutherland 1996). PFF functions primarily as an interface for the management of a landscape that is an amalgamation of independently operating stands simulated by FORECAST. At present, the habitat suitability component of F O R E C A S T has yet to be fully integrated into PFF. SIMFOR, however, is by its very nature a habitat suitability projection program. It utilizes GIS data layers, customizable estimates of vegetation dynamics, customizable habitat suitability relationships, and a separate harvest scheduler for simulating the effects of forest harvest activities. SIMFOR is also capable of simple spatial calculations, such as representing the impacts of patch size or edge effects.  4.4 Conclusion Management for timber and management for wildlife are generally compatible as long as wildlife habitat requirements are understood, and incorporated in planning along with requirements for timber management (Thomas 1979). Habitat suitability modeling is clearly an important decision support tool that helps to integrate wildlife habitat management into the overall management of forest ecosystems. The statistical methods for building mathematical models from empirical data that have been outlined and exemplified. within this thesis are reasonable. However, it seems likely that estimating parameters using a mixed modeling approach would have taken better advantage of the information contained in the available data. Nonetheless, the model selection approach of an a priori model building emphasis, using A I C for model ranking and producing weighted average parameter estimates, is sound. In all predictive models, assumptions and error are unavoidable, and therefore application requires careful consideration and common sense. Model producers have the responsibility to explain the processes, assumptions, strengths and weaknesses, of their models (Thomas 1986),. and I have attempted to do this to the best of my ability. Users of any model, at the very least, should understand the inherent model assumptions, the variability in the data used, and the best estimates of the uncertainty that is present in all model predictions (Thomas 1986).  76  References Cited Akaike, H . 1973. Information theory as an extension of the maximum likelihood principle. Pp. 267-281 in Petrov, B.N., F. Csaki, (eds.) Second International Symposium on Information Theory. Akademiai Kiado, Budapest. Akaike, H . 1981. Likelihood of a model and information criteria. Journal of Econometrics 16:3-14. Alberta Forestry, Lands and Wildlife, Fish and Wildlife Division. 1989. Management plan for mule deer in Alberta. Wildlife Management Planning Series Number 1. Altendorf, K . B . , J.W. Laundre, C A . Lopez Gonzalez, J.S. Brown. 2001. Assessing effects of predation risk on foraging behavior of mule deer. Journal of Mammology 82(2):430-439. Ammann, A.P., R.L. Cowan, C L . Mothershead, B.R. Baumgardt. 1973. Dry matter and energy intake in relation to digestibility in white-tailed deer. Journal of Wildlife Management 37(2): 195-201. Anderson, A . E . , D.E. Medin, D.C. Bowden. 1972. Mule deer numbers and shrub yieldutilization on winter range. Journal of Wildlife Management 36(2):571-578. Anderson, A.E., D.E. Medin, D.C. Bowden. 1972. Mule deer fecal group counts related to site factors on winter range. Journal of Range Management 25(l):66-68 Anderson, A.E., D.C. Bowden, D.E. Medin. 1990. Indexing the annual fat cycle in a mule deer population. Journal of Wildlife Management 54(4):550-556. Anderson, D.R., K . P . Burnham, W.L. Thompson. 2000. Null hypothesis testing: problems, prevalence, and an alternative. Journal of Wildlife Management 64(4):912-923. Anderson, D.R., K.P. Burnham, W.R. Gould, S. Cherry. 2001. Concerns about finding effects that are actually spurious. Wildlife Society Bulletin 29(1):311-316. Anderson-Sprecher, R. 1994. Model comparisons and R . The American Statistician 48(2): 1132  117. • Apps, C D . , B . N . McLellan, T.A. Kinley, J.P. Flaa. 2001. Scale-dependent habitat selection by mountain caribou, Columbia Mountains, British Columbia. Journal of Wildlife Management 65(l):65-77. Armleder, H . M . , R.J. Dawson. 1992. Logging on mule deer winter range: an integrated management approach. Forestry Chronicle 68(1): 132-137.  77  Armledder, H . M . , M J . Warerhouse, D.G. Keisker, R J . Dawson. 1994. Winter habitat use by mule deer in the central interior of British Columbia. Canadian Journal of Zoology 72: 17211725. Armledder, H . M . , M J . Waterhouse, R J . Dawson, K . E . Iverson. 1998. Mule deer response to low-volume partial cutting on winter ranges in central interior of British Columbia. British Columbia Ministry of Forests Research Report #16. i Arnold, J.F., H.G. Reynolds. 1943. Droppings of Arizona and antelope jack rabbits and the "pellet census". Journal of Wildlife Management 7(3):323-327. Atilgan, T. 1996. Selection of dimension and basis for density estimation and selection of dimension, basis and error distribution for regression. Communications in Statistics: Theory and Methods A25(l): 1-28. Austin, M.P., P.C. Heyligers. 1989. Vegetation survey design for conservation: gradsect sampling of forests in north-east New South Wales. Biological Conservation 50:13-32. Austin, D.D. and P J . Urness. 1985. Values of four communities for mule deer on ranges with limited summer habitat. Journal of Range Management 38(2):167-171. Ballard, W.B., D. Lutz, T.W. Keegan, L . H . Carpenter, J.C. deVos, Jr. 2001. Deer-predator relationships: a review of recent Noth American studies with emphasis on mule and blacktailed deer. Wildlife Society Bulletin 29(1):99-115. Bartmann, R . M . 1984. Estimating mule deer winter mortality in Colorado. Journal of Wildlife Management 48(l):262-267. Battin, J. 2004. When good animals love bad habitats: ecological traps and the conservations of animal populations. Conservation Biology 18(6): 1482-1491. Bedrick, E.J., C-L. Tsai. 1994. Model selection for multivariate selection in small samples. Biometrics 50:226-231. Bennett, L.J., P.F. English, R. McCain. 1940. A study of deer populations by use of pelletgroup counts. Journal of Wildlife Management 4(4):398-403. Best, L.B., D.F. Stauffer. 1986. Factors confounding evaluation of bird-habitat relationships. Pp. 209-216 in Verner, J., M . L . Morrison, C J . Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Box, E.P. 1979. Some problems of statistics and everyday life. Journal of the American Statistical Association 74(365): 1-4. Box, G.E.P., D.R. Cox. 1964. A n analysis of transformations. Journal of the Royal Statistical Society, Series B 26(2):211-252. 78  Box, G.E.P., D.R. Cox. 1982. A n analysis of transformations revisited, rebutted. Journal of the American Statistical Association B 77(377):209-210. Box, G.E.P., G.C. Tiao. 1962. A further look at robustness via Bayes' theorem. Biometrika 49(3 /4):419-432. Boyce, M.S. 1992. Population viability analysis. AnnualReview of Ecology and Systematics 23:481-506. Boyce, S.G., and McNab, W.H. 1994. Management of forested landscapes: simulations of three alternatives. Journal of Forestry 92(l):27-32. Brand, G.J., S.R. Shifley, L.F. Ohmann. 1986. Linking wildlife and vegetation models to forecast the effects of management. Pp. 383-387 in Verner, J., M . L . Morrison, C J . Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Brieman, L . 1995. Better subset regression using the nonnegative garrotte. Technometrics 37(4):373-384. British Columbia Ministry of Environment, Land and Parks. 1998! Ground-based inventory methods for selected ungulates: moose, elk, and deer. Standards for Components of British Columbia's Biodiversity No.33. British Columbia Ministry of Forests. 1994. Catalogue of curves for curve fitting. Biometrics Information Handbook No.4. 110pp. http://www.for.gov.bc.ca/hfd/pubs/Docs/Bio/Bio04.htm Bronowski, J. 1978. The origins of knowledge and imagination. Yale University Press, New Haven. Buckland, S.T., K . P . Burnham, N . H . Augustin. 1997. Model selection: an integral part of inference. Biometrics 53:603-618. Bunnell, F.L., R.S. McNay, and C S . Shank. 1985. Trees and snow: the deposition of snow on the ground: a review and quantitative synthesis. Research, Ministries of Environment and Forests, IWLFR-17. Victoria, B.C. Bunnell, F.L., K . L . Parker, L.L. Kremsater, F.W. Hovey. 1986. Thermoregulation and thermal cover of deer and elk on Vancouver Island: problem analysis. Research, B . C . Minsitry of Environment and Parks and B.C. Minsitry of Forests and Lands, rWTFR-28. Victoria, B.C. Bunnell, F.L., F.W. Hovey, R.S. McNay, K . L . Parker. 1990. Forest cover, snow conditions, and black-tailed deer sinking depths. Canadian Journal of Zoology 68:2403-2408. Bunnell, F.L., L . L . Kremsater, E. Wind. 1999. Managing to sustain vertebrate richness in forests of the Pacific Northwest: relationships within stands. Environmental Reviews 7:97-146. 79  Burnham, K.P., D.R. Anderson. 2002. Model Selection and Multimodel Inference. Springer, New York. Campbell, E.G. and R.L. Johnson. 1983. Food habits of mountain goats, mule deer, and cattle on Chopaka Mountain, Washington, 1977-1980. Journal of Range Management 36(4):488-491. Carpenter, L.H., O.C. Wallmo, and R.B. Gill. 1979. Forage diversity and dietary selection by wintering mule deer. Journal of Range Management 32(3):226-229. Carson, R . G . , and J . M . Peek. 1987. Mule deer habitat selection patterns in northcentral Washington. Journal of Wildlife Management 51(1) :46-51. Cavanaugh, J.E. 1999. A large-sample model selection criterion based on Kullback's symmetric divergence. Statistics and Probability Lectures 42:333-343. Cavanaugh, J.E. 2004. Criteria for linear model selection based on Kullback's symmetric divergence. Australian and New Zealand Journal of Statistics 46(2):257-274. Chamberlin, T.C. 1890. The method of multiple working hypotheses. Science 15(366):92-96. Chatfield, C. 1995. Model uncertainty, data mining and statistical inference. Journal of the Royal Statistical Society A 158(3):419-466. Christensen, N.L., A . M . Bartuska, J.H. Brown, S. Carpenter, C. D'Antonio, R. Francis, J.F. Franklin, J.A. MacMahon, R.F. Noss, D.J. Parsons, C H . Peterson, M . G . Turner, R . G . Woodmansee.  1996. The report of the Ecological Society of America Committee on the  scientific basis for ecosystem management. Ecological Applications 6(3): 665-691. Clark, J.D, J.E. Dunn, K . G . Smith. 1993. A multivariate model of female black bear habitat use for a Geographic Information System. Journal of Wildlife Management 57:519-526. Coe, P.K., B.K. Johnson, J.W. Kern, S.L. Findholt, J.G. Kie, M J . Wisdom. 2001. Responses of elk and mule deer to cattle in summer. Journal of Range Management 54:A51-A76. Collins, W.B. 1981. Habitat preferences of mule deer as rated by pellet-group distributions. Journal of Wildlife Management 45(4):969-972. Collins, W.B., and P J . Urness. 1983. Feeding behaviour and habitat selection of mule deer and elk on northern Utah summer range. Journal of Wildlife Management 47(3):646-663. Constan, K . J . 1972. Winter foods and range use of three species of ungulates. Journal of Wildlife Management 36(4): 1068-1076. Conroy, M J . 1993. The use of models in natural resource management: prediction, not prescription. Transactions of the North American Wildlife and Natural Resources Conference 58:509-519..  80  Conroy, M . J . , Y . Cohen, F.C. James, Y . G. Matsinos, B . A . Maurer. 1995. Parameter estimation, reliability, and model improvement for spatially explicit models of animal populations. Ecological Applications 5(1): 17-19. Cote, S.D., T.P. Rooney, J-P. Tremblay, C. Dussault, D . M . Waller. 2004. Ecological impacts of deer overabundance. Annual Review of Ecology, Evolution, and Systematics 35:113-147. Cox, D.R. 1958. Some problems connected with statistical significance. The Annals of Mathematical Statistics 29(2):357-372. Dale, V . H . , A . W . King, L . K . Mann, R.A. Washington-Allen, R.A. McCord. 1998. Assessing land-use impacts on natural resources. Environmental Management 22(2): 203-211. Daust, D.K., and Sutherland, G.D. 1996. SIMFOR: Software for Simulating Forest Management and Assessing Biodiversity. Pp. 15-29 in: The status of forestry/wildlife decision support systems in Canada: symposium proceedings. Thomson, J.D. (ed). 1994. Natural Resources Canada, Canadian Forest Service, Great Lakes Forestry Center, Sault Ste. Marie, Ontario. Dawson, R.J., H . M . Armleder, M.J. Waterhouse. 1990. Preferences of mule deer for Douglasfir foliage from different sized trees. Journal of Wildlife Management 54(3):378-382. D'Eon, R.G. 2001. Using snow-track surveys to determine deer winter distribution and habitat. Wildlife Society Bulletin 29(3):879-887. D'Eon, R . G . 2004. Snow depth as a function of canopy cover and other site attributes in forested ungulate winter range in southeast British Columbia. B.C. Journal of Ecosystems and Management 3(2): 9 pp. D'Eon, R.G., S.M. Glenn, I. Parfait, M.-J. Fortin. 2002. Landscape connectivity as a function of scale and organism vagility in a real forested landscape. Conservation Ecology 6(2): 10. http://www.consecol.org/vol6/iss2/artl0. D'Eon, R.G., R. Serrouya. 2005. Mule deer seasonal movements and multi-scale resource selection using GPS radiotelemetry. Journal of Mammalogy, in press; 32 pp. Deschamp, J.A., P.J. Urness, D.D. Austin. 1979. Summer diets of mule deer from lodgepole pine habitats. Journal of Wildlife Management 43(1):154-161. Donovan, T.M., F.R. Thompson III. 2001. Modeling the ecological trap hypothesis: a habitat and demographic analysis for migrant songbirds. Ecological Applications 11(3):871-882. Draper, D. 1995. Assessment and propagation of model uncertainty. Journal of the Royal Statistical Society, Series B 57(l):54-97.  81  Dutilleul, P., P. Legendre. 1992. Lack of robustness in two tests of normality against autocorrelation in sample data. Journal of Statistical Computation and Simulation 42:79-91. Dwernchuk, L.W., D.A. Boag. 1972. Ducks nesting in association with gulls: an ecological trap? Canadian Journal of Zoology 50:559-563. Edwards, R . Y . 1956. Snow depths and ungulate abundance in the mountains of western Canada. Journal of Wildlife Management 20(2): 159-168. Edge, W.D., S.L. Olson-Edge, L.L. Irwin. 1990. Planning for wildlife in national forests: elk and mule deer habitats as an example. Wildlife Society Bulletin 18:87-98. Fahrig, L . , G. Merriam. 1985. Habitat patch connectivity and population survival. Ecology 66(6):1762-1768. Fahrig, L., J. Paloheimo. 1988. Effect of spatial arrangement of habitat patches on local population size. Ecology 69(2): 468-475. Farmer, A . H . , M.J. Armbruster, J.W. Terrell,. R.L. Schroeder. 1982. Habitat models for landuse planning: assumptions and strategies for development. Transactions of the North American Wildlife and Natural Resources Conference 47:47-56.Hilborn, R, M . Mangel. 1997. The ecological detective. Princeton University Press, Princeton, N.J. Fielding, A . H . , J.F. Bell. 1997. A review of methods for the assessment of prediction errors in conservation presence/ absence models. Environmental Conservation 24(l):38-49 Fielding, A . H . , P.F. Haworth. 1995. Testing the generality of bird-habitat models. Conservation Biology 9(6):1466-1481. Flack, V.F., P.C. Chang. 1987. Frequency of selecting noise variables in subset regression analysis: a simulation study. The American Statistician 41(l):84-86. Freeland, W.J., D . H . Janzen. 1974. Strategies in herbivory by mammals: the role of plant secondary compounds. American Naturalist 108:269-289. Fretwell, S.D., H . L . Lucas, Jr. 1970. On territorial behaviour and other factors influencing habitat distribution in birds. I. Theoretical development. Acta Biotheoretica 19:16-36. Fryxell, J:M., A.R.E. Sinclair. 1988. Causes and consequences  of migration by large  herbivores. Trends in Evolution and Ecology 3(9):237-241. Fujikoshi, Y , K . Satoh. 1997. Modified A I C and C  p  in multivariate linear regression.  Biometrika 84(3):707-716. Gilbert, P.F., O.C. Wallmo, B . Gill. 1970. Effect of snow depth on mule deer in Middle Park, Colorado. Journal of Wildlife Management 34(l):15-23.  82  Gorman, J.W., R J . Toman. 1966. Selection of variables for fitting equations to data. Technometrics 3(1):27-51. Guisan, A . , N . E . Zimmerman. 2000. Predictive habitat distribution models in ecology. . Ecological Modelling 135:147-186. Hamel, P.B., N . D . Cost, R . M . Sheffield. 1986. The consistent characteristics of habitats: a question of scale. Pp. 121-138 in Verner, J., M . L . Morrison, C . J . Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Hamlin, K . L . , S J . Riley, D. Pyrah, A.R. Dood, R J . Mackie. 1984. Relationships among mule deer fawn mortality, coyotes, and alternate prey during summer. Journal of Wildlife' Management 48(2):489-499. Hanley, T.A. and J.D. McKendrick. 1985. Potential nutritional limitations for black-tailed deer in a spruce-hemlock forest, southeastern Alaska. Journal of Wildlife Management 49(1): 103114. Hansen, A J . , S.L. Garman, B. Marks, D.L. Urban. 1993. A n approach for managing diversity across multiple-use landscapes. Ecological Applications 3: 481-496. Hansen, A.J., S.L. Garman, J.F. Weigand, D.L. Urban, W.C. McComb, M . G . Raphael. 1995. Alternative silvicultural regimes in the Pacific Northwest: simulations of ecological and economic effects. Ecological Applications 5(3): 535-554. Hansen, R . M . , L.D. Reid. 1975. Diet overlap of deer, elk, and cattle in southern Colorado. Journal of Range Management 28(l):43-47. Hassell, M.P., R . M . May. 1974. Aggregation of predators and insect parasites and its effect on stability. TheJournal of Animal Ecology 43(2):567-594. Hicks, C.R., K . V . Turner. 1999. Fundamental concepts in the design of experiments. Oxford University Press, New York. Hilborn, R, M . Mangel. 1997. The ecological detective. Princeton University Press, Princeton, NJ. Hicks, C.R., K . V . Turner, Jr. 1999. Funadamental comcepts in the design of experiments. Oxford University Press, New York, N Y . Hirzel, A . A . Guisan. 2002. Which is the optimal sampling strategy for habitat suitability modelling. Ecological Modelling 157:331-341. Hobbs, N.T., D . L . Baker, R.B. Gill. 1983. Comparative nutritional ecology of montane ungulates during winter. Journal of Wildlife Management 47(1): 1-16. 83  Hobbs, N.T. 1989. Linking energy balance to survival in mule deer: development and test of a simulation model. Wildlife Monographs 101:1-39. Hobbs, N.T. 2003. Challenges and opportunities in integrating ecological knowledge across scales. Forest Ecology and Management 181:223-238. Hobbs, N.T., T.A. Hanley. 1990. Habitat evaluation: do use/ availability data reflect carrying capacity? Journal of Wildlife Management 54(4):515-522 Hocking, R.R. 1976. The analysis and selection of variables in linear regression. Biometrics 32(l):l-49. Hodges, J.S. 1987. Uncertainty, policy analysis and statistics. Statistical Science 2(3):259-275. Holling, C S . 1992. Cross-scale morphology, geometry, and dynamics of ecosystems. Ecological Monographs 62(4):447-502. Holt, R.D., S.W. Pacala, T.W. Smith, J. Liu. 1995. Linking contemporary vegetation models with spatially explicit animal population models. Ecological Applications 5(1): 20-27. Homyack,  J.A., D.J. Harrison, W.B. Krohn.  2004.  Structural differences  between  precommercially thinned and unthinned conifer stands. Forest Ecology and Management 194:131-143. Huggard, D.J. 1999. Developing animal-habitat models for management of high-elevation forests. Ph.D. Thesis. University of British Columbia, Vancouver, B.C. Huntly, N . 1991. Herbivores and the dynamics of communities and ecosystems. Annual Review of Ecology and Systematics 22:477-503. Hurlbert, S.H. 1984. Pseudoreplication and the design of ecological field experiments. Ecological Monographs 54(2): 187-211. Hurvich, C M . , R. Shumway, C - L . Tsai. 1990. Improved estimators of Kullback-Leibler information for autoregressive model selection in small samples. Biometrika 77(4):709-719. Hurvich, C M . , C-L. Tsai. 1989. Regression and time series model selection in small samples. Biometrika 76(2):297-307. Hurvich, C M . , C-L. Tsai. 1990. The impact of model selection on inference in linear regression. The American Statistician 44(3):214-217. Irby, L.R., J.E. Norland, J.A. Westfall Jr., M . A . Sullivan. 2002. Evaluation of a forage allocation model for Theodore Roosevelt National Park. Journal of Environmental Management 64:153-169. Jackson, L . J., A . S . Trebitz, K . L . Cottingham. 2000. A n introduction to the practice of ecological modeling. Bioscience 50(8):694-706. 84  James, F.C., C E . McCulloch. 1990. Multivariate analysis in ecology and systematics: panacea or Pandora 's Box? Annual Review of Ecology and Systematics 21:129-166. Jameson, D.A. 1967. The relationship of tree overstory and herbaceous understory vegetation. Journal of Range Management 20:247-249. Jefferys, W.H., J.O. Berger. 1992. Ockham's razor and Bayesian analysis. American Scientist 80:64-72. Johnson, C.R., B.A. Reiling, P. Mislevy, M . B . Hall. 2001. Effects of nitrogen fertilization and harvest date on yield, digestibility, fiber, and protein fractions of tropical grasses. Journal of Animal Science 79(9):2439-2448. Johnson, D.H. 1999. The insignificance of statistical significance testing. Journal of Wildlife Management 63(3):763-772. Johnson, J.B., K . S . Omland. 2004. Model selection in ecology and evolution. Trends in Ecology and Evolution 19(2): 101-108. Johnson, C.R., B . A . Reiling, P. Mislevy, M . B . Hall. 2001. Effects of nitrogen fertilization and harvest date on yield, digestibility, fiber, and protein fractions of tropical grasses. Journal of Animal Science 79(9):2439-2448. Julander, O. 1955. Deer and cattle range relations in Utah. Forest Science 1(2): 130-139. Julander, O., W.L. Robinette, D.A. Jones. 1961. Relation of summer range condition to mule deer herd productivity. Journal of Wildlife Management 25(l):54-60. Keay, J.A. and J.M. Peek. 1980. Relationship between fires and winter habitat of deer in Idaho. Journal of Wildlife Management 44(2):372-380. Kennard, R.W. 1971. A note on the C statistic. Technometrics 13(4):899-900. p  Kessler, W.B., H . Salwasser, C.W. Cartwright, Jr., J.A. Caplan. 1992. New perspectives for sustainable natural resource management. Ecological Applications 2(3):221-225. Kie, J.G., C J . Evans, E.R. Loft, J.W. Menke. 1991. Foraging behavior by mule deer: the influence of cattle grazing. Journal of Wildlife Management 55(4):665-674. Kie, J.G. 1996. The effects of cattle grazing on optimal foraging in mule deer. Forest Ecology and Management 88:131-138. Kim, H-J., J.E. Cavanaugh. 2004. Model selection criteria based on Kullback information measures for non-linear regression. Journal of Statistical Planning and Inference In Press, Corrected Proof: 18pp.  Kimmins, J.P., D. Mailly, B. Seely. 1999. Modelling forest ecosystem net primary production: the hybrid simulation approach used in FORECAST. Ecological Modelling 122: 195-224. 85  Kimmins, J.P. 2004. Forest Ecology: a foundation for sustainable forest management and environmental ethics in forestry. Prentice Hall, Upper Saddle River, N J . Kirchhoff, M . D . and J.W. Schoen. 1987. Forest cover and snow: implications for deer habitat in southeast Alaska. Journal of Wildlife Management 51(l):28-33. Kittredge, J. 1953. Influences of forests on snow in the ponderosa-sugar pine-fir zone of the Sierra Nevada. Hilgardia 22(l):l-96. Kliskey, A.D., E.C. Lofroth, W.A. Thompson, S. Brown, H.Schrier. 1999. Simulating and evaluating alternative resource-use strategies using GIS-based habitat suitability indices. Landscape and Urban Planning 45:163-175. Krebs, C J . 1994. Ecology. HarperCollins College Publishers, New York, New York. Krebs, C J . 1999.  Ecological Methodology. Addison Wesley Education Publishers, Menlo  Park,CA. Kristov, V . 2004. Investigations of indirect relationships in ecology and environmental sciences: a review and the implications for comparative theoretical ecosystem analysis. Ecological Modelling 174:37-54. Kufeld, R.C., O.C. Wallmo, and C. Feddema. 1973. Foods of the Rocky Mountain mule deer. U S D A Forest Service Research Paper RM-111, Rocky Mountain Forest and Range Experimental Station, Colorado Division of Wildlife, Department of Natural Resources. Kullback, S., R . A . Leibler. 1951. On information and sufficiency. The Animals of Mathematical Statistics 22(l):79-86. Kuusipalo, J., J. Kangas. 1994. Managing biodiversity in a forestry environment. Conservation Biology 8(2):450-460. Laird, N . M . , J.H. Ware. 1982. Random-effects models for longitudinal data. Biometrics 38:963-974. Lancia, R . A . , S.D. Miller, D . A . Adams, D.W. Hazel. 1982. Validating habitat quality assessment: an example. Transactions of the North American Wildlife and Natural Resources Conference 47:96-110. Law, R., U . Dieckmann, J . A J . Metz. 2000. Introduction. Pp. 1-6 in Dieckmann, U , R. Law, J.A.J. Metz (eds.). 2000. The geometry of ecological interactions: simplifying spatial complexity. Cambridge University Press, Cambridge. Legendre, P. 1993. Spatial autocorrelation: trouble or new paradigm? Ecology 74(6):16591673. Levin, S.A. 1992. The problem of pattern and scale in ecology. Ecology 73(6):1943-1967. 86  Levin, S.A. 1998. Ecosystems and the biosphere as complex adaptive systems. Ecosystems 1:431-436. Levin, S.A., S. Barrett, S. Aniyar, W. Baumol, C. Bliss, B . Bolin, P. Dasgupta, P. Ehrlich,.C. Folke, I-M. Gren, C.S. Holling, A . Jansson, B-O. Jansson, K - G . Maler, D. Martin, C. Perrings, E. Sheshinski. 1998. Resilience in natural and socioeconomic systems. Environment and Development Economics 3:222-235. Lindgren, P.M.F., T.P. Sullivan. 2001. Influence of alternative vegetation management treatments on conifer plantation attributes: abundance, species diversity, and structural diversity. Forest Ecology and Management 142:163-182. Lippke, B.R., J.T. Bishop. 1999. The economic perspective. Pp. 597-638 in Hunter, M . L . Jr. (ed.). Maintaining Biodiversity in Forest Ecosystems. Cambridge University Press, Cambridge. Liu, J. 1993. E C O L E C O N : an ECOLogical-ECONomic model for species conservation in complex forest landscapes. Ecological Modelling 70: 63-87. Lohr, S.L. 1999. Sampling: design and analysis. Duxbury Press, Pacific Grove, C A . Longhurst, W . M . 1954. The fecal pellet group deposition rate of domestic sheep. Journal of Wildlife Management 18(3):418-419. Loveless, C M . 1964. Some relationships between wintering mule deer and the physical environment. Transactions of the North American Wildlife Conference 29:415-431. Lyon, L.J. and C E . Jensen. 1980. Management implications of elk and deer use of clear-cuts in Montana. Journal of Wildlife Management 44(2):352-362. Mac Arthur, R. 1955. Fluctuations of animal population and a measure of community stability. Ecology 36(3):533-536. MacNally, R. 2000. Regression and model-building in conservation biology, biogeography and ecology: the distinction between - and reconciliation of - 'predictive' and 'explanatory' models. Biodiversity and Conservation 9:655-671. Mallows, C L . 1973. Some comments on C . Technometrics 15(4): 661-675. p  Marcot, B.G., M . G . Raphael, K . H . Berry. 1983. Monitoring wildlife habitat and validation of wildlife-habitat relationships models. Transactions of the North American Wildlife and Natural Resources Conference 48:315-329. Martinka, C J . 1968. Habitat relationships of white-tailed and mule deer in northern Montana. Journal of Wildlife Management 32(3):558-656.  87  Marzluff, J.M., J.J. Millspaugh, K.R. Cedar, C D . Oliver, J. Withey, J.B. McCarter, C L . Mason, J. Comnick. 2002. Modeling changes in wildlife habitat and timber revenues in response to forest management. Forest Science 48(2): 191 -202. Maurer, B . 1986. Predicting habitat quality for grassland birds using density-habitat correlations. Journal of Wildlife Management 50(4):556-566. Mclnnes, P.F., R J . Naiman, J. Pastor, Y . Cohen. 1992. Effects of moose browsing on vegetation and litter of the boreal forest, Isle Royale, Michigan, USA. Ecology 73(6):20592075. McKelvey, K.S., McDaniel, G.W., Mills, L.S., and Griffin, P . C 2002. Effects of plot size and shape on pellet density estimates for snowshoe hares. Wildl. Soc. Bull. 30: 751-755. McNay, R.S. 1985. Forest crowns, snow interception, and management of black-tailed deer winter habitat. Research, Ministries of Environment and Forests, IWIFR-19, Victoria, B.C. Meidinger, D., and Pojar, J. 1991.  Ecosystems of British Columbia.  Research Branch,  Ministry of Forests, Victoria, B.C. Special Report Series No. 6. Ministry of Sustainable Resource Management, Government of British Columbia. CaribooChilcotin Land Use Plan. 1995. http://srmwww.gov.bc.ca/md/regionaycarichil/90dayrep.htm Ministry of Sustainable Resource Management, Government of British Columbia. OkanaganShuswap Land and Resource Management Plan. 2001. http://srmwww.gov.bc.ca/rmd/lirnp/okanshus/docs/osliTnpfiill.pdf Morris, M.S., J.E. Schwartz. 1957. Mule deer and elk food habits on the National Bison Range. Journal of Wildlife Management 21(2):189-193. Morris, D.W., D.L. Davidson. 2000. Optimally foraging mice match patch use with habitat differences in fitness. Ecology 81(8):2061-2066. Murray, D.L., Roth, J.D., Ellsworth, E., Wirsing, A J . , and Steury, T.D. 2002. Estimating lowdensity snowshoe hare populations using fecal pellet counts. Can. J. Zool. 80: 771-781. Mysterud, A., E. Ostbye. 1999. Cover as a habitat element for temperate ungulates: effects on habitat selection and demography. Wildlife Society Bulletin 27(2):385-394. Neff, D J . 1968. The pellet-group count technique for big game trend, census, and distribution: a review. Journal of Wildlife Management 32(3):597-614. Neter, J, M . H . Kutner, C J . Nachtsheim, W. Wasserman. 1996. Applied linear statistical models. McGraw-Hill, Boston, Ma. Nicholson, M . C , R.T. Bowyer, J.G. Kie. 1997. Habitat selection and survival of mule deer: tradeoffs associated with migration. Journal of Mammalogy 78(2):483-504. 88  Nyberg, J.B., F.L. Bunnell, D.W. Janz, R . M . Ellis. 1986. Managing young forests as blacktailed deer winter ranges. British Columbia Ministry of Forests, Land Management Report #37. 58 pp. Oh, J.H., M . B . Jones, W . M . Longhurst, G.E. Connoly. 1970. Deer browsing and rumen microbial fermentation of Douglas-fir as affected by fertilization and growth stage. Forest Science 16(l):21-27. Oliver, C D . , B.C. Larson. 1990. Forest stand dynamics. McGraw-Hill, New York. O'Neil, L.J., A . B . Carey. 1986. Introduction: when habitats fail as predictors. Pp. 207-208 in Verner, J., M . L . Morrison, C..J. Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Oreskes, N . , K . Shrader-Frechette, K . Belitz. 1994. Verification, validation, and confirmation of numerical models in the Earth Sciences. Science 263:641-646. Parish, R., R. Coupe, D. Lloyd. 1996. Plants of southern interior. British Columbia and the Inland Northwest. B.C. Ministry of Forests and Lone Pine Publishing, Vancouver. Palmer, M . W . 1990. The estimation of species richness by extrapolation. Ecology 71(3):11951198. Parker, K . L . , C T . Robbins, T.A. Hanley. 1984. Energy expenditures for locomotion by mule deer and elk. Journal of Wildlife Management 48(2):474-488. Pearce, J.L., K . Cherry, M . Drielsma, S. Ferrier, G. Whish. 2001. Incorporating expert opinion and fine-scale vegetation mapping into statistical models of faunal distribution. Journal of Applied Ecology 38:412-424. Pearce, J., S. Ferrier. 2000. A n evaluation of alternative algorithms for fitting species distribution models using logistic regression. Ecological Modelling 128:127-147. Peek, J.M., J J . Korol, D. Gay, T. Hershey. 2001. Overstory-understory biomass changes over a 35-year period in southcentral Oregon. Forest Ecology and Management 150:267-277. Peek, J.M., M . D . Scott, L.J. Nelson, D J . Pierce. 1982. Role of cover in habitat management for big game in northwestern United States. Transactions of the North American Wildlife Conference 47:363-373. Parker, K . L . , C T . Robbins. 1984. Thermoregulation in mule deer and elk. Canadian Journal of Zoology 62:1409-1422. Parker, K . L . , C T . Robbins, T.A. Hanley. 1984. Energy expenditures for locomotion by mule • deer and elk. Journal of Wildlife Management 48(2):474-488.  89  Parker, K . L . , M.P. Gillingham. 1990. Estimates of critical thermal environments for mule deer. Journal of Range Management 43(1):73-81. Pederson, J.C., K.T. Harper. 1978. Factors influencing productivity of two mule deer herds in Utah. Journal of Range Management 31(2): 105-110. Peek, J.M., B . Dennis, T. Hershey. 2002. Predicting population trends of mule deer. Journal of Wildlife Management 66(3):729-736. Picton, H.D. 1984. Climate and the prediction of reproduction of three ungulate species. Journal of Applied Ecology 21:869-879. Plante, M . , K . Lowell, F. Potvin, B . Boots, M . Fortin. 2004. Studying deer habitat on Anticosti Island, Quebec: relating animal occurrences and forest map information. Ecological Modelling 174:387-399 Preece, D.A. 1990. R.A. Fisher and experimental design: a review. Biometrics 46:925-935. Puoli, J.R., G.A. Jung, R . L . Reid. 1991. Effects of nitrogen and sulfur on digestion and nutritive quality of warm season grass hays for cattle and sheep. Journal of Animal Science 69(2):843-852. Pullin, A.S., T . M . Knight, D . A . Stone, K . Charman. 2004. Do conservation managers use scientific evidence to support their decision-making? Biological Conservation 119:245-252. Raudenbush, S.W., A.S. Bryk. 2002. Hierachical linear models: applications and data analysis methods. Sage Publications, London. Rencher, A.C., F.C. Pun. 1980. Inflation of R in best subset regression. Technometrics 2  22(l):49-53. Reynolds, C S . 2002. Ecological pattern and ecosystem theory. Ecological Modelling 158:181200. Robinette, W.L., O. Julander, J.S. Gashwiler, J.G. Smith. 1952. Winter mortality of mule deer in Utah in relation to range condition. Journal of Wildlife Management 16(3):289-299. Robinson, G.K. 1991. That B L U P is a good thing: the estimation of random effects. Statistical Science 6(1):15-51 Rodgers, G., O. Julander, W.L. Robinette. 1958. Pellet-group counts for deer census and range use index. Journal of Wildlife Management 22(2): 193-199. Rogers, L . L . 1987. Seasonal changes in defecation rates of free-ranging white-tailed deer. Journal of Wildlife Management 51(2):330-333. Rotenberfy, J.T. 1986. Habitat relationships of shrubsteppe birds: even "good" models cannot predict the future.  Pages 217-221 in J. Verner, M . L . Morrison, and C J . Ralph, editors. 90  Wildlife 2000: Modelling Habitat Relationships of Terrestrial Vertebrates. The University of Wisconsin Press, Madison, Wisconsin. Rykiel, E J . Jr. 1996. Testing ecological models: the meaning of validation. Ecological Modelling 90:229-244. SAS Systems. 2001. SAS / STAT Release 8.02. SAS Institute, Cary, N C . Salwasser, H . 1992. From new perspectives to ecosystem management: response to Frissell et al. and Lawrence and Murphy. Conservation Biology 6(3):469-472. Sargeant, G.A., L.E. Eberhardt, J.M. Peek. 1994. Thermoregulation by mule deer (Odocoileus hemionus) in arid rangelands of southcentral Washington. Journal of Mammalogy 75(2):536544. Saunders, D.A., R.J. Hobbs, C.R. Margules. 1991. Biological consequences of ecosystem fragmentation: a review. Conservation Biology 5(l):18-32. Schabenberger, O., F.J. Pierce. 2001. Contemporary statistical models for the plant and soil sciences. C R C Press, Boca Raton, FL. Schamberger, M . L . , L.J. O'Neil. 1986. Concepts and constraints of habitat-model testing. Pp. 5-10 in Verner, J., M X . Morrison, C.J. Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Schoen, J.W. and M . D . Kirchhoff. 1985. Seasonal distribution and home range patterns of sitka black-tailed deer on Admiralty Island, southeast Alaska. Journal of wildlife Management 49(1):96-103 Schwarz, G. 1978. Estimating the dimension of a model. The Annals of Statistics 6(2):461-464. Seymour, R.S., M . L . Hunter Jr. 1999. Principles of ecological forestry. Pp. 22-61 in Hunter, M . L . Jr. (ed.). Maintaining biodiversity in forested landscapes. Cambridge University Press, Cambridge. Shackleton, D . M . 1999. Hoofed mammals of British Columbia. U.B.C. Press, Vancouver. Shannon, C E . 1948. A mathematical theory of communication. Bell System Technical Journal 27:379-423, 623-656. Short, H.L., D.R. Dietz, E.E. Remmenga. 1966. Selected nutrients in mule deer browse plants. Ecology 47(2):222-229. Silen, R.R., W . K . Randall, N . L . Mandel. 1986. Estimates of genetic parameters for deer browsing of Douglas-fir. Forest Science 32(1):178-184. Smith, A . D . 1953. Consumption of native forage species by captive mule deer during summer. Journal of Range Management: 30-37. 91  Smith, A . D . 1964. Defecation rates for mule deer. Journal of Wildlife Management 28(3):435444. Smith, T . M . 1986. Habitat-simulation models: integrating habitat-classification and forestsimulation models. Pp. 389-393 in Verner, J., M . L . Morrison, C.J. Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Sokal, R.R., F.J. Rohlf. 1995. Biometry. W.H. Freeman and Company, New York. Spalinger, D.E., C T . Robbins, T.A. Hanley. 1986. The assessment of handling time in ruminants: the effect of plant chemical and physical structure on the rate of breakdown of plant particles in the rumen of mule deer and elk. Canadian Journal of Zoology 64:312-321. Spies, T.A., M . G . Turner. 1999. Dynamic forest mosaics. Pp. 95-160 in Hunter, M . L . Jr. (ed.). Maintaining biodiversity in forested landscapes. Cambridge University Press, Cambridge. Stauffer, D.F., L . B . Best. 1986. Effects of habitat type and sample size on Habitat Suitability Index models. Pp. 71-77 in Verner, J., M . L . Morrison, C J . Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Stewart, K . M . , R.T. Bowyer, J.G. Kie, N J . Cimon, B . K . Johnson. 2002. Temporospatial distributions of elk, mule deer, and cattle: resource partitioning and competitive displacement. Journal of Mammology 83(l):229-244. Stickney, P.F. 1980. Data base for post-fire succession, first 6 to 9 years, in Montana larch-fir forests. U S D A Forest Service General Technical Report INT-62. Stoica, P., Y . Selen, J. L i . 2004. Multi-model approach to model selection. Digital Signal Processing 14:399-412. Store, R., J. Jokimaki. 2003. A GIS-based multi-scale approach to habitat suitability modeling. Ecological Modelling 169: 1-15. Sullivan, T.P. 2005. Influence of repeated fertilization on forest ecosystems: I. Abundance and diversity of understory vegetation. In preparation. Sullivan, T.P., D.S. Sullivan, P.M.F. Lindgren, J.O. Boateng. 2002. Influence of conventional and chemical thinning on stand structure and diversity of plant and mammal communities in young lodgepole pine forest. Forest Ecology and Management 170:173-187. Sugiura, N . 1978. Further analysis of the data by Akaike's Information Criterion and the finite corrections. Communications in Statistics: Theory and Methods A7(l): 13-26.  92  Swift, R.W. 1948. Deer select most nutritious forages. Journal of Wildlife Management 12(1):109-110. Tabachnick, B.G., L.S. Fidell. 2001. Using multivariate statistics. Allyn and Bacon, Boston. Telfer, E.S. 1978. Cervid distribution, browse, and snow cover in Alberta. Journal of Wildlife Management 42(2):352-361. Ter Braak, C.J.F. 1986. Canonical correspondence analysis: a new eigenvector technique for multivariate direct gradient analysis. Ecology 67(5): 1167-1179. Thomas, J.W., H . Black, Jr., R J . Scherzinger, R J . Pederson. 1979. Deer and elk, pp. 104-127. In: J. W. Thomas (ed.). Wildlife habitats in managed forests-the Blue Mountains of Oregon and Washington. U.S.D.A Handbook No. 533. Washington, D.C. Thomas, J.W. 1986. Wildlife-habitat modeling - cheers, fears, and introspection. Pp. xix-xxv in Verner, J., M . L . Morrison, C J . Ralph (eds.). Wildlife 2000: modeling habitat relationships of terrestrial vertebrates. The University of Wisconsin Press, Madison, WI. Tomm, H.O., J.A. Beck, Jr., R J . Hudson. 1981. Responses of wild ungulates to logging practices in Alberta. Canadian Journal of Forest Research 11:606-614. Torbit, S . C , L . H . Carpenter, D . M . Swift, A.W. Alldredge. 1985. Differential loss of fat and protein by mule deer during winter. Journal of Wildlife Management 49(l):80-85. Tucker, R.E., W.Majak, P.D. Parkinson, A.McLean. 1976. Palatability of Douglas-fir foliage to mule deer in relation to chemical and spatial factors. Journal of Range Management 29(6):486489. Turner, M . G . 1989. Landscape ecology: the effect of pattern on process. Annual Review of Ecology and Systematics 20: 171-197. Turner, M . G . , G.J. Arthaud, R.T. Engstrom, S J . Hejl, J. Liu, S. Loeb, K . McKelvey. 1995. Usefulness  of spatially explicit population models in land management. Ecological  Applications 5: 12-16. Uresk, D.W., V . A . Uresk. 1982. Diets and habitats of mule deer in South-central Washington. Northwest Science 56(2): 138-147. Urness, P J . , D J . Neff, J.R. Vahle. 1975. Nutrient content of mule deer diets from ponderosa pine range. Journal of Wildlife Management 39(4):670-673. U.S. Fish and Wildlife Service. 1981. Standards for the Development of Habitat Suitability Index  Models  (103  ESM). U.S. Fish  and  http://policy.fws.gov/ESMindex.html.  93  Wildlife  Service,  Washington, D C .  Vanclay, J.K., J.P. Skovgaard. 1996. Evaluating forest growth models. Ecological Modelling 98:1-12. Vangilder, L.D., O. Torgerson, W.R. Porath. 1982. Factors influencing diet selection by whitetailed deer. Journal of Wildlife Management 46(3):711-718. Van Home, B . 1983. Density as a Misleading Indicator of Habitat Quality. Journal of Wildlife Management 47(4):893-901. Wallmo, O.C. 1981. Mule and black-tailed deer of North America. University of Nebraska Press, Lincoln. Wallmo, O . C , W.L. Reglin, D.W. Reichert. 1972. Forage use by mule deer relative to logging in Colorado. Journal of Wildlife Management 36(4): 1025-1033. Wallmo, O . C , L . H . Carpenter, W.L. Regelin, R.B. Gill, D.L. Baker. 1977. Evaluation of deer habitat on a nutritional basis. Journal of Range Management 30(2): 122-127. Walters, C. 1986. Adaptive management of renewable resources. Macmillan Publishing Company, New York. Walters, C.J., C S . Holling. 1990. Large-scale management experiments and learning by doing. Ecology 71(6):2060-2068 Wambolt, C L . and A.F. McNeal. 1987. Selection of winter foraging sites by elk and mule deer. Journal of Environmental Management 25:285-291. Waterhouse, M.J., H . M . Armleder, R.J. Dawson. 1991. Forage litterfall in Douglas fir forests of the central interior of British Columbia. British Columbia Ministry of Forests Research Note No.108. Waterhouse, M.J., H . M . Armleder, R.J. Dawson. 1994. Winter food habits of mule deer in the central interior of British Columbia. Ministry of Forests Research Note 113. Victoria, B.C. Wessels, K.J., A.S. Van Jaarsveld, J.D. Grimbeek, M.J. Van der Linde. 1998. A n evaluation of the gradsect biological survey method. Biodiversity and Conservation 7:1093-1121. White, K . L . 1960. Differential range use by mule deer in the spmce-fir zone. Northwest Science 34(4): 118-126. Whittaker, R.H., S.A. Levin, R.B. Root. 1973. Niche, habitat, and ecotype. The American Naturalist 107:321-338. Wickstrom, M . L . , C.T. Robbins, T.A. Hanley, D.E. Spalinger, S.M. Parish. 1984. Food intake and foraging energetics of elk and mule deer. Journal of Wildlife Management 48(4): 12851301.  94  Willms, W., A . McLean, R. Ritchy. 1976. Feeding habits of mule deer on fall, winter, and spring ranges near Kamloops, British Columbia. Canadian Journal of Animal Science 56:531542. Willms, W., A . McLean, R. Tucker, R. Ritcey. 1980. Deer and cattle diets on summer range in British Columbia. Journal O f Range Management 33(l):55-59. Woolf, A . C K . Nielson, T. Weber, T J . Gibbs-Keininger. 2002. Statewide modeling of bobcat, Lynx rufus, habitat in Illinois, U S A . Biological Conservation 104: 191-198. Wu, J, D. Marceau. 2002. Modeling complex ecological systems: an introduction. Ecological Modelling 153:1-6. Yang, Y . , R.A. Monserud, S. Huang. 2004. A n evaluation of diagnostic tests and their roles in validating forest biometric models. Canadian Journal of Forest Resources 34:619-629. Yarmoloy, C , M . Bayer, V . Geist. 1988. Behaviour responses and reproduction of mule deer, Odocoileus hemionus, does following experimental harassment with an all-terrain vehicle. The Canadian Field Naturalist 102:425-429.  95  Appendix 1: Extended Model and Data Analysis Results Table A l . l . Full A I C output for summer data, including AIC weights and evidence ratios calculated for the full model set. C  ;K 1 2 3 4 5 6  fertilizer, shrub volume, herb volume, herb species richness, shrub species richness, density, density x shrub volume, density x herb volume fertilizer, shrub volume, herb volume, herb species richness, shrub species richness, density, density x shrub volume fertilizer, shrub volume, herb volume, herb species richness, shrub species richness, density, density x herb volume fertilizer, shrub volume, herb volume, shrub species richness, density, density x shrub volume fertilizer, shrub volume, herb volume, herb species richness, density, density x herb volume fertilizer, shrub volume, herb volume, density, density x shrub volume, density x herb volume  Ih  10  33.64  28.26  0.00000  1370303.82  0.486  9  26.57  21.19  0.00001  39841.24  0.486  9  26.57  21.18  0.00001  39819.24  0.486  8  20.57  15.18  0.00019  1981.89  0.482  8  20.53  15.15  0.00019  1946.96  0.485  8  20.71  15.33 . 0.00017  2129.32  0.476  7  fertilizer, shrub volume, herb volume, density, density x shrub volume  7  15.46  10.07  0,00241  153.83  0.476  8  fertilizer, shrub volume, herb volume, density, density x herb volume  7  15.45  .10.07  0.00241  153.70  0.476  9  fertilizer, shrub volume, density, density x shrub volume  6  10.94  5.55  0.02304  16.07  0.471  10  fertilizer, herb volume, density, density x herb volume  6  12.95  7.57  0.00841  44.03  0.359  11  fertilizer, shrub volume, herb volume, shrub species richness, herb species richness, density  8  20.50  15.12  0.00019  1919.36  0.486  12  fertilizer, shrub volume, herb volume, shrub species richness, herb species richness  7  17.02  11.64  0.00110  336.53  0.391  13  fertilizer, shrub volume, herb volume, herb species richness  6  12.62  7.24  0.00993  37.30  0.379  14  fertilizer, shrub volume, herb volume, shrub species richness  6  12.47  7.08  0.01073  34.52  0.388  15  fertilizer, shrub volume, herb volume, density  6  10.86  5.47  0.02401  15.43  0.476  16  fertilizer, shrub volume, herb volume  5  8.90  3.51  0.06398  5.79  0.359  17  fertilizer, shrub volume  4  5.38  0.00  0.37032  1.00  0.353  18  fertilizer, herb volume  4  7.38  2.00  0.13644  2.71  0.216  19  fertilizer, density  4  5.52  0.13  0.34647  1.07  0.344  96  Table A1.2. Full AIC output for winter data, including AIC weights and evidence ratios calculated for the full model set. C  ©  (M@0  CS  i  m§®  m •  "'DEOO®  7  fertilizer, shrub volume, herb volume, herb species richness, shrub species richness, density, density x shrub volume, density x herb volume fertilizer, shrub volume, herb volume, herb species richness, shrub species richness, density, density x shrub volume fertilizer, shrub volume, herb volume, herb species richness, shrub species richness, density, density x herb volume fertilizer, shrub volume, herb volume, shrub species richness, density, density x shrub volume fertilizer, shrub volume, herb volume, herb species richness, density, density x herb volume fertilizer, shrub volume, herb volume, density, density x shrub volume, density x herb volume fertilizer, shrub volume, herb volume, density, density x shrub volume  8  fertilizer, shrub volume, herb volume, density, density x herb volume  7  18.18  10.91  0.00130  234.40  0.328  9  fertilizer, shrub volume, density, density x shrub volume  6  14.41  7.15  0.00852  35.65  0.272  10 fertilizer, herb volume, density, density x herb volume fertilizer.shrub volume, herb volume, shrub species richness, herb species 11 richness, density fertilizer, shrub volume, herb volume, shrub species richness, herb species 12 richness 13 fertilizer, shrub volume, herb volume, herb species richness  6  14.37  7.11  0.00869  34.95  0.275  8  21.72  14.46  0.00022  1378.09  0.430  7  17.59  10.32  0.00174  174.26  0.365  6  13.11  5.84  0.01638  18.55  0.358  14 fertilizer, shrub volume, herb volume, shrub species richness  6  14.29  7.02  0.00906  33.53  0.280  15 fertilizer, shrub volume, herb volume, density  6  13.59  6.32  0.01288  23.60  0.328  16 17 18 19  5 4 4 4  10.31 7.31 7.48 7.27  3.05 0.05 0.22 0.00  0.06627 0.29672 0.27260 0.30385  4.58 1.02 1.11 1.00  0.275 0.232 0.219 0.235  1 2 3 4 5 6  fertilizer, fertilizer, fertilizer, fertilizer,  shrub volume, herb volume shrub volume herb volume density  97  10 34.78 27.52  0.00000  944823.18  0.435  9  27.71  20.44  0.00001  27498.11  0.434  9  27.76 20.49  0.00001  28155.90  0.432  8  23.06  15.79  0.00011  2683.99  0.352  8  21.71  14.44  0.00022  1365.53  0.431  8  23.37  16.11  0.00010  3145.75  0.332  7  18.15  10.89  0.00131  231.21 0.330  Table A1.3. Final models. Listed values are weighted averages of parameter estimates for multivariate linear regression models of summer and winter data; for the full models created using all a priori selected variables (models cut off at a Aj value of 10), as well for as a reduced (management) model list composed of models with a Aj value less than 10 that do not include shrub or herb species richness. Standard errors (SE) of the model averaged parameter estimates are present in parentheses. Recall that models arose using pellet-group densities that were standardized by region, and expressed on a per ha basis. Stand densities were also on a per ha basis, and volume measurements are actually crown volume indices (calculated as described in the Methods section of Chapter 3) expressed as m /0.01 ha. 3  ftwrih  mm, Full Model  Summer Winter  Management Summer Model Winter  s  i  s  0.3957  -0.1935  3.06E-03  3.72E-04  -8.07E-05  3.27E-04  1.04E-04  -3.48E-09  7.94E-10  (SE=0.4504)  (SE=0.3417)  (SE=0.0603)  (SE=0.0311)  (SE=0.0105)  (SE=0.0282)  (SE=0.0165)  (SE=3.47E-04)  (SE=1.67E-04)  0.6471  -0.2175  -1.74E-03  9.55E-04  4.46E-05  1.20E-04  3.76E-04  6.40E-09  -6.38E-09  (SE=0.4045)  (SE=0.3545)  (SE=0.0541)  (SE=0.0420)  (SE=0.0087)  (SE=0.0213)  (SE=0.0285)  (SE=2.59E-04)  (SE=1.90E-04)  NA  NA  0.4025  -0.1939  3.02E-03  3.52E-04  -8.24E-05  (SE=0.4428)  (SE=0.3377)  (SE=0.0596)  (SE=0.0301)  (SE=0.0105)  0.6631  -0.2217  -1.60E-03  8.77E-04  2.90E-05  (SE=0.3919)  (SE=0.3515)  (SE=0.0521)  (SE=0.0407)  (SE=0.0077)  98  NA  NA  -3.55E-09  8.11E-10  (SE=3.47E-04)  (SE=1.88E-04)  6.57E-09  -6.55E-09  (SE=2.67E-04)  (SE=2.24E-04)  A1.4. Summary of canonical correlation analysis results for fertilized and unfertilized plots analyzed separately, by site. Correlations listed are between available independent variables and their canonical variable. Site 1 = Summerland, Site 2 = Kelowna, Site 3 = Gavin Lake (Cariboo).  Table  S  fl  2  0.30  0.70  -0.71  0.2,1  0.27  -0.01  -0.34  0.32  0.12  0.23  0.15  -0.17  0.71  -0.03  0.20  0.40  -0.36  0.49  0.95  -0.06  -0.36  0.39  0.68  -0.61  -0.02  0.91  0.58  0.50  -0.48  -0.04  0.26  0.15  0.03  0.49  0.28  -0.04  0.23  0.15  0.08  -0.56  -0.28  -0.14  Stand Density  -0.43  0.31  0.47  -0.29  -0.28  0.17  -0.37  0.35  0.28  0.41  0.41  0.25  Density X Shrub Volume  -0.14  0.29  0.56  -0.19  -0.29  0.48  -0.20  0.52  0.56  0.14  0.69  0.78  Density X Herb Volume  -0.16  -0.06  -0.21  0.22  0.14  -0.19  -0.13  0.39  0.24  -0.09  0.46  0.34  Explained Standardized Variance of PelletGroup Counts  0.41  0.60  0.53  0.72  0.83  0.89  0.42  0.47  0.68  0.46  0.20  0.74  Herb Volume  -0.34  -0.69  Shrub Volume  0.12  0.06  Herb Species Richness  0.44  Shrub Species Richness  99  fl  .  2  ®  g  •a 0.42  •  Appendix 2 Habitat Suitability Component Design for FORECAST  A2.1 General Structure The wildlife habitat suitability component described here has been designed to integrate into the ecosystem simulation program FORECAST, and includes three main components; 1) The wildlife suitability interface ('Wildlife Habitat Suitability' window), a main page for defining the number of species to model, naming the species, giving access to the variable definition and equation building tools, and for displaying constructed equations; 2) The 'Wildlife Variable Selection' window, a tool for defining variables for use in habitat suitability equations; 3) The 'Build Habitat Suitability Equations' window, an equation building and editing tool.  A2.2 Wildlife Habitat Suitability Interface The general structure of the wildlife habitat suitability interface is fairly simple (see Fig. A2.1). The functions to be performed by the main interface page are; 1) To identify the number of wildlife species to be represented; 2) To name the species; 3) To provide access to the "Wildlife Variable Selection" window; 4) To provide access to the "Build Habitat Suitability Equations" window; 5) To display the defined variables and relationships for easy reference. The first section of the wildlife interface involves setting the number of wildlife species to be represented. The only other viable option would have been to pre-determine how many wildlife species could be represented in a F O R E C A S T run. However, this would have compromised flexibility in application. Next, the interface includes a location for naming the species. A l l F O R E C A S T output referring to the species uses the text in this location for identification. The sections to the right of the main interface include one button each for providing access to the "Wildlife Variable Selection" (the "Define Variables" button) and "Build Habitat Suitability Equations" (the "Build Equations" button) windows. The "Variable Definitions" box displays all defined variables for easy reference. Similarly, the "Habitat Suitability Equations" box lists the defined habitat suitability relationships. Pn this way, the user can quickly review defined relationships at a glance without having to open the editor window and go over the defined relationships one at a time. It would likely be an appealing additional feature for the interface to have the capability to display the simple curves definable in the "Manual Input" box 100  of the "Relationship Workshop" (described later). However, the displaying of equations and curves within the main interface window has not been included in the interface design at this point.  • Wildlife Habitat Suitability  Represented Species Wildlife Species  Loaded Data Set  C:\Horizori\DataSets\Lodgepolelfds  Detete  Variable Definitions  Number of Wildlife Species  1  (fa highlighted species)  Define Variables  shrub volume shrub species richness herb species richness  Habitat Suitability Equations (for highlighted species)  Save Lists  Restore Lists  Help  Build Equations  Close  Figure A2.1. Wildlife interface main page.  A2.3 Wildlife Variable Selection Perhaps the most complex component of the F O R E C A S T wildlife interface involves defining variables out of F O R E C A S T data outputs for later use in constructing mathematical relationships. The structure for the "Wildlife Variable Selection" window is shown in Figure A2.2. It is composed of 6 parts; a "Variable Identity" box, a "Feature" box, a "Species" box, a "Trait" box, basic mathematical buttons, and a "Variable Definition" box.  101  I Wildlife Variable Selection  Variable Identity 1  shrub volume  2  herb volume  3  shrub species richness  4  r  Trait Modifiers Number of Variables Units Conversion X Add Variable  herb species richness stem density  Feature  Variable Definitions  |  Delete Variable j  Species  I  Range Limits  > I  Trait  Variable Number  Close  Figure A2.2. The Wildlife Variable Selection component.  The first step in defining a variable will be to select a variable and to name it. This is done in the "Variable Identity" box in the top left corner. To add a new variable, the user simply clicks into an empty space next to a variable identification number, and type. Next, the user will specify the feature, species (one or more), and trait for the variable. To specify the feature, the user will click in the "Feature" box and select from the available options of "Trees", "Plants", "Snags" and " C W D " . Only one feature can be included per variable. Options available in the "Species" and "Trait" boxes will be dependent on the option selected in the "Feature" box. Next, in the "Species" box, all species to be referred to by the variable must be selected.  Options  available will depend on the option selected in the "Feature" box, and species names listed will come from species names specified in the "Header" tab of FORECAST'S Ecodata file. Species are selected from the list by simply clicking on them. Next, in the trait box, the particular trait of the "Feature" and "Species" in question will be specified. In the final version, the user will likely select from pre-defined "Traits", rather than be free to define them independently. Options to be available are as follows;  102  F o r "Feature" "Trees"; SPH (stems per hectare) @ D B H (diameter at breast height), SPH @ height, SPH @ age, B A (basal area) @ D B H , B A @ height, B A @ age, average D B H , average height, average age, canopy closure (and/or light reaching ground), and canopy depth. F o r "Feature" "Plants"; height, volume, fruit, new volume. F o r "Feature" "Snags";  decay state, SPH @ D B H , SPH @ height.  F o r "Feature" " C W D " ;  decay state, SPH @ D B H , SPH @ length.  There is an interesting point made by K i m Scoular (personal communication) that should be noted here regarding defining edible biomass of plants. Here I simply defined the trait as new volume (growth since the previous winter), however this may be an over simplification. As K i m Scoular stated "perhaps moose eat more twigs than deer". At the time of writing, no solution had been found regarding this problem. Once the variable (in the "Identity" box) "Feature", "Species", and "Trait" are all selected, the user will click the "Define Variable" button below the "Trait" box, and an expression describing the variable thus far will appear in the "Variable Definitions" box at the bottom of the "Variable Selector" window. The appropriate variable identification number will simultaneously appear in the "Variable Number" box. Finally, the user will specify the final details of the variable with the mathematical buttons in the top right section of the "Variable Selector" window and the numerical pad on a keyboard. Currently, only lesser-than and greater-than symbols are available, but lesser-than-orequal-to and greater-than-or-equal-to functions will likely be available in the final version. For an example of application, the species "Douglas-fir" and trait "SPH @ D B H " followed by the sign ">" and the number "40" would specify that the variable refers to the number of stems per ha of Douglas fir above 40 cm in dbh. If, D B H was unimportant, no range limit button would be selected, specifying that Douglas fir stems per ha at all D B H are referred to (e.g. large Fd = Douglas fir SPH @ DBH).  If it simplifies matters, "Traits" of just " S P H " and " B A " may be  included in the final version to provide the same feature. Variable definitions entered will be saved automatically upon closing the "Wildlife Variable Selection" window. Also, variable definitions entered can be reviewed and edited by again clicking on the appropriate variable names in the "Identity" box.  103  A2.4. Build Habitat Suitability Equations I Build Habitat Suitability Equations  Equation Identity:  Manual Input: Manual  r Functions:  VariableSelect:  shrub volume shrub species richness herb species richness stem density  Relationship: Mule Deer = 0.4414 + (-1 .60E-03)(shrub volume) + (877E-04)(herb volume) + (2.90E-05)(stem density) + (6.57E-09p(sh  Clear  Vallidate  Close  Figure A2.3. The Build Habitat Suitability Equations component.  The "Build Habitat Suitability" component of the wildlife habitat suitability interface (see Figure A2.3) provides the user an opportunity, through a number of possible methods, to express relationships between variables defined in the "Wildlife Variable Selection" component and the habitat suitability of species of interest. Habitat suitability relationships can be defined by one equation, or by more than one component equation combined by a total habitat suitability equation.  Although incomplete at the time of writing, relationships between variables and  habitat suitability will be able to be defined mathematically or by manually entering values in the final version. This component has 7 main parts; a "Equation Identity" box, a "Variable Select" box, a "Functions" box, a "Manual Input" box, a "Feedback" check box, "Validate" and "Clear" buttons, and the "Relationship" box. The first step in defining relationships will be to select and name the curve to be defined. The user simply clicks the appropriate space in the "Equation Identity" box, and types in the desired name. Unless the "Feedback" or "Manual" boxes are checked, after selecting a curve the name of the curve followed by an equal sign will appear in the "Relationship" box.  104  Variables in the "Variable Select" box will be those variables defined in the "Wildlife Variable Selection" window (see section A2.3), as well as each curve defined. Equations already defined will also appear in the "Variable Select" box. This will give the user the option of defining habitat suitability of a species through the combination of multiple habitat suitability curves that each define habitat suitability in relation to a particular component of the habitat. Although incomplete at the time of writing, habitat suitability equations defined for other defined species may also appear here, clearly marked as such, so that the user may define how the habitat suitability of one species will impact on the habitat suitability of another. Variables will be selected by clicking on the appropriate space in the "Variable Select" box. If the "Feedback" and "Manual" boxes are not checked, the name of the selected curve will appear as a predictor variable wherever the cursor is placed in the "Relationship" box. If the "Feedback" box is checked, the dependent variable will be the first variable selected in the "Variable Select" box, and should automatically be followed by an equal sign for clarity. Subsequently selected variables will appear, as normal, on the right side of the equation, wherever the cursor is placed. If the "Manual" box is checked, the name of variable selected will appear in the left column of the "Manual Input" box as the independent variable. Equations can be built in the "Relationship" box through a combination of selecting desired variables, using the numerical keypad on a keyboard (including brackets), and adding mathematical functions from the "Functions" box. Functions will be selected by simply clicking upon them with a computer mouse. Once clicked, the function will enter into the "Relationship" box wherever the cursor is placed. The functions to be included are all standard mathematical functions; exponent, logarithm, natural log (LN), and e, with the possible inclusion of trigonometric functions. If the "Feedback" box is checked, as mentioned above, the first variable selected will appear as the dependent variable followed by an equal sign, and each subsequent variable selected will appear as a predictor variable on the opposite side of the equation. The "Feedback" box will give the user the opportunity to define impacts on defined variables (or perhaps even the habitat suitability of other defined species) in relation to calculated habitat suitability of a species. This would only be of use for animals whose populations are tightly tied to habitat components, and no such animals that I currently know of would achieve abundances that would appreciably impact habitat components. Therefore, this feature would not likely be useful without the addition of population demographic details, and may be excluded from the final version 105  Whenever the "Manual" box is checked, the identified equation will refer to the values manually entered into the right column of the "Manual Input" box, and the name of the selected variable will appear in the left column. This option will allow the user to more easily take a traditional habitat suitability index (HSI) approach, where responses to variables are generally defined as linear changes of habitat suitability in relation to the supply of a particular habitat variable. This is the general approach adopted by SIMFOR, and the widespread use of habitat suitability index modelling in wildlife management necessitates it's inclusion into the FORECAST wildlife interface as well. model is shown in Figure A2.4.  A n example of one component of a traditional HSI  The "Build Habitat Suitability" component will be required to  automatically form an equation of a line between specified low and high values (calculating slope and intercept behind the scenes) in order to be able to assign a particular habitat suitability index value for an intermediate level of supply for the variable in question. The habitat variable of interest will be selected from the "Variable Select" box, and will appear in the independent variable column (left column) of the "Manual Input" box. It would be helpful at this stage i f a simple graph, similar to the structure of Fig. 4, could appear in the "Relationship" window to visually represent curves entered into the "Manual Input" box, however this is not an essential function.  % herb cover  Figure A2.4. Example of a traditional HSI component model. The user specifies at what level in the independent variable the habitat suitability index reaches it's high and low points, and simply draws a straight line between the two. 106  Appendix 3 Univariate Responses of Model Variables to Standardized Pellet-Group Densities A3.1. Graphs of Univariate Responses in Summer 1.2  0.8 A unfertilized • fertilized  0.6  0.4  0.2  A  ^  20  40  60  80  shrub volume index (m /0.01 ha) 3  Figure A3.1. Standardized pellet-group density vs. shrub volume index for fertilized and unfertilized plots in summer.  1.2  0.8  0.6 A unfertilized • fertilized  ! — 0.4 A  0.2  A* A  20  40  60  80  100  120  herb volume (m /0.07 ha)) 3  Figure A3.2. Standardized pellet-group density vs. herb volume index for fertilized and unfertilized plots in summer.  107  1 9  0)  I §1 £  «  •o  2  £  -  0.8  0-6  QQ. 3 O  A unfertilized • fertilized 0.4  A  0.2  A  •  10  12  14  shrub species richness (# species per stand)  Figure A 3 . 3 . Standardized pellet-group density vs. shrub species richness for fertilized and unfertilized plots in summer.  1 9  •a  0.8  o ^ O) o ». Q«> in  — a.  0.6  A unfertilized  o £  • fertilized 0.4  S  0.2  10  15  20  25  herb species richness (# species per stand)  Figure A 3 . 4 . Standardized pellet-group density vs. herb species richness for fertilized and unfertilized plots in summer.  108  1.2  a. ^  ,  0.8  P  A unfertilized • fertilized  O) Q. Q. 3  -a  2  |  S  o.4  0.2  ^  A A  • I  A  500  1000  1500  2000  2500  stand density (stems per ha over 3 m in height)  Figure A3.5. Standardized pellet-group density vs. stand density for fertilized and unfertilized plots in summer.  A3.2. Graphs of Univariate Responses in Winter  1.2 1 Q. O  f  0.8  0.6 CL  = (A 0)  a. a  A'^ • •  A unfertilized  •  • fertilized  • A A  0.4  •  AA A  0.2  •  •A  • A  A  A •  A  0 20  40  60  80  shrub volume index (m /0.07 ha) 3  Figure A3.6. Standardized pellet-group density vs. shrub volume index for fertilized and unfertilized plots in winter.  109  1.2  1 u 0.8  Q.  A unfertilized = (A O Q. O. 3 TJ  • fertilized  0.6  o  S « 0.4  A A  0.2  0 20  40  60  80  100  120  herb volume index  (m /0.01 ha) 3  Figure A3.7. Standardized pellet-group density vs. herb volume index for fertilized and unfertilized plots in winter.  1.2  1  >.  'co  c o  •a  0.8  A unfertilized  i ! O)  o  °<D 01  • fertilized  0.6  1 =  °. o "g a 0.4 N TJ kQ TJ C  tt  0.2  (fl 5  10  15  shrub species richness (# species per stand)  Figure A3.8. Standardized pellet-group density vs. shrub species richness for fertilized and unfertilized plots in winter.  110  1.2 1 0.8  Q.  3 o =  10  a. s •a £ N  —  T3  A unfertilized  0.6  W  • fertilized  0.4 0.2 0 10  15  20  25  herb species richness (# species per stand)  Figure A3.9. Standardized pellet-group density vs. herb species richness for fertilized and unfertilized plots in winter.  1.2 1  2 "> = (A 0) Q . Q. 3  •o 2 2  0.8  •  m A  m •  •  •  •  AA  0.4  •  H  A  A  •  A  • fertilized A  A A  A  0.2  A unfertilized  •  0.6  •  0 500  1000  1500  2000  2500  stand density (stems per ha over 3 m in height)  Figure A3.10. Standardized pellet-group density vs. stand density for fertilized and unfertilized plots in winter.  Ill  A3.3. Simple Correlations Table A3.1. Simple Pearson's correlation coefficients between available variables for summer data. eftrafe  ssreftnrass  agiilitgamrgiil  spaa®  gtefe  8§8dl  .  gfarfo I  > m$KMi  ' sireftnrass  o  t£S§§  1.00  0.44  0.39  0.53  0.17  0.10  -0.35  -0.18  0.55  -0.36  -0.30  0.04  0.41  0.16  0.44 0.39  1.00 0.65  0.65 1.00  0.35 0.34  -0.22 -0.34  -0.45 -0.20  0.08 0.01  0.22 0.08  0.15 0.29  -0.64 -0.71  -0.64 -0.40  -0.27 0.05  0.01 0.07  -0.35 -0.30  0.53  0.35  0.34  1.00  0.34  0.22  -0.07  0.10  0.45  -0.36  -0.64  -0.17  0.16  0.29  0.17  -0.22  -0.34  0.34  1.00  0.63  -0.16  0.00  0.18  0.61  0.20  -0.13  0.07  0.93  0.10  -0.45  -0.20  0.22  0.63  1.00  -0.29  -0.06  0.34  0.40  0.53  0.51  0.40  0.87  -0.35  0.08  0.01  -0.07  -0.16  -0.29  1.00  0.89  -0.74  -0.05  -0.21  -0.15  -0.49  -0.24  -0.18  0.22  0.08  0.10  0.00  -0.06  0.89  1.00  -0.65  -0.12  -0.18  -0.01  -0.35  -0.02  0.55  0.15  0.29  0.45  0.18  0.34  -0.74  -0.65  1.00  -0.27  -0.10  0.16  0.64  0.27  -0.36  -0.64  -0.71  -0.36  0.61  0.40  -0.05  -0.12  -0.27  1.00  0.58  -0.09  -0.20  0.57  -0.30  -0.64  -0.40  -0.64  0.20  0.53  -0.21  -0.18  -0.10  0.58  1.00  0.62  0.26  0.42  0.04  -0.27  0.05  -0.17  -0.13  0.51  -0.15  -0.01  0.16  -0.09  0.62  1.00  0.66  0.23  0.41  0.01  0.07  0.16  0.07  0.40  -0.49  -0.35  0.64  -0.20  0.26  0.66  1.00  0.28  0.16  -0.35  -0.30  0.29  0.93  0.87  -0.24  -0.02  0.27  0.57  0.42  0.23  0.28  1.00  112  Table A3.2. Simple Pearson's correlation coefficients between available variables for winter data. stalk)  M b * . :  Iterib :  , IWb  mm  0  o  fSSSg..  ^^^^^^^^ jsreifeToissa  fMnrasa waferars g3CT<SnW5ifl  mclSirarass  1.00  0.44  0.41  -0.05  -0.01  -0.18  0.21  0.22  -0.18  -0.06  -0.06  -0.14  -0.14  -0.09  0.44 0.41  1.00 0.65  0.65 1.00  0.32 0.36  -0.17 -0.30  -0.37 -0.14  0.00 -0.10  0.21 0.08  0.11 0.26  -0.60 -0.69  -0.58 -0.37  -0.29 0.02  0.01 0.08  -0.30 -0.27  -0.05  0.32  0.36  1.00  0.37  0.30  -0.04  0.14  0.39  -0.38  -0.62  -0.13  0.15  0.36  -0.01  -0.17  -0.30  0.37  1.00  0.59  -0.07  0.04  0.18  0.51  0.11  -0.15  0.02  0.91  -0.18  -0.37  -0.14  0.30  0.59  1.00  -0.25  -0.02  0.37  0.29  0.49  0.55  0.43  0.86  0.21  0.00  -0.10  -0.04  -0.07  -0.25  1.00  0.88  -0.79  -0.01  -0.18  -0.13  -0.48  -0.16  0.22  0.21  0.08  0.14  0.04  -0.02  0.88  1.00  -0.69  -0.14  -0.19  0.01  -0.35  0.03  -0.18  0.11  0.26  0.39  0.18  0.37  -0.79  -0.69  1.00  -0.21  -0.03  0.14  0.60  0.28  -0.06  -0.60  -0.69  -0.38  0.51  0.29  -0.01  -0.14  -0.21  1.00  0.50  -0.12  -0.22  0.45  -0.06  -0.58  -0.37  -0.62  0.11  0.49 .  -0.18  -0.19  -0.03  0.50  1.00  0.64  0.30  0.35  -0.14  -0.29  0.02  -0.13  -0.15  0.55  -0.13  0.01  0.14  -0.12  0.64  1.00  6.69  0.24  -0.14  0.01  0.08  0.15  0.02  0.43  -0.48  -0.35  0.60  -0.22  0.30  0.69  1.00  0.28  -0.09  -0.30  -0.27  0.36  0.91  0.86  -0.16  0.03  0.28  0.45  0.35  0.24  0.28  1.00  113  Appendix 4 Graphs of PRESS Residuals vs. Standardized Pellet-Group Densities  0.80  < 0.60  < •  0.40  tn CO  •  J!  •  0.20  '(A «  •  *  or 0.00  •  DO.  0.20 •  1  0.40  •  •  0.60  0.80  1.  00  •  •  -0.20  *  •  •  •  -0.40  Observed Standardized Pellet-Group Density (# groups per ha)  Figure A4.1. Residuals vs. standardized observations for full summer model.  0.80  0.60  0.40  (0 !5  0.20  •  ^  •  0.00  o.bo  •  0.20  •  0.40  0.60  0.80  •  P0  -0.20  -0.40  Observed standardized Pellet-Group Density (# g r o u p s per ha)  Figure A4.2. Residuals vs. standardized observations for the summer management (reduced) model. 114  0.60 0.40 0.20 w  03 3  Ti  0.00  o.po  0.20  '0.40  0.60  0.80  Do  d) -0.20  -0.40 -0.60 -0.80  Observed standardized Pellet-Group Density (# groups per ha)  Figure A4.3. Residuals vs. standardized observations for the full winter model.  0.60 0.40 0.20 0.00 0. 00 TS '55 CD -0.20 « 3  oi  -0.40  0.20  '0.40  •  •  0.60  0.80  DO  •  -0.60 -0.80  Observed Standardized Pellet-Group Density (# groups per ha)  Figure A4.4. Residuals vs. standardized observations for the winter management (reduced) model.  115  

Cite

Citation Scheme:

    

Usage Statistics

Country Views Downloads
United States 30 0
Canada 6 9
India 4 0
China 3 1
Sweden 2 0
Saudi Arabia 1 0
Australia 1 0
France 1 0
New Zealand 1 0
City Views Downloads
Washington 22 0
Unknown 8 3
Ashburn 4 0
Shenzhen 3 1
Stockholm 2 0
Calgary 2 0
Kamloops 2 0
Federal Way 1 0
Whitehorse 1 0
Los Angeles 1 0
Jeddah 1 0
Melbourne 1 0
Etobicoke 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}
Download Stats

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0075008/manifest

Comment

Related Items