Study area

El Palmar National Park (31°55′S, 58°16′W), in Entre Ríos, north-eastern Argentina, covers ∼8500 ha of savannas, grasslands, scrublands and gallery forests, including high-density stands of the palm Butia yatay, limited by the Uruguay River in the east and a fast highway (Route 14) in the west (Gürtler et al. 2017). The park is mainly surrounded by dense eucalyptus and pine plantations, sometimes combined with cattle for firebreak maintenance, some cattle ranching (on both park borders by the Uruguay River), and a small patch for crops (on permanent rotation including maize, soybean, sorghum and wheat, adjacent to the southern border). Native scrubland occupies a very small fraction of lands adjacent to the park. Intensified agriculture is restricted to the upper watershed of Arroyo Palmar, some 30 km to the west of the park. A permanent water course (Arroyo del Palmar) divides the park into two zones; another three minor water courses cross the park from west to east. The northern zone is open for public use. The southern zone is a restricted area closed to tourists and considered of greater conservation value. Both zones are covered by palm trees, and are poached for deer, boar and capybaras. The increasing expansion of several exotic trees and shrubs, closely linked to the suppression of spontaneous fires and exclusion of cattle in 1970, led to the current lignification of both zones (Ruiz Selmo et al. 2007; Batista et al. 2018). Prescribed fires conducted irregularly since 2002 have not been able to reduce the coverage of woody plants and restore the savanna landscape.

The mean annual temperature over 2006–2015 was 19°C, and the annual mean rainfall was 1389 mm at the closest weather station in Concordia. Long droughts (≥3 months) occurred from the 2005 spring–summer to the 2006 winter; over the same period in 2008–2009; from summer 2011 to winter 2012; over the 2013 fall–winter, and from summer 2014 to fall 2015.

Field data

The data come from the ongoing exotic-mammal management program at EPNP initiated in 2006. The program recruited local sport and subsistence hunters to cull wild boar and exotic deer under a regulated framework. The goals were to reduce the negative impacts of wild boar (on ground rooting and the recruitment of yatay palm trees) and to reduce the abundance of axis deer and blackbucks. The details of the program and a map of the study area have been described elsewhere (Gürtler et al. 2017).

Wild boar and axis deer were culled by controlled shooting, hunting with dogs, and spotlight hunting with spotlights and firearms from pick-up trucks. Hunting with trained dogs targeted wild boar only (although the axis deer incidentally found were dispatched) mainly over 2006–2007, and was accomplished by hunters using several packs of mixed-breed dogs (Table 1 in Gürtler et al. 2018). Spotlight hunting (which had a marginal role over 2006–2015) and dog-hunting are excluded from further consideration here.

Table 1.

Annual crude hunting effort and catch of wild boar and axis deer killed or wounded by controlled shooting with rifles from elevated blinds according to year post-intervention (YPI) at El Palmar National Park, 2006–2015

Controlled shooting using rifles with designated calibres was conducted from elevated blinds (watchtowers) distributed rather uniformly across the park (see fig. 2b in Gürtler et al. 2017). The number of active blinds (i.e. blinds that returned at least one hunting form over any year, regardless of whether they fired any shot) increased from 29 to 40 between year 0 and 3, then decreased to 36 by year 6, and increased to 46 in year 9; 47 identified blinds ever existed over the decade. Each blind had 5–20 shooting lanes in which the vegetation had been cleared. Each lane was from 200 m to 300 m long. Most blinds remained at a fixed position over the decade; a small fraction of them was relocated in the vicinities of the former position for improved accessibility, safety and hunting success. Typically, each hunting party had one or two authorised shooters, and one person in charge of illuminating the target and assisting in retrieving the quarry. However, during some overnight sessions conducted in the initial years and rarely thereafter, two self-identified groups of hunters occupied the same blind over a given session, although not at the same time; for simplicity, we considered only the first of them. Overall, this occurred in less than 1% of all active hunting parties. Hunters and individual blinds were not in strict correspondence at all times: some hunters used the same blind for all or most of their affiliation to the program; some changed blinds on a few occasions after an extended trial period; some shifted between two blinds over several years, and some drop-outs were replaced by newcomers as time went by. Hunters were required to bait the surroundings of their assigned blind with corn or ground pet food and blocks of salt once or twice a week before hunting sessions; in practice, the type, amount and frequency of baiting varied widely.

Hunting sessions were generally conducted between 1600 hours and 1800 hours and 1000 hours and 1200 hours (diurnal or short sessions) every 1 or 2 weeks over 10–11 months a year, but the separation between successive sessions varied widely across the decade. Overnight (or long) sessions frequently started in the evening (between 1600 hours and 1900 hours up to April 2013, and from 2200 hours on since September 2014) and extended up to 0600–0800 hours. All sessions conducted between January 2006 and end of March 2007 were overnight, and thereafter occurred 0–2 times a year. The mean duration of hunting sessions decreased from 11.7 h in 2006 (overnight) to ∼5 h (diurnal) from 2008 onward (Table 1), depending on park decisions, weather conditions, and the choice of individual hunting parties. No catch quotas or limits to the number of shots allowed during a hunting session were set. Each hunter was allowed to take back home up to one annual trophy to minimise selective hunting.

Every hunting session had a park ranger who supervised the activities in situ or at the central operating post (since 2008) through a VHF radio, with which he or she communicated with all hunting parties and other assisting park staff to announce the start and termination of a hunting session and grant permission for ground movement. Hunters communicated to the supervisor any relevant event, reported every shot, and requested permission to descend to finish off or retrieve the quarry or search for wounded specimens within 200 m from the hunting blind. The radio operator recorded in a dedicated logbook the timing and number of shots that each self-identified hunting party made, whether they killed a deer or boar, and any other relevant events (e.g. poaching).

