Evaluating predator control using two non-invasive population metrics: a camera trap activity index and density estimation from scat genotyping

Study sites

Our broad study design was replicated across the Wimmera and Mallee regions of Victoria, Australia. Within each region, we conducted a before–after control–impact (BACI) experiment, with repeated surveys through time at one non-treatment and one treatment site (Fig. 1). Baiting had previously occurred at both treatment (or ‘impact’) sites prior to initial surveying. Treatment sites in our study were managed for foxes by Parks Victoria and the Mallee Catchment Management Authority after initial surveying, respectively, with nearby non-treatment sites left unmanaged as experimental controls.

Fig. 1.

Location of the four study sites in (a) the Wimmera region and (b) the Mallee region in Victoria, Australia, used for camera trap surveys and scat collection. Baiting was conducted at treatment sites (black transect), whereas the non-treatment sites (grey transect) remained unbaited. Approximate camera sites are spread within the areas represented by black stippling, close to transects; precise locations are not shown for security reasons.

In the Wimmera region, the paired sites comprised the non-treatment site, Nurcoung Flora Reserve (36°41′39.8″S, 141°42′42.0″E), and the treatment site approximately 10 km away, the Eastern block of Little Desert National Park (36°34′50.0″S, 141°48′55.0″E) (Fig. 1a). The area of native vegetation at the non-treatment site is 5.9 km², and the treatment site has an area of 476.8 km². Both sites are bordered by pastoral and arable agricultural land. The Wimmera bioregion experiences hot, dry summers and cool winters with an average rainfall of 414 mm, falling primarily from June to August (data for Nhill, the nearest meteorological station; Bureau of Meteorology 2020a). The landscape typically supports Lowan Sands Mallee (Ecological Vegetation Class 87, Eucalyptus spp, Leptospermum spp.) (Victorian Government Department of Sustainability and Environment 2005a).

The Mallee region paired sites comprised the non-treatment site at Wandown Flora and Fauna Reserve (34°48′24.9″S, 142°59′20.9″E), and the treatment site approximately 30 km away at Annuello Flora and Fauna Reserve (34°53′38.7″S, 142°33′33.7″E) (Fig. 1b). The area of native vegetation at the non-treatment site has an area of 25.3 km², and the treatment site an area of 361.0 km². Both sites are bordered by irrigated agriculture on the northern side, and pastoral and arable agricultural land elsewhere. The Mallee bioregion also has a hot climate, with dry summers and cool winters. Rainfall averages 328.3 mm, with a strong bias across August–October (Data for Ouyen, the nearest meteorological station; Bureau of Meteorology 2020b). Vegetation is dominated by Woorinen Sands Mallee (Ecological Vegetation Class 86, Eucalyptus spp, Triodia spp.) (Victorian Government Department of Sustainability and Environment 2005b).

Fox-baiting regime

The non-treatment sites had no history of fox management and were not managed during our study. Pre-existing management at treatment sites ceased a month prior to our study commencing in both regions: June 2019 in the Wimmera, and May 2019 in the Mallee (Fig. 2). We consider this window to be unbaited in our study design but note that fox densities may be lower in treatment sites as a result of prior management.

Fig. 2.

Timeline of study events at the (a) Wimmera and (b) Mallee region. Poison baits were laid (black) at the treatment sites. For each region, the months before the first baiting event were the unbaited period (light grey), and those after baiting were the baited period (dark grey). Note, baiting commenced at the Wimmera treatment site only on the last days of September, therefore we treated this month as unbaited.

Both treatment sites were ground baited (Animal Control Technologies, FOXOFF®) following our unbaited scat surveys. Although we consulted with the land managers, the timing, extent and intensity of baiting at treatment sites was largely out of our control and dependent on the objectives and constraints of the organisations.

At the Wimmera treatment site, the Eastern and Central blocks have been pulse baited since 2010. This involved three 9-week pulses per year with fortnightly bait replacement. The unbaited period of our study was from June 2019 to September 2019. An intensified baiting regime commenced in late September 2019, and our baited survey period followed in October 2019 (Fig. 2). While pulse baiting continued in the central block, baits were laid monthly at the eastern block. Across the two blocks, 345 bait stations were operated, positioned approximately 750 m apart (Supplementary Fig. S1).

