Current and predicted future distributions of wallabies in mainland New Zealand

ABSTRACT

Bennett's (Notamacropus rufogriseus) and dama (N. eugenii) wallabies have been increasing in numbers and distribution in mainland New Zealand. Here, we update their current distributions; estimate current rates of spread to predict their future distributions; and describe the extent of suitable habitat for each species. Current distributions based on breeding populations and probable distributions based on outlier confirmed sightings resulting from natural dispersal and illegal liberations suggest that Bennett's and dama wallabies currently occupy between 5322–15,229 km² (532,200–1,522,900 ha) and 1865–4126 km² (186,500–412,600 ha), respectively. In 50 years, best- and worst-case estimates predict that they could occupy as much as 5883–44,226 km² (588,300–4,422,600 ha) and 1912–40,579 km² (191,200–4,057,900 ha), respectively. Habitat suitability was predicted to be high in the North and South islands, except for areas of high elevation, intensive agriculture with limited woody vegetative cover and large urban centres. In the absence of widespread intensive control, the ranges currently occupied by Bennett's and dama wallabies are predicted to increase by up to 7- and 20-fold (respectively) in 50 years. As their distributions continue to expand, they will become more difficult to control and their impacts more widespread. We suggest that progressively containing wallabies within increasingly smaller distributions and reducing their numbers to minimise damage within those distributions should be considered the top management priorities. If progressive containment is feasible it could lead to eventual eradication.

KEYWORDS:

• Bennett’s wallaby
• dama wallaby
• habitat suitability
• invasive species
• Notamacropus eugenii
• Notamacropus rufogriseus
• pest control
• range expansion
• unwanted impacts

Introduction

Introduced wild mammalian herbivores pose a significant threat to native ecosystems and production landscapes in many invaded areas (e.g. Pimental et al. 2000; Vázquez 2002; van Wieren & Bakker 2008; Latham et al. 2017). Their feeding and activity habits impact native vegetation, plantation forests, and agricultural and horticultural crops (Leader-Williams 1988; Williams et al. 1995; Vázquez 2002), as well as intensify soil erosion (Bayne et al. 2004; Gangoso et al. 2006). Their impacts are often most severe on island ecosystems (Courchamp et al. 2003) and this is evident in New Zealand where the flora evolved in the absence of mammalian herbivory (Greenwood & Atkinson 1977).

Excluding moose (Alces alces), 22 mammalian herbivores have extant wild populations in New Zealand following introductions that occurred mostly during the mid to late 1800s and early 1900s (King 2005). Many of these species have been shown to cause significant damage to agriculture and the environment (see Cowan 2005; Forsyth & Tustin 2005; Norbury & Reddiex 2005; Nugent & Fraser 2005; Parkes 2005; and references therein). Among these are five species of wallabies that have extant wild populations in New Zealand following liberations that occurred in the late 1800s for private collections and sport hunting (King 2005). Two species, Bennett's wallaby (Notamacropus rufogriseus) and dama (or tammar as it is commonly called in Australia) wallaby (N. eugenii), have well-established populations in mainland New Zealand (Warburton 2005a, 2005b) and a further three species occur only on Kawau Island in the Hauraki Gulf, North Island (Warburton 1986). One of the species now confined to Kawau Island was eradicated from two islands (Rangitoto and Motutapu) in the Hauraki Gulf in the 1990s (Mowbray 2002).

Bennett's and dama wallabies have increased in numbers and distribution since their initial releases, and like some species of macropods (kangaroos and wallabies) in their native Australia (Caughley et al. 1987; Di Stefano 2004), they have had significant unwanted impacts in some invaded areas. The primary impacts that Bennett's and dama wallabies have in New Zealand are damage to agriculture, seedlings within exotic plantation forests (mostly Pinus radiata) and native vegetation. Bennett's wallaby has been recognised as a pasture pest since the 1940s (Warburton 2005b). The plant species that form the bulk of their diet are the same as those grazed by sheep (McLeod 1986). Accordingly, if wallabies are allowed to reach moderate-to-high densities (estimated at c. ≥ 2 per ha, or level 3 on the Guilford Scale; an index based on wallaby faecal pellet density and sign and sightings), they may reduce stocking rates of livestock (Warburton 2005b; NPCA 2015); although this relationship needs to be quantified. Dama wallabies are also primarily grazers, with up to c. 70% of their diet comprising pasture species (Williamson 1986). Although they are better adapted to living in forest interior than Bennett's wallabies, forest-dwelling dama wallabies tend to be smaller and have lower kidney fat reserves than individuals with ready access to grassed clearings or pasture (Williamson 1986).