The culled specimens were brought to the central operating post and measured. For each hunting party, park staff filled in a numbered form regardless of whether any wild boar or deer was dispatched, and the reporting hunter signed it. Each form included date, hunters’ names, hunting method, initiation and termination times of the session, firearms and calibres, catch time, whether any boar or deer was wounded and escaped, each specimen’s morphometric measurements, and other details.

Definitions and metrics of relative abundance

For present purposes, we adopt the following definitions: a ‘hunting party’ represents a group of hunters hunting together (using the same blind) on a particular occasion (session), regardless of whether any exotic ungulate was dispatched or any shot was fired. For the few cases in which two parties successively occupied the same blind on a particular session (which occurred in less than 1% of all hunting parties), for simplicity, we only considered the first of them. A ‘hunting session’ is the activity of a hunting party on a given day over one continuous time interval. If the same hunting party hunts on the same session running overnight, we count this as one hunting session. Active hunting parties are those who participated in a hunting session, regardless of whether they were effective in culling or not and of how long they hunted.

Catch (C_j, or harvest) is the number of specimens culled by a hunting party (j) during a fraction or the whole of a hunting session using a given hunting method. C_j includes a few culled specimens in a poor body condition which were discarded and not measured. For each ungulate species separately (and for both pooled) culled by a defined method, aggregate catch (C._t) is the sum, over all hunting parties j, of the catch of hunting party j at time t (i.e. session, trimester or month, and year). Wounded specimens that escaped (as reported by hunters) were considered separately as putative crippling losses; these data were widely variable over time and thus were considered less reliable. Therefore, we alternatively excluded them or pooled them with killed boar, deer and undefined ungulate species (i.e. prey).

For a given hunting party, session and method, hunting effort (E_j) is the sum of hunting-party hours regardless of whether the party caught any boar or deer, and regardless of whether there was one shooter or there were two shooters on a hunting blind in a given session. E_j was computed as the effective duration of individual sessions (π_j,t, the difference between termination and initiation times of each hunting party at time t). Standardised hunting effort is constrained to hunting activity executed over 17–23 h (i.e. the most typical time band and duration across the decade) to allow for hypothetically greater ungulate activity over this period and the much longer duration of sessions during 2006 and 2007 (see above). Aggregate hunting effort (E_.t) was the total number of party hours invested by all hunting parties in a given hunting session or time period t by using a defined method.

The CPUE index (C/E, catch per unit effort) is linearly related to density D (abundance N relative to area A, N/A) through the catchability coefficient (q) under the assumptions of equal D across the area, units of effort (here hunting parties) acting independently, low catch probabilities, and constant q and A (Clark 1985; Hilborn and Walters 1992; Skalski et al. 2005). Explicitly, C/E = qD.

We used four partially related metrics to represent the relative abundance of wild boar and axis deer, and of both species combined (including killed and wounded specimens of any type, i.e. prey, as if hunters were not able to distinguish between exotic ungulate species) over a given time period (Table 1). The underlying rationale was to allow for the limitations of crude indices as surrogates of population size, and because to the best of our knowledge, whether crude or standardised measures of CPUE are better metrics for wildlife has not been addressed before. For crude catch per hunting-party session (C, metric 1), the time unit is the entire hunting session regardless of its duration. Standardised catch per hunting party (metric 2) is the sum of specimens culled over 17–23 h by each hunting party, regardless of the total duration of each session. Crude CPUE (metric 3) is C/E (crude catch per crude effort). Standardised CPUE (metric 4) restricted both the catch and hunting effort to a fixed time window (17–23 h), because hunting sessions differed in total duration and not all hunting parties spent the same amount of time on a given session. For each metric, we calculated the mean and variance across active hunting parties for each trimester or month (spatial TL) or for each hunting blind (temporal TL) across the decade.

Datasets

The management program kept hunting records and compiled a digital database for interim assessment and decision-making. Our current analysis is based on a recent complete review of the digital database against the original hunting forms (in paper format), with the goals of verifying data, filling in missing information, and using revised criteria (e.g. alternative exclusion of wounded specimens). These revisions led to slight changes in catch and hunting effort numbers relative to the data used in previous analyses (Gürtler et al. 2017, 2018). The revised database for controlled shooting includes 6104 records, 1712 wild boar and 2388 deer specimens culled with exact date and time of session, hunting method and hunting party form. Two relevant variables had missing data, namely, hunting blind identity (28 of 6104 records, 0.5%), and time of catch (668 of 4100 quarry, 16.3%). Specimens that were reportedly wounded and escaped were considered separately from killed specimens. For 158 specimens with complete morphometric data (excluding tail length and weight) but lacking species identification, we assigned them to species on the basis of a discriminant analysis of log₁₀-transformed measurements as described in Gürtler et al. (2017). Here, we used the available hunting sessions’ logbooks to trace back the reported species identity of specimens with missing identification in the hunting forms; the concordance between morphometric assignments and logbook records among 55 specimens with data was 100%. Data for the 10-year period (  Table S1 (WR20119_AC.zip)) and the code used for analysis (  Appendix S1 (WR20119_AC.zip)) are included as supplementary files.

Fitting and testing Taylor’s law