The Mallee treatment site was first baited in 2019. Baits were deployed across two pulses, one in April and May. The unbaited period of our study was from June 2019 to November 2019. The baited period of our study commenced with pulse baiting in December 2019 and continued into the next baiting event in April 2020. Baits were laid at 223 stations, spaced approximately 260 m apart. Each baiting pulse lasted for 2 weeks, and any removed baits were replaced daily.

Camera traps

All four sites were monitored with cameras traps (Wimmera: Scout Guard SG560k and Ltl Acorn Ltl-5310; Mallee: Scout Guard SG560k) to estimate fox activity. The AMPE continuously operates 10 solar-powered cameras at all sites (N = 40; Fig. 1). These cameras have all been operating since 2015, except for the Mallee treatment site (deployed in 2019).

We obtained 12 months (June 2019–May 2020) of raw camera data from the Victorian Malleefowl Recovery Group for each site. We grouped these camera data into unbaited and baited periods, specific to each region (Fig. 2). Note, baiting commenced at the Wimmera treatment site only on the last day of September, therefore we treated this month as unbaited.

The camera setup followed the approach described by van Hespen (2015) and van Hespen et al. (2019). Cameras were 1 m above ground in an off-track, grid formation, and spaced 0.8 km apart (Fig. 1). They were programmed to capture one photograph with a 5-min interval between sequential triggers to minimise repeat captures, and were not fitted with a lure to avoid altering predator behaviour (Benshemesh et al. 2014). Although camera spacing was not independent relative to fox movements, the design captured the ‘background’ rate of fox activity (therefore predation opportunities) at sites and accounted for potential spatial and temporal heterogeneities in fox activity (van Hespen et al. 2019).

Scat collection

We established a 14-km scat collection transect at each site, with a minimum of 8-km displacement between transect start and end (Fig. 1). Transects were based along roads adjacent to the camera traps at each site. This placement ensured that we maximised the number of scats collected (because roads are regions of high fox activity in mixed-agricultural landscapes; Hradsky et al. 2017) and that we surveyed a similar area as the camera traps. Transect length was approximately double the length of average fox home range estimates (95% kernel) from the nearest analogous bioregion (Carter et al. 2012). This ensured that the sampling area was likely to be greater than the range of any individual, an important model assumption for spatially explicit capture–recapture analysis (Efford and Boulanger 2019). Where feasible, transects were continuous for logistical reasons; otherwise, we split transects to preferentially survey along connected native vegetation (Fig. 1a).

We surveyed each transect across two primary survey sessions per region, one in the unbaited period and another in the baited period (Fig. 2). We cleared each transect of fox scats 2 days prior to a survey session so that only fresh scats would be collected. We aimed to collect scats across four secondary survey occasions (survey days) per survey session, alternating between the non-treatment and treatment transects over eight consecutive days. This ensured that scats were less than 48 h old and reduced the likelihood of DNA degradation (Piggott 2004). In total, there were eight survey sessions: two regions × two sites × two periods.

The survey sessions in the unbaited periods were conducted in full in both regions: Wimmera (September 2019) and Mallee (November 2019). For the baited period in the Wimmera region, one survey occasion at the non-treatment site was lost due to a Total Fire Ban day (December 2019). In the Mallee, the first survey occasion in the baited period was delayed until 7 days after scats were cleared due to a Bushfire State of Emergency, and we were only able to complete two and three survey occasions at the non-treatment and treatment site, respectively (January 2020).

Transects were walked simultaneously by two observers to collect fox scats from the middle of the road and up to 2 m either side. At all times, one observer was experienced in scat identification and the other had been trained for at least 1 day prior. We identified the species based on visual and odour specifications (Triggs 2004) and associated signs (Moseby et al. 2012). One disposable glove was used per scat to avoid cross-contamination of DNA. Scats were placed into individual zip-lock bags and marked with date, site, unique identifier and collector according to Piggott and Taylor (2003). We preserved scats frozen at −18°C immediately upon returning from collection (Piggott and Taylor 2003).

Analytic approach

We used a model comparison approach to evaluate support for an effect of baiting on fox activity (derived from camera records) and fox density (derived from scat genotypes).

