Planning for wolf-livestock coexistence: landscape context predicts livestock depredation risk in agricultural landscapes

Extensive pastoral livestock systems in Central Europe provide multiple ecosystem services and support biodiversity in agricultural landscapes but their viability is challenged by livestock depredation ( LD ) associated with the recovery of wolf populations. Variation in the spatial distribution of LD depends on a suite of factors, most of which are unavailable at the appropriate scales. To assess if LD patterns can be predicted sufficiently with land use data alone at the scale of one federal state in Germany, we employed a machine learning supported resource selection approach. The model used LD monitoring data, and publicly available land use data to describe the landscape configuration at LD and control sites (resolution 4km*4km). We used SHapley Additive exPlanations to assess the importance and effects of landscape configuration; and cross-validation to evaluate the model performance. Our model predicted the spatial distribution of LD events with a mean accuracy of 74%. The most influential land use features included: grassland, farmland and forest. The risk of livestock depredation was high if these three landscape features co-occurred with a specific proportion. A high share of grassland, combined with a moderate proportion of forest and farmland, increased LD risk. We then used the model to predict the LD risk in five regions; the resulting risk maps showed high congruence with observed LD events. While of correlative nature and lacking specific information on wolf and livestock distribution and husbandry practices, our pragmatic


Introduction
Pastoral livestock systems in Europe are valued for providing a suite of ecosystem services and supporting biodiversity (Dean et al., 2021), but their persistence and viability is challenged by the recolonization and expansion of legally protected carnivores . This is particularly the case for wolves (Canis lupus), the most widespread large carnivore in Central Europe (Chapron et al., 2014). Livestock depredation (LD) by wolves includes mostly sheep but also other ruminants, such as goat, cattle and farmed deer (Pimenta et al., 2017) and causes animal stress (Janczarek et al., 2021) as well as serious concerns, economic costs, and negative attitudes among livestock herders and parts of the general public (Arbieu et al., 2019).
How sustainable coexistence of wolves and pastoral livestock system can be achieved is part of an ongoing societal learning and negotiation process (Kiffner et al., 2019). Reducing LD is arguably one of the most important tasks for facilitating a more sustainable coexistence of people, livestock and wolves (Van Eeden et al., 2018). To achieve this overarching goal, identifying the effect of land use characteristics on LD could yield valuable information for livestock herders, politics and administration (Marucco and McIntire, 2010). Identifying such spatial associations could support a more focussed planning and implementation of prevention, support and compensatory measures, and thus ensure a more cost-effective use of public funds.
Predicting the LD risk by large carnivores would ideally be based on mechanistic models, which focus on the causality of input-output relationships. However, their use is constrained by the scarcity of relevant data at the appropriate scale. On the one hand, such models require information on the spatial distribution of wolves (Marucco and McIntire, 2010), which could itself be influenced by the distribution and density of wild prey (Janeiro-Otero et al., 2020).
Wolves typically occupy large home range sizes (200 km² territory size on average in Germany, www.dbbw.de) and have a wide habitat niche (Mancinelli et al., 2018). For example, habitat suitability models for wolves in Europe suggest that wolves mostly prefer forests, meadows and wetlands, and typically avoid anthropogenic landscape features such as roads, settlements and cropland (Jędrzejewski et al., 2008). Even though wolves seem to avoid anthropogenic structures, they clearly use these landscape features occasionally. Several studies reported that wolves use bridges, roads, pipelines, and other linear features as travel corridors (e.g. Plaschke et al., 2021). On the other hand, LD is conditional on the availability of different livestock species and this precondition is ideally incorporated in risk models (Kuiper et al., 2021). Additional factors to be considered in LD risk modelling include the presence and type of practices to prevent depredation, as well as landscape features that mediate hunting success by large carnivores (Gable et al., 2021).
However, the full suite of these data is typically not available at scales useful for planning livestock prevention methods, e.g. for an entire state or landscape of several thousand square kilometres. For example, ruminant livestock management (especially of sheep) often combines keeping animals in fixed paddocks during the winter time, and rotational grazing with movable fencing during the vegetation period. Particularly in the latter case, the distribution of small livestock is dynamic and thus difficult to monitor and predict. Similarly, while wolves are monitored in Germany, movement data of wolves are only available for a few individuals. Currently, the most widely available data sources originate from monitoring programs that focus on LD events (Khorozyan and Heurich, 2022). This monitoring typically includes systematic recording of the location of the lethal interactions between wolves and livestock. Although the data collected through these efforts are not covering all the above described information, they are standardized, regularly updated and have a large spatiotemporal coverage.
In this paper, we make use of these readily available LD monitoring data, link them to open source spatial datasets, and employ state-of-the art machine learning models to predict the LD risk. Our intention is not to improve or reveal a clear mechanistic understanding of this risk (Kuiper et al., 2021), as it will be unclear whether the identified associations between LD and specific landscape features are primarily due to wolf-related or livestock-related variables.
However, our data-driven approach can be used to plough through the available data of inputs and outputs and predict the LD risk without focusing on the underlying causalities. This allows us to explore which land use configuration is currently facilitating LD by wolves across wider areas. Finally, we use our model to develop risk maps, which in turn can inform decisions on appropriate livestock management actions, and thus can contribute to a more sustainable coexistence between wolves and ruminant livestock husbandry in central Europe.

