Changing patterns of conflict between humans, carnivores and crop-raiding prey as large carnivores recolonize human-dominated landscapes

wolf population in the Italian Apennines to study how this affected both types of human-wildlife conflic. We used large datasets of wolf occurrences ( n = 351), livestock depredation events ( n = 165), and crop damage events by wild boar ( n = 3442) to independently model the determinants of livestock depredation and crop damage distribution in relation to wolf habitat suitability over a ten-year period of increasing wolf numbers. These analyses yielded two major insights. First, livestock depre-dations were mainly related to insufficient prevention measures (e.g. lacking fencing) rather than landscape context, providing a clear pathway to conflict mitigation. Second, crop damage decreased in areas of higher wolf habitat suitability and became more likely in areas of lower wolf habitat suitability, closer to settlements. This suggests increasing predation pressure forces wild boars to avoid the most suitable wolf habitat, leading to a redistribution of crop damage in the landscape. More generally, our study highlights complex human-wildlife interactions as large carnivores recover in human-dominated landscapes, suggesting that multiple, co-occurring conflicts need to be assessed jointly and adaptively in order to foster coexistence between humans and wildlife.


Introduction
Large carnivores play a key role in ecosystem dynamics, but are disappearing across the world (Ripple et al., 2014). These species require extensive areas of habitat and are particularly vulnerable to habitat loss and fragmentation (Romero-Muñoz et al., 2019;Wolf and Ripple, 2017). Furthermore, large carnivores are often considered to be a threat, real or perceived, to humans and livestock. This often leads to the legal or illegal killing of carnivores to prevent, or in response to, conflicts (Morehouse et al., 2018). As a result, conflict is currently the main reason why many large carnivore populations are on the brink of collapse (Loveridge et al., 2017). Understanding why and where humancarnivore conflict occurs is therefore critically important for fostering coexistence between wildlife and humans (Miller, 2015).
In some regions, such as Europe and North America, large carnivores are currently re-expanding their ranges after centuries of persecution and local extirpation (Chapron et al., 2014;Smith et al., 2016). Several factors contribute to these trends, including revised hunting regulations, better protection (Boitani and Linnell, 2015), and structural changes in agriculture that lead to farmland abandonment and outmigration from rural areas (Kuemmerle et al., 2016;Navarro and Pereira, 2015). The recovery of large carnivore populations is encouraging, given that most of these species are vulnerable or endangered (Chapron et al., 2014). Additionally, the return of large carnivores can restore lost ecological functions, including predatory pressure on herbivores (Berti et al., 2019;Ripple et al., 2014). However, a major concern surrounding the comeback of large carnivores is the potential escalation of conflicts with humans (Broekhuis et al., 2017), which may lead to declining public acceptance of carnivores and undermine conservation success in the long run (Behr et al., 2017;Kuijper et al., 2019).
A key concern in this context is that traditional husbandry methods to prevent conflict with large carnivores have declined or were abandoned as landscapes became predator-free (Ciucci and Boitani, 1998). As a result, livestock today often graze freely and unsupervised, making them vulnerable to attacks (Hanley et al., 2018). A wide range of factors have been found to influence livestock depredation risk (Miller, 2015). For instance, preventive measures such as proper fencing and surveillance of livestock strongly lower the risk associated with carnivores presence (Dondina et al., 2015). Conversely, landscape patterns, and particularly a high interface between pastures and carnivore habitat, can increase the risk for livestock (Boitani and Linnell, 2015). Protected areas provide refuge to large carnivores in human-dominated landscapes (Grilo et al., 2019), and can potentially increase livestock vulnerability in their surroundings. On the other hand, protected areas typically have a higher natural prey base and can thus also serve to mitigate conflict (Wolf and Ripple, 2018). The relative importance of these spatial factors remains poorly understood, particularly for rebounding carnivore populations. Closing this knowledge gap is important for proactively implementing conflict mitigation strategies and for the spatial prioritization of interventions Wolf and Ripple, 2017).
Conflicts associated with rebounding large carnivore populations, however, potentially extend beyond direct livestock depredation (Muhly et al., 2011). A central, yet frequently overlooked aspect in humanwildlife conflict mitigation is that recovering carnivores might influence the behavior of their wild prey, which can alter local densities and activity patterns (Cunningham et al., 2019;Laundré et al., 2014). Many prey species are themselves associated with conflict with humans in the form of crop damage, whose costs typically outweigh economic losses due to livestock depredation. For example, economic damages caused by wild ungulates in Europe reach €80 million/year (Valente et al., 2020).
In Germany alone, crop damage by wild boar (Sus scrofa) is estimated at more than €17 million/year, compared to about €1 million/year of damages caused by wolves (Canis lupus; Welt.de, T-online.de, accessed June 2019). As a consequence, even if returning large carnivores have a relatively small effect on the behavior of their wild prey and associated prevalence of conflict, the economic impacts of such effects could be very substantial in the long run.
Rebounding gray wolf populations in Europe provide a unique opportunity to better understand dynamics in human-wildlife conflict. Wolves have been heavily persecuted in Europe for centuries, leading to their widespread extirpation in the 19th and early 20th centuries (Meriggi et al., 2011). Wolves persisted in Eastern Europe and the Balkans, as well as in several small and isolated populations in Italy, Spain, and Finland. From there, wolves have recolonized much of their historical range since the mid-20th century , which has led to increasing conflict with humans. However, to what extent a rebounding wolf population might impact agriculture through changes in prey density or behavior remains unassessed. Here, we make use of the natural experiment of a rebounding wolf population in the northern Apennines, Italy, to investigate the changing dynamics of livestock depredation and of crop damage by wild boar, the primary prey of wolves in this region (Torretta et al., 2017). We used a dataset of wolf presences, livestock depredation events, and crop damage claims to address two main research questions: 1. Which factors, particularly those related to husbandry systems and landscape patterns, determine livestock depredation risk? 2. Did the increasing wolf population in the study area affect the distribution of crop damage over the period 2009-2018?