We hypothesised that three main effects could influence fox activity and density: treatment (non-treatment or treatment site); period (unbaited or baited); and region (Wimmera or Mallee). We derived a set of candidate models for each response variable (activity rate and density), based on single, additive and interactive combinations of these effects, and a null model. In each case, support for models that included an interaction effect between treatment and period would be considered evidence that baiting affected fox activity or density. The region covariate was included to represent differences in the landscape context and baiting regime between the regions.

We ranked candidate models using Akaike’s information criterion adjusted for small sample size (AICc) (Burnham and Anderson 2002). We regarded the ‘best’ model to be the one with the lowest AICc value but considered those within two AICc as equally plausible. To measure precision around our estimates, we calculated the percentage relative standard error (RSE) as the standard error of the coefficient divided by the estimate; an RSE > 20% was assumed to indicate low precision and, therefore, low statistical power to detect a change in the population size (Williams et al. 2002; Palmer et al. 2018; Efford and Boulanger 2019). All statistical analyses were performed in in R v4.2.1 (R Core Team 2022).

Activity analysis

We identified fox captures from camera traps using FastStone software (Benshemesh et al. 2018; FastStone 2020). To minimise non-independent captures of the same individual, we classed repeat captures within a 30-min time frame as a single capture event using the ‘CamTrapR’ package (Niedballa et al. 2016). We then extracted the total number of independent fox captures per site for each camera-month. We assigned months into unbaited and baited periods (as detailed above; Fig. 2).

We modelled fox activity as the mean fox capture rate (number of captures per month) using generalised linear mixed effect models (GLMMs). To account for excess zeros in monthly fox capture rates (Blasco-Moreno et al. 2019), we assessed the fit of four alternative distributions (Poisson, negative binomial, zero-inflated Poisson and zero-inflated negative binomial) using simulation for dispersion and zero-inflation with ‘DHARMa’ (Hartig 2020). P-values suggested a zero-inflated negative binomial provided the best fit to the data, given by:

where P_c,m is the number of captures taken by camera c in month m, η_c,m is the mean capture rate and α is the variance of the distribution.

We fitted models using maximum likelihood methods, using the ‘glmmTMB’ package (Magnusson et al. 2021). We derived a set of 12 candidate models, based on combinations of the fixed effects and a null model (Table 1). To decide whether we needed to account for resampling at the same camera trap points each month, we evaluated our most complex fixed-effects model (treatment × period × region), with or without inclusion of camera trap point identity (1–10) as a random effect. We then used the structure of the most strongly supported model to fit the 12 candidate fixed effects models. We assessed competitive models with a combination of model weights, 95% confidence intervals of interaction terms and R² values (calculated as per Nakagawa and Schielzeth (2013)).

DNA extraction and species identification

Scats were removed from the freezer prior to DNA extraction and allowed to thaw before swabbing the outer surface of the scat repeatedly with a Copan dry swab (155CIS, Copan). Swabs when then placed in lysis buffer (500 μL ATL and 40 μL Proteinase K) in a 1.5-mL tube, with the swab tip snapped so that the tube lid could close. The tube was vortexed thoroughly, incubated at 56°C for 2 h before the sample was centrifuged at 6000g to pellet any debris. We then removed 200 μL of supernatant for extraction following the DNeasy spin column protocol with final elution in 100 μL buffer AE. Initially, we screened extracted DNA to confirm presence of fox DNA using the vv-CytB fox-specific primers in a qPCR assay described in Berry et al. (2012) on a Roche LightCycler 480 II system. A TaqMan Exogenous Internal Positive Control (VIC probe) was run for each sample to test for the presence of PCR inhibitors. Samples that were inhibited were then diluted 1/10 and rescreened with the qPCR assay.

All samples that were deemed positive for fox DNA (including diluted samples) were then genotyped by microsatellite loci and a single sex-specific marker to identify individuals and their sex. Microsatellite loci genotyped were V142, V374, V402, V602 (Wandeler and Funk 2006), C2054 (Francisco et al. 1996), DB1, DB3 and DB4 (Holmes et al. 1995), and the sex-determining marker was vvSRY (Berry et al. 2007). All markers were tagged with a unique fluorescent label during PCR using the method outlined in Blacket et al. (2012). Reactions included 5-μL Qiagen multiplex mix (Qiagen, Chadstone, Victoria, Australia), 4 μL of primer mix (0.2 μM of each primer) and 2 μL of template DNA. PCR conditions consisted of an initial 15-min denaturing step at 94°C, followed by 40 cycles of 94°C for 30 s, 59°C for 90 s and 72°C for 60 s, with a final extension step of 60°C for 30 min. All PCRs were undertaken on an Eppendorf Mastercycler S gradient (Eppendorf). All samples were genotyped in duplicate to avoid issues with allelic dropout.

