Study area

The study was conducted on the southeast facing slope of mount Chasseral, a prominent rise in the central Jura chain in western Switzerland (47°07′N, 07°03′E). The study area covered 25 km² of forest situated at altitudes of 940–1,560 m a.s.l., and included the locations where Zbinden (1979) studied habitat requirements and diet of the hazel grouse. East, southeast and west of the study area, the forest is more or less contiguous, while it borders on cattle pastures on the upper side and on gradually open grasslands along its southwestern margin. Forest structure and tree species composition have been shaped by forestry and cattle grazing for a long time. Forests consisting of stands of either pure beech Fagus sylvatica or mixed beech-spruce Picea abies make up 45% of the forested area. Another 35% are highly structured mixed stands of trees and shrubs of sycamore Acer pseudoplatanus, whitebeam Sorbus aria, rowan S. aucuparia, hazel Corylus avellana, willow Salix sp., buckthorn Rhamnus cathartica, and red-berried elder Sambucus racemosa, which indicate former cattle grazing. The remaining 20% are windthrow areas (5%) and other non-forested patches (15%), e.g. rocky outcrops.

Mapping the dependent variable: hazel grouse distribution

Hazel grouse were mapped in the field in October 1999. The grouse show territorial activity in autumn and are then usually found close to the centre of their home ranges (Kämpfer-Lauenstein 1995). Following the method of Huboux et al. (1994), a whistle (Bergmann et al. 1996: Fig. 11d) was blown to imitate hazel grouse song every 100 m along the 100 m contour lines across the entire study area. Distances on the ground between contour lines ranged within 300–400 m, which resulted in an average number of six lines along the 15–30 degree slopes of the study area. The position of the responding birds was then recorded with a handheld GPS receiver (GARMIN GPS 45). Additionally, all signs of grouse presence (e.g. droppings and feathers) were mapped by thoroughly searching a 100-m wide band along the transects (i.e. within 50 m on each side of the contour lines). This produced a coverage of more than a third of the study area. Altogether, hazel grouse presence was recorded at 143 sites (hereafter referred to as ‘grouse record’). These consisted of 131 sites with droppings, three with feathers, three song responses, and six sight records. Singing activity seems to be generally low in Swiss hazel grouse territories when compared to grouse in boreal forests (W. Suter, pers. obs.), and low willingness to react to imitated song was recorded earlier (Zbinden 1979). Because of the pronounced territorial activity of hazel grouse in October, it is unlikely that many droppings were found at feeding sites that were otherwise unsuitable habitat.

Sampling habitat variables from aerial photographs

We used infrared aerial photographs (scale 1:10,000–1:25,000 mostly from 1998 and 1999, plus a few from 1995 to complete coverage) to assess habitat variables at the landscape level (Table 1). All interpretations of aerial photographs were done manually using an analytical stereoscope (Leica Aviopret APT 2). A bird-centred analysis, BCA (Larson & Bock 1986), was performed to identify hazel grouse habitat variables on the photographs. In doing so, each grouse record was located on the photographs, and a circle with a radius of 50 metres on the ground was drawn around it, resulting in a circular interpretation area of 0.8 hectares (hereafter referred to as ‘grouse plot’). The radius used was a compromise since it had to be large enough to reflect the spatial requirements of hazel grouse yet small enough to be homogenous with respect to the main habitat variables. Therefore, grouse records were admitted to further analysis only if the centres were at least 100 m apart. The selection process was random and ensured that spatial autocorrelation was minimised by keeping the records spatially separate. The resulting data set contained N_gp = 101 grouse plots (out of the 143 original grouse records). We then generated a corresponding pseudo-absence data set with N_ngp = 101 data points for locations without hazel grouse records (hereafter referred to as ‘non-grouse plot’). To do so, we divided the study area into N_gridcell = 101 equal-sized grid cells, each of which contained one non-grouse plot. Its exact position in the cell was determined with random numbers representing the coordinates, but it had to be at least 100 m away from any grouse plot centre. We chose this combination of regular and random design to make sure that the ‘non-grouse’ habitat variables were sampled in their full range within the relatively diverse study area.

Table 1.

