Inferring Species Interactions from Long-Term Monitoring Programs: Carnivores in a Protected Area from Southern Patagonia

Protected areas recently created in Argentina often include previously degraded lands, such as sheep ranches in the Patagonian deserts. We show the results of a 14-year monitoring program of three formerly persecuted carnivores, the culpeo fox (Lycalopex culpaeus), the South American grey fox (Lycalopex griseus) and the puma (Puma concolor), in two abandoned sheep ranches that were incorporated into a Patagonian national park approximately 25 years ago. The culpeo fox population underwent an average annual decline of 10–23%, whereas the grey fox and puma populations increased at an average annual rate of 7% and 19%, respectively. The grey fox’s increasing trends were strongly correlated with the decline of the culpeo fox, whereas the correlations between the fox and puma trends were weaker. Culpeo fox decline was stronger in the ranch where sheep and predator controls had been removed earlier. These relationships between species trends support the competitive release hypothesis, assuming that puma competition with the culpeo fox for trophic resources is stronger than competition with the grey fox, and that the puma can exclude culpeo foxes through interference. Species trends suggest a competitive hierarchy between fox species, with grey fox being the inferior competitor. However, mechanisms other than competition could not be discounted. Our study illustrates how long-term monitoring of interacting species allows a better understanding of ecological processes and wildlife ecology.


Introduction
Few areas across the planet remain where natural systems have not been damaged, directly or indirectly, by human activities, which in many cases pose a threat for biodiversity conservation [1][2][3]. The establishment of protected areas is a common conservation strategy to reduce the intensity of prey on lambs and young sheep [22,31], intensive culpeo fox control has been a common activity in sheep ranches [25], reaching annual reductions of over 75% in some local populations [32]. Although the grey fox does not pose a risk to livestock, it has also been culled, either for the value of its fur or indirectly due to the extensive use of nonselective predator control methods [25]. Patagonian sheep stock peaked in 1952 at 25 million head and decreased continuously thereafter, as a consequence of overgrazing effects, desertification and wool price devaluation [18,33]. During the second half of the 20th century, livestock activity was completely abandoned in large regions of Patagonia [34], which has allowed the recovery of some native species such as the puma and the guanaco [23,35]. However, predator persecution of foxes and pumas has persisted in some areas [29,36]. Accordingly, the persecution exerted on carnivores over time might have altered the competitive interactions between these species, as well as their relative abundances.
Changes in the relative abundance of sympatric carnivores can have far-reaching ecological consequences, including the precipitation of trophic cascades and species declines [37,38]. Carnivores play an important role in preserving ecological function and dynamics of ecosystems, driving top-down effects such as controlling populations of native or introduced species [39,40]. Larger predators may regulate populations of smaller predators through exploitative competition upon shared resources, or through interference competition by harassment and intraguild predation [41]. One consequence of losing larger predators is an increase of mesopredator abundance [38,42,43]. As a rule, larger species tend to dominate competitive interactions over smaller species [44]. The competitive release hypothesis states that a reduction in the abundance of a dominant competitor triggers an increase in the numbers of subordinate competitors [45,46], whereas the recovery of the superior competitor would produce the opposite effect.
The puma plays a dominant role over both fox species, including high potential to kill them [47]. Going down in the competitive hierarchy, a dominant role of the larger culpeo fox (weighing~8 kg) over the smaller grey fox (weighing~2.6 kg) [48] would be expected. Although there is evidence of the recovery of the puma throughout many regions of Patagonia after the reduction of predator control associated with the abandonment of sheep ranches [35,36], it is unknown how puma resurgence has affected other carnivores such as the culpeo fox and the grey fox and their competitive interactions. The diets of these three species overlap. Pumas mainly prey on larger prey such as guanacos, but European hares are an important diet item [24,49]. Grey foxes mainly feed on small mammals and invertebrates, although European hares are also included in their diet [48,50]. Culpeo foxes feed mainly on European hares and small mammals [48,50].
Twenty years ago, our research team started a long-term monitoring programme of the populations of the two fox species present in Patagonia [30], whose main objective was to document population trends in order to guide conservation management [51]. This monitoring programme was carried out in a protected area which had recently incorporated abandoned sheep ranches. Adopting an adaptive monitoring framework [11], in 2002, we included the puma as a monitored species in order to better understand changes in fox populations. In this work, we present population trends of the three species and explore whether the observed trends support the predictions of dominance hierarchy and the competitor release effects. We also highlight the relevance of long-term monitoring programs in recently created protected areas for a better understanding of ecological processes triggered after the cessation of human activities.