After PCR, genotyping was undertaken using an Applied Biosystems 3730 capillary analyser, with product lengths determined relative to a GS500LIZ_3730 size standard. GeneMapper ver. 4.0 (Applied Biosystems) was used to score allele sizes, with all profiles examined manually. Individual fox haplotypes (multi-locus genotypes) were then identified with GenAlEx v6.5 (Peakall and Smouse 2012), as in Le Pla et al. (2022). If replicates matched at a locus, we assumed this to be the correct genotype, whereas if replicate genotypes differed at a locus, we then manually re-examined each histogram output to determine if there was an error in scoring, and if still conflicting, we excluded the locus from analysis. Individuals with six or more microsatellite loci scored across both replicate samples were included in all subsequent analyses to limit errors when identifying individuals (Waits et al. 2001).

Spatially explicit capture–recapture analysis

We estimated fox density using spatially explicit capture–recapture (SECR) analysis (Borchers and Efford 2008). SECR models the probability of detecting an animal as a declining function of distance from its range centre at an unknown location, and distance to the nearest trap location (Royle and Young 2008). Resulting models produce estimates of population density unbiased by edge effects and incomplete detection (Borchers and Efford 2008).

We fitted multisession SECR models to the data from the eight survey sessions. The multisession formulation was adopted because data aggregation can improve the precision of parameter estimations by sharing detection information across sessions (Efford et al. 2009) and regions (Morin et al. 2018; Proffitt et al. 2020).

Each 14-km transect was discretised into 200-m segments to create conceptual traps with scats assigned to the nearest segment (Efford 2011; Lonsinger et al. 2018; Le Pla et al. 2022). Each trap was considered a proximity detector to allow multiple detections per individual per occasion (Efford 2011). We defined a buffer distance of 3002 m around each transect as outer limit of the area of integration, where foxes beyond this estimate were considered unlikely to be detected in the survey (Borchers and Efford 2008). This value was based on double the diameter of the average home range (95% kernel) of GPS-collared foxes in the nearest similar climate zone (Carter et al. 2012). The buffer distance exceeded the minimum suggested buffer size of 2664 m, based on an estimate of dispersion, the root-pooled spatial variance (Calhoun and Casby 1958; Efford et al. 2009).

We modelled density with the ‘secr’ package (Efford 2020) by fitting the full-likelihood model assuming a half-normal detection function:

where D is density, g₀ is the probability of detecting an animal at a trap located at its activity centre (per occasion) and σ is a home range spatial scale parameter (Efford et al. 2009).

As for the activity analysis, we had a set of 12 density models. Further, we considered that fox detectability might vary with period due to seasonal differences in behaviour, or with sex if females and males had differences in movement behaviour (e.g. male foxes generally disperse further than females) (Trewhella et al. 1988). Therefore, for each density model, we varied the detectability parameter (g₀) and spatial scale parameter (σ) as constant, or as a function of period or sex. This resulted in 22 candidate models (Table S1).

Fox activity

Across 11 and 12 survey months, we recorded 106 fox detections (27 unbaited, 79 baited) from 4503 trap nights in the Wimmera, and 483 detections (267 unbaited, 216 baited) from 6782 trap nights in the Mallee, respectively. This corresponded to a higher camera-trapping rate in the Mallee at 7.12 foxes per 100 trap nights, as compared with 2.35 foxes per 100 trap nights in the Wimmera. In the Wimmera region, foxes were detected at a total of 44% and 57% of the functioning camera traps at the non-treatment (n = 9) and treatment (n = 7) sites, respectively. In the Mallee region, foxes were detected at 100% of the functioning camera traps at the non-treatment (n = 9) and treatment (n = 10) sites.