Habitat variables interpreted on infrared aerial photographs within a 50-m radius around grouse and non-grouse records (NFI = Swiss National Forest Inventory; Keller 2001)

We determined the habitat variables to be used in the models in two steps. First, we used the literature on habitat requirements of the hazel grouse to establish a comprehensive set of forestry and biology related variables. From this list, we chose the variables that were retrievable from aerial photo interpretation, and built an interpretation key. This key was then tested on the ground (ground truthing) and improved accordingly. Two variables (shrub and herb layer) were difficult to assess under fully closed canopy. Using field data, we defined cover of shrub and herb layer to be 0% in closed stands, while herb layer was set at 5% when there were small canopy gaps present. Finally, 18 variables were interpreted at the 202 grouse and non-grouse plots on the infrared aerial photographs, and 12 of these were used in the models (see Table 1).

Statistical analyses

We applied univariate descriptive statistics to explore differences in relevant habitat variables between grouse and non-grouse plots. We then fitted habitat suitability models using the statistical R environment 2.0.0 (R Development Core Team 2005). To do so, the number of variables was first reduced on the basis of a Bonferroni-corrected Pearson correlation analysis (see Cabin & Mitchell 2000). Where two variables correlated significantly (P < 0.05), we chose the one better explaining bird presence/absence. Explanatory power of individual variables was determined based on the Akaike information criterion, AIC (Chambers & Hastie 1992), by applying simple logistic regressions including only one independent variable at a time. If two correlated independent variables had a similar predictive power, the one that was more management oriented was used in the further modelling procedure. This approach is compliant with criterion 5 of Mosher et al. (1986) since all structural variables are assumed to be relevant to the distribution of hazel grouse at the landscape scale. A logistic regression model (link function = logit) of the GLM family (McCullagh & Nelder 1989) was then calibrated because bird presence/absence values are expected to follow a binomial distribution. All powers and interactions up to the second order were included for model calibration. Standardisation through a z-transformation of the variables and conversion of ordinal into continuous variables (see Table 1, where unit = cat.) were studied, but did not improve the modelling output. Finally we optimised a full model including all independent variables through a stepwise variable selection procedure based on the AIC.

The performance of the resulting habitat suitability model was evaluated through a k-fold cross-validation (Fielding & Bell 1997, Guisan & Zimmermann 2000), a procedure we chose because of the limited sample size. Grouse and non-grouse plots were each randomly partitioned into five equal-sized groups (k = 5). Model parameters were then fitted five times on ⅘ of the data and evaluated on the remaining ⅕ so that each point was left out once for testing. Model performance was then evaluated based on the AUC, the threshold-independent area under the receiver operating characteristic curve (ROC; Hanley & McNeil 1982). To estimate the robustness of the 5-fold data partitioning we repeated the procedure 10 times with individual random sample splits each. Kappa (Cohen 1960) was used to determine the appropriate cut-off level of the habitat suitability predictions to discriminate between suitable and unsuitable habitats (Fielding & Bell 1997).

A classification tree (Breimann et al. 1984) was additionally calibrated with all independent modelling variables as predictors. We selected this non-parametric modelling approach (in contrast to the parametric GLM) to evaluate the robustness and generality of the GLM-based habitat suitability model, as well as the hierarchical interactions of the model variables. First, a full tree was grown with a minimum of three observations per terminal node. Second, a 20-fold cross-validation was repeated 20 times in order to find the optimised number of terminal nodes for pruning the tree, and ultimately to avoid over-fitting of the classification tree model (Chambers & Hastie 1992). Evaluation of the tree model was performed based on a repeated 5-fold cross-validation, similar to the GLM. This latter cross-validation was not intended to prune the tree any further, but to evaluate the accuracy and robustness of the model.

Spatial model application