Study area
We focused on a 1950 km 2 study area in the Italian northern Apennines, corresponding to 34 municipalities in the lowlands and mountains south of the city of Bologna (Fig. 1). Elevation ranges from 50 to 2000 m above sea level. Climate is temperate-continental, with annual average temperature of 11.5 • C and rainfall of 900 mm (Cervi et al., 2018). Land cover consists of mixed forests and other semi natural areas (56%), as well as croplands (33%) and urban areas (9%) . Protected areas cover 11% of the study region. These areas are Regional Parks where historical sites and cultural land-use practices, including livestock farming, are protected alongside local flora and fauna, with strict regulations on hunting and timber harvesting. Livestock farming is widespread in the area, primarily of sheep, goats, and cattle, and to a lesser extent horses and donkeys. During the grazing season (April to September), livestock roam freely in fenced pastures (Dondina et al., 2015). Transhumance is very rare. Wolves were historically widespread in the region but were extirpated after World War II. Starting in the late 1990s, wolves began returning to the area and their population has been increasing since (Galaverni et al., 2016;Loy et al., 2019). The preferred wolf prey species in the area is wild boar, but also roe deer (Capreolus capreolus) and red deer (Cervus elaphus) (Meriggi et al., 2020).

Datasets used
We used data on livestock depredation events caused by wolves during the period 2011-2016 from the official, verified records of regional authorities (http://regione.emilia-romagna.it; accessed May 2018). These reports contain farm names, date of depredation events, livestock species, and number of lost animals. Using these data, we contacted, visited and interviewed 66 farmers who claimed depredation events, as well as 19 neighboring farmers who did not suffer livestock depredation during 2011-2016. Neighboring farmers were selected within a distance of 3 km from pastures with depredation events (following Mech et al., 2000). We considered only pastures with ≥8 goats/sheep heads for at least five years during 2011-2016. For each farm, we collected information on the husbandry systems (Table 1), including type of reared livestock species, presence of multiple livestock species on the same pasture, presence and quality of fences around the pasture, surveillance by human shepherds or guardian dogs, and presence of night-time corrals for livestock.
Our wolf occurrence data included 351 georeferenced presence locations, which were genetically-verified from samples that were noninvasively collected (e.g., hair, feces, urine, blood traces and tissues from carcasses) between 2008 and 2012 across the northern Apennines (Caniglia et al., 2014). Data on crop damage by wild boar consisted of georeferenced and dated events (n = 3442) that were compiled by regional authorities for the period 2009-2018 (Centro Servizi ATC Bologna; contacted April 2020). These data were collected by expert surveyors as part of the official damage compensation process provided by the regional government to farmers.
We gathered a dataset of geospatial data characterizing landscape patterns (Table 1) of the study area from the regional government geodatabase (http://geoportale.regione.emilia-romagna.it; accessed December 2018). We considered features that were previously found to explain wolf habitat suitability in Italy (e.g., Falcucci et al., 2013;Milanesi et al., 2015), including variables related to land cover, human disturbance, topography, and distance from refuge-like areas (i.e., protected areas where wildlife hunting is forbidden and human activities are limited), which could affect patterns of livestock depredation risk and, more generally, predators' space use. Notably, land-cover features were derived from a regionally-adapted CORINE land-cover map with 90 classes a spatial resolution of 40 m, and a minimum mapping unit of 0.16 ha. We resampled all features to an equal-area coordinate system (ETRS89/UTM Zone 32 N) using bilinear interpolation, and to spatial resolutions of 100 m, 200 m, 300 m, 400 m, 500 m, and 600 m (see next section for rationale) using a moving window approach. Note that each of these resolutions is much larger than the original spatial resolution of the land-cover map.

Analytical framework
Our overall approach consisted of three main steps (Fig. 2). First, to understand the main drivers of livestock depredation, we modelled the number of depredation events per pasture (considering farms surveyed in the field; Fig. 1) as a function of landscape patterns and husbandry systems. Second, we modelled and mapped wolf habitat suitability across the study area as a function of landscape patterns. Third, we modelled the relationship between wolf habitat suitability and the distribution of crop damage by wild boar as function of time in 2009-2018, a period over which the wolf population size has increased significantly (Galaverni et al., 2016).
To understand the determinants of livestock depredation events, we considered only depredation events on sheep and goats, which accounted for about 90% of all events registered during 2011-2016. We did this to reduce data heterogeneity due to different husbandry systems for cattle and equids (Pimenta et al., 2018) and the overall small number of depredation events for these species. As response variable, we used the number of depredation events per pasture (maximum = 6). To assess the importance of each potential determinant of livestock depredation (Table 1 -all variables), we used generalized linear models (GLM). First, we tested for multi-collinearity among variables using the variance inflation factor (VIF), by fitting an exploratory model including all variables, with threshold set at VIF > 3, but no collinearity was found. We used a Poisson error distribution, and the Kolmogorov-Smirnov test to assess data dispersion (using the R package DHARMa; Hartig and Hartig, 2017). This confirmed a Poisson distribution of errors as fitting our response variable. We also tested for spatial autocorrelation in our response variable using a test based on Moran's I, as implemented in the R package ape (Paradis and Schliep, 2019), and found no indication of spatial autocorrelation. We then built three GLMs: one with all variables (n = 9), a second with husbandry system variables only (n = 5), and a third with the landscape pattern variables only (n = 4), the resolution of which was correspondent to the wolf habitat suitability model with highest AUC score (see next section). We compared the performance of these three models using Nagelkerke's R 2 calculated through the R package Performance (Nakagawa et al., 2017). For each of the three models, we compared all possible combinations of variables using the R package MuMIn (Barton, 2016), and ranked each individual (sub-)model based on the Bayesian Information Criterion (BIC). We then averaged the best-performing models (i.e., models within a ΔBIC ≤4 from the best-performing model; Anderson and Burnham, 2002;Grueber et al., 2011) and assessed final univariate variable importance by hierarchical partitioning (Mac Nally, 2002).
To predict wolf habitat suitability across our study area, we used Maximum Entropy modelling (MaxEnt; Phillips et al., 2006). This algorithm performs well with small sample sizes, ranks consistently among the highest performing algorithms, and can consider spatial bias in occurrence data by using similarly biased background data, as needed in our analysis (Elith et al., 2011). We used the wolf occurrence points as presence locations, and as they were spatially clustered towards areas where wolf packs occurred during the time of data collection, we created a Gaussian kernel density map and sampled ten times as many background points as wolf presence records (Harris et al., 2014;Merow et al., 2013), stratified according to this density pattern (Kramer-Schadt et al., As predictor variables, we tested the landscape pattern variables (Table 1 -Landscapes patterns), summarized at different resolutions (i. e., 100 m, 200 m, 300 m, 400 m, 500 m, 600 m). We systematically compared MaxEnt models for each set of predictors (i.e., each resolution) to find the best-performing set given our dataset. For every model, we tested for potential multi-collinearity among predictors using VIF > 3 as a threshold, but found no indication of multi-collinearity. To identify the best MaxEnt parameterization for each model, we used the R package ENMeval (Muscarella et al., 2014), aiming at a balance between goodness-of-fit and model complexity (Merow et al., 2013). We tested a wide range of penalties for model complexity (i.e., 0.5 to 8, at 0.5 steps) and combinations of predictor transformations (i.e., linear, hinge, quadratic, linear-quadratic). We then selected the best-performing parametrization based on the Akaike Information Criterion (AIC), calculated using five-fold cross-validation (James et al., 2013). Finally, we compared among the best-performing models per resolution to identify the single-best model using the averaged area under the curve (AUC) statistic, also calculated based on five-fold cross validation.
To estimate whether the rebounding wolf population affected the distribution of crop damage across the 10-year period, we built two, complementary regression models. First, we modelled wolf habitat suitability values at crop damage sites as a function of time (model one), with years as the time step. Second, we modelled damage occurrence at damage sites as function of wolf habitat suitability (model two), where our response variable consisted of damage occurrence vs. absence at these sites. We note that modelling crop damage occurrence vs. a random sample of (pseudo-)absences is not feasible in our case as the occurrence of crop damage depends on wild boar distribution and farmers' willingness to claim for damages. Not controlling for these variables would introduce considerable bias, but spatial data for these factors is unavailable. We therefore focused on damage sites only (i.e., where we know that wild boar occur and that farmers are reporting damages). Both of our models assume wolf habitat suitability patterns to be stable during our observation period, which appears justified because wolves had recolonized the entire area before that period and generally showed a stable distribution in the Apennines, also at fine spatial scales (Caniglia et al., 2014;Mancinelli et al., 2018). However, the wolf population increased substantially during the observation period (Galaverni et al., 2016), suggesting a potentially increasing effect on wild prey.
As the response variable in model one (wolf habitat suitability ~ time) was scaled between 0 and 1, we used quasi-binomial regression modelling (Consul, 1990). We extracted the wolf habitat suitability and year of damage occurrence for event and modelled wolf habitat suitability as a function of time (years). A positive correlation between these two variables indicates that damage increasingly occurs in better wolf habitat, for example as wild boars increasingly use crops as refuges against increasing predation pressure. In contrast, a negative correlation indicates decreasing damage prevalence in areas most suitable for wolves, for example because wild boars seek to avoid such areas due to a higher risk of predation. The response variable in model two (damage occurrence ~ wolf habitat suitability * time) was binary (presence vs. absence of damage) so we used logistic regression modelling. We identified the presence and absence of damage for every year at sites where at least one damage was reported in the study period. We then modelled damage likelihood as a function of wolf habitat suitability and time (years). A positive coefficient of the interaction of suitability * time indicates an increasing likelihood that damages occur in better wolf habitat, whereas a negative correlation indicates decreasing damage recurrence in areas most suitable for wolves.
After detecting spatial autocorrelation in model residuals for both models using Moran's I test (R package spdep; Bivand and Piras, 2015), we decided to use Generalized Linear Mixed Models using Penalized Quasi-likelihood (R package MASS; Ripley et al., 2013), which perform well for both logistic and quasi-binomial distributions (Zuur et al., 2009). To limit effects of spatial and temporal autocorrelation, we included spatial correlation structure between sites of crop damage  Pinheiro et al., 2017), and we used the municipalities of the study area as random intercept factor.

Determinants of livestock depredation
Our data included a total of 528 livestock heads that were killed in 165 individual depredation events between 2011 and 2016. The mean annual number of animals killed was 88 (standard deviation SD = 38.5) in 29 events (SD = 11.0). Sheep were killed in the majority of depredation events (69.4% of events), followed by goats (18.8%), cattle (10%), and equids (1.8%). Our averaged model of livestock depredation risk considering all variables, based on depredation events of sheep and goats only, had a Nagelkerke's R 2 of 0.44. Eleven models were equally well-supported  (BIC ≤ 4) as the best-performing model. We thus averaged them, and our final average model included six variables (presence of equids/cattle, quality of fences, pasture area, degree of surveillance, number of buildings, presence of corrals; Fig. 3). All eleven models included the variable quality of fences, eight included the variable presence of equids/cattle, and five the variable pasture area. Only one out of the four landscape pattern variables (i.e., number of buildings) was selected in the final model.
In terms of relative variable influence (Fig. 3-A), livestock depredation risk increased with lower fence quality (coefficient ± standard error = 1.20 ± 0.39, p = 0.002), the presence of equids/cattle on the same pasture with sheep and goats (0.94 ± 0.37, p = 0.012), and larger pasture area (0.02 ± 0.01, p = 0.054). Lacking herd surveillance, the absence of corrals, and lower number of buildings were also all positively related to increasing depredation risk, but these relationships were not significant (p > 0.05; Table 2). Comparing the model with all variables, based on both husbandry systems and landscape patterns, to models using only one of these variable groups showed that husbandry systems model had only a marginally lower Nagelkerke's R 2 (0.41) than the full model, whereas the model based on only landscape pattern variables performed poorly (R 2 = 0.08).

Wolf habitat suitability
The sensitivity analyses of habitat suitability models calculated for predictor variables of different spatial resolutions showed the highest goodness-of-fit for a resolution of 500 m. We report here only results for this resolution and refer to the Supplementary Material for the full results. Our wolf habitat suitability model overall had a high goodness-offit, with an AUC of 0.81. Our Maxent model-selection procedure suggested a model with low complexity penalty (i.e., 1) and only hinge features as best-performing. The most important explanatory predictor in the model was the share of woodlands (with a modelling contribution of 69%), followed by distance from protected areas (21%) and terrain ruggedness (10%). The predictor number of buildings did not contribute appreciably to the model (< 1%), although this variable was sometimes more important in coarser resolution models (Fig. S1). The response curves of our best-performing model (Fig. 4-B) showed generally plausible relationships, highlighting an increase in wolf habitat suitability with an increasing share of woodlands and a decrease of wolf habitat to matrix/unsuitable levels distant from protected areas and in rugged terrains.
Predicting wolf habitat suitability showed that patches of highquality habitat occurred across the study area and at all elevations, generally in areas where human disturbance is lower in this landscape ( Fig. 4-A). Matrix habitat was particularly widespread in the norther study area, where intensive human pressure prevails, coinciding with the suburban area of Bologna and the margin of the Po plain. In the southern study area, widespread suitable habitat was predicted along the northern Apennines ridge.

Crop damage risk and wolf habitat suitability
Our crop damage dataset included a total of 3442 events occurring in all 34 municipalities (min = 3 events/municipality, max = 810 events/ municipality, average = 101 ± 142.2 SD). Crop damage by wild boar was distributed fairly evenly over time, with average 344 damage events per year (SD = 64.0). The highest annual number of damages occurred in 2012 (497 events). Crop damage was generally widespread in the study area, but more frequent in the southern part (i.e., in more mountainous and densely forests areas; Fig. 5-B). On average, 3.4 (SD = 3.4) crop-raiding events occurred at the same site.
Our first model (model one), which assessed wolf habitat suitability at crop damage localities as a function of time (with years as a time step), revealed a significant correlation between these variables (− 0.0191 ± 0.0077, p = 0.014; Table S2). This correlation was negative, suggesting a tendency of crop damage to occur in areas of lower habitat suitability for wolves over time and of less crop damage events to occur in highly suitable wolf habitat (Fig. 5-A). While the coefficient for this model was relatively small, due to the temporal unit of one year in our model, the overall effect is considerable when aggregated over time. For example, our results suggest a decrease in habitat suitability values at crop damage sites of about 1.9% per year, or 19.1% over the entire time period covered by crop damage dataset (2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018). Our second model (model two), which assessed the likelihood of crop-raiding recurrence at same sites of damage as a function of wolf habitat suitability and time (again, with annual time step), also revealed significant correlation between these variables (− 0.1089 ± 0.0412, p = 0.008; Table S3). The negative coefficient here suggests a tendency of crop damage likelihood to decrease in sites located in areas of higher wolf habitat suitability (Fig. 5-A). This result is consistent with those from our first model of crop damage.

Discussion
Human-carnivore conflicts translate into major challenges for the conservation of large carnivores, particularly where these species return to human-dominated landscapes after historical extirpation (Ghoddousi et al., 2021;Wolf and Ripple, 2017). Understanding what determines patterns of conflict in such situations is key for fostering long-term coexistence of humans and large predators (Chapron and López-Bao, 2016;Lamb et al., 2020). Efforts have so far mainly focused on conflicts directly attributable to carnivores, such as what determines livestock depredation risk and where depredation risk is high. However, this neglects the additional and potentially widespread indirect effects that returning top predators can have on the space use of their natural prey, which themselves are associated with conflict, such as in the form of crop damage.
Focusing on the case of a rebounding wolf population in the Italian Apennines, we carried out, to the best of our knowledge, the first exploration of both types of conflict, using spatial models of wolf habitat suitability, livestock depredation risk, and crop damage distribution by wild boar. Two main findings emerge from our work. First, our results indicate a clear association of livestock depredation events with inadequate livestock protection, independent of landscape patterns. This suggests that re-adopting livestock protection measures, which have been largely abandoned due to the long absence of wolves in the area, can substantially minimize direct conflict between wolves and people. Second, we found crop damage to decreasingly occur in areas most suitable for wolves, according to our habitat mapping. We suggest this indicates a redistribution of crop damage in our study area as predation pressure by wolves grew during the time when the wolf population increased, and wild boar frequented areas associated with higher human disturbance, where they can avoid wolves, more often. On a more general level, our study highlights that managing for human-carnivore coexistence must consider the geography of multiple, co-occurring conflicts, as carnivores are recolonizing their historical ranges.
Inadequate livestock protection, particularly weak fences, increased the risk of wolf attacks to livestock in our case (Fig. 3). Fences are costly and difficult to maintain, particularly for large pastures and in rugged terrains, often resulting in inefficient prevention (Dondina et al., 2015). Large pastures are frequently associated with mixed livestock groups, including equids and cattle, and such husbandry systems were related to higher depredation risk in our study (Fig. 3-A) and elsewhere (e.g., northern Portugal; Pimenta et al., 2018). Indeed, cattle and equids require large pastures, which are harder to keep fenced effectively (e.g., maintenance of holes in fences). Moreover, larger pastures are usually more isolated than smaller pastures, thus having a greater interface with woodlands. Similarly, smaller pastures are often closer to settlements. Together, these factors likely explain why sheep and goats reared alongside horses and cattle suffer from higher depredation risk. In the Italian northern Apennines, pasture fencing for preventing wolf attacks has only recently been adopted (Gazzola et al., 2008) and our results highlight the potential for lowering depredation risk if fences were upgraded and more widely used. Additionally, our results corroborate views that sheltering livestock during the night and surveilling livestock on pastures, both of which has characterized livestock husbandry systems historically (Boitani et al., 2010), could further decrease the depredation risk (Dondina et al., 2015;Ogada et al., 2003;Pimenta et al., 2017). Husbandry system variables were generally more important than landscape patterns in explaining depredation risk in our case. This is in line with recent work by Pimenta et al. (2018), which stresses the importance of preventive measures on pastures located inside the home ranges of wolf packs in northern Portugal, regardless of landscape context. In contrast, prior work in our region highlighted the importance of specific landscape variables in determining depredation risk, including elevation and the share of grasslands (Milanesi et al., 2019). However, this study did not consider variables related to the husbandry system, which is not surprising, because such variables are typically hard to obtain and require extensive fieldwork. Our work thus highlights the importance of accounting for characteristics of the livestock husbandry systems, besides landscape patterns, and cautions against causal interpretations of models solely based on landscape variables. While  landscape variables can indeed be useful indirect proxies for assessing depredation risk, allowing for upscaling and predictions across larger areas, variables that more directly describe livestock husbandry features, such as the existence of fencing, the composition of livestock herds, or the presence of herd surveillance and overnight structures, are likely more powerful predictors of depredation risk, and therefore more useful to inform conflict management.
Interestingly, we did not find a direct relationship between the presence of guardian dogs and livestock depredation risk, although guarding dogs are usually recognized as an effective prevention measure (Mosley et al., 2020). However, we note that guarding dogs, particularly dogs staying with livestock herds day and night, are not common in the northern Apennines. Indeed, in many of the farms we visited, the guardian dogs were juvenile or lacking proper training, which reduces their effectiveness in protecting livestock (Khorozyan et al., 2017). As a result, our findings should not necessarily be interpreted as discouraging the use of guarding dogs or questioning their usefulness in the Apennine context or elsewhere.
A key finding of our study was that crop damage by wild boar, the principal prey of wolves in Italy (Mori et al., 2017), decreased in areas most suitable for wolves and increased elsewhere. A possible explanation for this finding is the phenomena of the 'landscape of fear' that predators create and that prey respond to (Laundré et al., 2014). There is increasing evidence for top-down effects of carnivores on the behavior of prey species, although these effects have predominantly been shown for undisturbed, natural systems (Ripple et al., 2014). In human-dominated landscapes, these effects are typically weaker, and the landscape of fear is driven by human disturbance, to which both predators and their prey respond by avoiding human-dominated areas (Ciuti et al., 2012;Dorresteijn et al., 2014;Gaynor et al., 2019). In our study area, populations of wolves and their wild prey mainly cluster in patches of semi-natural vegetation surrounded by a matrix of human-dominated areas, including agricultural areas and settlements (Mancinelli et al., 2018). These semi-natural patches thus function as refuges for both wolves and prey. With wolf populations increasing, prey must trade-off between avoiding humans and avoiding predation-risk by carnivores (Gordon, 2009). Prey species often manage this trade-off by becoming more tolerant to human disturbance and sometimes even seeking out such areas, as this reduces predation risk (Berger, 2007;Palmer et al., 2021;Schmitz et al., 2004). In our situation, with wolf population size and predation pressure increasing, we would expect wild boar to increasingly avoid the most suitable habitat patches for wolves. This, in turn, would increase wild boar presence in more human-disturbed areaswhich is what we found, as crop damages were more likely to occur in less suitable areas for wolves over time (Fig. 5-A).
Key assumptions of this interpretation are that wolves had recolonized available habitat in our study area prior to our observation period, and that despite this, the number of wolves increased in our observation period. Both assumptions are supported by available observations and data, with wolves numbers doubling across the Apennines since 2008 (Galaverni et al., 2016), when wolves were already found all across our study area (Caniglia et al., 2014). A plausible explanation is thus that an increasing wolf population led to an expanding landscape of fear through increasing predation pressure, thereby redistributing wild boars and their crop damage in the landscape. If the wolf population would not have been in equilibrium with its environment in our study area, wolves likely have continued to fill in less suitable, secondary habitat after 2012. This cannot be ruled out, although we consider it to be an unlikely scenario, since wolves already occurred in the north of our study area close to the Po plain, one of the most human-modified areas in Europe, in 2012 (Caniglia et al., 2014). Regardless, the observed response we found could still be explained by wild boars increasingly seeking shelter in the more human-dominated areas frequented by wolves in the last years of the period covered by our study, such as agricultural areas within a natural matrix. Unfortunately, we were unable to quantitatively test whether wolf habitat suitability was stable during our study period due to missing robust wolf occurrence data after 2012.
To our knowledge, this is the first study providing indication for returning large carnivores having an effect on prey-related humanwildlife conflicts (i.e., the redistribution of crop damage by wild boar in our case). Our findings are broadly in accordance with the predator-prey space race theory (Muhly et al., 2011;Sih et al., 2005) and in line with studies suggesting large predators exert some level of top-down control in human-dominated landscapes as well (Dorresteijn et al., 2015). From a conflict perspective, changing prey behavior in landscapes recolonized by carnivores may have negative effects (e.g., higher intensity or redistribution of crop damage) or positive effects (e.g., lower risk of car collisions; Gilbert et al., 2017). Quasi-experimental studies, for example comparing recolonization areas with appropriate counterfactual areas where carnivores are still absent, could quantify these effects, thereby informing land-use and conservation planning. Still, we caution against over-generalizing our results. The outcomes of large carnivore recovery in human-dominated landscapes are likely dependent on the socialecological context in which recovery happens. For instance, if carnivores are culled when exceeding tolerated population densities, which are usually very low, their ecological effects on wild prey's behavior are small (Kuijper et al., 2016). Likewise, intensive poaching of carnivores may result in high rate of pack turnover, which can increase livestock depredation (Lovari et al., 2007).
Many species of ungulates are considered overabundant in Europe, leading to widespread damage to agriculture and forestry (Apollonio et al., 2010). Given the continued range expansion and population growth of large carnivores in Europe, and particularly of wolves, the effects of re-distributing patterns of damage could be widespread, but remains so far largely unexplored. This is unfortunate, given that damages by ungulates are economically much more significant than the direct damages by carnivores through livestock depredation (Carpio et al., 2021). Furthermore, damages by vehicle collisions with ungulates are also economically relevant (Raynor et al., 2021), and frequencies and spatial patterns of such collisions could be affected by the return of large carnivores. From a conservation planning and wildlife management perspective, this means that the prediction of carnivore range expansion should be accompanied by assessments of potential effects by carnivores on wild prey behavior. More generally, holistic assessments of large carnivore recoveries should consider trophic interactions as large carnivores return to human-dominated landscapes, as well as the services and disservices that large carnivores deliver (Rode et al., 2021). Designing effective conservation strategies to achieve coexistence thus requires that we simultaneously address multiple types of humanwildlife conflict, often with stakeholders with diverging interests (König et al., 2020), to achieve long-term coexistence of humans and wildlife in shared landscapes (Carter and Linnell, 2016;Kuijper et al., 2019Kuijper et al., , 2016.

Author statement
M.D. and T.K. conceived the ideas. M.D., A.G., F.M.S. and T.K. designed the study. M.D., E.F. and R.C. provided key datasets. M.D., A.G. and F.M.S. implemented the analyses. M.D. and T.K. led the writing. All co-authors contributed critically to interpretation and manuscript drafts.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. particular A. Romero-Muñoz and J. Oeser for their help in processing input data of the livestock depredation risk modelling and the model of wolf habitat suitability. We are grateful to D. Zanin (ATC Bologna) and M.L. Zanni (Regione Emilia-Romagna) for providing data and technical support. We also thanks Rhys Lemoine of the BIOCHANGE Group at Aarhus University and Casey Lawrence (Department of English, Trinity College Dublin) for English language editing on the manuscript. M.D. was supported by the German Academic Exchange Service (DAAD) via a Short-Term Grant. A.G. appreciates the financial support of the German Research Foundation (DFG; PArCS project #409732304). We thank two reviewers and Editor G. Cummings for encouraging and very helpful comments that improved this manuscript in major ways.