Study Area
The study was conducted in the Bosques Petrificados de Jaramillo National Park (hereafter denoted as BPJNP) located in the central north of the Santa Cruz province, Argentine Patagonia (47.66 • S 65.89 • W; Figure 1). Currently, the BPJNP extends over 78,543 ha of Patagonian shrub-steppe, an ecoregion of special interest for global conservation [52]. In 1997, two neighbouring sheep ranches were incorporated into the national park, increasing its size by about 30,000 ha. Prior to its incorporation, Cerro Horqueta ranch (CH) had~1000 sheep, while El Cuadro ranch (EC) had~500, but the livestock load for both Diversity 2020, 12, 319 4 of 17 fields together exceeded 5000 sheep in the 1950s. Livestock raising and predator persecution went on at least until 1999 in EC, while these activities ceased in CH around 1992. In both ranches, the vegetation is dominated by tussock grasses and low, dome-shaped spiny shrubs [53], whose cover ranges from less than 10% in the most arid sites to 60%, suggesting a marked degradation that could be attributed to many decades of overgrazing, as observed in other Patagonian steppes [54,55]. Summer temperatures average 17 • C, and winters are relatively mild with an average of~12 freezing days/year. Annual rainfall ranges between 100 and 300 mm, and snow is rare. denoted as BPJNP) located in the central north of the Santa Cruz province, Argentine Patagonia (47.66° S 65.89° W; Figure 1). Currently, the BPJNP extends over 78,543 ha of Patagonian shrubsteppe, an ecoregion of special interest for global conservation [52]. In 1997, two neighbouring sheep ranches were incorporated into the national park, increasing its size by about 30,000 ha. Prior to its incorporation, Cerro Horqueta ranch (CH) had ~1000 sheep, while El Cuadro ranch (EC) had ~500, but the livestock load for both fields together exceeded 5000 sheep in the 1950s. Livestock raising and predator persecution went on at least until 1999 in EC, while these activities ceased in CH around 1992. In both ranches, the vegetation is dominated by tussock grasses and low, dome-shaped spiny shrubs [53], whose cover ranges from less than 10% in the most arid sites to 60%, suggesting a marked degradation that could be attributed to many decades of overgrazing, as observed in other Patagonian steppes [54,55]. Summer temperatures average 17 °C, and winters are relatively mild with an average of ~12 freezing days/year. Annual rainfall ranges between 100 and 300 mm, and snow is rare.

Field Surveys
The monitoring program was designed and launched in 2000 [30,51]. To assess the distribution and relative abundance of foxes, we used bait stations, a noninvasive method to survey carnivores at large scales [56]. A bait station consisted of a 1-m-diameter circular surface of sifted, smoothed earth, with a 30 g bait buried at its centre [57]. Bait composition consisted of a mixture of minced meat, hydrogenous oil, cornstarch and the commercial lure Cat Passion (O´Gorman Enterprise Inc., Broadus, Montana, USA) [58].

