Using occupancy models to assess the direct and indirect impacts of agricultural expansion on species’ populations

Land-use change is a global threat to biodiversity, but how land-use change affects species beyond the direct effect of habitat loss remains poorly understood. We developed an approach to isolate and map the direct and indirect effects of agricultural expansion on species of conservation concern, using the threatened giant anteater (Myrmecophaga tridactyla) in the Gran Chaco as an example. We reconstructed anteater occupancy change between 1985 and 2015 by fitting single-season occupancy models with contemporary camera-trap data and backcasting the models to 1985 and 2000 land-cover/use maps. Based on this, we compared the area of forest loss (direct effect of agricultural expansion) with the area where forests remained but occupancy still declined (indirect effect of agricultural expansion). Anteater occupancy decreased substantially since 1985, particularly after 2000 when agriculture expanded rapidly. Between 1985 and 2015, ~ 64,000 km2 of forest disappeared, yet occupancy declined across a larger area (~ 102,000 km2), extending far into seemingly untransformed habitat. This suggests that widespread sink habitat has emerged due to agricultural land-use change, and that species may lose their habitat through direct and indirect effects of agricultural expansion, highlighting the urgent need for broad-scale conservation planning in the Chaco. Appropriate management responses could proactively protect more habitat where populations are stable, and restore habitat or address causes of mortality in areas where declines occur. Our work also highlights how occupancy modelling combined with remote sensing can help to detect the direct and indirect effects of agricultural expansion, providing guidance for spatially targeting conservation strategies to halt extinctions.