Bennett's wallabies browse a number of native plant species (McLeod 1986), and can contribute significantly to the depletion of forest understorey and prevent regeneration of palatable species (Warburton 2005b). Although dama wallabies prefer edge habitats they are able to live solely within native podocarp/mixed hardwood forests (Williamson 1986). In forested habitats with few grassy clearings, the diet of dama wallabies comprises a high percentage of native plants. For example, Williamson (1986) reported that only 3% of the diet of forest interior dama wallabies was composed of grasses, whereas native kāmahi (Weinmannia racemosa) and māhoe (Melicytus ramiflorus) together made up c. 50% of their diet. Knowlton & Panapa (1982) erected three permanent 20 m × 20 m plots in the Lake Okataina Scenic Reserve, North Island. One plot excluded only dama wallaby, one plot excluded red deer (Cervus elaphus) and dama wallaby, and one plot was left open to both mammalian herbivores as a control. Results from these plots showed that browsing of seedlings by dama wallabies prevented the regeneration of the most palatable species and that plant species diversity was 57% higher where dama wallabies were excluded (Knowlton & Panapa 1982; Wallace & Wallace 1995; Wright 2017). Benes (2001) and Wright (2017) concluded that the combined effects of dama wallabies and red deer noticeably altered forest understorey composition and structure, and that the regeneration of the most palatable plant species will require the removal or reduction to low densities of both taxa. At high densities, wallabies can also impact soils and accelerate erosion rates (Warburton et al. 1995).

The management of wallabies in New Zealand is the responsibility of affected regional councils (Warburton 2005a, 2005b) and all affected regions include wallabies in their regional pest management plans (RPMP). These plans have objectives to keep wallabies at low abundances, prevent spread outside of containment areas delineated by councils and eliminate isolated populations outside of containment areas (Figure 1). Importantly, the responsibility for wallaby control on private land lies with the landowner and regional councils can issue ‘notices of compliance’ to landowners when wallaby numbers on their properties are assessed and found to be too high (e.g. if Bennett's wallaby density is ≥ level 4 on the Guilford Scale; NPCA 2015; G. Sullivan, Environment Canterbury, pers. comm., 15 January 2016). The most common methods currently used for reducing wallaby numbers include shooting (primarily from the ground) and poisoning using sodium fluoroacetate (compound 1080) or potassium cyanide in Feratox pellets deployed in bait stations (Choquenot & Warburton 2006; Shapiro et al. 2011; NPCA 2015). If landowners fail to manage wallabies in accordance with RPMP rules they can be subject to legal action under the Biosecurity Act 1993.

Figure 1. Historical and current distribution and confirmed sightings and kills (since 2000) of Bennett's wallabies in South Island and dama wallabies in North Island. The black dots represent the sighting and kill locations used to estimate probable distributions using kernel density estimators (see Methods for details); the white dots (Bennett's wallaby only) represent outlying observations suspected to be illegally-liberated individuals. Containment areas were delineated by regional councils. They were designed to contain all known breeding populations of wallabies and have geographic boundaries, such as large rivers, urban areas and mountains that provide defensible boundaries against range expansion by wallabies. The details associated with confirmed wallaby sightings and kills were provided by regional councils and the New Zealand Department of Conservation.

As of 2016, Environment Canterbury's annual budget for Bennett's wallaby management included c. $66,000 for compliance inspections, c. $16,000 for faecal pellet transect monitoring and c. $20,000 for incursion response (G. Sullivan, Environment Canterbury, pers. comm., 15 January 2016). Similarly, Bay of Plenty and Waikato Regional Councils had a joint dama wallaby management programme with a combined annual budget of c. $210,000, plus an additional contribution of c. $10,000 from the New Zealand Department of Conservation (DOC) (G. Corbett, Bay of Plenty Regional Council, pers. comm., 15 January 2016). This expenditure did not include control work paid for by landowners affected by wallabies, for which we do not have estimates. The combined council and agency budget of c. $320,000 per year for control and surveillance of wallabies in mainland New Zealand, and additional landowner expenditure for control, have probably been effective at slowing their rate of spread, but it has failed to prevent them spreading via natural dispersal from their designated containment boundaries (Figure 1). For example, Bennett's wallabies have used bridges and dams to cross the Waitaki River at the southern edge of their containment area and established a breeding population that straddles southern Canterbury and northern Otago (Warburton et al. 2014). Moreover, the number of enforcement notices issued to properties with higher than permissible numbers of Bennett's wallabies has on average increased since 2004 (none prior to 2004; 5.2 enforcement notices from 2004–2017; max. = 16 in 2015), despite comparable effort being put into compliance inspections during that period (B. Glentworth, Environment Canterbury, unpubl. data).

