Recognition of Animals’ Affiliation with Habitat Patches Using Machine Learning Techniques

1. Introduction

Wisents/European bison (Bison bonasus L.), after their extirpation from the wild over 100 years ago, were maintained in captivity, and, beginning in 1952, have been gradually reintroduced to certain designated sites that they once naturally occupied. The present world population of the species originates from only 12 founders who survived in captivity, which makes it one of the most inbred mammals in the world [1].

There has been only a few successful reintroductions of this species, all in Central and Eastern Europe, and the newly-established herds have practically remain in isolation from each other. Therefore, the species, due to its very low genetic diversity, is potentially highly threatened by inbreeding depression. This may affect the fitness of individuals, their ability to adapt to changes in environment, and, finally, the chances for long-term survival of the species [2,3,4].

Nevertheless, because of socio-economic problems resulting from damages in agriculture and cultivated forest stands caused by some large wisent herds, it is suggested for the future to maintain rather small populations of wisent within any given area [5,6].

There are two choices regarding the management of small, isolated populations. One is an attempt to establish gene flow, which may lead to outbreeding depression, and the other is to maintain subpopulation isolation , risking inbreeding depression and potential genetic drift [7,8].

Natural exchange of individuals among particular herds is not an easy process be-cause of the general lack of suitable migration corridors for large wildlife species in con-temporary Europe, as well as the social organization of wisent herds, where females with young remain mostly within stationary ranges and only single males are are migratory[9,10].

The reintroduction of wisents to the wild was performed in sites selected according to their known habitat preferences. In some case, however, animals did not remain around the site of the release but moved to other habitat patches or split into separate groups occupying different home range territories [11,12].

Monitoring of spontaneous movements of animals among neighboring herds is very difficult and, in most cases, impossible. Hence, in order to calculate the degree of potential gene exchange, it is important to be able to estimate the intensity of such migrations. That can be achieved by identifying the affinity of particular individuals to their maternal groups.

An earlier attempt to distinguish wisents occupying different home ranges (i.e. be-longing to different herds) during both the vegetative and winter seasons was done through kernel density estimation and spatial analysis. The results obtained using the different methods show that there is a possibility of generalization [13,14].

The advent of modern data collection techniques and geospatial technology including real-time location tools and high-resolution satellite imagery has led to the generation of vast amounts of data. In this regard, machine learning algorithms facilitate their deeper analysis and bring about better understanding of the actual situation. In contrast, classical statistical methods are constrained by numerous limitations and assumptions regarding the type of data and the relationships between variables.

Therefore, in this paper we tested the machine learning approach for spatial data classification, seeing this as a new method for recognizing the affiliation of particular individuals to their maternal home ranges. This approach, we think, could be, useful for future management of free-ranging populations. In this regard, the eXtreme Gradient Boosting (XGBoost) algorithm was applied for classification tasks due to its high efficiency and superior performance in processing large datasets. The present study was primarily concerned with the following: (1) introducing a novel ML approach for subpopulation identification based on animal presence data, (2) developing the XGBoost method for classification by using CLC data, and (3) generating comparative analysis of the classification efficiency of the vegetative and winter seasons.

2. Materials and Methods

2.1. Study Area