Case study region
As case study area for our research (Figure 1) we chose the federal state of Brandenburg, Germany with a total area of 29°654 km 2 , half of which is used for agriculture (48.6%), followed by cover of forests (34.8%), settlements, roads, water bodies and other vegetation (combined 16.6%) (Statistik, 2022).
<<insert Figure 1 here>> Brandenburg is located in the centre of Europe, and until the Middle Ages, the wolf population was relatively large and stable. With the expansion of deforestation and livestock husbandry, LD by wolves was considered a major problem resulting in organised wolf eradication programmes, which eventually led to the local extinction of wolves by the end of the 18 th century (Ludwig, 2017). During the 19 th and 20 th centuries, reports suggest occasional occurrences of wolves immigrating to Germany from Eastern Europe. However, these wolves were persecuted and typically killed. Since Germany's reunification in 1990 the wolf is under legal protection (EEC, 1992). In 2000, wolves started reproducing in Saxony, and in 2007, the first territorial wolves in Brandenburg were documented Schade, 2010). Since their recolonization, Brandenburg's wolf population grew quickly and expanded along a South-East -North-West gradient, while some parts of Brandenburg are not yet colonized (DBBW, 2020). Concurrent with the expansion of the wolf population, the frequency of LD events increased (LfU, 2021), and state authorities initiated systematic monitoring of LD events, provided funding for livestock protection (e.g. for fencing, livestock guarding dogs and operating costs) and damage compensation. In 2020, annual expenditures in Germany for herd protection amounted to € 9.5 million, while compensation payments for livestock losses amounted to € 0.8 million (DBBW, 2020).
The increase in LD events caused by wolves (Khorozyan and Heurich, 2022) Maintaining and restoring grassland biodiversity via supporting and expanding extensive grazing-based livestock systems is thus a key political priority in Brandenburg. Additional benefits are new farm income sources through meat production including novel value chains, such as direct marketing and improving sustainability through better linking meat production in rural and meat consumption in urban areas, and EU-wide promotion of the "Farm to Fork strategy". From a socioeconomic and circular economy perspective, the low supply of locally produced grazing-fed, and thus more sustainable, livestock is not desirable. For example, only ten percent of the organic meat consumed by residents of Germany's capital Berlin, which is located in the center of Brandenburg, is actually produced in Brandenburg, while the rest originates from other regions (Scholl, 2007). However, livestock owners, particularly shepherds, are questioning the future of their farming system given the extent of livestock losses due to wolves and a perceived lack of sufficient political action to address this issue (NABU, 2022).
As a consequence, the recovery of the wolf population and the desired expansion of grazinglivestock systems are conflicting. In addition to other socio-ecological stressors, impacts of LD may challenge the sustainability of this type of land use. Improved knowledge on where LD events are more likely to occur could help formulate more targeted measures and thus contribute to resolving this conflict and contribute to more sustainable wolf-livestock coexistence.