A zero-inflated negative binomial distribution with a random effect for camera traps provided the best fit to the fox activity data. The top-ranked model was strongly supported over the other candidate models (ΔAICc to next model 8.16; model weight = 0.94; Table 1). The top model included a three-way interaction among treatment, period and region. However, this was primarily driven by a change in fox activity at the non-treatment site in the Mallee.

There was no evidence that fox baiting reduced fox activity in either region. In the Wimmera region, fox activity was low and similar at both sites across the unbaited and baited periods (Fig. 3a). In the Mallee region, fox activity was higher at the non-treatment than treatment site during the unbaited period; estimates overlapped during the baited period due to a large but uncertain decline in fox activity at the non-treatment site (Fig. 3b). The coefficient for the three-way interaction term was moderate and positive, but uncertain (1.19; 95% CI: −0.20–2.57). This uncertainty likely reflects the high level of variation in capture rates among camera traps within sites, which resulted in poor precision around most of the activity estimates (Table S2).

Fig. 3.

Predicted fox activity estimates (±95% CI) from the top-ranked model. Activity was derived from the monthly fox capture rate from 10 camera traps at each site. Estimates for non-treatment sites are light grey and treatment sites are black. Panels show estimates for the unbaited and baited periods in (a) the Wimmera and (b) the Mallee region.

There was no support for the simpler model that could also indicate an effect of baiting on fox activity (treatment × period); Table 1. The fixed effects for this model had little explanatory power (R²m = 0.04; Table 1).

Scat collection

In total, we collected 638 scats across the four sites. The number of scats collected per survey session varied between regions and periods. The greatest number of scats was detected from the Mallee region (Table 2). Scat numbers were lower in the baited periods, with approximately 60–70% fewer scats detected at all sites (Table 2), likely reflecting a seasonal difference.

Genotyping success and recapture rates

Of the 638 scats, 285 (45%) yielded fox DNA and 266 (42%) amplified enough microsatellite loci (n ≥ 6) to assign individual identities; 372 samples were unable to be genotyped (Table 2). We identified 108 individual foxes (Table 3). No individual foxes were detected moving between sites or regions, suggesting spatial independence of study sites.

In the Wimmera region, we detected 51 unique individuals; 67% of these were at the non-treatment site (Table 3). At the non-treatment site, no fox contributed to more than 8% of all detections (n = 78 scats; Table 3). However, at the treatment site, one fox accounted for 26% of all detections (n = 42 scats; Table 3). In the Mallee, we detected a similar number of individuals to the Wimmera (57), but the number of individuals was approximately evenly split between the non-treatment and treatment sites (Table 3). No fox accounted for more than 13% of detections at the non-treatment site (n = 67 scats: Table 3). At the treatment site, recaptures were heavily skewed towards a few individuals, with three foxes contributing 41% of all detections (n = 79 scats; Table 3).

Across all sites, many individuals detected during the baited periods had previously been detected in the unbaited periods. Recapture rates between the periods ranged between 50% and 82% (Table 3) and were similar or higher in the treatment sites than the non-treatment sites, suggesting that these individuals survived baiting, and that baiting did not increase population turnover (Table 3).

Spatially explicit capture–recapture density modelling

Nine models of fox density received substantial support (Table 4). Model uncertainty was high and individual model weight was low, partly due to uncertainty around the detectability parameters: four density models with an identical density model but with or without an effect of sex on the detectability parameters received very similar support (Table 4).

The top-ranked model contained an effect of sex on the detectability parameters, an interaction between treatment and period and an additive effect of region, indicating some evidence for an effect of baiting on fox density (Table 4; Fig. 4). However, uncertainty around the interaction term was wide (0.55; 95% CI: −0.17–1.26), and two other models without this interaction term received similar support (Table 4). The top-ranked model suggested that g₀ was lower for female foxes (0.05; 95% CI: 0.04–0.08) than male foxes (0.08; 95% CI: 0.06–0.12). However, the σ value was larger for female foxes (726.78 m; 95% CI: 606.52–870.88), as compared with male foxes (569.01 m; 95% CI: 481.51–672.41).

Fig. 4.

Predicted fox density estimates (±95% CI). Density was derived from scat genetic analysis and spatially explicit capture–recapture analysis. Estimates for non-treatment sites are light grey and treatment sites are black. Panels show estimates for the top-ranked model in the (a) Wimmera and (b) Mallee regions.

