A combination of discriminant function analysis and multiple regression was used to develop a linear model of American marten Martes americana winter habitat use within second-growth boreal forest in northeastern Ontario, Canada. Four structural variables significantly discriminated between sites that were used or not used by martens: the percentage of spruce or fir trees, tree height, the number of downed logs, and canopy closure. The model was tested against a second data set and was not invalidated. The results demonstrated that martens were using second-growth forests in Ontario, and that their response to structural characteristics was similar to responses described previously in uncut forests.
The American marten Martes americana has received attention in recent years from forest managers because its habitat preferences appear to be in conflict with the interests of the timber industry. Martens generally prefer mature and overmature coniferous forests (Marshall 1951, Koehler & Homocker 1977, Soutiere 1979, Buskirk & Powell 1994) and these same forests are being heavily logged in many areas of North America. Forests in northeastern Ontario, Canada, have been commercially harvested for the past 100–200 years, and many are currently ‘second- growth’, ranging from 50 to 100 years in age; the oldest seral stages are becoming increasingly rare in the region, and it is apparent that martens are making some use of the widely available second-growth forests (Ontario Ministry of Natural Resources, unpubl. data). The habitat relationships of martens in these second-growth forests are unclear although regeneration is believed to be suboptimal due to a reduction of necessary structure and diversity compared to mature and decadent forests (Thompson & Harestad 1994). As loss of habitat has been linked to the extirpation of martens from parts of their historic range in North America (Gibilisco 1994), an understanding of what second-growth forest types are used for by martens is important, so that these forest types may be managed for in timber harvest plans. As severe conditions limit martens to a narrower range of habitats during winter (Buskirk & Powell 1994), winter habitat relationships are of particular interest.
This study was an effort to develop a model of marten winter habitat use in the second-growth forests of northeastern Ontario. We report on a data- based model developed with multivariate statistics and on a test of that model performed on a second data set.
Methods
Study area
The 1,100-km2 study area was located on Crown Land south of Timmins, Ontario (48°30′N 81°40′W) and was typical of the Northern Clay Belt forest region (Rowe 1972) which is predominantly a black spruce Picea mariana lowland. The topography was nearly level and Precambrian granite was abundant; drainage in the region was poor. Black spruce dominated the overstory although Balsam fir Abies balsamea was abundant as well and often associated with trembling aspen Populus tremuloides and white birch Betula papyrifera. Eastern white cedar Thuja occidentalis occurred in swampy patches and stands of jack pine Pinus banksiana were common in drier, sandy sections. White spruce Picea glauca was scattered along the few upland sections, but not dominant. The understory was characterised by thick patches of speckled alder Alnus rugosa, while other species such as mountain maple Acer spicatum and beaked hazel Corylus cornuta were distributed throughout the area. Common herb layer species included Sphagnum spp., Clintonia borealis, Trientalis borealis, Streptopus roseus, and Aralia nudicaulis. (For a more detailed description of vegetation in the study area, see Bowman, Robitaille & Watt 1996). The region has a history (100–200 years) of clearcut logging, and stands were predominantly mid- to late-successional (75% of stands were classed between 50 and 89 years; Ontario Ministry of Natural Resources, unpubl. data). Although registered traplines existed in the study area, we attempted to minimise potential impacts that removing martens might have on the study. Transects were not established within 1 km of traplines, and no traplines were active for martens during the snow-tracking period. Mean temperature during snow-tracking in February and March of 1994 was -13.5°C. In that period, 41 cm of precipitation fell, 100% of which was snow.
Snow-tracking
We snow-tracked to monitor martens, which limited us to making inferences about habitat use, rather than preference (Peek 1986) as we could obtain no information about fitness with our method. If the limits of interpretation are recognised, snow-tracking can be an efficient means of monitoring martens (Raphael 1994). Line transects (1 to 3 km) were first established from Forest Resource Inventory (FRI) maps of the area. Transects were then identified in the stand itself and navigated on snowshoes using hip chains and compasses. All transects were sampled once for marten use, during February and March of 1994. The intersection of marten tracks with a transect was recorded and considered to indicate marten habitat use. All 100-m intervals and points of marten intersection were marked with flagging tape. We did not perform snow-tracking after more than two consecutive days without fresh snowfall (Thompson, Davidson, O'Donnell & Brazeau 1989). We snow-tracked 142 kilometres on 57 transects.
Habitat sampling
Snow-free habitat sampling occurred during June–October 1994, and June-September 1995. As one part of a provincial government program to study wildlife habitat (Watt 1990), we assessed marten habitat use at a 100-m scale of resolution, to be consistent with other provincial studies. Vegetation was sampled systematically every 100 m (‘non-use’), or, at a point of marten snow track intersection (‘use’). To minimise spatial autocorrelation, only one track intercept was sampled when tracks were clustered within the 100 m scale (the centre intercept, systematically). Sampled plots were 10 by 10-m quadrats.
Within each quadrat, several vegetation descriptors were measured (see Table 1 for descriptions of variables). Canopy closure was determined with a spherical densiometer (Lemon 1957). A Basal Area Factor (BAF) 2 prism was used to count trees from the centre of the plot. Species and diameter at breast height (dbh) were recorded for all trees that fit the prism. Tree heights were measured with a clinometer. The number and dbh of all snags within the plot was recorded. All downed logs that intersected the plot were measured and counted. The number of stumps, rock outcrops and other irregularities were recorded. Shrub closure was determined using ocular estimates. Over the two summers, 877 quadrats were sampled for habitat descriptors on 41 transects.
Model development
Twenty-five transects were randomly selected to develop the habitat model; the remaining subset of data was used for testing the model. One factor analysis of variance was used to test for differences in habitat descriptors between 151 (marten) use and 357 non-use quadrats. Quadrats with no trees (fitting BAF2 prism; less than 5%) were eliminated from all analyses. Correlation coefficients between pairs of habitat descriptors were determined and where r ≥ 0.50, the least significant descriptor of the pair was eliminated from further analysis. A combination of linear multiple regression and stepwise discriminant function analysis (using SPSS software, SPSS Inc.) was used to develop a habitat model using presence or absence of marten as the dependent variable. Mathematically, these two multivariate approaches are similar and the resulting coefficients are proportional (Tatsuoka 1971, Brennan, Block & Gutiérrez 1986). Discriminant function analysis was used to determine classification coefficients and group separation, while multiple regression was used to calculate standard errors for the coefficients. For discriminant analysis the method of Rao's (1970) V was used, and the probability-of-F-to-enter was set at P ≤ 0.05. Prior probabilities were based on sample size. The habitat model was based on the general linear model:
where Y = a calculated value; a = a constant; b = the unstandardised coefficients (from regression or discriminant); and X = habitat descriptors entering the model. This general linear model was included in the following posterior probability equation (Green 1978, Brennan et al. 1986). The posterior probability P(l/X) of an unclassified plot being used by marten was based on a vector of habitat measurements (equation 1) and was given by:
where ql = the prior probability that the habitat is used by marten; q2 = 1 - ql; x′ = a vector of habitat measurements - the X variables from (1); k = a vector of unstandardised discriminant function coefficients and constant; tl and t2 = the mean discriminant scores (group centroids) of the use and non-use quadrats. To ensure that a linear model described these relationships, classes were developed for all descriptors in the model and the bivariate relationship for marten usage of each class was tested for linearity.
Model testing
We assessed our habitat model by trying to invalidate it through comparison with a second data set. The test sample was 293 (69 use and 224 non-use) quadrats, and was used to assess individual variables and overall model output (sensu Schamberger & O'Neil 1986). The performance of variables in the habitat model was tested with a second discriminant function analysis, again using use and non-use quadrats as the dependent variable. Univariate-F-ratios were used to examine between group differences. Methodology of the first discriminant was repeated, however the test sample was classified twice, using functions from the model test discriminant analysis and functions from the original model development. Chance-corrected classification rates were compared using Cohen's (1960) Kappa-statistic, as suggested by Titus, Mosher & Williams (1984), Capen, Fenwick, Inkley & Boynton (1986) and Rexstad, Miller, Flather, Anderson, Hupp & Anderson (1990).
The performance of model output was assessed by generating posterior probability scores (equation 2). The distributions of these scores for use and non-use quadrats were compared using a Kolmogorov-Smimov (K-S) test. Further, we tested for a positive linear relationship between marten habitat use and the posterior probability of quadrats being used by marten. Ten classes of probabilities (0.0 – 0.1 to 0.9 – 1.0) were developed and compared to marten selection of each class (use/availability) with a Spearman rank correlation.
Results
Model development
After controlling for possible multicollinearity, 16 environmental variables were included in multivariate analyses - nine of these variables had significant between-group differences at the bivariate level (Table 2).
Marten-use and non-use quadrats could be discriminated significantly (Wilks' lambda = 0.83; df = 1, 4; P < 0.001; Rao's V) on the basis of a linear combination of four weighted variables: percentage of spruce or fir, tree height, number of downed logs, and canopy closure. Prior probabilities were 0.70 and 0.30 for use and non-use quadrats, respectively. Unstandardized discriminant function coefficients served as the basis for the model and were proportional to multiple regression coefficients and standard errors by a factor of 5.8 (Table 3). Mean discriminant scores were -0.29 for the non-use quadrats and 0.70 for the use quadrats. The discriminant function classified 72% of cases correctly (Kappa = 0.22 ± 0.11; z = 3.86; P < 0.001). Pooled within-group covariance matrices were not equal (F = 6.69; df = 1; P < 0.001; Box's M test). The posterior probability of a quadrat containing marten tracks was given by:
where x′k = [0.0259 x spruce or fir (%)] + [0.0814 × median tree height (m)] + [0.1231 × number of downed logs] + [0.0140 x canopy closure (%)] - 3.9071.
Table 1.
Description of variables sampled on 10 by 10-m quadrats in boreal forest near Timmins, Ontario, Canada, and subsequent categories developed for analyses.
Table 2.
Mean and standard error of variables sampled on 10 by 10-m quadrats in boreal forest near Timmins, Ontario, Canada. Quadrats had previously been sampled for use or non-use by martens.
Table 3.
Order of variable selection, coefficients, and error terms for multivariate analyses of habitat descriptors from 10 by 10-m marten-use and non-use quadrats sampled in boreal forest near Timmins, Ontario, Canada.
Analysis of use-availability data for the four model variables indicated that all of the bivariate relationships were significantly linear: percentage of spruce or fir (F = 50.76; df = 4, 503; P < 0.001), median tree height (F = 13.39; df = 18,489; P = < 0.001), downed logs (F = 10.01; df = 5, 502; P < 0.001), and canopy closure (F = 7.98, df = 4, 503; P < 0.01).
Model testing
Univariate-F-ratios demonstrated that all four variables in the model had significant differences between 69 use and 224 non-use quadrats in the test sample: percentage of spruce or fir (F = 4.23; df = 1, 292; P < 0.05), tree height (F = 5.15; df = 1, 292; P < 0.05), number of downed logs (F = 17.44; df = 1, 292; P < 0.001), and canopy closure (F = 4.18; df = 1, 292; P < 0.05). Use and non-use quadrats were discriminated (Wilks' lambda = 0.90; df = 1, 3; P < 0.001; Rao's V) on the basis of a linear combination of three variables: number of downed logs, percentage of spruce or fir trees, and canopy closure. Tree height did not significantly improve the fit (P = 0.10) and did not enter the function. Prior probabilities were 0.24 and 0.76 for use and non-use quadrats, respectively. Mean discriminant scores were 0.60 for use quadrats and -0.18 for non-use quadrats. Pooled within-group covariance matrices were not equal (F = 5.53; df = 1; P < 0.001; Box's M test). Correct classification was not significantly higher than chance (Kappa = 0.17 ± 0.17; z= 1.88; P> 0.05); 13.0% and 95.5% of use and non-use quadrats, respectively, were correctly classified. When this sample was classified using discriminant functions from model development, correct classification was significantly greater than according to chance (Kappa = 0.11 ± 0.10; z = 2.24; P < 0.05); 75% and 40% of use and non-use quadrats, respectively, were correctly classified.
Figure 1.
Selection indices (SI) of marten (based on use/availability) for 10 classes of posterior probability (PP) scores: 1 (0.0 - 0.10), 2 (0.11 - 0.20), 3 (0.21 - 0.30), 4 (0.31 - 0.40), 5 (0.41 - 0.50), 6 (0.51 - 0.60), 7 (0.61 - 0.70), 8 (0.71 - 0.80), 9 (0.81 - 0.90), and 10 (0.91 - 1.0). The Spearman correlation between SI and PP class was significant (rs = 0.96; N = 293; P < 0.001).
A Kolmogorov-Smimov test demonstrated that the distribution of posterior probability scores for use quadrats (mean = 0.40 ± 0.02) was significantly shifted to the right of the distribution of scores for non-use quadrats (mean = 0.28 ± 0.01) (z = 2.56; P < 0.001). Marten selection indices were positively correlated with 10 classes of posterior probability scores (rs = 0.97; N = 293; P < 0.001) (Fig. 1).
Discussion
During winter, cold temperatures and difficulties in foraging through deep snow constrain martens to habitats which offer thermal microenvironments, and sufficient resting, denning, and foraging opportunities (Buskirk 1984, Bateman 1986, Thompson 1986, Buskirk & Powell 1994). Previous winter habitat models for martens (models typically developed from literature reviews) have been based on the ‘martens require old-growth’ paradigm (Harrison, Phillips, Chapin, Katnik & Hodgeman, in press) and have consistently included variables representing tree species composition, seral stage, coarse woody debris, and canopy closure (Allen 1982, Ritter 1985, Larue 1993). The included variables are believed to be the important, measurable habitat characteristics which offer martens thermal microenvironments and sufficient resting, denning, and foraging opportunities in winter. Our data-based model was the first to describe winter habitat use in the second-growth forests of Ontario, and the inclusion of similar variables, compared to the previous ‘old-growth’ models, demonstrated that martens in second-growth were responding to similar structural characteristics as those described previously in mature and overmature forests. Martens in the study area were constrained to habitats with high values for four variables: percentage of spruce or fir trees, tree height, number of downed logs, and canopy closure. More generally, closed spruce-fir forest with ample downed logs.
The structure characterised by the habitat model reflected the suggested winter survival strategies of martens. These strategies include: staying under cover to avoid avian (Hawley & Newby 1957) and terrestrial (Thompson & Harestad 1994) predators; and accessing the subnivean layer to rest, den, and forage (Hargis & McCullough 1984, Bateman 1986, Com & Raphael 1992). Closed spruce-fir forests with ample downed logs provide more cover in winter than other forest types, and facilitate martens penetrating the subnivean layer by providing snow breaks caused by the woody debris (Buskirk & Powell 1994).
Larsen (1980) states that in Precambrian shield communities (such as the study area), closed spruce forest is the most advanced stage of forest succession (i.e. Clementsian succession; Clements 1916). An alternative suggestion is that forest ecosystems are dynamic, but contain stable (and different) ‘site type’ assemblages, often on the same landscape (West, Shugart & Botkin 1981). Larsen's interpretation suggests that martens in the study area were using the most advanced seres available, a conclusion supported by a majority of marten research (for reviews see Buskirk & Powell 1994, Buskirk & Ruggiero 1994). The more recent interpretation indicates that martens were using particular site type assemblages (i.e., closed spruce-fir forests; sensu Bowman et al. 1996). Regardless of interpretation, martens exhibited significant selectivity within the range of ‘second-growth’ cover types.
Performance and validity of habitat model
Testing of a second data set at both the variable and output levels (Schamberger & O'Neil 1986) did not invalidate the model developed here. Each of the four model variables had significant between-group differences, however, there was shared variance among the variables, as only three entered the discriminant function. Correct classification of the test sample with discriminant functions from model development was significantly greater than by chance, although this classification rate was low. Titus et al. (1984) suggest that correct classification rates are often low as a result of unequal prior probabilities. In the test sample, prior probabilities were 0.76 and 0.24 for non-use and use quadrats, respectively. A second potential cause of low correct classification rates was our inability to be certain that a non-use quadrat was in fact, not used at all by martens.
Overall model output was significantly higher in use plots than in non-use plots, which indicated that the model was successful at predicting when a sampled quadrat would contain marten tracks. The positive linear relationship between posterior probability scores and marten selection indices was further evidence of this success.
The scope of this model is somewhat limited by the methods used in its development. First, snow-tracking only reflects habitat use by martens which are travelling on the snow's surface. While this provides useful insight (Hargis & McCullough 1984, Thompson et al. 1989, Raphael 1994), it means that there is no direct consideration of arboreal- or subniveanspace use in the model, and therefore, not all habitat use by martens is described. Conservatively, our model describes third-order habitat selection (sensu Johnson 1980) of martens travelling on the snow's surface. Secondly, the model was developed to describe habitat use within a narrow range of conditions (Clay Belt second-growth forests); the study area was relatively homogeneous compared to those of other published studies (e.g. burned vs unburned sites, cut vs uncut sites; Soutiere 1979, Thompson 1994). This is one reason why martens exhibited low (although significant) selectivity in the study; 17% of variance between use and non-use quadrats was explained by the model. In addition to being limiting, the efficacy of the methods used to develop the model can be questioned. The scale of resolution chosen for this study (100 m) may have been inadequate; some martens range widely (up to 1,600 ha; Buskirk & Macdonald 1989) and sampling at a 100-m scale could result in spatially autocorrelated data. Regardless of sampling scale, recent thinking suggests that independent samples might be an unrealistic objective (Legendre 1993); autocorrelation is often present in ecological data. We are proceeding with a spatial analysis of these data to take advantage of potential autocorrelation (J.-F. Robitaille, unpubl. data). The lack of a spatial component in the present approach required ‘traditional’ parametric multivariate techniques, which involved assumptions of multivariate normality, homogeneity of dispersions, and known prior probabilities. Violation of these assumptions is common in ecology (for reviews see Williams 1983, Rexstad, Miller, Flather, Anderson, Hupp & Anderson 1988). Testing of univariate distributions indicated that several environmental variables departed from normality and pooled within-group covariance matrices were not equal during model development. Williams (1983) states that when assumptions are violated, analysis must be considered exploratory (sensu Tukey 1980). Capen et al. (1986), Rexstad et al. (1990), and Taylor (1990) agree that discriminant analysis models should be derived from large data sets, and applied to independent data to test validity. Our model was derived from a large data set (508 cases x 16 variables) and was not invalidated when tested against an independent sample; these findings suggest that our habitat model may be considered confirmatory (sensu Tukey 1980). We think that despite its limitations, the model provides useful insight into the habitat use of martens in second-growth boreal forests.
Models of ecological systems will never completely mimic the real world; they are simplifications of relationships. A successful model provides the essence of a relationship while removing extraneous information (Levin 1992). Our model indicates that closed spruce-fir forests with ample downed logs were the most used (winter) habitat of martens in the second-growth forests of northeastern Ontario. These same forests are eventually harvested for wood pulp - a conflict which must be recognised if martens (and associated species) are to have adequate resources in managed boreal forests.
Acknowledgements
financial support for this study was provided by Northeast Science and Technology, Timmins Region of the Ontario Ministry of Natural Resources, and we thank the Timmins District Office of the OMNR for logistic support and technical advice. We also thank winter and summer field crews for their enthusiastic effort. This material, in various stages of development, has benefited from the thoughtful comments of RJ. Beckett, J.A. Bissonette, J. Kierstead, H.D. MacPherson, F.F. Mallory, R. Otto, J.-R Ouellet, I.D. Thompson, W.R. Watt, and two anonymous reviewers.