The resulting predictive models were applied to the entire study area. Thereby, the significant model variables first had to be derived in a spatially explicit way from the aerial photographs. We delineated areas that were homogenous with respect to the significant model variables and used the same photo interpretation key as in the sampling approach. Of these areas, 70% were < 1 ha, 16% 1–2 ha, and 14% > 2 ha. Subsequent analytical photogrammetry (Leica BC2000S) included geo-referencing with spot heights from the Swiss topographic map (1: 25,000) and restitution of vector data. We imported the resulting polygon layers in a Geographic Information System (ArcInfo^® 8.3) and completed them by adding the pertaining variable information to each polygon (N = 802). Afterwards a grid layer was created for all but one of the variables and geo-coded with the Swiss digital elevation model (spatial resolution = 25 m). Only the forest edge boundaries were first calculated per hectare and then interpolated through the inverse distance weighting method (routine IDW in ArcInfo^® 8.3) to produce an analogous grid layer. By running the final predictive habitat models with the generated variable grids, final maps containing habitat suitability information in each grid cell were created.

Habitat variables

The predominant habitat on 84% of all hazel grouse plots at the Chasseral was classified as true forest. The surroundings (i.e. the matrix; see Table 1) of the grouse plots, however, consisted of only 60% forest; the remaining 40% were more or less open forest-pastures. Grouse records clustered near the upper and lower altitudinal limits of forests. Plotting the sample plots onto the forest map showed that grouse plots tended to be closer to outer forest edges than non-grouse plots. The mean forest edge density on grouse plots (155 m/ha) was more than twice that on non-grouse plots (64 m/ha; W = 8,407, P < 0.001).

Forest areas containing hazel grouse were vertically well structured. Stands were multi-storied on 75% of the grouse plots, but two- and single-storied on only 13% and 12%, respectively. Vertical forest structure on non-grouse plots was significantly more uniform (χ² = 57.9, df = 4, P < 0.001), with only 27% multi-storied, 13% two-storied, but 44% single-storied. Likewise, hazel grouse habitat was significantly more diverse in horizontal structuring (χ² = 38.8, df = 4, P < 0.001). Among all hazel grouse plots, 42% were in semi-open stands consisting of tree groups of different sizes, 28% on closed stand edges, 20% within homogeneously closed stands, and 10% on open glades or clearings. By contrast, 47% of the non-grouse plots were in homogeneously closed stands. Stands with hazel grouse had mainly intermediate canopy closure and were significantly different from non-grouse plots which tended to be either open or very closed (χ² = 62.7, based on 10,000 Monte Carlo simulations, P < 0.001).

A similar distribution pattern emerged if cover of the tree layer was expressed as a continuous variable (Fig. 1A). Cover of shrub (see Fig. 1B) and herb layer (see Fig. 1C) also differed between grouse and non-grouse plots. Hazel grouse habitats tended to have intermediate cover by shrub and a dense herb layer whereas non-grouse plots had little or no shrubs at all (see Fig. 1B). Stands with a herb cover below 20% were rarely used by hazel grouse (see Fig. 1C). On average, forest stands were younger on grouse plots than on non-grouse plots, moderately so for coniferous but strongly so for deciduous trees (Fig. 2). Only with respect to tree type were the two grouse categories rather similar, the ratio of deciduous to coniferous trees being 1.5 on grouse plots and 2.0 on non-grouse plots.

Figure 1.

Distribution of tree (A), shrub (B) and herb (C) cover at grouse plots (▪) and at non-grouse plots (□).

Figure 2.

Distribution of development stages (grey shades) of deciduous and coniferous trees at grouse and non-grouse plots. Each bar represents all 101 grouse or non-grouse plot samples. Hence, plots without the respective tree type present were assigned to the ‘absent class’ (no shade).

Predictive habitat suitability models

The final generalised linear model retained five of the 12 independent variables entered (Table 2). The most influential variable, forest edge density, was related to general forest structure. Then, two direct resource-related variables (Guisan & Zimmermann 2000) contributed significantly to the final model, i.e. cover of shrub and herb layer. Evaluation of the GLM through a 5-fold cross-validation resulted in a mean AUC of 0.90. Kappa statistic was highest (0.62) at a cut-off level for predicted habitat suitability of 0.50, where the correct classification rate was 0.81.

Table 2.

Final generalised linear model, GLM, to predict hazel grouse habitat suitability (0–1). ** indicates significance at < 0.01 and * indicates significance at < 0.05.