The top-ranked model suggested that fox densities were similar between regions (Fig. 4). Density estimates in the unbaited period ranged from 0.69 (95% CI: 0.47–1.0) to 1.06 (95% CI: 0.74–1.51) foxes km⁻² across regions. Precision was highest at each site and region in the unbaited period, with RSE ≤ 20% (Table S3). Fox densities in both regions declined more at treatment sites (35–40% decline) than non-treatment sites (12–14% decline) between the survey periods, although confidence intervals still overlapped pre-baiting estimates (Fig. 4).

Fox activity

Our top-ranked model of fox activity, derived from the AMPE camera trap data, provided no evidence that baiting decreased fox activity. Instead, fox activity differed greatly between the two regions. This variation was likely due to camera error and low sample size rather than a reflection of true difference. Three of the 20 camera traps in the Wimmera failed and did not record any photos. Others produced false zeros or were poorly aligned, leading to consistently low and imprecise fox capture rates.

Environmental factors may also have influenced fox behaviour and therefore detectability and activity. Fox activity was higher at the Mallee non-treatment site, perhaps because cameras at this smaller site tended to be closer to habitat edges and roads, where fox activity is known to concentrate (Šálek et al. 2010; Hradsky et al. 2017). Furthermore, fox breeding (winter) and cub births (spring) coincided with unbaited survey periods, and the baited survey periods included the dispersal season (autumn) (Pech et al. 1992). Seasonal changes in movement behaviour (such as females provisioning for cubs) may have also contributed to the higher initial fox activity at the Mallee non-treatment site. However, the expected increase in fox activity from cub independence and dispersal is not seen in the baited period. Although activity counts from camera traps are commonly employed to evaluate the outcomes of fox baiting (Towerton et al. 2011; van Hespen et al. 2019), these metrics may have limited value when population densities are low, or when environmental factors – or baiting – alter fox behaviour (Sollmann et al. 2013; Stephens et al. 2015). Our findings reiterate that under such circumstances, camera traps may not receive sufficient capture rates to detect population change.

Trapping rates on cameras in the Wimmera were similar to rates from two baited studies in similar bioregions: nearby in the Grampians, approximately 100 km to the south-east (2.35 foxes per 100 trap nights; Ramsey et al. 2015), and from semi-arid Western Australia (2.28 foxes per 100 trap nights; Thompson et al. 2019). However, the Mallee camera trapping rate was substantially smaller than the rate of approximately 24.67 foxes per 100 trap nights recorded at the same non-treatment site 4 years prior to our study (van Hespen et al. 2019). The discrepancy between these trapping rates is likely caused by differences in how consecutive captures on cameras were handled: we used a 30-min interval to class detections as ‘independent’, whereas van Hespen et al. (2019) used 5 min.

Spatially explicit capture–recapture

Our top model of fox density, derived from the scat genotyping data, provided weak evidence of an effect of baiting on fox density. Fox densities were similar across both regions prior to baiting, and there was some evidence that densities at the treatment sites declined more than at the non-treatment sites among the survey periods. It is important to note that fox densities at the treatment sites may have been impacted by prior baiting events. This may explain the tendency toward lower densities in treatment than non-treatment sites prior to baiting. The BACI design allowed us to quantify the effects of the next management pulse on the fox populations, while controlling for seasonal changes by monitoring changes in the non-impact landscapes. However, it was likely to be conservative in terms of quantifying an overall impact of baiting on unmanaged fox populations.

The identification of individual foxes using microsatellite genetic markers provided powerful information on fox survival rates in the baited period, revealing that more than half the foxes at each treatment site were unaffected by management actions. These foxes may not be encountering baits, not consuming baits or not consuming baits while the baits contain toxic levels of 1080. This suggests that managers need to explore options for improving population knockdown rates (rather than focus on reducing immigration). Such options could include modifying the baiting regime or method of bait delivery. For example, to increase the number of baits available to each fox, managers could increase bait density or the frequency of baiting (Saunders and McLeod 2007; Towerton et al. 2016). The likelihood of foxes consuming a toxic dose will decrease with bait avoidance and caching (Allsop et al. 2017). Trialling other bait matrices to improve palatability (Van Polanen Petel et al. 2001), trialling other poisons with greater efficacy (Southwell et al. 2011; Gentle et al. 2017) or combining baits with canid pest ejectors (Marks et al. 2002) may result in more rapid consumption of baits. Scat genetic analysis is evidently a valuable technique to refine baiting and monitoring regimes and for estimating fox densities.