Despite recent range expansions, no primary research has been done on either species of wallaby on mainland New Zealand since the 1980s, and no reviews have been undertaken since 2005 (Warburton 2005a, 2005b). The objectives of this work were: 1. to update the current distributions of Bennett's wallaby and dama wallaby in mainland New Zealand; 2. estimate current natural rates of spread and use these to predict future distributions over the next 50 years; and 3. describe the extent of suitable habitat for each species in mainland New Zealand.

Methods

Wallaby distributions and rates of spread

Historical and current distributions for Bennett's and dama wallabies have been delineated by affected regional councils (Environment Canterbury, Otago Regional Council, Bay of Plenty Regional Council and Waikato Regional Council) based on qualitative and quantitative surveys conducted by them and DOC. We define these as ‘current distributions’ of known breeding populations of wallabies. However, since the early 2000s regional councils and DOC have collated numerous confirmed sightings and kills of wallabies outside of the current distributions and these records are likely the result of illegal liberations and/or natural dispersal (Figure 1). As Lindström et al. (2013) point out, including outlying observations of individuals may be particularly important for estimating population rates of spread and predicting future distributions, because these estimates will include individuals dispersing at the invasion front, as opposed to residents located proximally to the front. In this sense, omitting dispersing individuals may underestimate actual rates of spread at the leading edge of the invasion front. Accordingly, we used 265 and 82 outlying sighting/kill locations (all of which were assumed to have originated from natural dispersal) for Bennett's wallaby and dama wallaby, respectively (see Figure 1, ‘Animals shot or sighted’), to estimate the current ‘probable distribution’. These boundaries were estimated as the 95% isopleth of kernel density estimators (KDE) around the outlying observations plus 200 random points drawn within the current distribution. The 200 points were used to depict the area identified as currently occupied by wallabies and were randomly drawn to avoid biasing the KDE towards specific areas within this distribution. The 2015 boundaries of the current and probable distributions represent best- and worst-case management scenarios, respectively, from which to predict future distributions of wallabies.