We fitted TL and tested its adequacy as a description of the data in three steps, following the detailed account in Cohen et al. (2017b). First, for each ungulate species separately and for both combined, we computed the sample mean (m_h) and the sample variance (v_h) of the catch and CPUE achieved by each active hunting party at a given blind and session over each trimester. For example, for each trimester, we pooled the numbers of boar killed by a hunting party in a hunting blind at a given hunting session (variable X, either catch or CPUE), then computed the mean and variance of X across active hunting parties in this trimester (regardless of the number of hunting sessions or active parties during the 3-month period), then log₁₀-transformed both the mean and the variance of X, and performed an ordinary least-squares linear regression of log₁₀ v_h on log₁₀ m_h, t = 1, …, 38, across all trimesters with data. We repeated these procedures for each ungulate species culled (taken separately and combined), hunting shift (diurnal versus overnight), and human disturbance (public vs restricted zones) using Stata 15.1 (Stata Corp 2018). We also included a least-squares linear regression for the combined data for ungulate species, including the number of killed or wounded specimens regardless of whether they had been identified to species level (i.e. prey). Analyses always complied with the requirements that the mean abundance was greater than 0, that at least 15 observations (here hunting parties for each 3-month or 1-month period) were used to calculate each mean and variance, and that the linear regression (here, for a given species, zone, and hunting shift) should include at least five paired data of v_h and m_h (Taylor et al. 1988, p. 721).

For this and all subsequent statistical tests of a null hypothesis, the critical P-value was set at 0.01. We did not attempt to correct the critical value of these tests for simultaneous statistical inference (Miller 1966) because of the complex dependence of the data used in these tests. For example, the widely used Bonferroni correction assumes that all tests are statistically independent, which is not true here because the four different metrics of abundance are based on overlapping data. Hence, we rely on a combination of statistical hypothesis testing and visual judgments.

Second, we tested for curvature in the relation of log₁₀ v_h to log₁₀ m_h, by fitting a quadratic regression log₁₀ v_h = a + b × log₁₀ m_h + c × (log₁₀ m_h)² by least-squares (Taylor et al. 1978, p. 388, their eqn 14). We refer to this quadratic generalisation of TL as QTL. If the confidence interval of the coefficient c did not include 0, the data rejected TL. We examined the residuals of the linear regression models for heteroskedasticity, normality, skewness and kurtosis by using the commands swilk, estat hottest and estat imtest using Stata 15.1. The results of these tests were not always consistent, even when these tests were supposed to evaluate the same characteristic of the residuals. In such cases, we examined the data graphically and used our judgment to decide on the overall pattern of the data, rather than relying on contradictory evaluations from different tests of statistical significance. Each residual of the TL linear regression measured the stability of population abundance following Döring et al. (2015).

Third, when the analyses in steps 1 and 2 did not reject TL, we used analysis of covariance (ANCOVA, implemented in the anova command of Stata 15.1) to test for differences in the slope and intercept of the species-specific TLs fitted to different subsets of the data (i.e. ungulate species, zones, and hunting shifts). For example, one ANCOVA treated ‘species’ as a categorical variable and asked whether ‘species’ or the interaction term ‘species × log sample mean’ significantly influenced log sample variance. If the interaction term influenced log sample variance, then ‘species’ affected the slope. If both ‘species’ and the interaction term influenced log sample variance, then both the intercept and the slope of TL depended on the ungulate species. We used Welch’s t-test for two quantities with unequal variances to compare slope estimates under two conditions. We present back-transformed means and 95% confidence intervals of indices of ungulate abundance for convenience of interpretation.

Finally, we asked whether the choice of the time scale over which relative abundances were computed would affect the parameters of Taylor’s law. Therefore, we repeated steps 1–3 described above for 1-month periods over the 10-year study.

Spatial Taylor’s law

Wild boar

The spatial TL described controlled shooting of wild boar reasonably well (Fig. 1, Table 2). TL, which corresponds to the null hypothesis that c = 0 in the QTL, was rejected (P < 0.01) only for crude CPUE. The low-variance outlier in trimester 2 and the high-variance outlier in trimester 22 (Fig. 1c) probably contributed to the rejection of TL for crude CPUE. We discuss outliers below. For the remaining three of the four metrics, the slope coefficients b of TL were insignificantly different from 1 (Table 2, Fig. 1a, b, d). However, the intercept a of TL was significantly positive for standardised catch and significantly negative for standardised CPUE, rejecting the Poisson model of random catches in these cases.

Fig. 1.

The spatial TL described the relationship between y = log₁₀ v and x = log₁₀ m of the relative abundance of wild boar, under conditions of controlled shooting from a fixed hunting blind. Each point represents the (a) spatial mean and spatial variance of crude catch per hunting party-session, (b) standardised catch per hunting party-session, (c) crude catch per unit effort, and (d) standardised catch per unit effort across El Palmar National Park at one trimester over 2006–2015. The solid straight lines are fitted by least-squares regression to the data from each trimester separately. The dashed belt above and below each solid regression line represents a 95% confidence interval for individual residuals. Numbers next to the circles represent trimesters over the 10-year period.

Table 2.

Linear regression estimates (by ordinary least-squares) of the parameters of the spatial Taylor’s law (TL) log₁₀ v = a + b × log₁₀ m, where m = sample mean, v = sample variance of various metrics of the relative abundance of wild boar, axis deer and both species combined (including killed and wounded specimens, or excluding wounded prey) computed for 3-month periods at El Palmar National Park, 2006–2015

Minimum log mean, minimum over all hunting blinds with at least one prey of the log₁₀ of the sample mean m_t catch or CPUE in time t. Maximum log mean, maximum over all hunting blinds with at least one prey of the log₁₀ of the sample mean m_t catch or CPUE. Range log mean, maximum log mean minus minimum log mean, to indicate the range of the abscissa in the scatterplot of TL. Intercept_a, least-square estimate of the intercept a of the log–log form of Taylor’s law (TL) log₁₀ v = a + b × log₁₀ m. s.e._a, standard error of the estimate of a. Slope_b, least-square estimate of the slope b of TL. s.e._b, standard error of the estimate of b. Curvature_c, least-square estimate of the curvature c of the quadratic Taylor’s law (QTL). s.e._c, standard error of the estimate of ∼adj R², adjusted R² (adjusted for the number of fitted parameters). d.f., error (residual) degrees of freedom (d.f.) equals number of observations minus the number of fitted parameters. Test_a_0, P-value of the test of the null hypothesis that the intercept a = 0. Test_c_0, P-value of the test of the null hypothesis that the coefficient c of the quadratic term equals 0. A P-value of <0.01 rejects the adequacy of TL because the relationship of log₁₀ v to log₁₀ m has statistically significant curvature. Swilk, Shapiro–Wilk test of the normality of residuals from the linear regression of TL. A low P-value rejects the assumption of normally distributed residuals. Heteroskedasticity, P-value of the test of the null hypothesis that the residuals from the linear regression are homoskedastic, that is, all have the same variance. Information matrix, P-value of the test of the null hypothesis for heteroskedasticity, skewness, and kurtosis. Catch, catch per hunting session; Stand., standardised; CPUE, catch per unit effort; Het., heteroskedasticity; IM, information matrix