Our study was able to establish the first empirical density estimates of fox populations for the semi-arid regions of Victoria. These estimates are broadly comparable with those from other semi-arid and arid grazing habitat in Australia (e.g. Marlow 1992; Marlow et al. 2000; Berry et al. 2012). They are higher than estimated densities from the Grampians, approximately 100 km south-east, which estimated 0.22 foxes km⁻² using presence–absence camera trap data of foxes without individual identification (Ramsey et al. 2015). Our density estimates are also substantially greater than those from the closest analogous bioregion (approximately 200 km north-west). Here, concurrent GPS collaring and camera trap spatially explicit mark–resight analysis provided estimates of 0.05–0.14 foxes km⁻² (Carter et al. 2019). In contrast, however, densities were substantially lower than those observed in Australian urban and temperate agricultural habitats (Coman et al. 1991; Thompson and Fleming 1994; Marks and Bloomfield 1999; Berghout 2001; Gentle 2005); this was expected because foxes have larger home ranges and occur at lower densities when productivity is low.

Comparison of metrics

Despite considerable uncertainty associated with all metrics, density estimates produced consistently lower RSEs than the activity indices, indicating that scat genetic analysis had a greater capacity to detect a change in fox populations, if one occurred. One potential reason for this is that our scat survey transects covered a larger area than the camera arrays and were therefore likely to produce less variation in sample size and more-precise estimates (Güthlin et al. 2014b).

Further, operating error and seasonal variation in fox behaviour led to low and variable fox capture rates on the camera traps. Activity indices can be skewed if detections are sparse or variable (Sollmann et al. 2013). Where capture rates are low, as exemplified in the Wimmera region, an activity index will be biased by detectability and thus have limited capacity to infer trends in density (Sollmann et al. 2013; Stephens et al. 2015). Nonetheless, our activity models indicated some promise for inferring population density: Mallee fox activity approximately corresponded with the suggested decline in density produced by SECR analysis. This suggests that camera trap activity indices may have the potential to parallel trends in density, given appropriate capture rates.

Survey design – camera activity index

Simple changes that optimised survey effort would likely improve the statistical power of the camera survey design. Although increasing the survey period would reduce uncertainty in our spatially replicated study, it is irrelevant for long-term monitoring programs such as the AMPE. Instead, survey effort would best be increased by ensuring that all cameras are functioning and not obscured. The detection rate of foxes at each camera could also be increased by placing cameras at areas of higher fox activity, such as on tracks (Mann et al. 2015; Carter et al. 2019; Geyle et al. 2020). However, camera placement on tracks has several limitations that need to be considered, including spatial variation in track networks and capture rates (e.g. Carter et al. 2019; Kolowski et al. 2021), and temporal variation in fox behaviour, such as dispersal patterns (Soulsbury et al. 2011). On-track cameras may also be logistically difficult for long-term monitoring programs due to concerns of theft and vandalism (e.g. Foster 2008; A. Watkins, pers. comm., May 2020). Lures and attractants may be an alternative to increase fox visitation rate and therefore capture rate (Meek et al. 2014) but again, may raise logistical concerns with maintenance.

Survey design – scat genotyping for density estimates

With spatially explicit capture–recapture analyses, the precision of the density estimate increases as the number of unique individuals and recaptures increases (albeit, non-linearly). That is, the relative standard error of the density estimate is inversely related to the square root of the minimum of either the number of unique individuals genotyped or the number of recapture events (Efford and Boulanger 2019). In our study, the number of individuals and recaptures were acceptably high (>20) for most of the unbaited surveys but were generally low (<10) during the baited surveys at both treatment and non-treatment sites. This likely reduced the statistical power of our study to detect an effect of fox baiting on fox densities, if one occurred (Efford and Boulanger 2019). Sample sizes could be boosted by improving genotyping success. Of the collected scat samples, only 42% of scats with extracted DNA were assigned individual identities, a rate comparable with similar nearby studies (Piggott et al. 2008; A. Weeks, unpubl. data). Degradation of scat DNA from environmental exposure can significantly increase microsatellite genotyping error (Murphy et al. 2007; Brinkman et al. 2010). Our survey periods overlapped extreme weather events, which exposed scat DNA to thunderstorms, dust storms and long intervals between collections. These conditions also prevented us conducting several surveys in the post-baiting sessions; more repeat surveys would have resulted in the collection of more scats for these sessions.