The final classification tree (Fig. 3) confirmed the findings from the GLM, having three out of six relevant variables in common (forest edge density, shrub cover and vertical stand structure). Shrub cover was the main variable to partition the grouse and non-grouse plots. Similar to the GLM structural variables, such as forest edge density, canopy closure and vertical stand structure, were relevant for the primary partitions of the TREE. The more terminal partitions additionally included development of deciduous trees. Evaluation of the final TREE through a 5-fold cross-validation resulted in a mean Kappa of 0.80, and a mean correct classification rate of 0.90.

Figure 3.

Final classification tree (TREE) to predict hazel grouse habitat as suitable or unsuitable.

The two models thus performed similarly. Both indicated that at the landscape level primarily structural and secondarily resource-related variables were relevant for hazel grouse habitat suitability mapping on aerial photographs, with shrub cover being a key variable. The spatially explicit application of the models to the study area identified suitable habitat mostly concentrated near the outer forest edges or near larger forest openings. The spatial patterns of suitable habitat were similarly predicted in both model applications. However, as was the case with the accuracy measures, the GLM was more restrictive than the TREE.

Hazel grouse habitat characteristics and local distribution

Habitat requirements of the hazel grouse have repeatedly been studied at various sites across its vast distribution area, and results range from meticulous descriptions (summarised e.g. by Swenson 1995, Bergmann et al. 1996, Zbinden & Blattner 1998) to more formal quantitative analyses (e.g. Swenson & Angelstam 1993, Åberg et al. 2000a, Åberg et al. 2003, Sachot et al. 2003). Habitats used span a broad spectrum, from ancient deciduous coppice forests in central Europe to untouched taiga forests in the Asian Far East. However, they show strong structural similarities, having in common multiple tree layers, dense shrub cover mixed with small openings, and the presence of certain softwood species providing winter food in the form of fruits, buds and catkins. Hazel grouse may thus be seen as a habitat specialist adapted to early succession on small-scale gaps along old-growth forests (Swenson 1995, Bergmann et al. 1996).

Our results of the habitat analysis fully comply with the previous knowledge. Hazel grouse at Mt. Chasseral selected multi-storied stands that were also well structured in the horizontal dimension, with unevenly distributed trees resulting in intermediate canopy closure, and a moderately dense shrub cover with young deciduous trees and a thick herb layer. These conditions were fulfilled in the study area almost exclusively along the upper and lower altitudinal limits of the forested area, primarily as a result of former cattle grazing (Zbinden 1979), which is now on the decline. In the forest interior, the few scattered records were mostly along road and windthrow openings, probably of birds feeding on berries in autumn (Zbinden 1979). Parts of the study area of Zbinden (1979: Fig. 1), situated in the westernmost quarter of our study area and containing a large former windthrow area inside the forest, are now largely abandoned, probably because the regeneration has grown up to homogeneous stands unsuitable for hazel grouse. Thus, contemporary distribution of hazel grouse depends strongly on past forest uses and occasional natural disturbances.

Habitat variable assessment from aerial photographs

The principal aim of our study was to evaluate whether suitable habitats could be identified on aerial photographs. Thus, we used habitat variables that were identifiable on colour infrared photographs up to the scale of 1:25,000. Since structural variables are particularly important for defining hazel grouse habitat, and since structural qualities, such as stand composition, developmental stage or cover of the different layers, generally involve height differences, they can easily be identified on aerial photographs by using stereoscopic analysis techniques. We also found that the spatial resolution on the photographs was fine enough to represent vertical and horizontal structural properties relevant to hazel grouse, although it was at its lower tolerable limit. An emerging remote sensing technique with strong characteristics in extracting structural forest properties over areas that are comparably large is LiDAR, Light Detection And Ranging, (e.g. Riaño et al. 2003, Mathys 2006), therefore representing a potential future source for habitat data. But the importance of resource related variables, e.g. shrub or herb layer, shows that spectral information is needed as well. We conclude that the stereoscopic aerial colour infrared photography is an appropriate tool to derive hazel grouse habitat properties from both structural and spectral features visible on the photographs.