The extreme data point labelled 2 (corresponding to the second trimester in 2006) in three of the metrics differed from other trimesters in having few hunting parties on each session, long overnight sessions, and mostly comprising blinds located in the restricted zone, and hunting procedures being in the adjustment stage (Fig. 1a–c). Trimester 1 suffered the same problems. Trimester 22 had fewer hunting sessions than did preceding and subsequent trimesters, mostly because of management-related issues. When we re-ran the analyses in Table 2 after deleting trimesters 2 and 22, of the four metrics for each species, none rejected TL in favour of QTL.

The Shapiro–Wilk test showed significant deviations from normality for crude CPUE only in both TL and QTL, but the other tests showed significant deviations of residuals from normality, homoskedasticity and normal kurtosis for crude catch and crude CPUE.

Axis deer

The data for axis deer were consistent with (i.e. failed to reject) the spatial TL for the four metrics, with no significant curvature. The log mean crude catch was highly significantly correlated with the log variance of crude catch across all active hunting blinds over both zones of the park (adj. R² = 0.779); the slope coefficient did not differ significantly from 1 (b = 1.175 point estimate ± 0.103 standard error) but the intercept differed highly significantly from 0 (a = 0.250 ± 0.051; Fig. 2a, Table 2). The first and second trimesters appeared repeatedly as outliers (Fig. 2a, b, d), along with trimesters 7 and 17 sporadically. The slope coefficient b was not significantly different from 1 for standardised catch (1.000 ± 0.076, adj. R² = 0.749) and standardised CPUE (b = 1.068 ± 0.106, adj. R² = 0.731), unlike for crude CPUE (b = 1.307 ± 0.120, adj. R² = 0.760; Fig. 2b–d, Table 2). Intercepts were significantly different from 0 for crude and standardised catch only. Residuals were not normally distributed in most cases, according to the Shapiro–Wilk test, and showed significant deviations from homoskedasticity and normal kurtosis only for standardised CPUE.

Fig. 2.

The spatial TL described the relationship between y = log₁₀ v and x = log₁₀ m of the relative abundance of axis deer, under conditions of controlled shooting from a fixed hunting blind. Each point represents the (a) spatial mean and spatial variance of crude catch per hunting party-session, (b) standardised catch per hunting party-session, (c) crude catch per unit effort and (d) standardised catch per unit effort across El Palmar National Park at one trimester over 2006–2015. Key as in Fig. 1.

ANCOVA showed no significant effects of species on coefficients a and b of TL for crude catch, standardised catch and standardised CPUE, suggesting that if species affected the slope at all, it did so weakly for deer crude catch only. We did not compare crude CPUE between species by ANCOVA because QTL described crude CPUE better for boar, and TL described CPUE better for deer.

Human disturbance

When we separated the data from the sites in the public-use zone from those from the restricted zone, TL was not rejected for any metric or either species (Table 3). Hence, it was valid to compare these zones using ANCOVA. ANCOVA yielded no significant differences in the parameters of the spatial TL between the public-use zone and the restricted zone for any of the four metrics for wild boar and axis deer (Table 3, Fig. 3), indicating no significant effects of anthropic disturbances and other undefined landscape variation (i.e. forest structure and adjacent crops) on TL parameters. The adj. R² was greater than 0.75 in 12 of 16 cases (16 cases = 4 metrics × 2 species × 2 zones; Table 3).

Table 3.

Linear regression estimates (by ordinary least-squares) of the parameters of the spatial Taylor’s law (TL) log₁₀ v = a + b × log₁₀ m, with m = sample mean, v = sample variance of various metrics of the relative abundance of wild boar, axis deer and both species combined (including killed and wounded specimens, or excluding wounded prey) at El Palmar National Park, 2006–2015, according to hunting zone

Key as in Table 2. Catch, catch per hunting-session; Stand., standardised; CPUE, catch per unit effort; Het., heteroskedasticity; IM, information matrix

Fig. 3.

The spatial TL described the relationship between y = log₁₀ v and x = log₁₀ m of the crude catch per unit effort of exotic mammals, under conditions of controlled shooting from a fixed hunting blind located in the public-use zone (triangles, red regression line) or restricted zone (diamonds, green regression line) of the park. Each point represents the spatial mean and spatial variance of crude catch per unit effort for (a) wild boar and (b) axis deer across El Palmar National Park at one trimester over 2006–2015. Key as in Fig. 1.

For wild boar, the slopes of TL for crude catch per hunting-party session varied little, from 0.972 to 1.016 in the restricted and public-use zones respectively (Table 3). Very slight differences between zones were recorded in the slopes for standardised catch (1.112 and 1.101), standardised CPUE (1.065 and 1.198), and crude CPUE (1.343 and 1.380). Intercepts were significantly different from 0 for standardised catch and standardised CPUE in the public-use zone only.

