IntroductionTop

Hunting has been part of the history of mankind. Initially had a utilitarian purpose, to bring down wild beasts and/or procure sustenance. Although subsistence hunting still remains in many regions of the world (Shanley et al., 2013; Fragoso et al., 2016), in developed countries, hunting has become an important traditional and recreational activity (Fischer et al., 2012). For example, in Europe there are 7 million hunters and the hunting activity worth 16 billion Euros (FACE, 2016). However, there are well-documented examples of undesirable consequences caused by hunting (Coltman et al., 2003; Jerozolimski & Peres, 2003; Torres-Porras et al., 2014; Perea et al., 2014). On the other hand, hunting can be a conservation tool, playing a role in biodiversity conservation and in the rehabilitation of wildlife areas, permitting income generation without compromise wildlife population growth (Fletcher et al., 2010; Mustin et al., 2011; Estrada et al., 2015). For this reason, hunting has become an important issue in conservation biology, therefore monitoring and control of this activity is primordial to guarantee its sustainability.

In Spain, hunting is an important socio-economic activity, there are more than 43 million hectares dedicated to hunting activity, divided in more than 30,000 hunting states (MAGRAMA, 2015). Each one governed by a mandatory technical document called Technical Hunting Plan (THP; see more information about THP in methods section). Hence, the governmental agencies responsible to manage hunting (we shall use Administration from here on), commonly are understaffed and underfunded and are unable to adequately do this job, because this high number of small private hunting states, associated exclusively by political-administrative criteria hinder the sound management of hunting territories.

For this reason, managing these thousands of hunting states is not a pragmatic solution; it is necessary to use a management unit larger than the hunting state. The game management needs to be understood at a landscape scale, and recently it has been claimed the need of analyse the management at the level of ecologically significant management units (Selier et al., 2014; Torres-Porras et al., 2014). In this sense, the development and implementation of innovative management techniques as hunting regionalization, is a new management tool that will provide improvements in planning and development sustainable hunting and biodiversity conservation. The hunting regionalization has its basis in ecoregionalization, that is the fusing of the ecological concept of ecosystems with the geographic concept of regions (Loveland & Merchant, 2004). Accordingly, a hunting region is an area that has, by environ­mental, ecological, socioeconomic and game species characteristics, a degree of internal homogeneity and significantly different from other areas (Vargas, 1997).

There are two approaches to defining ecoregions, one is the environmental approach, which involves classification on the basis of environmental factors such as climatic and topographic parameters; the other approach is using biological variables to define regions. The use of biological variables can be problematic because its long-term instability. For example, the use of vegetation for the establishment of hunting regions is inappropriate because depends directly on land use changes (Fernández Alés et al., 1992). The same applies to the use of distribution data (e.g. The Global Biodiversity Information Facility [GBIF] or Atlas of the Terrestrial Mammals of Spain), because are a temporal approximation to the chorology of the species and therefore are subject to continuous modifications, extensions and reinterpretations (Delibes & Palomo, 2007). Thus, we propose a method of defining an environmental classification scheme within which hunting units may then be described. This first step is required in order to define and delineate biogeographic communities and to produce consistently defined provincial or national maps (Vargas et al., 2006a; Lechner et al., 2016).

The main objective of this study was to develop a new game management tool through the establishment of a hunting regionalization system, in order to facilitate and improve the management of hunting in a large region. The regionalization system proposed require three partial objectives: (i) to determine the environmental units on a set of climatic and topographic variables; (ii) the characterization of the environmental units on a set of hunting variables; and (iii) to generate different hunting regionalization models to select the best of them. A secondary objective was to define the harvest rate and the main species harvested in the new hunting regions.

Material and methodsTop

Study area

We used Castilla-La Mancha region, in central Spain, as case study since it has important role in the Spanish hunting activity. There are seven big game species and 40 small game species (DOCM, 1996), around 7 million hectares with active hunting management (more than 84% of the surface of Castilla-La Mancha) and more than 5,000 hunting estates (MAGRAMA, 2015). The main small game species are red-legged partridge (Alectoris rufa), wild rabbit (Oryctolagus cuniculus) and Iberian hare (Lepus granatensis), while the principal big game species are red deer (Cervus elaphus) and wild boar (Sus scrofa) (Ríos-Saldaña, 2010). The number of hunting estates may vary from year to year, although during the study period 5,322 preserves were registered, of which 1,096 (20.5%) correspond to big game estates, 4,194 (78.4%) to small game and 32 (1.1%) to intensive estates (characterized by the legal release of farm-reared animals throughout the hunting season, in order to artificially increase game harvest). Physiographically, Castilla-La Mancha has a plain (La Mancha Plain), surrounded by mountain chains including the Central Range to the north, the Iberian Range to the northeast, and the Sierra Morena to the south (Fig. 1). The predominant altitude (600-1000 m) encompasses 66% of the total surface area. The climate is typically Mediterranean, with relatively mild winters and severe summer droughts (Pons, 2011).