Field Surveys
The monitoring program was designed and launched in 2000 [30,51]. To assess the distribution and relative abundance of foxes, we used bait stations, a noninvasive method to survey carnivores at large scales [56]. A bait station consisted of a 1-m-diameter circular surface of sifted, smoothed earth, with a 30 g bait buried at its centre [57]. Bait composition consisted of a mixture of minced meat, hydrogenous oil, cornstarch and the commercial lure Cat Passion (O´Gorman Enterprise Inc., Broadus, Montana, USA) [58].
We set 8 lines of 6 stations in each abandoned sheep ranch. Stations were set on alternate sides along secondary, unpaved roads. Separation of adjacent stations within each line was 0.5 km, and lines were at least 1 km apart. We assumed that the probability of any individual fox visiting more than one line per survey approached zero. Therefore, we considered each line as an independent sampling unit. Lines were surveyed twice per year, at six-month intervals, for 14 years (2000-2013) in autumn and 13 years (2000-2012) in spring. For each survey, we operated bait stations during three consecutive nights. Each morning, we recorded whether fox tracks were printed on the smoothed soil. Fox tracks were easily distinguishable between the species, since grey fox tracks are narrower and much shorter than culpeo fox tracks [58]. Tracks for which accurate species identification was not possible due to Diversity 2020, 12, 319 5 of 17 strong winds and/or poor printing were excluded from the analyses. In order to maximise accurate track identification, all surveys were conducted in dry weather.
Bait stations were designed to monitor fox populations, and they did not perform well in detecting puma trends in our study area (unpublished data). Therefore, between 2002 and 2009, we surveyed the puma population using fresh signs (scats, tracks and scrape marks) detected along foot transects as indices of relative abundance [59,60]. The study area was divided into 134 4-km 2 square cells (72 in CH and 62 in EC; Figure 1). With the aid of a GPS (Garmin 3Plus; Lenexa, Kansas, USA), a mean distance (± SD) of 3.01 ± 0.71 km was erratically travelled in each cell, starting at one of the cell edges and passing through the cell's centre. All fresh signs of pumas within~2 m of each side of the transect were recorded. Transects were mainly performed by the same observer (95%; D. Procopio), who was occasionally supported by A. Travaini and J.I. Zanón. Logistical difficulties from financial constraints and a lack of human resources prevented us from sampling all cells each year. Despite this, on average, more than 50% of cells from each ranch were visited annually (mean% ± SE for CH: 58% ± 8%; for EC: 60% ± 7%). Surveys were carried out in autumn, except for the first two years, when they were done in spring. By pooling spring and autumn data, we assumed that season had a small effect on the index of puma abundance.

Analyses
Visitation indices denote the possible ways of converting field data into a quantitative indicator of species activity at a bait station. Different visitation indices are built using different combinations of daily records throughout the three-day operation period. Travaini et al. [30] used power analysis to examine the effect of choosing each of seven possible visitation indices on the ability of the monitoring program to detect changes in fox abundance. They found that the proportion of zeroes (no visits) had a large effect and that recording fox visits within 72 h after activation showed the lowest variation and improved the expected power to detect population trends of both fox species monitored simultaneously. We used this visitation index, defined as the number of bait stations visited by each fox species per line during any of the three consecutive nights [30]. Thus, our index of relative abundance could take integer values between zero and six visited stations per line.
Population trends were estimated using the R package rtrim [61], a reimplementation of the original TRIM software [62]. TRIM was specifically developed for the analysis of wildlife monitoring data and has been widely used in the study of population trends of several vertebrate taxa, e.g., [63,64]. TRIM allows dealing with common issues in long-term time series: missing counts, overdispersion and serial correlation. Missing counts are replaced with predicted counts estimated by log-linear Poisson regressions models, allowing for analysing data sets with up to 50% of missing counts [62,65]. In our data set, missing values were between 6% and 12% for foxes, and 41% for puma. TRIM includes a correction parameter when overdispersion occurs and also controls for serial correlation among time series counts by using a Generalised Estimating Equations (GEE) approach [62]. TRIM computes two types of indices: (i) model indices, based entirely on the model; and (ii) imputed indices based on the observations plus, for missing counts, estimated values based on the model. These indices represent change between years in terms of relative variation in total population size. Imputed indices are employed to estimate a mean annual change rate, since they show a more realistic course in time and a trend category is assigned according to its statistical significance and the estimated magnitude of change [62].
For each fox species and season, we fitted a linear trend model for the visitation index with all years as change points (time-effect model), using the complete data set. We set 2000 as the base year, so the index for 2000 took the value 1. We used a stepwise procedure to identify significant changes in the slope using Wald tests to remove nonsignificant change points from the model (default significance threshold: 0.2) according to the parsimony principle [64]. This procedure provides a high confidence on true variations in the studied populations [63]. In order to check this, we built a model where all change points were included as mandatory, and we compared it with the model built under a stepwise Diversity 2020, 12, 319 6 of 17 procedure using the Akaike Information Criterion (AIC) and considering differences between models in ∆AIC > 2 to reject the less parsimonious model [66]. The same procedure was used to evaluate the puma population trend for the 2002-2009 period, using the number of fresh signs as count data and setting 2002 as the base year. All models were run with overdispersion and serial correlation taken into account (default TRIM thresholds: >3.0 and >0.4, respectively). We estimated overall trends using the multiplicative trend slope (β), which reflects the annual average percent change [62]. When β < 1, the population decreases; when β = 1, the population remains stable; and when β > 1, the population increases. In addition, we used β to estimate the overall population size change (PSC) for the whole time series using the formula: (PSC) = (β (n − 1) − 1) × 100, where n is the total number of years for the whole monitoring period, and change is expressed as a percentage [67].
Seasonal effects on population trends of foxes were tested using the goodness of fit of models and 95% confidence intervals for the PSC. Goodness of fit of models was tested by the likelihood ratio (LR) and chi-squared tests. Nonsignificant values (i.e., p > 0.05) for one or both tests indicate that the data fitted a Poisson distribution [62].
We also evaluated how differences between sheep ranches in the cessation time of livestock activities (including predator control) affected the slopes of both fox species and puma overall trends, assuming a similar management model in both ranches. To do this, we included ranch identity as a covariate in a simple linear trend model for each species and season. The data structure prevented the fit of a time-effect model with a stepwise procedure [62]. Differences on trend slopes were evaluated by the Wald test, using a significant level value of p < 0.05 [62].
We used simple regression models to evaluate how species population trends were related to each other according to the hypothesis of size-based dominance hierarchy [37,44]. We used the trend index of the putative dominant species at year "t" as the explanatory variable for the observed trend index for the putative subordinate species at year "t + 1", following a procedure to assess top-down regulation relationships [37,68].
Trend graphs were built using imputed trend indices, which were more realistic and more easily visualised than model indices [62]. Analyses were developed using "rtrim", "stats" and "graphics" R packages [69].