For axis deer, the slopes of TL for the public-use zone had a wider range (1.072–1.260) than did those for the restricted zone (1.202–1.341), but for each metric the difference between zones was slight except for standardised CPUE (b = 1.072 and 1.283; Table 3). When intercept estimates differed strongly from 0, they tended to do it in both zones. Residuals for deer metrics frequently departed from normality and homoscedasticity.

On average, across 36–38 trimesters with data (depending on zone and metric), the back-transformed mean crude catch (per hunting-party session) and crude CPUE of wild boar in the restricted zone were 43% (t = 2.5, 66.5 d.f., P = 0.015) and 34% (t = 2.3, 73.9 d.f., P = 0.023) greater than those in the public-use zone respectively (Fig. 4a). Conversely, the mean crude catch and crude CPUE of axis deer were 60% (t = 3.4, 75.9 df, P = 0.001) and 87% (t = 4.6, 74.0 df, P < 0.001) greater in the public-use zone respectively (Fig. 4b). Standardised catch and standardised CPUE displayed similar patterns and narrower differences between zones for each ungulate species taken separately.

Fig. 4.

Back-transformed spatial mean and 95% confidence intervals of the relative abundance of (a) wild boar and (b) axis deer by four metrics (×100) according to hunting zone (public use versus restricted use) across El Palmar National Park at one trimester over 2006–2015.

Hunting shift

When we separated the data from the diurnal hunting shifts from those from the overnight hunting shifts, TL was not rejected for any metric or either species, except crude CPUE of axis deer in diurnal shifts (Table 4), which demonstrated statistically significant convexity. Hence, it was valid to compare these shifts using ANCOVA in all cases except crude CPUE for axis deer. ANCOVA found no significant differences in the slopes and intercepts of the spatial TL between hunting shifts (diurnal versus overnight sessions) for crude or standardised catch and standardised CPUE of both ungulates taken separately (Table 4, Fig. 5). For wild boar, the slopes of TL slightly varied between diurnal and overnight shifts for most metrics (e.g. for crude catch, from 1.080 to 1.086 respectively; Table 4, Fig. 5a, c). The hunting shift weakly modified the slope of crude CPUE (P = 0.018) for wild boar. For axis deer, standardised catch and standardised CPUE coefficients for overnight shifts tended to exceed those for diurnal shifts (Table 4, Fig. 5b, d). Nearly all intercept estimates for crude or standardised catch of axis deer were significantly different from 0.

Table 4.

Linear regression estimates (by ordinary least-squares) of the parameters of the spatial Taylor’s law (TL) log₁₀ v = a + b × log₁₀ m, with m = sample mean, v = sample variance of various metrics of the relative abundance of wild boar, axis deer and both species combined (including killed and wounded specimens, or excluding wounded prey) at El Palmar National Park, 2006–2015, according to hunting shift

Key as in Table 2. Catch, catch per hunting-session; Stand., standardised; CPUE, catch per unit effort; Het., heteroskedasticity; IM, information matrix

Fig. 5.

The spatial TL described the relationship between y = log₁₀ v and x = log₁₀ m of the relative abundance of exotic mammals, under conditions of controlled shooting from a fixed hunting blind in diurnal (triangles) and overnight (diamonds) hunting shifts. Each point represents (a, c) the spatial mean and spatial variance of crude catch per hunting-party session and (b, d) crude catch per unit effort for (a, b) wild boar and (c, d) axis deer across El Palmar National Park at one trimester over 2006–2015. Key as in Fig. 1.

Taylor’s law closely fitted the data in 11 of 16 cases (adj. R² > 0.75), with looser fits for standardised CPUE of deer in overnight shifts (adj. R² = 0.557). As mentioned above, TL was rejected for crude CPUE of deer in diurnal shifts, and the adj. R² = 0.449 of TL was notably lower than was the adj. R² = 0.574 of QTL, as expected (Table 4, Fig. 5). As for wild boar, residuals for half of the models of axis deer showed significant deviations from normality according to the Shapiro–Wilk test, and rarely deviated significantly from homoskedasticity and normal kurtosis.

On average, across the 10-year period, the back-transformed mean crude catch of wild boar in overnight shifts was 39% greater than that in diurnal shifts (which were much shorter), but the difference was not statistically significant (t = –1.7, 25.05 d.f., P = 0.11; Fig. 6a). Similarly, the mean crude catch of axis deer per hunting-party session was 36% greater in overnight shifts, but again the difference between shifts was not significant (t = 1.7, 22.6 d.f., P = 0.10; Fig. 6b). However, in terms of crude CPUE, the mean indices were marginally greater in diurnal shifts both for wild boar (by 38%; t = 1.7, 29.16 d.f., P = 0.098) and axis deer (by 39%; t = 1.9, 22.17 d.f., P = 0.074). For deer, both standardised catch and standardised CPUE weakly favoured diurnal shifts (P = 0.025).

Fig. 6.

Back-transformed spatial mean and 95% confidence intervals of the relative abundance of (a) wild boar and (b) axis deer by four metrics (×100) according to hunting shift (diurnal versus overnight) across El Palmar National Park at one trimester over 2006–2015.

Effects of time scale on Taylor’s law parameters

  Tables S2 and S3 (WR20119_AC.zip) provide the same statistics as Tables 2–4, but computed for monthly periods rather than trimesters. The spatial TL was not rejected for either species and for any metric except wild boar using standardised catch or crude CPUE for all data pooled (  Table S2 (WR20119_AC.zip)), and for wild boar crude CPUE in the diurnal shift only ( Table S3 (WR20119_AC.zip)).

In general, the choice of time scale affected little the slope coefficients for each metric, but estimates based on monthly periods were nearly always greater for axis deer and for the restricted zone ( Table S2 (WR20119_AC.zip)). Standard errors were smaller, and the fit of the models improved slightly in most cases.