The free-ranging wisent herds currently inhabiting the Bieszczady Mountains originate from two reintroductions: the first in 1963 to the eastern part of the range, and the second in 1980 to its western part. All introduced individuals originated from captivity and belonged to the Lowland – Caucasian line, bearing genes of the only representative of Caucasian wisents (the Caucasus #100) that survived in the breeding center. The eastern subpopulation (known as the Tworylne herd) dwells within the forest districts of Lutowiska and Stuposiany and partially within the Bieszczadzki National Park [1,15,16]. The home range of the western subpopulation (known as the Baligród herd) extends over the forest districts of Baligród, Komańcza, Cisna, and Lesko (Figure 1).

During the period of time covered by this study (2002 and 2021), wisent numbers in Bieszczady grew from about 150 to about 700 individuals, which is over 40% more than the target number originally planned for this population [17,18]. Such significant increase of population numbers could be associated with intensified use of the most preferred habitats and overexploitation of natural food base as well as depredation of crops [19, 20}. In this study we analysed data on the occurrence of wisents in habitat patches within the home ranges inhabited in Bieszczady by two subpopulations of this species.

The Bieszczady Mountains, being a part of the Western Carpathians, are situated in southeastern Poland. Their maximal elevation is at 1346 m above sea level, and they are mostly forested (almost 85% of their total land cover). The tree stands are dominated by beech-fir (Fagus sylvatica-Abies alba) associations, with alder-willow (Alnus glutinosa-Salix spp.) woods along watercourses. A considerable part of the presently forested area is also overgrown with Scotch pine (Pinus sylvestris), replanted or originating from secondary succession on former cultivated fields [16]. The majority of home ranges of both subpopulations are forested: for the western subpopulation, the forested area is 593 km² and the open area is 237 km², and for the eastern subpopulation, the forested area is 678 km², and the open area is 183 km².

2.2. Data Collection

We used for this analysis, 31480 records of wisent presence of two subpopulations: western (the Baligród herd) and eastern (the Tworylne herd) obtained between 2002 and 2021) from 6 radiocollared individuals (3 in every herd). The raw data we employed included telemetric fixes collected by the staff of the Carpathian Wildlife Research Station of the Polish Academy of Sciences (under the framework of the routine monitoring of free-ranging wisent in the Bieszczady Mountains). This survey was undertaken in the vegetative season (between April and October) and the winter season (between November and March). Table 1 shows the wisent presence records of the two subpopulations as broken down to season and in total.

Spatial distribution of free ranging wisent presence during the vegetative and winter seasons is shown in Figure 2 (a), (b).

The available set of telemetric data for every record consisted of the date and hour, elevation above sea level (a.s.l. [m]), and the location of the animal in WGS84 projection. Geospatial data were added to the numerical map of the Bieszczady Mountains via ArcGIS software. Next, thematic layers were created for further analysis (forested area, open land, composition of tree stands, infrastructure, elevation above sea level, records of wisent presence), which allowed for identifying their mutual relationships. Habitat parameters within both wisent ranges were classified according to CLC regulations. This land cover and land use inventory consists of 44 thematic classes, including forests, agricultural land, and human-related infrastructure [18]. Analyses were performed for three levels of CLC (Table 2, Table 3, Table 4, Table 5, Table 6, Table 7 and Table 8).

At the time of data collection, the average elevation above sea level was different within the home ranges of both subpopulations (Table 2).

The western subpopulation was found to inhabit generally higher elevations, by about 90 m on average. During the vegetative season, this difference was much lower – about 61 m, but in winter – it was over 118 m. Seasonal differences were much higher (about 83 m) within the area occupied by the eastern subpopulation, while between seasonal home ranges of the western subpopulation, elevation differed by only about 26 m.

On average, land cover and land use within the home ranges of both subpopulations were quite similar, with forests covering above 90% of their area, with much lower percentages of cultivated land being occupied (Table 3 and Table 4). The pertinent percentages are highlighted in bold text.

During the vegetative season, the home range of the eastern subpopulation contained more open spaces (by little over 4%), while that of the western population demonstrated more forested area occupation (by about 5%). During the winter season, the proportion of forests and open land was almost identical.

Types of land use and land cover within the home ranges of the Baligród and Tworylne herds according to CLC level 2 classification are presented in Table 5 and 6.

A comparison of the home ranges of both subpopulations according to CLC level 2 classification reveals that the percentage of forest inhabitation differed by about 7% during the vegetative season for the western subpopulation, in contrast to the eastern population. In the home range of the eastern subpopulation, however, the percentage of pastures and arable land habitation was greater during the vegetative season by about 2%, in contrast to the western population (Table 5). During the winter season, the differences between home ranges of both subpopulations were very slight, not exceeding 1%.

The analysis undertaken for CLC level 3 (Table 7 and Table 8) with regard to types of forest stands found in each subpopulation home range shows significant differences. Within the home range of the western subpopulation, coniferous stands were much more frequently evident (over 24% difference), while mixed forests were significantly less notable (by almost 16%). Moreover, the share of broad-leaved stands, was lower by about 7% within the home range of the western subpopulation (Table 7).