Modelling methodology
To predict the LD risk by wolves, we used a resource selection approach (i.e. comparing areas with LD and areas where no depredation occurred) and employed machine learning algorithms (Wieland et al., 2021) to identify spatial correlates of the LD risk (for an overall illustration of the methodology, see Figure ). we only considered events were wolves were either confirmed or very likely to have caused the death of the livestock; additionally, we excluded six cases with missing or erroneous coordinates, resulting in 1°094 data points used for modelling.
We obtained land use features from the OpenStreetMap (OSM) project. These data are freely available, have a high resolution (10 m * 10 m), and a high level of detail in terms of represented land use features. We defined the influence area of the landscape for an LD event with a cell size of 16 km² (4km*4km). We chose this cell size to approximate the daily movement of wolves which average 16 km per day (Wrangham et al. (1993): 14.3 km/day; Bryce et al. (2022): 18 km/day). As our approach relates to an area extent, we took the square root of this average daily movement rate to define the cell size. The coordinates of the LD event determine the centre of the influence area. The definition of the influence area is a compromise: it should be large enough to include all important landscape elements, but not too large to include unnecessary elements.
The chosen modelling approach requires presence and pseudo-absence sites for the training.
The presence sites are the LD, the pseudo-absence sites are sites without reported LD events.
Approximately the same number of presence and pseudo-absence sites should be provided for a training session.
The area of the state of Brandenburg (without the largest cities) is 28°922 km², thus containing 1,808 cells (28922 km² divided by 16 km²). The cell size of 4km*4km represents a compromise between important landscape elements and the bias caused by overlaps. We loaded 1°094 presence sites as images (png) from the OSM server. Subsequently, we randomly loaded pseudo-absence sites with a minimum (empirical) site centre distance of 2°500 m from each presence site. We chose this distance-based rule to minimize the likelihood that pseudoabsence sites coincide with LD events. Since we cannot exclude the possibility that come LD events were not reported we refer to these sites as "pseudo-absence". For both, pseudoabsence and presence sites, we extracted the proportion of land use variables. Initially, we included all available land use features. However, during the course of modelling we omitted some features (e.g. building or mining areas) as they were not influential. We based our final models on the proportion of ten features: (1) cropland and fallow crops; (2) grassland, meadow and paddock; (3) forest; (4) scrubland; (5) residential, industrial and commercial areas; (6) gardens, parks, cemeteries, and sports fields; (7) farmyard; (8) lake; (9) river; and (10)  <<insert Figure 2 here >>