The choice of a monthly time scale modified the effects of hunting shift and zone on some spatial TL parameters. ANCOVA detected highly significant differences in the slope, but not in the intercept, of TL between hunting shifts for crude CPUE of wild boar, and in the slope and intercept of standardised catch for axis deer only. Regarding human disturbance effects, ANCOVA yielded no significant differences in the parameters of the spatial TL between public- and restricted-use zones for any of the four indices of wild boar. However, for axis deer, both zone and the interaction term zone × log mean weakly (0.01 < P < 0.05) increased log sample variance of crude and standardised catch only.

All indices of relative abundance computed for 1-month periods differed more strongly between hunting zones and shifts than did those computed for 3-month periods and maintained the direction of differences between categorical levels. Wild boar was highly significantly more abundant in the restricted-use (southern) zone by the four metrics, whereas axis deer was captured more frequently in the public-use zone according to all metrics, with differences being highly significant in crude CPUE only (P < 0.001). Regarding hunting shift effects, the crude (P = 0.002) and standardised catch (P = 0.017) of wild boar and its standardised CPUE (P = 0.002) in overnight sessions significantly exceeded the corresponding indices recorded in diurnal sessions. For deer, all indices were significantly (0.001 < P < 0.002) greater in diurnal sessions, except crude catch, which was greater in overnight sessions, although not significantly so.

Temporal Taylor’s law

All metrics for axis deer and wild boar culled at identified blinds over a variable number of trimesters were consistent with (failed to reject) TL (Table 5). Adjusted R² ranged from 0.744 to 0.909. Three of the eight slope coefficients (8 = 2 species × 4 metrics) were significantly greater than 1 (for standardised CPUE of deer and boar and crude CPUE of boar), whereas all eight intercepts but one differed significantly from 0 (Figs 7, 8, Table 5). ANCOVA of TL parameters yielded no significant differences between ungulate species for any of the four metrics (P > 0.2).

Table 5.

Linear regression estimates (by ordinary least-squares) of the parameters of the temporal Taylor’s law (TL) log₁₀ v = a + b × log₁₀ m, with m = sample mean, v = sample variance of various metrics of the relative abundance of wild boar, axis deer and both species combined (including killed and wounded specimens, or excluding wounded prey) at El Palmar National Park, 2006–2015

Key as in Table 2, but the minimum, maximum and range log mean and log variance are here calculated over time at each given hunting blind. Catch, catch per hunting-session; Stand., standardised; CPUE, catch per unit effort; Het., heteroskedasticity; IM, information matrix

Fig. 7.

The temporal TL described the relationship between y = log₁₀ v and x = log₁₀ m of the relative abundance of wild boar, under conditions of controlled shooting from a fixed hunting blind. Each point represents the temporal mean and temporal variance of (a) crude and (b) standardised catch per hunting party-session, (c) crude catch per unit effort and (d) standardised catch per unit effort across El Palmar National Park at one hunting blind over 2006–2015. Key as in Fig. 1.

Fig. 8.

The temporal TL described the relationship between y = log₁₀ v and x = log₁₀ m of the relative abundance of axis deer, under conditions of controlled shooting from a fixed hunting blind. Each point represents the temporal mean and temporal variance of (a) crude and (b) standardised catch per hunting party-session, (c) crude catch per unit effort and (d) standardised catch per unit effort across El Palmar National Park at one hunting blind over 2006–2015. Key as in Fig. 1.

The hunting blinds showing maximum performance (i.e. greater mean crude catch or crude CPUE among those participating in more than 80 sessions) for a given species were mostly invariant of the individual metric considered, but the exact blinds and park sections greatly differed between wild boar (Fig. 7) and deer (Fig. 8). For wild boar crude catch, the three top hunting blinds (3S, 1S, and 1N) and some immediate followers (0S and 4S) were concentrated at the park’s western extreme close to a permanent water course (Arroyo El Palmar; Fig. 7). For deer crude catch, five of the top six blinds (10N, 8N, 21N, 25N and 18N) were concentrated in a stretch of 1 or 2 km of the public-use zone covered by a dense forest of exotic trees or shrubs adjacent to the river coastline, on the eastern extreme of the park (Fig. 8).

Principal findings

We found that the spatial and temporal Taylor’s law (TL) adequately described the relationship of the log sample variance to the log sample mean of four indices of the relative abundance of wild boar and axis deer at multiple fixed hunting sites of a protected area over a decade, except for 4 of 89 cases tested for curvilinearity of the relation. In most cases where TL held for both species and the slopes b could be validly compared, they did not differ significantly between wild boar and axis deer in both the spatial TL and the temporal TL, and ∼1 < b <2, as widely verified (e.g. Taylor et al. 1978, 1988).

If captures of boar and deer had been Poisson distributed across hunting blinds and hunting sessions in each trimester (or month), with different mean values in different trimesters (or months), then the slope of the spatial TL would have been expected not to differ significantly from b = 1 and the intercept of the spatial TL would have been expected not to differ significantly from a = 0, because the Poisson distribution (a commonly used standard for pure ‘randomness’) has a variance equal to its mean, for every value of the mean.

In our data, the slopes b of TL were not significantly different from 1 by three of our four metrics, but one of the most frequently used metrics (crude CPUE) was consistently greater than 1 in all cases we examined (e.g. ranges 1.307–1.669 for the spatial TL in Table 2, 1.127–1.323 for the temporal TL in Table 5), regardless of whether each ungulate species was taken separately or both were combined, and regardless of whether killed and reportedly wounded specimens were pooled (prey). Thus, the slopes of TL depended on the specific metric considered.

When TL was not rejected, tests of the null hypothesis that a = 0 rejected the null hypothesis with P < 0.01 in 5 of the 11 cases in Table 2 for the spatial TL and 8 of the 11 values in Table 5 for the temporal TL. In general, the null hypothesis that a = 0 was not rejected in precisely those cases where the slope b > 1 substantially, as if the TL regressions (both spatial and temporal) had either zero intercept and steep positive slope or positive intercept and lower slope close to 1. These results were robust to the choice of time scale (trimesters or monthly periods).