The aforementioned data demonstrates that the only significant differences between the home ranges of both wisent subpopulations lies in the proportions of various types of tree stands, mainly in the vegetative season (by about 35% higher for coniferous forest for the western subpopulation and 23% for mixed forest for the eastern subpopulation). Proportions of forest cover and open lands are, however, quite similar.

2.3. Extreme Gradient Boosting Algorithm

The XGBoost algorithm was created by Tiangi Chen and Carlos Gusterin in 2016 [21]. XGBoost combines a number of weak classifiers into a single strong classifier using the gradient tree boosting technique [22,23,24].

Let $D=\left\{\left(x_i,y_i\right)\right\}$where $i=\mathrm{1,2},\dots ,n,x_i\in R^m$ and $y_i\in R$ be a given data set with $n$ examples and $m$ features. A tree ensemble model using $t$ base models to predict the result ${\widehat{y}}_i^{\left(t\right)}$ can then be written in the form [21]:

$${\widehat{y}}_i^{\left(t\right)}={\sum }_{j=1}^t{f}_j\left({\mathrm{x}}_i\right)={\widehat{y}}_i^{(t-1)}+f_t\left({\mathrm{x}}_i\right).$$

(1)

where ${\widehat{y}}_i^{\left(t\right)}$ is the final tree model, ${\widehat{y}}_i^{(t-1)}$ is the previous tree model and $f_t\left({\mathrm{x}}_i\right)$ is the newly generated tree model.

The process of learning the model estimates the model parameters by minimizing the objective function [21]:

$${\mathcal{L}}^{\left(t\right)}={\sum }_{i=1}^{n}l\left(y_i,{\widehat{y}}_i^{\left(t\right)}\right)+\mathsf{\Omega }\left(f_t\right),$$

(2)

where ${\widehat{y}}_i^{\left(t\right)}$ is the prediction of the target $y_i$, $l$ is the training loss and $\mathsf{\Omega }$ represents the regularization term.

By substituting the formula (1) into (2) and simplifying the formula (2) the objective function can be transformed into the following form:

$${\tilde{\mathcal{L}}}^{\left(t\right)}=\sum _{i=1}^n\left(g_i{f}_t\left({\mathrm{x}}_i\right)+{\displaystyle \frac{1}{2}}{h}_i{f}_t^2\left({\mathrm{x}}_i\right)\right)+\mathsf{\Omega }\left(f_t\right).,$$

(3)

where functions $g_i$ and $h_i$ are the first and second order gradient statistics of $l$. Moreover, the regularization term $\mathsf{\Omega }$, can be written in the form:

$$\mathsf{\Omega }\left(f_t\right)=\gamma T+{\displaystyle \frac{1}{2}}\mathsf{\lambda }\sum _{j=1}^T{w}_j^2,$$

(4)

where parameters $\mathsf{\gamma }$ and $\mathsf{\lambda }$ are penalty coefficients, and the number $T$ is the total number of tree leaves.

The solution of (3) gives the optimal weight of leaf $j$ [21]:

$$w_i^{\mathrm{*}}=-{\displaystyle \frac{{\sum }_{i\in I_j}{g}_i}{{\sum }_{i\in I_j}{h}_i+\mathsf{\lambda }}},$$

(5)

where $I_j=\left\{i\right|q\left(x_i\right)=j\}$ denotes the instance set of leaf $j$.

2.4. Logistic Regression

The logistic regression model for the dichotomous variable $Y$ with values 0 and 1 is defined by the formula:

$$P\left(Y=1|x_1,x_2,\dots ,x_n\right)={\displaystyle \frac{\exp \left(a_0+\sum _{i=1}^n{a}_i{x}_i\right)}{1+\exp \left(a_0+\sum _{i=1}^n{a}_i{x}_i\right)}},$$

(6)

where $a_i$ for $i=\mathrm{0,1},\dots ,n$ are the regression coefficients and $x_1,x_2,\dots ,x_n$ are the independent variables. To predict a class for $Y$, the probability is compared to a classification threshold, which is usually equal to 0.5.

The regression coefficients are estimated using maximum likelihood estimation and gradient descent [25,26].

2.5. Evaluation Metrics