Technical hunting plan (THP)

THPs have been mandatory since 1990 (Gálvez, 2004) and periodically submitted to the appropriate Regional Government (every 5 years in Castilla-La Mancha). In these documents, all game activities (management practices, hunting modalities and bags, etc.) are reported with the ultimate aim of protecting and fostering game richness (Vargas et al., 2006b). A total of 5,322 THPs (100%), in force in 2005, were available for analyses. We excluded for analysis the intensive estates (32 THPs) due to their peculiarities.

Environmental variables

We determined environmental units on a set of climatic and topographic variables that were available as digital coverage’s (Table S1 [suppl]) and the 919 municipalities as work units. Altitude is released as a digital coverage by the Land Processes Distributed Active Archive Center, located at the US Geological Survey’s EROS Data Center (http://LPDAAC.usgs.gov). Slope was calculated from altitude through the Idrisi SLOPE command (Eastman, 2001). Climatic data variables result from records of 30 years and are generally representative of present climatic conditions (Font, 1983) and were digitized using CartaLinx 1.2 software and processed using the Idrisi32 GIS software (Barbosa et al., 2003). Climatic and topographic variables have been normalized (Box-Cox transformation), as ­found necessary, for entry into a Principal Components Ana­lysis (PCA). We used the PCA to identify the main axes of environmental variation in study area. The positions (coordinates) of municipalities, with the axes selected by the PCA analysis were used as variables in the hierarchical clusters analysis (see below).

Clustering procedures

We employed a two-stage procedure where a hierarchical algorithm is used to define the number of clusters; these results then served as the starting point for subsequent nonhierarchical clustering (k-means method). Different studies have shown that this procedure increases validity of solutions (Milligan, 1980; Punj & Stewart, 1983).

In the first stage, we analyzed the agglomeration schedule, specifically the change in the agglomeration coefficient to decide how many clusters are needed to represent the data. The agglomeration coefficient is the sum of the within- group variance of the two clusters combined at each successive step. Thus, a small coefficient means that fairly homogenous clusters are being attached to each other, while large coefficients indicate that dissimilar clusters are combining. Consequently, when the increase in the coefficients between two adjacent steps is large (showing that heterogeneous clusters are being combined), indicates the turning point at which it makes sense to group data. The agglomeration coefficient showed an incre­mental change indicating the maximum number of environmental units is 22.

Afterwards, k-means partitioning was employed to cluster the municipalities in environmental units. In this approach, the range of potential k groups is established in advance; thus, we considered the number of clusters obtained from hierarchical cluster analysis (22 environmental units; previous stage) like starting seeds. Hence, the municipalities have been grouped according to their environmental and spatial characteristics, but always bearing in mind that the number of resulting environmental units could not exceed 22; therefore, this was done for all possible models (since 22 environmental units until one environmental unit). The k-means procedure appears to be more robust than any of the hierarchical methods with respect to the presence of outliers, error perturbations of the distance measures, and the choice of a distance metric (Punj & Stewart, 1983). Finally, we performed a filtering of environmental units obtained in the k-means using a criterion of continuity (i.e., that the same environmental unit cannot be divided).

Hunting characterization

After a previous detection and examination of outliers (using the interquartile range method; Swallow & Kianifard, 1996), we characterized environmental units using stepwise logistic regression (Hosmer & Lemeshow, 2000) on a set of hunting variables.

The hunting variables employed were hunting yields (HY) and hunting dominance (HD) considering 48 different game species (Table S2 [suppl]). The game species can vary year to year, therefore, we only consider the species included in the THPs in force in 2005. The average hunting yields were estimated from data from the game estates, reported in the THPs. As digital maps of the game estates are not available we ascribed each game estate to its corresponding municipality (to the one of largest surface when covered more than one municipality) and estimated the average hunting yields of individuals captured on each municipality according to the expression:

where HY is the hunting yield per municipality, expressed as the number of animals captured per 100 ha of game estate.

The hunting dominance (HD) of a municipality is the number of species that constitutes the main hunting use. It is calculated by the expression:

where e is the basis of the natural logarithm, and H' is a hunting diversity, calculated as diversity Shannon index (Shannon, 1949):

where n is the hunting yield for each game species captured in a municipality and N is the sum of the hunting yields from the different species captured in that municipality.

The criteria used to identify species that become part of HD in each municipality were to choose those with high HY in the municipality, according to the value of HD obtained.

Selection of the best hunting regionalization model

It has been assumed that the best hunting regio­nalization corresponds to the group of municipa­lities (environmental units) best characterized from the standpoint of hunting. The evaluation of performance measures for each family first required the derivation of matrices of confusion that identified true positive (a), false positive (b), false negative (c) and true negative (d) cases predicted by each model. From the values in the matrix of confusion, we therefore calculated alternative performance measures including sensitivity, specificity, false positive rate, false negative rate, positive predictive power (PPP), and negative predictive power (NPP) (Fielding & Bell, 1997). Also, Piédrola (2003) provided the following scale to interpret the PPP: >10 excellent, 5-10 good, 2-5 regular and 1-2 deficient. However, conventional statistics of association on these 2×2 tables are inappropriate for assessing model performance (Manel et al., 2001). For this reason, we included the threshold-independent measures: receiver operating characteristic (ROC) plots and area under the curve (AUC). We used this approach (ROC curves and AUC) to make an analysis of similarity of how well the model fits the data (Fielding & Bell, 1997; Manel et al., 2001; Brown & Davis, 2006; Stow et al., 2009). Methods involving ROC curves assess performance from model output at all possible probability thresholds. High performance models are indicated by large AUC. Usually AUC values of 0.5-0.7 are taken to indicate low accuracy, values of 0.7-0.9 indicate useful applications and values of > 0.9 indicate high accuracy (Swets, 1988). We report the sensibility and specificity so that the relative importance of commission and omission errors can be considered to assess the method performance (see Lobo et al., 2008). The use of AUC to measure the model performance was criticized by Lobo et al. (2008), but correctly applied and interpreted remain valid (Elith & Graham, 2009). Furthermore, Santika (2011) showed that AUC correctly ranked the performance of a model and provide useful information about the predictive performance of this.

In this research, we selected the hunting regions with higher AUC values. The municipalities that were not part of any of these regions elected were relocated to neighboring regions following a continuity criterion. Subsequently, the process was repeated.

When two hunting regions shared some municipalities five new hunting characterizations were performed to evaluate the following possibilities: i) region A + region B, ii) region A + shared municipalities, iii) region B - shared municipalities, iv) region A - shared municipalities, and v) region B + shared municipalities.