Either combination of TL parameters (zero intercept and slope greater than 1 or positive intercept and slope indistinguishable from 1) rejected Poisson distributions as a sufficient model for the variability in harvests. Although many factors could cause deviations from pure randomness, two mechanisms that have received much attention among ecologists are behavioural and environmental. Individual deer and boar are likely to be attracted or repelled by other individuals of the same or other species, and these attractions or repulsions are likely to change over time, for example, seasonally, resulting in spatially and temporally varying aggregations. Only foxes in the park may affect the offspring of either exotic ungulate, which lack an apex predator in the region. Axis and other deer species may escape from the immediate presence of wild boar (Tolleson et al. 1995; Latham 1999; Ferretti et al. 2008). Moreover, spatial and temporal environmental variations (e.g. availability of water, hunting pressure, food resources, and weather) are likely to induce or inhibit aggregations of deer and boar around the hunting blinds where controlled shooting was implemented.

We used four metrics to describe the relative abundance of wild boar and axis deer because crude catch or crude CPUE were sometimes considered poor indices of population density, especially in some fisheries (e.g. Walters 2003; Maunder and Piner 2015), but there is no broad consensus on this issue (Pauly et al. 2013). However, CPUE sometimes provided a valid index of wildlife density for management of deer and game birds (e.g. Novak et al. 1991; Lancia et al. 1996; Cattadori et al. 2003; Rist et al. 2010). The relative density of wild boar sighted correlated positively and significantly with other indices such as the fraction of plots with fresh feral pig dung (Hone 2012, p. 35). Here, all four metrics for log mean abundance (temporal TL) were highly significantly (P < 0.001) correlated for wild boar (range of r, 0.9090–0.988) and deer (range of r, 0.964–0.983) across all active hunting blinds. However, crude CPUE differed from all other metrics for boar and deer in that the slope coefficient was significantly greater than 1 in nearly every case, and tended to be much greater when pooled boar and deer or prey were considered in the spatial TL (b = 1.655 and 1.669). This was less obvious for the temporal TL (b = 1.114 and 1.310). These results imply that the variance in crude CPUE increased faster than in proportion to the mean CPUE. One possible interpretation is that when the local density of wild boar or deer increased or more specimens were attracted to the baits (or more bait was deployed regularly), more skilful hunters or those with better gear (or investment in bait) may have taken advantage of these opportunities, whereas less experienced hunters equipped with poorer gear (or with poorer baiting habits) continued to cull quarry at roughly the same rate. Similarly, the slopes b tended to depart more from 1 when both ungulate species were combined or wounded prey were included, again suggesting that hunters facing better chances (i.e. either species) would tend to depart more from their average performance. Standardised CPUE and standardised catch would be less affected by better conditions or chances because the restricted time window in which hunting efforts and outcomes were considered put a cap on both, and allowed fewer opportunities for improved performance; thus, as for other metrics, b was consistently around 1 or marginally different from 1 in both spatial and temporal TL.

Neither human disturbance (i.e. zone) nor the hunting shift significantly modified the parameters of the spatial TL, except in one case. However, Fig. 4 shows that, on average, across the 10-year period, the mean relative abundance of wild boar was consistently and significantly greater in the restricted zone, whereas that of axis deer was significantly greater in the public-use zone, and more specifically on the stretch of exotic trees or shrubs adjacent to the Uruguay River. Similarly, axis deer occupancy peaked in proximity to permanent water courses in Australia’s dry landscapes (Forsyth et al. 2019), as wild boar did in the lower Mississippi valley, USA (e.g. Paolini et al. 2019). Both zones in the park were embedded in forest plantations but differed in the intensity of disturbance (e.g. tourism-related vehicle use) and other features. Only the restricted zone had permanent crop fields adjacent to the border, and these were in proximity to better, more abundant refuge afforded by vegetation type, more permanent water courses and wetlands within the park. The between-zone differences in ungulate abundance were unexpected and have implications for dedicated research efforts and management (see below).

The effect of hunting shift on abundance depended on the metric of abundance. Mean crude catches of boar and deer were substantially greater in overnight shifts (which were much longer than diurnal shifts), but differences in mean crude CPUE were marginally significant and slightly favoured diurnal shifts. According to local hunters and rangers, both ungulate species apparently displayed a peak of activity right before sunset and at dawn. The available evidence does not allow us to conclude whether diurnal or overnight shifts are better for culling boar or deer, in part because in these unplanned comparisons, the exact timing of overnight shifts was quite variable over time, as was the duration of hunting sessions. These factors are likely to affect both the chance of hunter–prey encounter at the blinds and hunter performance over prolonged sessions (i.e. fatigue). Cost-effectiveness and safety issues may need to be factored in, and these depend on hunters’ motivations and setting. Notwithstanding these concerns, deer crude catch was the only metric that favoured overnight shifts substantially over the remaining metrics for deer (Fig. 6b). Camera-trap surveys conducted in the public-use section of the park over 2017–2018 unexpectedly showed that axis deer displayed major activity periods around 0100 hours and 0700 hours over late winter and early summer (Nicosia et al. 2019). On the basis of these results, we implemented a collaborative intervention trial with park management and program-affiliated hunters, in which a series of four or five diurnal hunting sessions was followed by another series of five overnight sessions by using standard hunting procedures. The trial confirmed that deer harvest measured both by crude catch and crude CPUE was three to four times greater in overnight sessions (Gürtler et al., unpubl. report to Administración de Parques Nacionales 2019).