In binary classification, only two classes are considered: one of which is usually referred to as 'positive' and the other as 'negative'. The number of correct classifications, represented by true positives (TP) and true negatives (TN), and the number of misclassifications, represented by false positives (FP) and false negatives (FN), are illustrated in Figure 3 [27]. The "Actual Class" refers to the true, correct label of a data point in a classification problem. It is the "ground truth" or the known correct classification.

The following indicators are used:

• Accuracy – the percentage of correct classifications:

$$Accuracy={\displaystyle \frac{TP+TN}{TP+TN+FP+FN}}.$$

(7)

2.

Precision – the number of cases, expressed as a percentage, in which the classifier yielded the correct result:

$$Precision={\displaystyle \frac{TP}{TP+FP}}.$$

(8)

3.

Recall (sensitivity, true positive rate) – the ability of the model to capture positive cases:

$$Recall={\displaystyle \frac{TP}{TP+FN}}.$$

(9)

4.

Specificity (true negative rate) – the ability of the model to capture negative cases:

$$Specificity={\displaystyle \frac{TN}{TN+FP}}.$$

(10)

5.

F1-score – the harmonic mean of precision and recall:

$$F1score={\displaystyle \frac{2·Precision·Recall}{Precision+Recall}}$$

(11)

6.

Receiver Operating Characteristic (ROC) curve – a visual representation of the relationship between the efficiency of classifying positive cases (sensitivity) and the inefficiency of classifying negative cases (1 – specificity) using all possible classification thresholds. The area under the perfect ROC curve (AUC) is equal to 1 [27].

2.6. Model Development

In our research, the XGBoost algorithm from the xgboost Python package [28] was employed to identify the affiliation of an individual with its proper subpopulation. Afterwards, a grid search method was applied to find the best XGBoost classification model. The entire machine learning modelling can be described thusly:

The geospatial data was processed via ArcGIS packages to create tabular data.

The tabular data was preprocessed in order to choose the proper variables and eliminate missing or incorrect values.

The final data set was divided into training and test subsets by applying the 80:20 rule.

A randomized search model was applied so as to reveal the most suitable hyperparameters. Here, 40 sets of random hyperparameters were selected from the grid search by applying cross-validation splitting strategies [28,29]. The optimal values of hyperparameters were then found based on the AUC ROC metric, through partitioning of the training data into 10 subsets. Finally, the best set of parameters found in the grid was returned so as to uncover the best estimator.

The optimized XGBoost model was employed to make predictions for training and test sets. Next, the evaluation metrics were calculated and compared.

3. Results

In the presented work, we applied the powerful machine learning algorithm XGBoost. The XGBoost model was compared to the logistic regression model, as the latter is a baseline binary classification algorithm.

3.1. Model Settings

The data used in the classification procedure included 17,737 records of wisent presence during the vegetative season and 13,743 records during the winter season.

The data contains the class variable describing the wisent subpopulation (Baligród herd – 0, Tworylne herd – 1), and four feature variables: the numeric variable elevation a.s.l. [m] and three categorical variables including thematic classes of the CLC level 1, CLC level 2 and CLC level 3 data. The conducted analysis was done for each season separately. The XGBoost model was optimized based on grid search parameter tuning with 10-fold cross-validation. The optimal hyperparameters values of the XGBoost model are presented in Table 9.

Figure 4 (a) and (b) illustrates the training and test curves of the developed XGBoost model for both vegetative and winter seasons.

The stable shape and descending trend of these curves clearly evident the correctness of the developed XGBoost behaviour. Accordingly, the training loss was reduced to less than 0.1. The test loss is, however, slightly greater than the training loss. Still, XGBoost avoided overfitting during the training process.

3.2. Model Performance

The evaluation metrics used to describe the performance of the optimized XGBoost model for the vegetative and winter seasons are presented in Table 10.

We can observed that the XGBoost model achieved very high performance metrics for both seasons, although they were higher for the winter season. Furthermore, all metrics achieved very high results for both the training and test sets. Here, the accuracy values were around 91% for the vegetative season and 96% for the winter season. It can also be noted that the F1-score values were around 85% and the ROC-AUC values were around 90%. Such results correspond to very high model classification efficiency. This is also confirmed by high values of precision, recall and specificity, the highest of which are specificity values. These values reveal that the model very successfully classified individuals belonging to the Tworylne herd. The confusion matrices of the XGBoost classifier for the vegetative and winter seasons based on the independent test sets are shown in Table 11 and 12. The best results are in bold font.