Model training, evaluation and validation
To compare land-use characteristics of sites with LD to pseudo-absence sites, we employed XGBoost (https://xgboost.readthedocs.io/en/stable/index.html), a modern and powerful machine learning algorithm suitable for solving classification tasks (Brownlee, 2016). Together with the scikit-learn library, which includes tools for data pre-processing tools and calculation of evaluation metrics, XGBoost can be used to train and validate a model. To illustrate the calculation of these parameters, assume following exemplary confusion matrix of one model run (Table 1): <<insert Similar to other machine learning algorithms, XGBoost models are complex models, which are different to interpret without additional tools. We use here the SHAP software (https://github.com/slundberg/shap), developed on the basis of cooperative game theory (Lundberg and Lee, 2017), which allows for the quantification of the feature importance in the training data.

Results
We used the waterfall visualization, as part of the SHAP library, to analyse the effect of each feature on the rating of the sites. The waterfall visualization is read from bottom to top. The blue arrows point towards no LD risk, while the red arrows point towards LD risk. If the result (top) is > 0 LD risk is assumed to be high; if it is < 0 the LD risk is assumed to be low. Figure 3 illustrates an example for a map with (Figure 3 A) and without LD (Figure 3 B).
<<insert Figure 3 here>> The trained model was applied to predict the LD risk in five 50 km x 50 km areas in Brandenburg (see Figure 1)  <<insert Figure 5 here>> We tested dispersion of the key metrics of the models using the feature set described in Figure   2 for 30 runs (Table 2). Accuracy values ranged from 0.70 to 0.79, with a mean value of 0.74.
Based on the statistics of key metrics as well as of the SHAP analysis, we chose one model for validation in the five sub-regions of Brandenburg (Figure 1). The validation scores of the trained model were: accuracy=75.8, precision=74.2 and recall=77.7.
<<insert Table 2 here>> Results of the fifth step (prediction), in which we applied the trained model to the five subregions in Brandenburg (Figure 1), are shown in Figure 6. The likelihood of LD is shown on the left side in a raster of 1km*1km ( Figure 6) on a scale from 0 (no risk) to 1 (highest risk).
The underlying OSM map of each area is shown in the central panel. On the right panel, this map is overlaid with a heatmap (bicubic interpolation from 1km*1km to 30m*30m). Light (dark) areas indicate a high (low) LD risk. LD events are marked as blue dots (middle and right panel). Simulation and observed data showed a high degree of correspondence, thus confirming the potential usefulness of our approach for the prediction of the LD risk.

Interpreting land use configuration -livestock depredation associations
Given the suite of factors that can potentially influence LD by wolves, the central question of the presented work was whether there are patterns in these events that can be explained exclusively by the surrounding land use features. Our study showed that the LD risk by wolves was strongly associated with the co-occurrence of distinct land use proportions of grassland, forest and farmland.
The identified strong influence of forest is in line with expectations, as wolves prefer forests as key habitat (e.g. Jędrzejewski et al., 2008); similarly, the presence of grasslands and farmsteads in a landscape patch implies a higher availability of livestock as potential prey.
Consequently, LD appears to primarily occur at the landscape interface, where suitable wolf habitat borders or overlaps with areas where livestock is kept. These results suggest that the occurrence of LD is driven by land uses that influence both resource selection decisions by wolves (i.e. positive association with forests) as well as the distribution of livestock (i.e. grassland and farmland). This aligns well with research on other human-wildlife conflict scenarios, suggesting that damage occurrence is fundamentally influenced by principles of resource selection decisions by both wildlife and humans (Bautista et al., 2021).
However, a more nuanced, functional interpretation of the identified co-occurrence is not possible given the interrelatedness of livestock presence and wolf habitat use, which cannot be disentangled with the currently used data. Considering that predictive ability does not imply causality, the importance of grasslands and farmland in our models does necessarily suggest that wolves have a tendency to hunt in the proximity of these land uses. In fact, farmland or human infrastructure such as buildings are generally avoided by wolves, as more detailed studies on wolf movement have clearly demonstrated (e.g. Mancinelli et al., 2019).
Simulation results for the five sub-areas confirmed the strong impact of the co-occurrence of these three key land use variables. In the region "Uckermark" (region 4), LD was less likely to occur. This region is dominated by large agricultural fields and features little forest cover and little grassland. Therefore, from its structural characteristics, it appears less suitable for animal husbandry. From the wolf's perspective, the landscape appears also less suitable due to the relative scarcity of forest cover (e.g. Jędrzejewski et al., 2008). In tandem, this landscape configuration results in a comparatively low LD risk. In contrast, region 3 ("Havelland") is a diverse, structural rich landscape with a large share of grassland and forests. Consequently, both observed and simulated data indicate high depredation risk ( Figure 6).
In sum, the modelling results suggest that LD is more likely to occur in diverse landscapes with a minimum of three co-occurring landscape elements (forest, grassland, and farmsteads), transboundary areas and smaller fields, while in cropland-dominated landscapes with little grassland and forest cover, depredation events are unlikely to occur.

Methodological aspects
Our approach used only few, yet readily available input data and achieved a moderate accuracy for predicting LD by wolves. However, the obtained accuracy also means that every fourth prediction is incorrect. Thus, landscape features are important, but other factors which could not be considered here due to the mentioned data availability constraints, are perhaps required for predicting depredation by wolves. The key advantages of the presented methodology are that it makes use of available data, and it is exclusively based on machine learning without requiring the involvement of experts or other potentially costly work steps.
It can be efficiently applied and transferred to other regions if training data (i.e. systematically collected LD data) are available.
Moreover, the results can be used in many ways. Not an obvious, but in our eyes important use is the analysis of waterfall diagrams. In Figure  A key prerequisite of a good model is a solid database of livestock predation events.
Monitoring of livestock predation events has considerably improved over the years and represents a robust database for the analyses presented in this paper. The available monitoring protocol, however, is currently wolf-centred, e.g. it does not, or not systematically, include additional variables that might affect LD, such as herd size and composition or the type and quality, of livestock protection and husbandry measures applied by herders. Data on finescaled wolf and livestock distribution are also not available at the scale of a federal state.
Despite the lack of such data, our model provided reasonable prediction accuracy. The achieved accuracy of 74 % is likely sufficient for the simulation shown, yet not ideal. Possibly, achieving a higher accuracy is affected by the choice of pseudo-absence data, which were generated by a random selection process. Especially since the wolf population in Brandenburg was still expanding during the time period of the study (DBBW, 2020) some of the pseudounaffected sites would in reality be risky. In addition to accuracy, precision and recall should also be considered. Precision also includes the false positives, recall the false negatives.
Looking at both precision (72 %) and recall (78%) indicates that there are more false positives than false negatives (Beattie et al., 2020). Thus, more risky areas are falsely identified as low risk areas than low risk areas are falsely predicted as high risk. Again, this can be explained by the expanding wolf population, and the fact that not all areas of the state had been occupied by wolves during the entire study period.
A possible source of bias is that land use may change over time, yet we used maps from a single time step, which may not adequately describe the land use at each time step. In particular, this applies to earlier recorded depredation events. Potentially, this could have contributed to a lower accuracy, yet the influence of this factor is likely to be small.
The inclusion of additional factors in our model, such as more detailed data regarding livestock densities (Kuiper et al., 2021), e.g. from the Integrated Administration and Control System (Uthes et al., 2020), data on wolf distribution or proxy data for key wolf habitat such as protected areas or military training areas , barriers limiting wolf expansion such as motorways without green bridges, or the type of forest and wildlife management, could further improve the accuracy of the model. For example, one could assume that the intensity of hunting activities differ between privately and state-owned forests, resulting in heterogeneous wild ungulate densities. In turn, changes in ungulate densities, the main prey for wolves in Germany, could influence LD patterns.
However, these additional data sources are either not available at all or more difficult to obtain compared to the OSM maps used in our approach. Including them in our model would cause higher costs and divert from our initial intention to base our modelling on publicly available data.

Use for decision support
Our model produced risk maps illustrating the LD risk based on available land use parameters, a useful tool for rapid classification of areas with high and low depredation risk. This information can assist zonation-planning projects and help authorities decide on the degree of preventative and compensatory measures required in different areas. For example, livestock herders with permanent grazing areas located in high risk landscapes (dark areas in Figure 6) could be encouraged to adopt more effective prevention measures (e.g. electric fence and livestock guarding dogs) while in low-risk areas (white areas in Figure 6) less costintensive prevention measures (e.g. electric fence only) could suffice. Similarly, livestock herders could be encouraged to move livestock, especially mobile small stock such as sheep, to suitable and available pastures in low-risk areas. Such evidence-based, spatio-temporal avoidance by livestock herders has been effective in reducing carnivore depredation on livestock in other parts of the world (Melzheimer et al., 2020). Coexisting with large carnivores requires substantial learning and adaptations among key stakeholders  and an integrated perspective rooted in a socio-ecological-systems understanding and aimed at win-win solutions (Gordon, 2018). To overcome the implementation gap between knowledge of effective prevention methods and relevant action, dedicated farm coaching workshops in areas with high LD risk may be an effective means to increase the adoption of effective livestock prevention methods.

Ethics approval
Not applicable.

Data and model availability statement
None of the data or models were deposited in an official repository. The data are not publicly available as they include sensitive information.

Declaration of interest
None.

Highlights
 We use a machine-learning model to predict livestock depredation risk by wolves.
 The model requires past depredation events and publicly available land use data.
 The trained model predicts livestock depredation with 74% accuracy.
 Thus, landscape context appears to influence livestock depredation.
 Results help practitioners on deciding where to carry out protection measures. Figure 2: Graphical overview of the methodology to model livestock depredation risk due to wolves in Brandenburg. The XGBoost algorithm compares the proportion of ten main land use features between livestock depredation areas (three example sites out of the 1°094 cases in the top part) and randomly selected areas without livestock-depredation cases (three example sites in the bottom part). Land use features are colour-coded as follows: yellow -farmland and fallow crops, dark green -forest, light green -meadows and grassland, grey -residential, industrial and commercial areas, blue -lake, lightblue -river, brown -farmyard.