We calculated natural rates of spread (m yr⁻¹) for Bennett's wallaby for each of the three time periods for which distributions had been delineated: 1916‒1956, 1956‒1975 and 1975‒2015. Similarly, we used four time periods for which distributions had been delineated for dama wallaby: 1947‒1954, 1954‒1979, 1979‒2000 and 2000‒2007. In addition, we also calculated rates of spread between the latest (pre-2015) council delineated distributions (1975 for Bennett's and 2007 for dama) and the probable distributions constructed using the additional confirmed outlying sightings and kills. Rates of spread for each period were estimated by calculating the straight-line distance between polygon boundaries depicting two consecutive time periods (e.g. 1916‒1956). Because range expansion does not occur at the same rate around the boundary, we sampled the distance between the two polygons every 100 m and then constructed a distribution of annual rates of spread for each period. For this analysis we excluded all established populations that have resulted from known illegal liberations (e.g. Bennett's wallaby population at the northern end of Lake Pukaki; Warburton et al. 2014).

We predicted wallaby range expansion under four natural spread scenarios: 1. low rate of spread (1st quartile annual rate of spread) from the current distribution; 2. high rate of spread (3rd quartile) from the current distribution; 3. low rate of spread (1st quartile) from the probable distribution; and 4. high rate of spread (3rd quartile) from the probable distribution (see values highlighted bold in Table 1). We used 1st and 3rd quartiles instead of the mean or median of the distribution of rates of spread because the former was inappropriate given that the data were right-skewed and the latter did not depict the large rates of spread reported from some parts of the invasion front of both species. We used the 2015 polygon boundaries (current or probable) as the starting distribution from which to incrementally model annual spread. To include biological realism to the modelling process, we constructed layers that biased rates of spread in ArcGIS 10.2 (ESRI 2015). These included a hard ‘barrier to movement’ layer that completely prevented the extension of wallaby range and included the following features: the South Island main divide (Bennett's wallaby only); major rivers; man-made canals in the Canterbury High Country (Bennett's wallaby only); large lakes; and large towns. However, we allowed spread to occur through this barrier at bridges and dams, i.e. wallabies could extend their range beyond barriers but this was biased towards and slowed by crossing features. In addition to hard barriers, we also included a semi-permeable barrier for mountain ranges (excluding the main divide) ≥ 1500 m above sea level (a.s.l.). This layer allowed wallaby range expansion to occur across these high elevations, but the rate of spread was set at half that of elevations < 1500 m a.s.l. Using this modelling approach, we ran each of the four scenarios previously described for 50 years to estimate Bennett's wallaby and dama wallaby distribution at 5, 10, 20 and 50 years into the future. Range expansion for each wallaby species was modelled separately. Finally, we ran a fifth scenario to estimate the impact of human-assisted spread on the future distribution of Bennett's wallaby. In this case, we included recent confirmed sightings and kills that occurred distally to the area included in the probable distribution assuming they represented illegally-liberated established populations (based on information from regional councils and DOC), and allowed these to expand following the rules described above and using a high rate of spread. We did not conduct a similar analysis for dama wallaby because the worst-case scenario (expansion from the probable distribution at a high rate of spread) resulted in a distribution that already included all extreme outlying observations.

Table 1. Bennett's wallaby and dama wallaby natural rates of spread (m yr⁻¹) in mainland New Zealand for three and four different time periods, respectively. Estimates were derived for the known historical distributions as well as the ‘current’ and ‘probable’ (defined as the 95% kernel density estimator around sightings and kill locations) distributions. Rates of spread for a given period were calculated as the distance between the two distributions, measured every 100 m along the distributional boundary. We estimate best-case (current distribution) and worst-case (probable distribution) rates of spread (see descriptions in Methods). Values used to model wallaby range expansion under four different scenarios are shown in bold: we refer to the 1st quartile as a low rate of spread, whereas the 3rd quartile as a high rate of spread.

Wallaby species	Distribution	Time period	No. of years	Min.	1st quartile	Median	3rd quartile	Max.
Bennett's	Current	1916–56	40	0.3	318.8	697.3	1067.0	1703.5
		1956–75	19	0.2	27.2	244.9	620.5	1555.3
		1975–2015	40	0.0	16 . 5	61.3	181 . 7	647.5
	Probable	1975–2015	40	178.0	353 . 5	515.0	827 . 8	1232.5
Dama	Current	1947–54	7	246.7	521.7	657.9	1057.6	1857.1
		1954–79	25	0.5	84.1	192.0	273.4	523.6
		1979–2000	21	0.8	160.8	300.2	450.9	773.3
		2000–07	7	0.0	3 . 1	3.5	92 . 9	976.9
	Probable	2007–15	8	22.3	592 . 9	1106.6	1981 . 3	3883.8

Wallaby habitat suitability

The range expansion models described above depict future distributions at a given point in time representing first-order selection (sensu Johnson 1980) by wallabies. These models do not describe habitat use within predicted future distributions. To identify suitable habitat for Bennett's and dama wallabies and predict potential distribution in New Zealand, we used Maxent 3.3.3 (Phillips et al. 2006). Maxent is a machine learning approach for modelling species distributions from presence-only records. Maxent has been shown to perform as well, or better than, other methods for modelling these kind of data (Elith et al. 2006). In addition to presence data, Maxent uses a number of randomly selected points (called ‘pseudo-absences’) that are then combined with biophysical covariates to construct an index of habitat suitability for each cell ranging from 0 (least suitable habitat) to 1 (most suitable habitat). Presence data for Bennett's and dama wallabies were obtained from three sources. First, locations of incidental observations or where animals had been killed outside known wallaby distributions were obtained from regional councils and DOC (Figure 1). Second, we randomly drew 300 and 267 locations for Bennett's wallaby and dama wallaby, respectively, from within their current distributions to represent occurrence within these areas. Finally, for Bennett's wallaby we also included the locations of 27 pellet transects surveyed since 1993 and a further 11 transects surveyed since 2008, all of which recorded wallaby presence (G. Sullivan and B. Glentworth, Environment Canterbury, unpubl. data). In total, we used 631 and 349 presence-only records for Bennett's wallaby and dama wallaby, respectively. Each wallaby species was modelled separately.

Based on wallaby biology in New Zealand (Warburton 2005a, 2005b), we identified 14 biophysical variables as potentially important predictors of habitat suitability (Table 2). These variables were generated from three different datasets: LENZ (Land Environments of New Zealand, Leathwick et al. 2003); LCDB (Land Cover Data Base v4.1, www.lcdb.scinfo.org.nz); and a digital elevation model (DEM) of New Zealand (www.lris.scinfo.org.nz). All these datasets were in raster format except for LCDB, which was converted to a 25 m resolution raster prior to analyses. Biophysical variables were collected within 1 km² (100 ha) cells; this cell size was chosen as a compromise between estimates of home range size for wallabies in New Zealand (only available for dama wallaby: c. 0.4 km² [40 ha]; Williamson 1986) and a practical unit size for conducting field surveys (Karanth et al. 2009).

Table 2. Biophysical variables used to model habitat suitability for Bennett's wallaby and dama wallaby in mainland New Zealand using incidental observations in Maxent. Variables were collected from the source layer within 1 km² (100 ha) cells. The last two columns show the relative contributions of each biophysical variable included in each model. Models were fitted separately for Bennett's wallaby and dama wallaby. Variables with a high relative contribution to model fit are shown in bold.

Variable	Source	Pixel size (m)	Percent contribution
			Bennett's wallaby	Dama wallaby
Average annual temperature (°C)	LENZ	100	3.2	12 . 6
Average elevation (m)	DEM	25	0.9	5.0
Average slope (°)	LENZ	100	16 . 9	21 . 0
Distance to cover (m)	LCDB	25	3.2	0
Distance to river (m)	LCDB	25	0.2	0.8
Min. temperature of coldest month (°C)	LENZ	100	48 . 9	33 . 0
% agriculture in flat areas^a	LCDB	25	0.1	0.1
% cover (forest, shrubs, tall tussock)	LCDB	25	5.7	0.3
% forest	LCDB	25	1.6	18 . 4
% pasture/grasses	LCDB	25	0.4	4.3
% pasture/grasses in hill and high country^a	LCDB	25	18 . 0	2.7
% shrubs	LCDB	25	0.5	0.2
% urban	LCDB	25	0.2	0.3
% water	LCDB	25	0.2	1.3

^aFlat areas are defined as those between 0 and 100 m a.s.l., hill and high country areas are defined as those above 100 m a.s.l.

We assessed the performance of the Maxent models by determining how well the model discriminated between suitable and unsuitable habitat over a range of thresholds. To do this, we estimated the area under the receiver operating characteristic curve (ROC–AUC), which represents the probability that a randomly chosen presence site will be ranked as more suitable than a randomly chosen pseudo-absence site (Fielding & Bell 1997). A model with AUC = 0.5 performs no better than random, a model with an AUC value between 0.7 and 0.8 has acceptable discrimination, and a model with an AUC value > 0.8 has excellent discrimination and indicates that the model has high predictive power (Fielding & Bell 1997). An additional measure of performance is the regularised training gain, which describes how much better the Maxent distribution fits the presence data compared to a uniform distribution, where all cells have the same probability and a training gain of 0. The exponential of the training gain is a measure of how many times higher the sample likelihood is compared to a random cell (Yost et al. 2008). For example, if the gain is 2, it means that the average sample likelihood is exp(2) ≈ 7.4 times higher than that of a random background pixel.

Results

Wallaby distributions and rates of spread

Chronological distributions for Bennett's wallaby and dama wallaby show significant geographical range expansion since early 1900s, and recent (post-2000) confirmed sightings and kills outside of current containment areas indicate continuing expansion (Figure 1). The predicted future distributions of Bennett's and dama wallabies are shown in Figures 2–3 and Table 3. Currently, Bennett's and dama wallabies occupy c. 5322 km² (532,200 ha) and c. 1865 km² (186,500 ha), respectively. Assuming low natural rates of spread (Table 1) from the current distributions, future distributions for Bennett's and dama wallabies are predicted to increase by c. 550 km² (55,000 ha) and c. 50 km² (5000 ha), respectively, over 50 years (Figures 2A, 3A; Table 3). This represents an increase of 10% for Bennett's wallaby, but only a trivial increase (2.5%) for dama wallaby. Assuming high rates of spread (Table 1) from the current distributions, distributions are predicted to increase to c. 9621 km² (962,100 ha) for Bennett's wallaby and c. 3265 km² (326,500 ha) for dama wallaby (Figures 2B, 3B; Table 3).

Figure 2. Predicted distributions of Bennett's wallaby at five time periods using four different estimates of rate of spread and two distributional polygons. A–B, Expand Bennett's wallaby 2015 ‘current distribution’ using either a low (16.5 m yr⁻¹) or a high annual rate of spread (181.7 m yr⁻¹), respectively; C–D, expand Bennett's wallaby 2015 ‘probable distribution’ using either a low (353.5 m yr⁻¹) or a high annual rate of spread (827.8 m yr⁻¹), respectively; E, predicted distribution of Bennett's wallaby including six illegally-liberated outlying populations. The 2015 ‘probable distribution’ was allowed to expand using a high annual rate of spread (827.8 m yr⁻¹).

Figure 3. Predicted distributions of dama wallaby at five time periods using four different estimates of rate of spread and two distributional polygons. A–B, Expand dama wallaby 2015 ‘current distribution’ using either a low (3.1 m yr⁻¹) or a high annual rate of spread (92.9 m yr⁻¹), respectively; C–D, expand dama wallaby 2015 ‘probable distribution’ using either a low (592.9 m yr⁻¹) or a high annual rate of spread (1981.3 m yr⁻¹), respectively.

Table 3. Size of predicted distributions (km²) for Bennett's wallaby and dama wallaby in mainland New Zealand. Distributions are predicted at five time periods using four different estimates of rate of spread for each species and two different distributions, ‘current’ and ‘probable’.

Wallaby species	Type of distribution	Rate of spread (m yr^−1)	Year
			2015	2020	2025	2035	2065
Bennett's	Current	16.5^b	5322	5395	5443	5553	5883
		181.7^c	5322	5947	6477	7434	9621
	Probable	353.5^b	14 135	14 925	15 703	17 018	20 631
		827.8^c	14 135	15 949	17 444	20 257	28 447
	Probable, IL^a	827.8^c	15 229	18 328	21 529	28 048	44 226
Dama	Current	3.1^b	1865	1873	1881	1887	1912
		92.9^c	2051	2165	2283	2508	3265
	Probable	592.9^b	4126	4789	5489	6780	11 070
		1981.3^c	4126	6335	8649	13 988	40 579

^aProbable, IL is predicted Bennett's wallaby range expansion from the 2015 probable distribution using a high rate of spread (827.8 m yr⁻¹) and assuming spread by illegally-liberated established populations in addition to spread by the main population within the current distribution.

^bLow rate of spread, corresponding to the 1st quartile of the distribution of rates of spread.

^cHigh rate of spread, corresponding to the 3rd quartile of the distribution of rates of spread.

The probable distribution of Bennett's wallaby is 14,135 km² (1,413,500 ha), which is almost three times larger than the current distribution (Figure 2; Table 3). Assuming a low or a high rate of spread from the probable distribution, Bennett's wallaby range is predicted to increase in 50 years to 20,631 km² (2,063,100 ha; Figure 2C) or 28,447 km² (2,844,700 ha; Figure 2D), respectively. The probable distribution of dama wallaby is 4126 km² (412,600 ha), which is more than double the size of the current distribution (Figure 2; Table 3). Assuming a low or a high rate of spread from the probable distribution, dama wallaby range is predicted to increase in 50 years to 11,070 km² (1,107,000 ha; Figure 3C) or 40,579 km² (4,057,900 ha; Figure 3D), respectively.

Including illegally-liberated established populations, the probable distribution of Bennett's wallaby is 15,229 km² (1,522,900 ha; Figure 2E). Predicting range expansion using these populations as the starting range and a high rate of spread resulted in large range expansions into northern Otago and mid-Canterbury (Figure 2E; Table 3). The distribution in 50 years under this scenario is predicted to be 44,226 km² (4,422,600 ha) or 55% larger than under the maximum size predicted under a solely natural spread scenario (28,447 km² [2,844,700 ha]).

Wallaby habitat suitability

The Maxent model for Bennett's wallaby predicted abundant suitable habitat within most of the eastern parts of the South Island and some suitable habitat in West Coast and Southland regions (Figure 4A). The primary lower elevation areas from which Bennett's wallabies are predicted to be absent are associated with high production exotic grassland (e.g. dairy farms) and larger urban centres. Similarly, the model predicted low habitat suitability at high elevations in the South Island's mountains, with the caveat that most valleys provide at least moderately suitable habitat for Bennett's wallaby. Three variables (minimum temperature of coldest month, percent pasture in hill and high country, and average slope) had a combined relative contribution of 84% (Table 2). The AUC (0.92) and Gain (1.31) values indicated that the Maxent model for Bennett's wallaby performed well.

Figure 4. Predicted habitat suitability derived using Maxent, estimated using incidental observations of wallabies (i.e. presence-only data). A, Bennett's wallaby in South Island; B, dama wallaby in North Island. Bright red and orange colours indicate areas predicted to be excellent to good wallaby habitat (HSI = 1), pale blue indicates moderately suitable habitat and dark blue represents poor habitat (HSI = 0).

The Maxent model for dama wallaby predicted abundant suitable habitat within most central parts of the North Island (Bay of Plenty, Waikato and Manawatu–Wanganui), much of Hawke's Bay and Gisborne regions, and smaller, but not insignificant parts of Northland, Auckland, Taranaki and Wellington regions (Figure 4B). High production exotic grassland (e.g. dairy farms) was the primary variable responsible for dama wallaby absence at lower elevations. The model also predicted low habitat suitability at high elevations and around major urban centres. Four variables (minimum temperature of coldest month, average slope, percent forest and average annual temperature) had a combined relative contribution of 85% (Table 2). The AUC (0.96) and Gain (1.72) values indicated that the Maxent model for dama wallaby performed well.

Discussion

Empirical estimates of rate of spread associated with the current distributions were highly variable across time periods for both species. Additionally, historical rates of spread were usually far greater than more recent estimates. This could arise from three non-mutually exclusive factors. First, rates of spread could have slowed recently as the wallaby invasion front met unsuitable habitat. Second, control of wallabies, particularly at the invasion front, may have been adequate to halt or minimise dispersal and thus population spread. This seems unlikely (at least) for Bennett's wallaby, because their distribution and numbers have increased since the disbandment of the wallaby control unit in South Canterbury in 1992 (Choquenot & Warburton 2006). Third, estimates of rate of spread associated with the current distributions may be representative of movement by individuals that reside on the proximal edge of the core breeding population, rather than of movement by dispersers (that generally occur at low densities) at the invasion front (e.g. Lindström et al. 2013). Interestingly, we found that recent rates of spread estimated using the probable distributions were comparable to the earliest historical estimates, suggesting that estimates using the current distributions do not accurately depict dispersal behaviour/movement at the invasion front. Greater rates of spread than those used here have been estimated for brushtail possums (Trichosurus vulpecula) in some areas (e.g. c. 8400 m yr⁻¹ in Northland [Te Kao to Te Paki], 1983–late 1986; Clout & Ericksen 2000). Similarly, the estimated rates of spread for several wild ungulate species introduced into New Zealand are either much larger (e.g. 8690 m yr⁻¹, chamois Rupicapra rupicapra) or comparable (e.g. 1770 m yr⁻¹, red deer; 805 m yr⁻¹, fallow deer) (Caughley 1963) to the spread values we used. As with possums in New Zealand (Cowan 2005), dispersal of wallabies in Australia is strongly biased towards young males attempting to avoid male–male competition (Eldridge & Coulson 2015), although the biological drivers of wallaby dispersal in New Zealand remains a knowledge gap.

High variation in estimated rates of spread created difficulty in determining appropriate values with which to predict future distributions. We consider that projections using low rates of spread (1st quartile) from the current distributions are unrealistic, as they do not account for the large number of recent confirmed sightings and kills outside of the current distributions. In the absence of widespread future control they are, at best, optimistically low. It is also unclear whether it is more appropriate to model future wallaby distributions using the current breeding distributions delineated by councils or (what we term) the probable distributions (i.e. including the large number of confirmed sightings and kills that occur outside council delineations). If these outlying locations represent breeding populations, even at low densities, then the use of the higher estimated rates of spread to model future distributions (Figures 2C or 2D, 3C or 3D) may provide more realistic projections. New surveillance methodology may help confirm the presence of low density breeding populations outside of containment areas, and thus refine estimates of rates of spread. In the absence of these data, we speculate that future wallaby distributions are likely best represented by Figures 2B or 2C, or 3B or 3C. Although not worst-case scenarios, these represent significant range increases for both species. For example, in 50 years the distribution of Bennett's wallaby is predicted to increase in area by 81%–288% with respect to the currently delineated range. For dama wallaby, the range increase would be between 75%–494%. If existing illegal liberations (assuming they are established breeding populations) are included into the estimate for Bennett's wallaby, then the distribution at 50 years is predicted to increase in area by c. 730%.

Further evidence for potential large increases (under current control levels, which appear insufficient to contain populations) in the distributions of both wallaby species comes from our habitat suitability modelling. Bennett's wallaby has been confirmed south of the Waitaki River (Warburton et al. 2014), and the availability of suitable habitat for wallabies invading southwards suggests that their range will increase substantially in the absence of effective control. There is also ample suitable habitat north of the current distribution of Bennett's wallaby, so we would also predict continuing expansion in that direction. The availability of suitable habitat for Bennett's wallabies means that they could occupy most of the South Island if allowed to expand. Similarly, the predicted continued expansion of dama wallabies into Mamaku and Horohoro, and invasion of places like Oropi and Te Urewera, seems reasonable given that forest types around these locations are similar to those already occupied by dama wallabies (Williamson 1986). Although the beech and podocarp dominated forests on the upper slopes of Te Urewera may be less suitable, lower altitude forests dominated by kāmahi (Weinmannia racemosa) and māhoe (Melicytus ramiflorus) are likely to be suitable habitat for dama wallabies (Williamson 1986). The habitat suitability model for dama wallaby supports these predictions, with this species likely to occupy most of the North Island (except for areas heavily influenced by intensive agriculture) if not contained.

Our rate of spread models predict future distributions of Bennett's and dama wallabies at five different time periods. This is important for management agencies because it incrementally illustrates the extent of each region likely to be impacted by wallabies if they are not sufficiently controlled. However, these models were not designed to assess habitat preferences, occupancy or densities within future predicted distributions. We acknowledge that the future predicted wallaby distributions shown here will probably overestimate the actual area occupied by wallabies. For example, dama wallabies are believed to avoid areas of intensive dairying with inadequate scrub/forest cover and thus those should be subtracted from their predicted range, that is 15% (5983 km² [598,300 ha] of 40,579 km² [4,057,900 ha]) of their predicted range in 50 years. Even accounting for this, dama wallaby are still predicted to increase by c. 1750% under a worst-case scenario. It would have been useful to constrain predicted future distributions by estimated habitat suitability to estimate actual area impacted by wallabies at various points in the future. However, we did not do this because of the subjectivity in choosing an appropriate cut-off value to translate habitat suitability into occurrence (Liu et al. 2005), and future changes in land use that could alter the availability and extent of suitable wallaby habitat.

We used environmental data surrounding confirmed wallaby locations to estimate suitable habitat for Bennett's and dama wallabies, and extrapolated this information to get an indication of areas that may be affected by wallabies in the future. Because these analyses used environmental data from wallaby presences (i.e. within the areas they currently occupy), they are probably biased towards predicting similar habitats elsewhere in each of the islands as being the most suitable. For areas that have not yet been invaded by wallabies it is more difficult to predict their suitability because of a lack of data. Thus, habitat suitability models should be considered to represent the minimum area likely to be occupied by wallabies in the absence of effective control. Moreover, because wallabies have not reached their equilibrium distribution in New Zealand, predictions from the Maxent model should be taken cautiously (Elith et al. 2010). One way of addressing these limitations would have been to fit the suitability models with wallaby presence data from their native Australia and use these models to extrapolate suitability in New Zealand. However, many uninvaded habitat types in New Zealand differ substantially from those in their native ranges, suggesting that such an approach may not be very informative; although some environmental variables, such as variation in rainfall in Tasmania, could be. Nevertheless, both models performed well (based on AUC), and thus we can be confident that areas where wallabies are likely to occur have been highlighted. Wallabies are likely to occur throughout the New Zealand mainland, with the exception of high elevations, urban areas and high production exotic grassland (e.g. dairy farms). It is possible that suitable habitat might have been underestimated in high production exotic grassland because the lack of fine-scale vegetation data meant that refugia (patches of scrub, hedgerows and riparian corridors) within these areas could not be easily identified. However, the spatial and temporal resolution of our analyses should be sufficient to guide wallaby management in both islands.

Populations of Bennett's and dama wallabies currently persist outside of containment areas and these, as well as confirmed new out-of-containment-area sightings and kills, appear to be increasing. Based on moderate empirical estimates of rates of spread, and in the absence of widespread sustained control or containment, the distributions of both species are predicted to increase substantially over the next 50 years, resulting in further impacts on agriculture, native vegetation and, potentially, contributing to soil erosion within newly invaded areas. If allowed to increase in distribution and numbers, they may have additional unwanted impacts through, for example, increased rates of collisions with motor vehicles as occurs in their native Australia (e.g. Bond & Jones 2014). The current levels of surveillance and control done by agencies and landowners is likely insufficient to prevent further geographical spread of wallabies. Significantly more funding and a well-coordinated strategy will be needed to achieve this. Effective management of these species will also require surveillance, detection and control methods that are fit for purpose, and these are seen as the highest priorities for the progressive containment of wallabies in New Zealand.

Acknowledgements

We thank personnel from Environment Canterbury, Bay of Plenty, Otago and Waikato regional councils, and the New Zealand Department of Conservation for providing data and valuable discussion. P. Cowan (Manaaki Whenua Landcare Research), E. van Eyndhoven (Ministry for Primary Industries) and two anonymous reviewers provided invaluable comments on earlier versions of the manuscript.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This study was supported by funding provided by the Ministry for Primary Industries, Manaaki Whenua Landcare Research and Strategic Science Investment Funding.