Main game species and harvest rate by hunting region

One of the main aims of hunting regionalization is to optimize the game resources based on their abundance, derived from hunting yields, and demand. In this regard, it is important to establish regional harvest rates and to know which are the main species hunted in each hunting region. For this, we used the interquartile range (IR), the median and the upper limit of hunting yields of each region.

In addition, as proposed by Vargas et al. (2006a), the game amplitude (GA) was used as an additional variable. The GA represents the number of game species more frequent harvested in hunting states belonging at the same hunting region. Moreover, the usefulness of this variable is to identify on which species is focusing the hunting in the hunting region, regardless of their hunting yields. The value of the GA is expressed by an integer number followed by a decimal, indicating the first digit of the number of preferred species (hunted in over 80% of the hunting states in the hunting region) and the second, the number of representative species (hunted in more than 50% of the hunting states of the region but without being preferred; Vargas et al., 2006a).

ResultsTop

The PCA result of environmental and spatial variables detected eight main axes (explaining 85% of the total variance). With the two-stages clustering procedure and the hunting characterization, we identi­fied 12 hunting regions as the best model (Fig. 2). Eight hunting regions had values typical for high accuracy models (AUC > 0.9), and the remaining four had values typical for useful applications (AUC > 0.8; Table 1). Moreover, ten hunting regions had values PPP >10 (Table 1).

Table 1. Area under the ROC curve (AUC) and values derived from matrix of confusion for assessing the hunting region’s model performance in Castilla-La Mancha.

In general terms, we found six big game regions (1, 2, 4, 8, 10 and 12) and six small game regions (3, 5, 6, 7, 9 and 11; indicated by the species and the sign of the coefficient in Table 2). We also identified, for example, that hunting regions 10 and 12, presents high HY of red deer (Cervus elaphus), and Spanish ibex (Capra pyrenaica), respectively (Table 2). The regions 1 and 8 present high HY of roe deer (Capreolus capreolus), while regions 2 and 4 have high HY of fallow deer (Dama dama). On the other hand, the regions 5 and 11 have high HY of rabbit (Oryctolagus cuniculus), while region 7 has for common pheasant (Phasianus colchicus), and Iberian hare (Lepus granatensis).

Table 2. Hunting characteristics of the 12 hunting regions in Castilla-La Mancha.

Furthermore, the main game species, identified by the value of the GA, comprise two species of big game and eleven of small game (Table 3). The region 1 had the highest game amplitude whereas the region 10 the lowest. Likewise, only three species are preferred in every region, wild rabbit, red-legged partridge (Alectoris rufa), and Iberian hare; contrarily, roe deer is preferred only in region 1 (Table 3).

Table 3. Regional game amplitude (GA) of the most characteristic big and small game species of the hunting regions.