Results
Stepwise procedure models (SWM) were the most parsimonious for every species (Supplementary  Material Table S1). All models fitted a log-linear distribution and were not affected by overdispersion and serial correlation, except the spring model for the culpeo fox, for which both the chi-squared and the likelihood ratio tests indicated a lower model fit (Supplementary Material Table S1), which, however, should not affect the reliability of this trend [62].
According to the best model, the overall population trend for the culpeo fox was significantly decreasing in both seasons ( Figure 2; Table 1), with an estimated average annual population decline of 23% in autumn and 10% in spring, which were classified by TRIM as a "steep decline" and "moderate decline", respectively (Table 1). In terms of population size, this means a population reduction of 97% (95% CI: 92-99%) in autumn and of 75% (95% CI: 42-89%) in spring at the end of the study period. These models included five significant change points, both for autumn and spring, which were similar between seasons (Supplementary Material Table S2; Figure 2). Diversity 2020, 12, 319 7 of 17 means a population reduction of 97% (95% CI: 92-99%) in autumn and of 75% (95% CI: 42-89%) in spring at the end of the study period. These models included five significant change points, both for autumn and spring, which were similar between seasons (Supplementary Material Table S2; Figure 2).   On the contrary, the overall population trend for grey foxes was significantly positive and similar between seasons, with an estimated average annual population increase of 7%, which was classified by TRIM as a "moderate increase" (Table 1; Figure 2). The grey fox population size doubled during the study period (a 125% increase; 95% CI: 80-180% in autumn and 44-247% in spring). These models included four significant change points, both for autumn and for spring, which were similar between seasons (Supplementary Material Table S2; Figure 2).
The puma population trend significantly increased for the studied time period (2002-2009), with an estimated annual population increase of 19%, classified by TRIM as a "steep increase" (Table 1; Figure 3). This means a tripling of its population size in the eight years of monitoring, with an estimate increase of 238% (95% CI: 43-625%). The model included three significant change points with two periods of significant increase and one with a significant decrease (Supplementary Material Table S2; Figure 3). * (p < 0.05); *** (p < 0.001).
On the contrary, the overall population trend for grey foxes was significantly positive and similar between seasons, with an estimated average annual population increase of 7%, which was classified by TRIM as a "moderate increase" (Table 1; Figure 2). The grey fox population size doubled during the study period (a 125% increase; 95% CI: 80-180% in autumn and 44-247% in spring). These models included four significant change points, both for autumn and for spring, which were similar between seasons (Supplementary Material Table S2; Figure 2).
The puma population trend significantly increased for the studied time period (2002-2009), with an estimated annual population increase of 19%, classified by TRIM as a "steep increase" (Table 1; Figure 3). This means a tripling of its population size in the eight years of monitoring, with an estimate increase of 238% (95% CI: 43-625%). The model included three significant change points with two periods of significant increase and one with a significant decrease (Supplementary Material Table S2; Figure 3).  Observed trends were similar between abandoned ranches, following the same pattern of overall trends previously described, i.e., decreasing trends for the culpeo fox and increasing trends for the grey fox and puma. Despite this, trend slopes varied significantly between ranches for the culpeo fox in autumn (Wald test = 13.2, df = 1, p < 0.001), showing a decrease rate for CH (model slope ± SE: −0.49 ± 0.09; average annual population decline: 39%) higher than that of EC (model slope ± SE: −0.16 ± 0.08; average annual population decline: 15%; Figure 4). On the contrary, no significant differences were found during spring (Wald test = 0.059, df = 1, p > 0.1) (Supplementary Material Figure S1). Similarly, no significant differences were found for the grey fox slopes between ranches both in autumn (Wald test = 0.02, df = 1, p > 0.1; CH: 0.06 ± 0.01; EC: 0.06 ± 0.02) and spring (Wald test = 0.24, df = 1, p > 0.1; CH: 0.05 ± 0.02; EC: 0.07 ± 0.03), showing for both ranches similar average annual population increases of~7% (Supplementary Material Figure S1). Finally, no differences were found for puma populations at the ranch level (Wald test = 9.39, df = 1, p > 0.1), showing a similar slope (CH: 0.12 ± 0.08; EC: 0.13 ± 0.08) and similar average annual population increases of~13% (Supplementary Material Figure S1).