In addition, the performance of the XGBoost model was compared with the performance of a baseline model for binary classification tasks (logistic regression) using the same training and test datasets. For comparison, the evaluation metrics for logistic regression are shown in Table 13.

The ROC-AUC values for the vegetative season based on the training and test sets were very low, about 55%. In the case of the winter season, these values were slightly better (about 65%). The results showed that the quality of the classification was very low, especially during the vegetative season. We can also observe that XGBoost achieved significantly higher values for all evaluation metrics except specificity, demonstrating the superiority of this method.

4. Discussion

In contemporary Europe, wildlife species mostly occupy relatively small, isolated patches of suitable habitats. Moreover, ecological connectivity among such sites often does not exist due to considerable fragmentation of the environment and numerous man-made barriers such as industrial or transportation infrastructure, settlements, etc. The size of the remaining fragments of natural habitats is often too small to support viable populations. Therefore, subpopulations of wildlife species are frequently threatened by inbred and loss of genetic diversity, which finally may lead to their extinction [31,32,33,34].

A solution to this problem is the restoration of ecological corridors, and doing so will allow the spontaneous exchange of individuals through migrations. This will enable their occupying suitable habitat patches outside of maternal home ranges [35]. Unfortunately, the assessment of the intensity of such natural processes is very difficult and often impossible unless there is a sufficient number of animals being radio-collared or permanently marked, thus enabling identification of their origin. This, however, requires complicated and expensive procedures connected with capturing, immobilization, and further monitoring of a sizable group of individuals [36].

Therefore, other methods that make it possible to distinguish individuals of different derivations provide valuable assistance in the conservation and management of wildlife species. The application of the XGBoost algorithm, in enabling the classification of particular animals as native to a given home range of a herd indicating its actual habitat preferences, provides such an opportunity.

Earlier, such attempts were based upon the use of kernel density estimation and spatial analysis. Unfortunately, the obtained products did not prove the universality of those approaches. Indeed, the results of applying such models to the available data pertaining to different herds varied considerably in their quality [13,14]. Although ranges inhabited by both herds had similar percentage of forested area, their stand compositions were quite different. In both cases, the most intensively used parts of their home ranges were situated more than 10 km from the respective sites of release. Such situation reflects the ability of wisents to select optimal habitat patches [17}.

Regarding the limitations of our study, we investigated only two herds of wisents. This is because datasets on the spatial distribution and habitat use of these two herds are quite extensive. In contrast, the monitoring of neighbouring herds in Slovakia and Ukraine has been much less intensive and has been carried out over shorter periods of time.

Moreover, differences in habitat selection could result not only from local availability of habitat patches, but to a certain extent, from variation in genetic or physiological traits, since the animals used for introductions were obtained from various breeding centres of diverse conditions.

In this study, we compared the classification results to one traditional model – logistic regression. However, the main aim of our study was to develop the XGBoost algorithm for the classification task. Further research will focus on using a larger dataset and on applying different machine learning models.

Nevertheless, the obtained results confirm the usefulness of the applied tool for analyses of the spatial structure of wildlife populations and for the facilitation of selection of optimal habitat patches for introductions of wildlife species.

5. Conclusions

In our study, the powerful machine learning XGBoost model was developed and applied for wisent herd classification as based on CLC data. The proposed optimized XGBoost model was shown to be an effective method for recognizing individual affiliation to habitats within a home range of a given subpopulation. Our experiments, based on grid search parameter tuning revealed that the examined model provided very high performance metrics for both the vegetative season (accuracy of 91.63%, F1-score of 85.81% and ROC-AUC of 89.05% for the test set) and winter season (accuracy of 96.01%, F1-score of 84.64% and ROC-AUC of 89.34% for the test set). These demonstrate that if there is available a considerable set of telemetric data, the ML approach may provide very accurate identification of habitat patches preferred by animals within their home ranges, so therefore it can be truly useful in the decision-making process for the management of wisents, as well as for other wildlife populations.