The bird-centred analysis approach (BCA) used on the aerial photographs proved to be an efficient sampling technique, since explanatory habitat variables were derived from sites where birds were present (‘grouse plots’) or absent (‘non-grouse plots’), thus reducing the total area to be mapped for model generation (Karl et al. 1999). Having distinct reference areas, i.e. a defined interpretation circle, also facilitated the interpretation of variables with variability in space (e.g. canopy closure representing tree cover per area).

Predictive habitat suitability modelling

Both models, GLM and TREE, fulfil the six criteria identified by Mosher et al. (1986) as a prerequisite of predictive models. The selected interpretation and modelling variables had to meet the requirements set by existing knowledge on the ecology of the hazel grouse (criterion 5), the characteristics of the aerial photographs, as well as applicability to forestry. Therefore mostly forestry-related variables, compatible with the Swiss national forest inventory, NFI (Keller 2001), were chosen (criterion 4). This ensured that the final model was easy to apply on aerial photographs and in the field (criterion 2) by natural resource managers (criterion 3).

Despite their different statistical approaches (parametric vs non-parametric), the two models had similar prediction accuracies, and achieved correct classification rates of 81 and 90% which were above the suggested minimum of 80% (criterion 1). We believe that it is mostly the technique used, i.e. the interpretation of infrared aerial photographs in combination with the BCA, which accounts for the high prediction accuracy. Perhaps some site aspects, such as stand composition within study area, may also be responsible for the outcome. At the Chasseral there were zones of unfavourable, mature commercial forests rather distinctive from suitable hazel grouse habitat, and this contrast might have had a facilitating effect on model calibration and evaluation. However, the grid-based sampling of non-grouse plots should have reduced this potential bias. Nevertheless, we recommend further investigations when transferring our models to topographically different areas. As Åberg et al. (2000a) experienced in Sweden, parameterisation and subsequent model building can be much more difficult in areas with subtle habitat variation. But then, forests of the Chasseral and their composition are typical for the general area of focus, i.e. the Jura mountains. We therefore suggest that these models be applicable for the Jura mountain range and similar pre-Alpine forests shaped by forestry in Switzerland and elsewhere (criterion 6).

Additionally, the robustness of the approach is supported by the fact that the variables retained by both models were either the same or represented the same structural aspects. TREE is a non-parametric, rather data-driven method. This method seems to better explain the variability of our data set than the parametric GLM method. However, another relative performance of the respective models at different spatial scales cannot be excluded (Thuiller et al. 2003).

Habitat suitability model application

When applied to the whole study area, the two predictive models produced similar results that are mostly consistent with known hazel grouse territories. Both models tend to overpredict habitat suitability due to the hectare-based variable ‘forest edge density’, but to a different degree because of the different weights in the models. Areas modelled as suitable but actually without grouse records consisted of scattered forest fragments within a relatively open pasture. These areas tend to be too small or too isolated from each other, given the poor dispersing propensity of hazel grouse (Åberg et al. 1995, Saari et al. 1998, Aberg et al. 2000b, but see Montadert & Leonard 2006). Since our model is a habitat suitability model sensu stricto, it does not include variables relating to spatial configuration of habitat patches such as connectivity or isolation. One way of including this effect into a model would be to generate and include variables expressing landscape metrics. However, we do not consider ignoring landscape metrics a drawback but rather an advantage for the practical usefulness of the model, as it will identify suitable habitat strictly from resource properties. The model will thus not classify a suitably structured habitat patch as unsuitable simply because the patch is too isolated. Managers will be interested in knowing the full extent of potential habitats. Planning measures to improve unsuitable ‘matrix’ for linking suitable patches will increasingly become important in central Europe (Bergmann & Klaus 1994), and knowledge on dispersal ability of the hazel grouse can be easily applied in the process (e.g. Åberg et al. 1995, Saari et al. 1998, Aberg et al. 2000b, Montadert & Leonard 2006) to decide where measures are warranted. Combining structural and spatial aspects would however be important in a model that focuses on understanding a given hazel grouse distribution pattern (see Sun et al. 2003). We conclude that habitat suitability for hazel grouse can be efficiently modelled from data combining a bird-centred analysis with information obtained from aerial colour-infrared photographs. Such models will then be a powerful tool in grouse-oriented forest management at the landscape level, particularly in the planning phase that involves the production of species action plans.