The hunting blinds with top mean crude catch or CPUE for a given species over a decade exceeded the rest by at least an order of magnitude, and differed between wild boar (Fig. 7) and axis deer (Fig. 8). Both the top blinds and their immediate followers were concentrated in exactly opposite stretches of the park: axis deer closer to the river coastline on the eastern extreme, and wild boar on the western extreme by the entry of a permanent water course (believed to be a corridor for boar from the upper watershed) with a dense riverine forest in proximity to high-density stands of palm trees. The latter shed an enormous amount of fruit from mid-summer to mid-fall, which is readily consumed by wild boar (Ballari et al. 2015). These substantial differences are even more noteworthy given the 10-year-long period of observations, during which the holders of local blinds experienced variable rates of turnover, improvement of hunting gear, and long-term landscape change (i.e. lignification). As explained above, various sources of heterogeneity in the local environment may explain the differences in relative abundance between zones, some of them of anthropic origin.

We tested TL by using a large amount of relative abundance data that have some limitations. Several hunting-related factors (e.g. bait, night-vision gear, hunters’ experience, weather) are expected to affect every relative abundance metric via their influence on hunting effectiveness. Hence, slight deviations from linearity in TL should not be given too much weight. Unlike in another multi-site study in which the data analysed for TL had been collected systematically according to defined sampling protocols (e.g. Cohen et al. 2017b), here, how each hunting party operated varied widely in several key aspects. Wild boar was hunted with dogs over the first 2 years after program onset, and this may have affected (positively or negatively) the effectiveness of controlled shooting. Virtually all outlier trimesters corresponded to fall and summer, and none occurred in winter. Fog is more frequent over the fall season, mostly affecting wetland sections of the park and permanent streams, and preventing hunters from sighting the prey; this would explain widely variable catches at some locations. Summer trimesters had reduced or nil coverage of hunting sessions from year 7 to 9.

Implications for management of invasive species

Taylor’s law identifies key zones or sections with a high mean and high variance of boar or deer relative density. We are not aware that TL has been used to show the most productive hunting spots of any wildlife species. If the goal is wildlife control or elimination, these ‘hot spots’ and associated habitats should be preferentially targeted with adequate methods. Locations or habitats with exceptionally low mean density and low variance (i.e. in the lower extreme of the prediction interval) entail stable, low-density local populations, whereas spots of exceptionally high mean density and high variance may represent either sources or sinks that demand adequate prioritisation. Even for locations that lie within the prediction interval of the log–log regression, high mean catches or CPUE with a low or high variance over time may signal different processes, such as focal immigration of exotic ungulates to the protected area, inconsistent baiting practices or poaching (i.e. local depletion of the target species by unauthorised third parties). For example, local hunters and rangers recurrently blamed poachers for pre-emptive culling of wild boar at the western extreme of the restricted zone, which translated into large swings in the harvest achieved by program hunters. According to these sources, poaching activity has decreased over time in direct correspondence to the decline in wild boar abundance caused by program activities, increased with recurrent economic crises, and likely shifted to the increasingly abundant axis deer and capybaras. In the past, poachers used dogs, horses and knives to hunt wild boar when it was abundant. Alternatively, successful program hunters may locally deplete the target species so that a period of time is required for local repopulation.

The greater abundance of wild boar in the restricted zone (with greater conservation value) entails greater risks of predation of yatay palm tree saplings and would strongly reduce its recruitment, one of the main goals of the management program and main reasons why it was created and has been sustained. Identifying the functional relationship between wild boar abundance and palm tree recruitment at local levels would provide valuable information and improved metrics for performance monitoring of the program, especially in the restricted zone.

Taylor’s law can be used to design sequential sampling to attain a fixed-precision estimate of wildlife abundance, as was recently undertaken for fish assemblages in a lake in Africa (Xu et al. 2019). This framework and TL can be used to sample wild boar or deer more efficiently. To apply this method to the boar and deer data, the mean and the variance have to be calculated cumulatively, starting with the first 5 or 10 observations, and then incrementing with each additional observation. Here an ‘observation’ could mean a single hunting session, or a trimester (including all the hunting sessions in that trimester), for example. A cumulative sample size curve is plotted, where x = the cumulative number of observations, and y = the cumulative number of boar (or deer) harvested by the x observations. When a certain line constructed from the parameters of the cumulative TL and from the desired fixed precision (see eqn 2 in Xu et al. 2019) intersects the cumulative sample size curve, the desired fixed precision has been attained.

In addition to fixed-precision sequential sampling, TL can be used to establish an upper and a lower ‘economic injury limit’ (EIL) that determines whether to start control measures (above the upper EIL), stop control measures and sampling (below the lower EIL), or continue sampling to improve the estimated abundance (between the lower EIL and the upper EIL). For wild boar, the maximum extent of ground rooting (one of the typical damages) allowed may represent the EIL. Given the positive relationship between two indices of wild boar relative abundance and ground rooting (e.g. Hone 2002, 2012, pp. 35, 59), the threshold CPUE at which control measures should be started is relevant and needs to be identified with the required precision.

Questions for future research

As with hunting blind locations, TL may also help identify key habitats where intensified control measures may be more cost-effective. The underlying drivers for heterogeneities in wild boar and deer abundance remain to be identified and have implications for improved management.

The exotic ungulate management program was based on the hidden assumption that the catches of wild boar and axis deer are independent, in the sense that shooting specimens from one species would not affect the chance of sighting, shooting and dispatching specimens from the second species. However, shooting is known to affect the behaviour of ungulates, and more specifically, wild boar displays elaborate evasive tactics when confronting hunting (Ohashi et al. 2013; Thurfjell et al. 2013), as do various deer species when confronting wild boar (Latham 1999). Whether the actual catch of hunting parties is spatially dependent cannot be inferred from our current analyses. These issues may be relevant for TL and improved management.