Discussion
Our results for foxes highlight the importance of designing monitoring programs on the basis of pilot studies and power analysis for the early detection of significant changes in population size [9,10]. Well-designed monitoring programs may be able to detect population change over a short time, but others may require a decade or two to show, for example, a 25% population change with 80% power at α = 0.10 [9]. Our monitoring program was designed to detect significant changes in the trends of both fox species in the first 5 years of monitoring (50% decrease and even a 30% increase in three years) with an estimated power of ~76% [30]. Thus, our monitoring program was effective because it was able to detect quick significant

Discussion
Our results for foxes highlight the importance of designing monitoring programs on the basis of pilot studies and power analysis for the early detection of significant changes in population size [9,10]. Well-designed monitoring programs may be able to detect population change over a short time, but others may require a decade or two to show, for example, a 25% population change with 80% power at α = 0.10 [9]. Our monitoring program was designed to detect significant changes in the trends of both fox species in the first 5 years of monitoring (50% decrease and even a 30% increase in three years) with an estimated power of~76% [30]. Thus, our monitoring program was effective because it was able to detect quick significant changes in population trends for both fox species during the first 5 years, in autumn as well as in spring. Further, the sensitivity of the monitoring program was greater than expected in detecting both significant increases and decreases [30,51].
The estimated trends were concurrent between seasons for both fox species, showing a significant decreasing trend for the culpeo fox population and a significant increasing trend for the grey fox population. The similarity between seasonal trends was maximal for the grey fox, obtaining the same mean annual change rates (7% of annual increase). Conversely, the mean annual decrease rate estimated for the culpeo fox in autumn was twice the corresponding spring value (a 23% vs 10% annual decrease, respectively). Overall, autumn trends produced more accurate estimates with narrower confidence intervals for population change, thus reducing uncertainty. Red fox (Vulpes vulpes) populations reach maximum abundance in autumn during juvenile dispersal [70,71], increasing the probability of detection and improving the accuracy of population parameter estimates during this season as compared with other seasons, including spring [72,73]. Patagonian foxes have a similar phenology and annual peaks of abundance usually occur in autumn, during juvenile dispersal [74], resulting in higher visitation rates with lower variability and more precise estimates. On the contrary, minimum population abundances are found in spring, after winter mortality, coupled with a reduction in fox activity during the breeding season [74]. Using annual troughs in fox abundance might have increased variability in visitation rate as well as uncertainty in estimates of population change. This effect was especially pronounced for the culpeo fox population, which was apparently involved in a drastic population decline. Decreasing numbers and reduced activity could explain seasonal differences in the estimated magnitude and rate of population collapse and the lack of fit in spring models [75]. Conversely, the higher abundance of the grey fox during the study period as a result of steady population increase might have diluted seasonal differences in population trend. Likewise, the sensitivity for detecting population increases is usually greater than the sensitivity for detecting declines when the sampling effort is kept constant [76].
Once fox trends were detected, the next step was inferring possible processes causing them. Population trends of culpeo and grey foxes showed a specular antagonistic pattern, where the continued increase of the grey fox population size was significantly related to the decline of the culpeo fox. Culpeo and grey foxes are opportunistic predators with similar ecological requirements and compete for the same resources [77]. With a larger body size (more accentuated in our study area; body size ratio > 3) [48], the culpeo fox could play a dominant role in competitive interactions similar to those described in other carnivore guilds [44]. Culpeo foxes may exclude grey foxes from high-quality habitats with abundant prey [78,79]. Under unfavourable conditions (e.g., reduced prey diversity and availability), competitive interactions could be more intense [50,77] and culpeo foxes could even kill grey foxes [47]. According to the competitor release effect, a reduction in the abundance of a dominant competitor allows subordinate competitors to take advantage of released resources, increasing their abundance [45,46]. An example of this ecological process has been shown with mesocarnivores in the UK, where the intense culling of the Eurasian badger Meles meles, a dominant competitor, was associated with an increase in densities of a subordinate competitor, the red fox [80]. Thus, the reduction in the population size of the dominant culpeo fox would have favoured grey fox numbers once competitive stress was relaxed.
In general, both fox species show considerable resilience to extraction from predator control, maintaining fairly stable populations over time [32,81,82]. Furthermore, when the intensity of predator persecution is reduced, culpeo foxes can respond by numerically increasing their density [81]. Hence, the deep decline observed in the culpeo fox population was an unexpected result in a protected area such as BPJNP where predator control was forbidden years before monitoring began and poaching rarely occurs. The first hypothesis that we considered to explain culpeo fox trends was the recovery of the puma population in BPJNP triggering top-down regulation processes.
On the one hand, the significant increase in puma population observed during our study suggests a population recovery similar to that described in other Patagonian regions where livestock activity and predator persecution has ceased [35]. In the absence of human persecution, large carnivores self-regulate their populations through social mechanisms and density dependence [83]. Accordingly, a significant increase in abundance would not be expected in a medium-high-density puma population, in which case most of the territories would be occupied [84,85]. Therefore, it is likely that pumas in BPJNP occurred at very low densities. Unlike foxes, predator persecution has decimated puma populations throughout Argentina, which were extirpated in many regions of Patagonia [28]. Thus, the status of the puma population at the creation of the BPJNP could have been at low density or even absent from some areas of the BPJNP, which could have favoured the release of subordinate competitors [38,42], such as the culpeo fox.
On the other hand, population trends of culpeo foxes (both autumn and spring) were negatively related to puma population trends, supporting the recovery of the regulation exerted by pumas on the culpeo fox population. Conversely, grey fox population trends were positively related to the increasing population trend of pumas, showing no regulation of pumas over grey foxes. Although these relationships were not statistically significant, the pattern observed in our study was similar to that described for the interaction between three canid species in North America [37]. This study shows how the recovery of the dominant predator, the grey wolf (Canis lupus) suppresses coyote (Canis latrans) populations, which in turn releases red foxes from regulation by coyotes. Levi and Willmers [37] point out how these interactions are determined by the trophic overlap between species, in turn related to their body size. Thus, wolves find a greater potential competitor in the larger coyote than in the smaller red fox. Similarly, pumas are more likely to focus their efforts on displacing or even killing culpeo foxes (we have direct evidences of that in the BPJNP; authors' unpublished data) than grey foxes because they might perceive culpeo foxes as stronger competitors due to their larger overlap on food resources. Trophic niche overlap between both species can be high. In a protected area located in the high mountains of central Argentina, where native prey increased after livestock cessation, pumas and culpeo foxes showed a high diet overlap which was maintained over the 6 years of study (Pianka's overlap index values range: 0.92-0.99) [86]. In southwestern Patagonia, pumas and culpeo foxes select the exotic European hare even in protected areas where native large-bodied guanacos occur [22,24]. In our study area, pumas mainly consume guanacos (60% of consumed biomass) and European hare (30% of consumed biomass) [49]. In contrast, European hare (36% of consumed biomass) and small rodents and armadillos (30% of consumed biomass) make up the bulk of the culpeo fox diet in BPJNP [50]. Despite the fact that the guanaco is the main prey of pumas, the European hare is an important prey for both species in BPJNP. So, it is expected that the increase in puma abundance intensified competitive interactions between both species.
Although our results would support these hypotheses, other factors not considered here could have been involved in the observed ecological processes. For example, our results partially support the hypothesis of a regulation process of the culpeo fox population by pumas. However, accepting this as the sole mechanism provoking the culpeo fox decline would imply that puma interference would be more efficient in depressing fox numbers than human persecution, which is questionable. Furthermore, marked changes in the availability of food resources advantageously exploited by each species could have benefited the grey fox or harmed the culpeo fox [48,50], explaining the observed trends in both fox species populations. On the other hand, our results in autumn support how the differences in time since the abandonment of livestock activity in ranches have significantly affected the population reduction rate of culpeo foxes. The regression was much more accentuated in the ranch where sheep were previously removed, i.e., CH, than in the ranch where about 500 sheep along with predator persecution remained in the area seven years longer, i.e., EC. Therefore, the removal of sheep from the BPJNP would be a limiting factor for the culpeo fox. Sheep have been described as an important food resource for culpeo foxes in active Patagonian ranches, both as live prey (mainly lambs and young sheep) and carrion [77]. Further, an important prey for the culpeo fox in Patagonia, the European hare [22,77], could have been harmed by the disappearance of sheep grazing that favours the proliferation of introduced lagomorphs by reducing the height of the grass and shrub cover [87,88]. However, we do not know whether management in both ranches was different, and how this could have influenced the results. Unfortunately, we do not have enough data to explain the dramatic population reduction observed for the culpeo fox in BPJNP, so we can only present our discussion as hypotheses that could have acted: puma competition, reduction in key food resources for culpeo foxes or a synergy between them. All these alternative explanations deserve further investigation to achieve a better understanding of the competitive interactions between these species.
Our long-term study highlights the importance of wildlife monitoring programs not only as a tool for an effective management and conservation of biodiversity [8], but also for proposing hypotheses about ecological processes and interactions [12,13,15]. The implementation of monitoring programs for carnivores in newly protected areas created on human-altered landscapes can be a key tool to assess the relevance of predators to restore degraded ecosystems and confer resilience against globally threatening processes [40]. Our long-term monitoring may provide insight into the poorly understood competitive interactions among the three most common carnivore species of Patagonian ecosystems, and further specific studies are required to assess the hypotheses proposed here.
Supplementary Materials: The following materials are available online at http://www.mdpi.com/1424-2818/12/ 9/319/s1; Figure S1: Non-significant ranch effect (EC: El Cuadro; CH: Cerro Horqueta) on population trends of culpeo fox, grey fox and puma for each season, Table S1: Trend models fitted for each species and season, and ordered according to their AIC values, Table S2: Model slope and standard error (SE) for change points selected in the stepwise modelling process for the three species and seasons.