The regional hunting yields for each of the 12 hunting regions are shown in Table 4. The median represents the most common values in the hunting states of each region for a given species, while the maximum limit values determine the range of normal harvest for each region and species, therefore, any hunting state with hunting yields higher of this limit, can be considered an atypical hunting state.

Table 4. Harvest rate (of the main species), expressed as individual hunted/100 ha of land bounded, for each one of the 12 hunting regions of Castilla-La Mancha.

DiscussionTop

Our results show that the establishment of hunting regions represents a new system of administrative game management. Although this is the first model of hunting regions proposed for this Community in Spain, several regional governments have established hunting regions to improve and facilitate management actions; however, each has followed a different methodology to identify them (GLR, 2003; Vargas et al., 2006a). In this paper we used a cluster analysis in two phases to minimize the arbitrariness with which the number of resulting clusters is decided. Thus, a maximum number of regions was established from which the evaluation of different region models began to define the most appropriate, as the regional model should not be too complex but not too simple, so it really is a useful administrative tool that can be adopted by those responsible for hunting.

The THPs were employed in this analysis becau­se they are a well-documented information source. Although, the use of this type of information has its limitations (Blanco-Aguiar et al., 2012; Caro et al., 2014); this indirect sampling method provides homogenous data in a large temporal series, with evident savings in material and human resources (Vargas et al., 2006b). Also, indirect measures of population abundance, such as hunting harvest data used in this study, are useful for many scientific purposes (Cattadori et al., 2003; Vargas et al., 2006b; Farfán et al., 2008; Delibes-Mateos et al., 2009; Ríos-Saldaña et al., 2013). On the other hand, the application of clustering methods to classification problems, in a similar fashion to the method used in this work, is well established in the literature (Rueda et al., 2010; Vasconcelos et al., 2011). One of the novel features here presented is the use of a two-stage procedure to establish the number of clusters. This played a significant role in clarifying the correct number of environmental units, and finally, hunting regions (as mentioned above). Likewise, logistic regression is a widely used tool for modelling species distribution (e.g. Vargas et al., 2007; Caprio et al., 2011; Santangeli et al., 2013).

Otherwise, outliers were detected in the HY by muni­cipalities of all game species of Castilla-La Man­­cha, prior to hunting characterization of the municipalities, something that has not been done so far in hunting regionalization models in other territories, despite the great importance of these outliers. In the case of Castilla-La Mancha it was necessary to exclude some values of HY. As an example of the influence of these values we can remark the case of the rabbit, when including outliers, the average of the HY in Castilla-La Mancha is 42.9 rabbits per 100 ha; instead, when excluding outliers the average of the HY is 15.9 rabbits per 100 ha. This allows another important result, the establishment of extraction rates of each hunting region excluding outliers. In other studies has been calculated the harvest rate of regions as the average of the yields of the municipalities that make up each region, as was done in Andalusia (Vargas et al., 2006a). However, it is considered that this is not the best way to calculate, since outliers have a great influence on the mean, leading to the establishment of rates of inadequate extraction, which could lead to overexploitation of resources, or a use thereof below its potential. In this paper the interquartile range method was applied for to analyse the variability in rates, identifying those that might be considered atypical by its marked difference from the rest. Moreover, an aspect susceptible to discuss is the use of hunting data (contained in the THPs) for establishing extraction rates, because of that other specific factors such as the release of farm-reared game (Caro et al., 2014), the speed of the reproductive rate (Price & Gittleman, 2007) or well, how to care for the offspring, the group size or the presence of infanticide (Caro et al., 2009), can dramatically affect­ the esta­blishment of sustainable harvest rates. Har­vest rates proposed in this paper represents the current hunting demand in each region, but it is not enough to guarantee the sustainability of hunting. Hence, monitoring the distribution and abundance of hunted wildlife is essential to avoid overexploitation of game species (i.e. Caro et al., 2015). Likewise, any administration that use data for development of hunting game management strategies should evaluate the data collection method used (checkpoints, hunting preserve memories, hunters personal surveys, etc.), as the accuracy of these data may depend on the method used (Kilpatrick et al., 2005).

Throughout hunting regionalization it is pretended to establish basic guidelines to optimize the game management adjusting it to the reality of each region. The use of hunting regions facilitates the tasks of assessment, management and monitoring of game by reducing the number of units. In addition, hunting regionalization will optimize, among other things, the utilization of hunting resources in each region according to their wealth and demand. Finally, this represents the development of new game management strategies, important not only to satisfy the new demands of society, but to ensure the stability of the resource and its potential long-term productive. Furthermore, the recovery of some game species such as wild rabbits and red legged-partridges or the effective control of generalist predators, will be some of the advantages of the hunting regionalization, because these actions are necessarily tied to combine the efforts of neighbouring hunting states that have a common problem.