Three key assumptions may also have impacted the accuracy of our density estimates, primarily through bias in fox detectability. First, we assumed equal scat collection survey effort among observers. If observer experience affected detection of individual foxes, the detectability component of the SECR model may have underestimated fox density (Borchers and Efford 2008). Second, we assumed that fox detectability followed a half-normal function, and therefore that large detection distances were improbable. Given that a transect survey design only informs fox detection within two directions of its activity centre, a half-normal detection function is a conservative estimate of detection probability (Efford 2011; Efford et al. 2016; Efford 2019). Third, our primary survey sessions were paired and relatively short, so we assumed demographic closure within a survey session. The ‘secr’ closure test function confirmed all sites and periods were closed to gains and losses, except the Mallee non-treatment site in the baited period where data were insufficient to test for violation detect (Efford et al. 2009). However, baited surveys may have captured transient foxes due to timing of juvenile independence and subadult dispersal (Pech et al. 1992). Where possible, we advise future studies to time scat collection surveys during periods of stable fox population dynamics to minimise any confounding effects of season, or to repeat surveying to lower this uncertainty.

Management implications

Our study highlights the many challenges that face land managers when monitoring and managing a cryptic, wide-ranging invasive predator. We reiterate the importance of rigorous study design and setup in detecting population change. Camera traps can be a long-term and cost-effective option for detecting elusive, low-density species (Cutler and Swann 1999; Moruzzi et al. 2002). However, without an informed power analysis or pilot study, poor detection rates and limited numbers of cameras can reduce their power to detect changes in predator activity, so careful consideration must be given to the make and number of cameras, their set-up, positioning and maintenance (Meek et al. 2015; Newey et al. 2015; Findlay et al. 2020; Kolowski et al. 2021). Additionally, activity indices are prone to bias if behavioural changes from baiting, season or habitat affect detectability, obscuring changes in density (Sollmann et al. 2013; Güthlin et al. 2014b; Stephens et al. 2015).

Non-invasive scat genetic analysis is a powerful alternative for detecting predators and monitoring population change (e.g. Gopalaswamy et al. 2012; Anile et al. 2014; Rodgers et al. 2014; Wegge et al. 2019). This methodology can be used to inform management actions, providing useful information on relatedness, survival and population turnover, as well as direct density estimates to refine baiting and monitoring programs, or as an initial reference value for individual-based population simulations (e.g. Hradsky et al. 2019). Technological advances are expected to reduce costs and increase genotyping success rates, improving precision in density estimates and making scat collection more cost-efficient (Morin and McCarthy 2007; A. Weeks, unpubl. data). However, scat genetic analysis, paired with SECR, provides only ‘snapshots’ of density in time, and is unlikely to be financially suitable for continuous monitoring of invasive predator populations.

Camera trap surveys and scat genetic analysis are potentially most useful in tandem, because they can provide multiple lines of evidence about the effectiveness of invasive predator management. Where density estimates achieve a desirable level of precision, they could also be used to calibrate activity indices (e.g. Hopkins and Kennedy 2004; Siddig et al. 2015). Exploring these two metrics collectively has the potential to improve overall monitoring; for example, density estimates from SECR could be used to design appropriate camera trap arrays, and camera traps could locate hotspot areas of fox activity to target scat collection.

Our study is one of the few existing comparisons of relative and absolute predator population density metrics under lethal management actions. Our results suggested some evidence of an effect of baiting on fox density but also high survival rates, a finding supported by a similar study in the wet forests of southern Victoria, Australia (Le Pla et al. 2022). We therefore add to the growing body of literature emphasising the importance of evidence-based management actions. Above all, we demonstrate the many challenges associated with the monitoring and management of invasive mammalian predators, and the need for cautious interpretation of relative abundance indices as standalone metrics.