Introduction
We are currently experiencing the highest extinction rates in Earth's history, mainly due to habitat destruction (Ceballos et al. 2015). This is particularly the case in the tropics and subtropics, where agricultural expansion and intensification are causing widespread destruction of natural vegetation (Laurance 2010). Understanding the effects of land-use change on biodiversity in these regions is therefore essential (Newbold et al. 2015), especially where agriculture expands into forests (Kehoe et al. 2017).
The most direct impact of agricultural expansion on biodiversity occurs via habitat loss and degradation. Fragmentation can also affect species directly by decreasing dispersal among patches (Hanski 2015). Additionally, indirect effects of agricultural land-use change can be observed, for instance through edge effects (Smith et al. 2018). This is because the accessibility of forest patches to humans increases in fragmented landscapes, and thus hunting pressure or forest fires can exert increasing pressure on species (Corlett 2007;Tabarelli et al. 2004). For example, human presence and hunting pressure are typically higher along edges of forests and protected areas (Brashares et al. 2001;Wittemyer et al. 2008). Likewise, hunting commonly occurs in or around areas where agriculture expands (Peres 2001;Romero-Muñoz et al. 2019, 2020bSemper-Pascual et al. 2019). Moreover, road expansion due to agricultural expansion not only increases accessibility, but also facilitates movements by animals, thus increasing the risk of road kills (Laurance et al. 2006(Laurance et al. , 2009. Finally, habitat fragmentation also alters biophysical conditions which in turn affect species distributions (Brook et al. 2008). Such effects go beyond the direct conversion of habitat through agricultural expansion, but are a consequence of this landuse change. We therefore refer to such effects as indirect effects of agricultural expansion. Indirect effects are typically hard to detect and thus often overlooked, despite their importance (Brook et al. 2008). Approaches for identifying where and how indirect effects impact species of conservation concern are therefore needed.
Much of the work assessing the impacts of agricultural expansion on biodiversity has relied on models based on presence-only data (Gaston 2009;Guillera-Arroita et al. 2015), either for individual species [e.g., habitat, Romero-Muñoz et al. (2019)] or entire communities [e.g., richness, Ferrier and Guisan (2006)]. Models based on this type of data can capture where species potentially occur, but are insensitive to detecting periods of population decline that precede local extirpation (Hobbs and Mooney 1998). Yet local extinctions can be preceded by many years of population decline (Hylander and Ehrlén 2013;Norris 2016), thus providing a window of opportunity to prevent those losses. Unfortunately, identifying these windows of opportunity is hard, as parametrizing models that can capture population trends, such as capture-recapture models (Chandler and Royle 2013), are unfeasable or too costly for many species. Operational tools for detecting windows of opportunity preceding local extinctions are therefore often missing.
Estimating species occupancy can be a viable alternative. Occupancy is defined as the probability of a site being occupied by a species (MacKenzie et al. 2002), and changes in occupancy provide a robust proxy for population declines, especially useful for elusive species which are difficult to detect (Beaudrot et al. 2016). Yet, even though occupancy models are increasingly used in conservation (Guillera-Arroita 2017), studies typically focus on a single time period (Lindenmayer et al. 2012). Developing approaches which allow us to track occupancy in space and time would therefore be highly useful for mapping where population declines might be ongoing.
Furthermore, most occupancy studies assessing the impacts of agricultural expansion on biodiversity typically do not distinguish between direct and indirect effects (Barlow et al. 2016;Gibson et al. 2011;Raiter et al. 2014). However, distinguishing between direct and indirect effects is crucial to spatially targeting conservation actions, in order to halt population declines of species of conservation concern. For example, proactively protecting remaining habitat from being converted should focus on areas where populations are not already declining due to anthropogenic factors (Gilroy and Edwards 2017). Conversely, other reactive conservation measures might be warranted to reduce mortality in landscapes where threatened species are in decline (e.g., anti-poaching measures, habitat restoration, corridor planning), especially if declines occur inside protected areas.
Here, we develop a new approach that backcasts occupancy models to separate the direct and indirect effects of agricultural expansion, using the Giant Anteater (Myrmecophaga tridactyla) (hereafter: anteater) as an example. Anteaters occur widely in Central and South America and the species is listed as vulnerable by the International Union for Conservation of Nature (IUCN) due to a * 30% population loss during the last decades (Miranda et al. 2014). Forest is the most important vegetation type for anteaters (Di Blanco et al. 2015), however, other studies have shown that they can also be associated with open habitats [e.g., grasslands or savannas; Di Blanco et al. (2015), Medri and Mourão (2005)]. Anteaters can also occur in anthropogenic landscapes when these landscapes provide food resources (Pardo et al. 2019). Although habitat loss has been proposed as a main threat to the species (Miranda et al. 2014), the effects of agricultural expansion on anteaters are poorly understood. To our knowledge, no study has evaluated how agricultural expansion affects anteater populations directly (e.g., through habitat conversion) and indirectly (e.g., through edge effects) anywhere in their range.
The overarching goal of our study was to estimate giant anteater occupancy in the Argentine Dry Chaco, based on contemporary camera-trap data and contemporary landcover/use maps. We then projected this occupancy model to past land-cover/use maps from 1985 and 2000 to reconstruct anteater occupancy changes in relation to agricultural expansion. Specifically, we addressed three research questions: 1. What factors influence contemporary anteater occupancy in the Argentine Dry Chaco? 2. How did agricultural expansion influence anteater occupancy between 1985 and 2015? 3. What is the relative importance of direct versus indirect effects of agricultural expansion on anteater occupancy changes?

Study area
The Gran Chaco is the largest tropical/subtropical dry forest worldwide (1.1 million km 2 ), extending into Argentina, Paraguay, Bolivia, and Brazil. The ecoregion is characterized by heterogeneous landscapes, with forests as the dominant land-cover (* 60% of the study region), followed by open habitats such as pastures (* 15%) and natural grasslands and savannas (* 10%) . Although the Chaco harbors high biodiversity (TNC et al. 2005), growing demand for beef and soybean have led to rapid agricultural expansion (Fehlenberg et al. 2017), especially after 2000. In fact, between 1985 of all forest in the Chaco was lost due to agriculture expansion for cropland (38.9% of all forest losses) or pastures (61.1% of all forest losses) , strongly threatening biodiversity (Periago et al. 2015;Semper-Pascual 2018).
We assessed giant anteater occupancy across the Argentine Dry Chaco. By focusing on the Dry Chaco only (i.e., excluding the Wet Chaco), we ensured that environmental conditions were fairly homogeneous across our entire study area. To delineate our specific study area, we overlaid the boundary of the Argentine Dry Chaco ecoregion with the anteater range according to the IUCN (Fig. 1, * 33 million ha). Our study area has a highly seasonal climate, with precipitation ranging from 450 to 900 mm and annual temperatures ranging between ? 42°C in summer to -7°C in winter (Boletta et al. 2006;Minetti et al. 1999).

Anteater data
We combined camera-trap data from three surveys which were all established to assess the mammal community structure in the study area: the first survey (survey1) was carried out in 2013 (Decarre 2015), the second survey (survey2) in 2014-2015 (Gómez-Valencia 2017) and the last survey (survey3) in 2016-2017 (Proyecto Quimilero). Camera-trap sites were distributed across the extant range of the species (northern part of the study area, Fig. 1). We did not sample the southern part of the study area, as according to the current IUCN range map, anteaters are absent there (Miranda et al. 2014). Sites were selected randomly while avoiding trails, to sample the full diversity of vegetation types in the region, and to avoid bias due to high capture rates along trails (Kolowski and Forrester 2017;Wearn et al. 2013). In addition, such a sampling design has been recommended for estimating occupancy (Burton et al. 2015;Wearn et al. 2013). Camera-trap sites were separated by at least 1500 m from one another, a distance approximately representing the radius of an average circular anteater home range (Table S1). Further details on the camera placement are provided in Online Appendix S1.

Predictor variables
We used four variables for modelling anteater detection probability: camera-trap survey, sampling effort, temperature, and precipitation (Tables 1, S2). For occupancy probability, we selected four landscape variables known to influence anteater habitat selection: percentage of forest cover, percentage of open habitats, percentage of edge between forest and open habitats, and a diversity index of suitable habitats (i.e., forest and open habitats) (Tables 1, S2). Importantly, these variables by themselves are not measures of direct or indirect effects of agricultural expansion.
The entire range of each landscape variable was represented by our camera-traps (Fig. S1). Given that camera-trap data were collected during different years, we used two high-resolution (30 m) land-cover/use maps to relate each camera-trap site to the corresponding land-cover/use map. We therefore related the camera-trap sites from survey1 to a land-cover/use map from 2013 ) and the camera-trap sites from survey2 and survey3 to an updated land-cover/use map from 2015 consistent with Baumann et al. 2017 (Fig. S2). We thus parametrized one single model. Such a modelling approach (i.e., time-calibrated model) is widely used and has several advantages compared to models based on data from a single time period (Kuemmerle et al. 2012;Nogués-Bravo 2009;Romero-Muñoz et al. 2019;Sieber et al. 2015). Time-calibrated models use all available data to parameterize a single model not linked to a specific time period, to then project this model into different time periods (Nogués-Bravo 2009). Consequently, timecalibrated models describe species' niches more fully, are less prone to sampling bias and ensure that observed changes are solely due to changes in predictor variables (and not due to different model parameterizations) (Kuemmerle et al. 2012;Nogués-Bravo 2009). The latter is particularly important in the context of our goal to separate the direct from the indirect effects of agricultural expansion.
To derive landscape indices, we used Morphological Spatial Pattern Analysis (Vogt et al. 2007) and Shannon's diversity index (McGarigal et al. 2012) and summarized them using a 1500 m radius circular moving window. We also tested a larger radius, however, as models based on variables for the 1500 m circular radius fitted best and align with home range estimates, we only show these results. In addition to the landscape variables, we included three variables reflecting human disturbance: distance to settlements, distance to roads, and distance to paved roads (Tables 1, S2). Importantly, all these human-disturbance indicators were time-invariant (i.e., they were calculated using data from * 2015), as we assumed that expansion of main roads or settlements did not greatly affect the region during the study period (Romero-Muñoz et al. 2020a). Table 1 Variables for estimating anteater detection probability and occupancy probability based on Desbiez and Medri (2010), Di Blanco et al. (2015Blanco et al. ( , 2017 and Medri and Mourão (2005)  Roads are an indicator of disturbance (e.g., providing access for hunters) PavRoads b Euclidean distance to the closest paved road Paved roads are accessible after heavy rains contrary to non-paved roads and human disturbance is expected to be more persistent a Correlated variables (r C 0.6) that were not included in the same detection model b Correlated variables (r C 0.6) removed from the final dataset for modelling occupancy. Since forest is considered essential habitat for the species (Di Blanco et al. 2015), when two variables were correlated we retained the one including forest (i.e., forest and edge) We aggregated all variables to a 300 m pixel resolution, using an equal-area coordinate system. This pixel resolution is much smaller than the average anteater home range (* 9 km 2 , Table S1). We standardized all variables (0 = mean, 1 = standard deviation) and checked for collinearity. When two variables were correlated (Spearman correlation coefficient C 0.6), we retained the ecologically more meaningful variable (Table 1).

Modeling anteater occupancy
We used a likelihood-based, single-season occupancy model (MacKenzie et al. 2017) to estimate anteater occupancy (W) in relation to landscape and human disturbance variables, while accounting for detectability (p). During wildlife surveys, individuals are often unobserved, either because they are truly absent at the sampling site, or because they are present, but undetected. Visiting each sampling site on multiple occasions enables the estimation of detection probability, which can be incorporated into the occupancy modelling to correct for imperfect detection (MacKenzie et al. 2002).
We generated sampling occasions by pooling daily detection/non-detection records for each sampling site into consecutive camera-days. We carried out sensitivity analyses to check for differences in occupancy and detection estimates when pooling daily detection/ non-detection records into different number of camera-days. We defined an occasion as an interval of seven camera-days (see Fig. S3). We only included sites that were actively surveyed for more than 21 consecutive days to obtain a minimum of three occasions per site, and a maximum of 84 days (i.e., 12 occasions) to account for the closure assumption of occupancy models. Therefore, we estimated true occupancy as we assumed that changes in occupancy between occasions did not occur. The number of occasions per site ranged between 3 and 12 (mean = 6.1) with a sampling effort per site ranging from 21 to 84 days (mean = 43). The final sample included a total of 106 camera-trap sites.
Occupancy models link a state model determining occupancy at each site with an observation model for detection which is conditional on occupancy. We first modelled detection probability and then occupancy probability. To model detection probability, we only included detection variables, using a null occupancy model. We considered all possible combinations of variables, while excluding correlated variables in the same model, and tested a quadratic effect of temperature as we hypothesized that anteater detection is lower during extreme temperatures. We ranked the resulting models using Akaike's Information Criterion, corrected for small sample size [AICc; Burnham and Anderson (2002)], where the detection model with the lowest AICc is best-fitting.
To model occupancy probability, we kept the best-fitting detection model constant and included variables that could affect anteater occupancy. We excluded correlated variables (see Table 1) and built six candidate models ( Table 2). Given that habitat loss is the main threat to anteaters (Miranda et al. 2014), we built a first set of candidate models testing all possible combinations of landscape structure variables. In addition, to test whether human disturbance may affect anteater occupancy, we built a second set of candidate models in which we included human disturbance variables to the models from the first set (Table 2). We compared all resulting models, ranked them according to AICc and considered models with a DAICc \ 2 to be competing. We used model averaging by calculating weighted average estimates based on models with DAICc \ 2 and following Burnham and Anderson (2002). To assess model fit, we ran a goodness-of-fit test on the global model and calculated an overdispersion parameter (cˆ) (MacKenzie and Bailey 2004), using 10,000 bootstrap replicates and the Pearson v 2 statistics. We fitted our models using the R package unmarked (Fiske and Chandler 2011).

Mapping changes in anteater occupancy over time
To map changes in anteater occupancy over time, we used our single-season, time-calibrated model to predict anteater occupancy for the years 2015, 2000, and 1985. The year 2015 represents the average year of our sampling period and the most recent land-cover/ use map for our study area. The year 2000 marks the beginning of the period when land conversion drastically increased. Finally, 1985 marks a baseline before this wave of conversion and fine-scale satellite images are available for land-use classifications.
Our approach relied on several assumptions. First, we assumed that the relationship between anteater occupancy and landscape variables remained unchanged over the study period, meaning that factors describing contemporary anteater occupancy would also describe past anteater occupancy. Second, land-cover/use variables were the only timevariant variables in our models as we assumed that no major changes in terms of roads/ settlements expansion occurred during the study period. Finally, agricultural expansion is the dominant process leading to change in land-cover. Therefore, if our models find differences in anteater occupancy across time, we assume that these differences can only arise from effects of agricultural expansion on our predictor variables.
To predict current anteater occupancy across our study area, we used model-averaged estimates derived from our time-calibrated occupancy model and the land-cover/use map from 2015. We then backcasted anteater occupancy by projecting the same model to the past land-cover/use maps from 1985 to 2000 (Fig. S2). Habitat models have been widely The first part of the table lists the detection models (only top models and the null model are presented; bestfitting model in bold), the second part shows the occupancy models (models with DAICc \ 2 in bold). Akaike weights (AICc w ) quantifies the probability of each model being the best model p probability of detection, W probability of occupancy used to predict habitat across larger regions, as well as to reconstruct past habitat (Tapia et al. 2018;Varela et al. 2009). Predictions are considered to be reliable when they are made within the range of environmental conditions used for model building, in our case the contemporary data (Elith and Leathwick 2009). To ensure we predicted to similar environmental conditions as those represented by our data, we compared the distributions of predictor variables at the camera-trap sites with the overall distribution of the predictor variables for the study area and for each time period we assessed (Figs. S4, S5). This showed that even though we projected our models across large areas in space and back in time, we did not predict anteater occupancy to areas with environmental conditions not represented by the data used to calibrate our models. In other words, while we projected our models to large areas, we did not extrapolate.
To estimate changes in anteater occupancy over time, we calculated the difference in occupancy predictions for each 300 m pixel between (1)  We assessed the uncertainty of the resulting predictions using a parametric bootstrap procedure (Davison and Hinkley 1997) with 10,000 replications. Further details on the uncertainty assessment are provided in Fig. S6.
To assess the relative importance of direct versus indirect effects of agricultural expansion on occupancy, we compared the area in which occupancy declined with the area directly affected by forest loss due to agricultural expansion (i.e., conversion from forest to cropland and pasture). We only considered conversions from forest to cropland and pasture, not conversions of other habitats such as grasslands and savannas, since forest is considered essential habitat for the species (Di Blanco et al. 2015;Mourão and Medri 2007), and as other land uses are not widespread in the study area. We refer to the direct effect of agricultural expansion as those areas where forest was converted. An indirect effect of agricultural expansion occurs in those areas where occupancy declined yet forests were stable during the observation period (Fig. 2). In these areas, a decrease in anteater Fig. 2 Direct and indirect effects of agricultural expansion on anteater occupancy. Agricultural expansion affects anteater population declines directly where forest is converted to croplands or pastures. However, anteater population declines can additionally occur in nearby forests via indirect effects of agricultural expansion, for example due to the expansion of roads, heavier traffic, easier access for hunters, or fires escaping from agricultural land. In our model, such effects collectively lead to declining anteater occupancy in otherwise stable forests occupancy can only be attributed to indirect effects of agricultural expansion (e.g., increase of human pressure along forest edges or expansion of logging-trails).

Results
Our dataset contained 71 independent anteater captures over 4508 trap-days (mean capture rate = 1.57 captures/100 trap-days). The naïve occupancy (i.e., proportion of sites with at least one detection) was 42%. Our best-fitting detection model contained the variables survey and temperature (Table 2), indicating that detection probability varied among the three camera-trap surveys and that it increased slightly with increasing temperature. The estimated detection probability across all sites was 0.14 (SE = 0.03).
The landscape variable percentage of forest had the strongest influence on occupancy probability and was selected in all three top models (DAICc \ 2). The best-fitting model included only percentage of forest ( Table 2). The second-best model also included percentage of edge, and the third-best model additionally included distance to roads and settlements. Model-averaged coefficients also corroborated that percentage of forest was the most important variable, with a positive effect on anteater occupancy (Table 3 and Fig. S7). Percentage of edge had also a positive effect on occupancy, but weaker than percentage of forest (Table 3). Regarding the human disturbance variables, occupancy probability increased further away from settlements, but decreased with distance to roads ( Table 3). The goodness-of-fit test on the global model indicated that our model fitted the field data well (p value = 0.95) and that the data were not overdispersed (ĉ = 0.05).
Occupancy maps highlighted that areas with a high percentage of forest cover had the highest occupancy probabilities, while agricultural areas had the lowest occupancy probabilities (Figs. 3, S2). Areas with high occupancy values were more extensive in 1985 than in 2000, and experienced an even greater decline in area between 2000 and 2015 (Fig. 3). The mean estimated occupancy across 10,000 pixels randomly selected from the study area was 84.59% (SE = 13.23) in 1985 and 72.90% (SE = 16.18) in 2015.
Estimated occupancy changes for individual pixels ranged from -84 to ? 35%. We classified occupancy changes into four categories based on quartiles (Fig. S8): (a) increase/ stable occupancy where occupancy increased (increase of 6-35%) or stayed stable (increase or decrease up to 6%), (b) low decrease (6-19% decrease), (c) medium decrease (19-46% decrease), and (d) high decrease (46-84% decrease). High decreases in occupancy mainly occurred in agricultural areas (e.g., agricultural frontiers in the Chaco and Santiago del Estero provinces). Conversions to cropland typically resulted in highest occupancy decreases, whereas lower decreases were seen for conversions to pasture (Fig. S9). Areas where occupancy stayed stable or increased were widespread in 1985-2000, but became scarce in 2000-2015. The same was true for the 1985-2015 period, suggesting that anteater occupancy decreased especially after 2000. Only areas characterized by a high share of forest cover (e.g., the El Impenetrable region in the north of Chaco province or the Copo National Park in the north-east of Santiago del Estero province) had increasing or stable occupancy at the end of our study period (Fig. 3). The uncertainty assessment of our occupancy change map based on parametric bootstrapping indicated that the standard deviation of our predictions of occupancy change was ± 8.39% (Fig. S6).
To separate the direct and indirect effects of agricultural expansion, we compared the area of forest loss with the area where anteater occupancy declined in our study region. Anteater occupancy decreased over wide areas and these areas were much more extensive All variables are standardized to 0 mean and 1 standard deviation for comparison p probability of detection, W probability of occupancy than the area directly affected by forest loss. Decrease in anteater occupancy occurred over approximately 102,000 km 2 (31% of the study region, Fig. 4). In contrast, only 64,000 km 2 (20% of the study region) experienced forest conversion to agriculture, with 34,000 km 2 being converted to pastures and 30,000 km 2 to croplands (Fig. 4).

Discussion
Agricultural expansion in the tropics is a major cause of local extinctions. Focusing on anteaters in the Argentine Dry Chaco, we demonstrate how projecting a time-calibrated occupancy model back in time can reconstruct occupancy change and separate the direct and indirect effects of agricultural expansion on species of conservation concern. Forest cover was the main driver of anteater occupancy, and consequently, the widespread deforestation driven by agricultural expansion in the Chaco has led to drastic declines in anteater occupancy, with only a relatively small proportion of the study area characterized The maps of changes in occupancy show four categories: a increase/stable where occupancy has increased (between 6 and 35%) or stayed stable (change B 6%), b low decrease (between 6 and 19%), c medium decrease (between 19 and 46%), and d high decrease (between 46 and 84%). Permanent water bodies are depicted as light gray by stable or increasing occupancy in 2015. Most areas turning from stable to declining occupancy did so after 2000, when landscape transformation increased. Importantly, anteater occupancy declined over much wider areas than those directly affected by forest loss, and occupancy declines extended far into seemingly untransformed habitat. This suggests that agricultural expansion has substantial and widespread indirect effects on species of conservation concern -beyond the spatial footprint of agricultural expansion itself. Stopping the decline of anteaters and other large mammals in the Chaco will thus require swift conservation action and planning, and our study can help to target both proactive and reactive management interventions to tackle ongoing defaunation. More broadly, we show how relatively simple and broadly applicable tools can help to provide a spatial template for conservation planning to mitigate the direct and indirect impacts of agricultural expansion. Anteater occupancy was primarily affected by landscape variables, especially forest cover, and occupancy increased with forest cover. This aligns with prior research emphasizing the value of forest for anteaters (Di Blanco et al. 2015), and is reasonable given that forests provide refuge from predation and protect from extreme temperatures. Edges between forest and open habitats were also positively related to anteater occupancy, however, edges had a lower effect on occupancy than did forest. One possible explanation is that anteaters usually rest in forested areas, and use open habitats for feeding and for moving between forest patches (Desbiez and Medri 2010; Medri and Mourão 2005).
We found that anteaters avoided areas close to human settlements, but selected areas close to roads, a finding in line with other studies (Di Blanco et al. 2015;Quiroga et al. 2016;Vynne et al. 2011). Although vehicle collisions have been proposed as a main cause of anteater mortality (Miranda et al. 2014), anteaters are known to use roads to move Fig. 4 Decrease in anteater occupancy (low, medium, and high) versus agricultural expansion (forest conversion to cropland and to pastures) through agricultural landscapes (Vynne et al. 2011), and our results suggest that they do not perceive roads as a threat. Additionally, anteaters in the Chaco are rarely hunted by people (Altrichter 2006;Camino et al. 2018) and this may explain why roads, a proxy for accessibility for hunters, did not negatively affect anteater occupancy in our case. The avoidance of human settlements can be explained by the widespread presence of dogs, which can attack and kill anteaters (Lacerda et al. 2009). However, the effects of presence of dogs and humans on anteaters in the Argentine Dry Chaco were still less important than the effects of forest cover (Table S3).
Our predictions showed that much of the Argentine Dry Chaco still has a high probability of anteater occupancy, especially in the north. In contrast, the south has less potential habitat, and the species is locally extinct there. The high occupancy estimates in the north of our study area can be explained by a combination of three factors: relatively large remaining forest patches, low hunting pressure on anteaters (Altrichter 2006;Camino et al. 2018), and the absence or extremely low density of jaguars, the anteaters' main predator (Astete et al. 2008;McBride et al. 2010;Sollmann et al. 2013), in the study area (Quiroga et al. 2014;Romero-Muñoz et al. 2019).
Although our models predicted relatively high occupancy in areas with high forest cover, our occupancy change maps first and foremost showed a drastic decrease in occupancy due to forest loss driven by agricultural expansion. In areas where agriculture expanded, occupancy decreased by up to 84% between 1985 and 2015, especially after 2000, when deforestation rates in the Chaco soared due to technological innovations in agriculture (e.g., introduction of genetically modified soybean), rising global soybean prices, and government incentives to foster agricultural expansion Fehlenberg et al. 2017;Goldfarb and Zoomers 2013). Areas with increasing or stable occupancy were relatively widespread during the first period, but became scarce after 2000 and are now limited to remote areas with high forest cover. Few such areas remain in the Chaco. Some of these consist of larger protected areas (Fig. S10), yet others remain beyond the current agricultural frontier. Factors explaining why agriculture has not expanded in these areas include the presence of indigenous communities, the national zoning plan that prohibits deforestation in some areas, and low market access from such regions (Piquer-Rodríguez et al. 2018).
In the Chaco, anteater occupancy decreased over an area almost twice as large as the area affected by forest loss (i.e., anteater occupancy decreased across and area of * 102,000 km 2 and forest loss occurred on an area of * 64,000 km 2 ). This suggests that agricultural expansion not only directly affects anteater occupancy via habitat conversion, but also has strong indirect effects that contribute to population declines beyond the footprint of forest conversion. For example, agricultural expansion typically leads to improvements of infrastructure and therefore, an increase in traffic, which in turn increases road mortality (Cáceres et al. 2010). Similarly, where agriculture expands, an inflow of people can occur, and pets and humans have easier access to previously remote areas. Additionally, agricultural expansion and associated human accessibility is increasing fire frequency in the Chaco (Argañaraz et al. 2015), and such fires are considered to be a key threat for the species (Miranda et al. 2014). Finally, forest loss alters local and regional climate (Brook et al. 2008), leading to drier and hotter conditions, which in turn amplify fire risk (Argañaraz et al. 2015). These potentially strong indirect effects that exceed far beyond the actual agricultural footprint are in line with findings from other regions though few studies have assessed this. For instance, extinction rates reported from West Africa were higher than those predicted by species-area models, likely due to hunting pressure along reserve edges (Brashares et al. 2001). Similarly, accessibility and hunting pressure increased due to forest fragmentation in tropical Asian forests (Corlett 2007). Our study provides further empirical evidence for the strong synergistic effects of land-use change and other extinction drivers.
Land-cover/use variables were the only time-variant predictors included in our models. This allowed us to reveal that the indirect impact of agricultural expansion on occupancy can extend far into seemingly untransformed habitats. Three alternative explanations for occupancy declines in untransformed habitat are possible, but not plausible for anteaters in the Chaco: (1) other environmental changes, (2) hunting pressure, and (3) time-delayed effects. Regarding other environmental changes, climate change might impact anteater occupancy. However, the Argentine Dry Chaco has neither experienced marked changes in climate during our study period (Fig. S11), nor is the Chaco a climatically marginal area for anteaters, as anteaters occupy both drier and wetter areas today (Miranda et al. 2014). Hunting pressure has increased in the Chaco over the last decades (Romero-Muñoz et al. 2019), likely due to improved road infrastructure, population increase, and the decline of forested areas. However, as mentioned above, anteaters are not targeted by hunters. We are not aware of any variable not related to, and driven by, agricultural expansion that might plausibly explain the broad-scale declines in anteater occupancy.
We used a large camera-trap dataset and our occupancy framework resulted in robust models and plausible occupancy patterns. Still, some uncertainties remain. First, we used a time-calibrated model that ensures observed changes in occupancy are solely due to changes in predictor variables, and that makes full use of all available camera-trap data (Nogués-Bravo 2009). Estimating occupancy based on past anteater data would have been interesting, but no dataset that would allow this analysis exists. Second, some predictor variables did not have a strong effect on anteater occupancy, and there was considerable uncertainty surrounding the effect of these predictors (Table 3). Yet, the goodness-of-fit test indicated that our models fitted the field data adequately, and the estimated errors of our spatial predictions were small, together suggesting that our inferences were robust. This suggests robust model projections, which form the basis for our separation of direct and indirect effects of agricultural expansion (rather than directly interpreting model coefficients). Third, we predicted occupancy over a large geographic extent and back in time. These predictions were made to areas within the range of environments sampled by our camera-trap sites, and thus we did not extrapolate ( Fig. S4 and Fig. S5), but we cannot fully exclude non-stationarity in the relationship between anteater occupancy and variables. Fourth, we interpreted our occupancy estimates as a proxy for population change. A high correlation between occupancy and abundance estimates has been previously demonstrated for several species (Clare et al. 2015;Linden et al. 2017;Parsons et al. 2017), however, this correlation remains untested for anteaters in the Chaco. Finally, detection probability differed among the three camera-trap surveys. We controlled for this by (1) including the variable survey in our detection model, and (2) relating each observation to the site-specific environmental variables from the sampling period. Still, a more homogenized survey design (e.g., same camera-trap brand) would have been preferable.
Overall, our results provide further evidence of the impacts of agricultural expansion on biodiversity. The fact that anteaters, a species which can tolerate certain levels of anthropogenic disturbance (Quiroga et al. 2016), may be greatly affected by forest loss, suggests that the consequences of agricultural expansion may be even more drastic for strictly forest-dependent species. Identifying the effects of forest conversion on biodiversity is therefore crucial for conserving anteaters and many other species. Our results showed that besides the direct effects of agricultural expansion, indirect effects are widespread. Identifying these often-overlooked indirect effects is therefore crucial but also challenging, as they can go unnoticed in forests that seem untransformed. In such forests, however, other threats such as hunting and forest fires can exert a major pressure on wildlife (Benítez-López et al. 2017;Brook et al. 2008;Redford 1992). As the footprint of agriculture expands in many tropical areas, paying closer attention to the additional and larger footprint of indirect negative impacts on wildlife becomes crucial.
Methodologically, we demonstrated how projecting time-calibrated occupancy models over time may help to identify the indirect effects of agricultural expansion, and to separate areas affected by direct versus indirect effects. Such a spatially-explicit approach can be used to predict occupancy declines, and thus can potentially reveal population declines before critical thresholds are reached. More generally, understanding and mapping occupancy declines due to direct and indirect effects of agricultural expansion provides starting points for conservation planning and for targeting conservation action. Proactive conservation strategies (e.g., protecting remaining core habitat) should be targeted in areas where occupancy increases or remains stable, and reactive strategies (e.g., restoring sink habitat) in areas where declines are occurring due to direct forest conversion. Additionally, to counteract the indirect effects of agricultural expansion, conservation actions should focus on implementing measures to lessen such indirect effects, for example through fire regulation, especially in areas that seem untransformed but where species' declines are observed. regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.