virToad: simulating the spatiotemporal population dynamics and management of a global invader

The cane toad (Rhinella marina) is one of the most globally significant and well-studied invasive alien species, and the detrimental impacts of its invasions warrant the design and application of decision support tools. While many models have been developed for guiding policies addressing cane toad invasions, none reliably predict the species’ population dynamics at scales relevant to on-the-ground management. We describe virToad—an individual-based life-history simulator of the cane toad. We then illustrate virToad’s ability to forecast the cane toad’s spatiotemporal population dynamics at local- to landscape-scales, and its potential for improving management responses to cane toad invasions. We designed virToad to make population dynamics an emergent consequence of the cane toad’s fitness-maximising behavioural responses to mechanistic constraints (e.g., water availability, kin selection), and to management actions. We used virToad to simulate cane toad population dynamics in the absence of management, and under alternative management strategies implemented across a spectrum of effort: hand-capturing and trapping of juveniles and adults, fencing waterbodies, and trapping and chemically suppressing tadpoles. virToad produced plausible predictions of cane toad population densities, detection probabilities, distributions, and spatial segregation. Simulation experiments indicated that the efficacy of competing management actions varied significantly, and that only moderate to high effort hand-capturing and trapping of juveniles and adults had the potential to suppress invasions. virToad is an open-source, rigorous, and extensible decision support platform that will enable researchers and practitioners to defensibly forecast local- to landscape-scale cane toad spatiotemporal population dynamics and management outcomes.


Introduction
The global spread of invasive alien species (IAS) is a legacy of recent centuries' rapid globalisation and industrialisation (Seebens et al. 2018;Hulme 2021), and is a hallmark of the Anthropocene (Crutzen 2002). Invasive alien species are non-native species introduced by humans into areas beyond their natural ranges, where they subsequently establish and spread (Blackburn et al. 2011). These species can cause profound and irreversible changes across scales, from genes to ecosystems, which consequently can compromise the conservation of biodiversity and ecosystem integrity (Blackburn et al. 2014). The resulting ecological changes can impact cultural and socio-economic activities of invaded jurisdictions (Bacher et al. 2018;Shackleton et al. 2019), and result in significant economic losses (global economic loss of ~ 160 billion US$ in 2017; Diagne et al. 2021).
One of the most globally significant alien invaders is the cane toad (Rhinella marina; Cuthbert et al. 2021). The cane toad is a neotropical species wellknown for its widespread and ecologically impactful invasions (Shine 2010;Shine et al. 2020). In the earlytwentieth century, it was deliberately introduced in many countries across the globe as a biological control agent to address sugarcane pest infestations (Easteal 1981;Pili et al. 2019). In Australia, since its first release in the 1930s, the cane toad has invaded ~ 1.2 million km 2 of land area across northern Australia (Urban et al. 2008). Cane toads impact Australian ecosystems in part by poisoning frog-eating predators, and by competing for food and habitat resources with native species, ultimately resulting in functionally irreversible community-level changes (Doody et al. 2017;Shine 2010;Shine et al. 2020). Decades of sustained investment in cane toad research and management have made the species' one of the most extensively studied and monitored IAS globally, and thus an ideal study system for invasion science and other scientific disciplines (van Wilgen et al. 2018). To date, several approaches to managing cane toad invasions have been developed and implemented, and many more have been proposed ; however, evaluating the effectiveness of these management strategies has proven challenging.
Ecological models are indispensable tools for managing biological invasions (Venette 2015), and numerous models have been developed to understand and predict cane toad invasions (e.g., Perkins et al. 2013;Tingley et al. 2014;Pili et al. 2020). Albeit instrumental for answering theoretical questions, most cane toad models lack the spatiotemporal components and biological realism necessary to generate defensible forecasts of invasions under plausible competing management scenarios (Cuddington et al. 2013;Buchadas et al. 2017;Schuwirth et al. 2019). Of the few models that have been embraced for decision support (e.g., Southwell et al. 2017), none reliably extrapolate from local to landscape scales -the latter being that at which alien species invasions and impacts are effectively controlled and mitigated, and thus where management planning should ideally take place (Glen et al. 2013;Lurgi et al. 2016).
Building on the existing literature and models addressing cane toad invasions, we developed a novel spatially-explicit, temporally dynamic, mechanistic, individual-based life history simulator called virToad. Individual-Based Models (IBMs) can account for spatiotemporally dynamic and fine-scale life-history processes (e.g., interactions with the environment or other organisms) that are impossible to represent in traditional models; and this additional biological nuance is essential when extrapolating localscale mechanisms and drivers up to populations and landscapes (Cuddington et al. 2013;Buchadas et al. 2017;DeAngelis and Yurek 2016). Likewise, IBMs have the potential to mechanistically link life-history processes to landscape structures (e.g., habitat fragmentation, landscape connectivity), dynamic environmental conditions (e.g., daily rainfall), and simulated management actions, making it possible to generate forecasts of population dynamics that are responsive to environmental change and to competing management scenarios (McLane et al. 2011;DeAngelis and Yurek 2016). Together, these attributes make IBMs essential decision support tools for IAS management (Briscoe et al. 2019;Schuwirth et al. 2019).
We developed virToad to understand better cane toad spatiotemporal population dynamics (e.g., density, distribution, spatial segregation) realised at local to landscape scales. We anticipate that virToad will principally be used as a decision support tool, enabling researchers and practitioners to explore, optimise, and prioritise competing management responses to cane toad invasions. To this end, we have designed virToad to simulate five currently used and proposed cane toad management strategies: (i) hand-capturing post-metamorphs (i.e., juveniles and adults), (ii) trapping post-metamorphs, (iii) fencing waterbodies, (iv) trapping tadpoles, and (v) chemically suppressing tadpoles.

Methods
We developed virToad using HexSim (v.4.0.18;www. hexsim. net;Schumaker and Brookes 2018). HexSim is a flexible model development environment where users can simulate a broad range of ecological processes. The virToad computer code and the project workspace (described below) can be freely accessed at https:// github. com/ arman pili/ virTo ad. We focus below on model design and implementation, the simulation experiments we conducted, and the techniques we developed for analysing emergent population dynamics. In Online Resource S1, we provide a TRACE (transparent and comprehensive model evaluation) document detailing our model design, implementation, testing, and applications (Grimm et al. 2014). In Online Resource S2, we provide a detailed description of our experimental simulation of competing management scenarios.

Simulation landscape and environment
virToad has a daily time step and is run in a simulation landscape made up of a 2-dimensional array of 86.6 m 2 hexagonal cells, with the distance between hexagon centres being 10 m. We ran our simulation experiments for five years on a 3 × 3 km simulation landscape (103,800 hexagon cells: 300 wide × 346 tall), characterising the landscape structure and dynamic environmental conditions in Beatrice Hill Farm in Northern Territory, Australia (12.621794' S, 131.305016' E; Fig. 1). Beatrice Hill Farm is where many of the country's long-term studies on the cane toad's biology and ecology have been conducted (Brown and Shine 2016;. Input data required by the virToad model include a landcover map depicting the spatial arrangement of habitat types (following the classification of Brown et al. 2006 andGonzález-Bernal et al. 2015), and a table containing daily and cumulative rainfall records. Users may also provide maps indicating the locations of management areas and fences, but these inputs are only required when running management scenarios. Tables are provided in CSV format. Meanwhile, all map data are imported into HexSim in ASCII raster format, except for fences data, which are provided as an ESRI shapefile. Section 2.5 of Online Resource S1 describes in detail the preparation of input data.

Life history
In virToad, as is the case in the real world, cane toads have a semi-aquatic amphibian life cycle, where terrestrial adults lay aquatic eggs, and aquatic larvae metamorphose into terrestrial juveniles. Toads interact with each other, especially when breeding and in the aquatic phase, during which tadpoles cannibalise and chemically suppress eggs to selectively promote the development and survival of kin. Moreover, toads can perceive environmental moisture and distant habitats (foraging, rehydration, and breeding areas). Social interactions and environmental constraints affect the toads' internal physiology and fitness, and in response, toads modify their behaviour, with the goal of maximising fitness. Below, we outline model mechanisms that are key to simulating these behaviours and ecological processes. Section 2 of Online Resource S1 describes the model in detail following the Overview, Design Concepts, and Details (ODD) protocol for describing agent-based models and IBMs (Grimm et al. 2020). Unless noted otherwise, aspects of the cane toad's biology and ecology discussed below have all been replicated in virToad.

Growth and development
Real-world and simulated cane toads have a semiaquatic amphibian life cycle (Crump 2009), with life stages classified into aquatic phase egg and tadpole classes, and terrestrial phase metamorph, juvenile, and adult classes (see Zug and Zug 1979;Hearnden 1991;Cohen and Alford 1993;Cohen 1995). Toads start as eggs and, after three days from birth, hatch into tadpoles (Hearnden 1991). Cane toad tadpoles produce chemical pheromones that suppress the growth and development of eggs.
Simulated tadpoles not exposed to suppression pheromones as an egg metamorphose at about 34 days, whereas those exposed metamorphose at about 55 days (Hearnden 1991;. The growth rate of simulated terrestrial phase toads follows a von Bertalanffy growth model (Kaufmann 1981), fitted using data obtained from field experiments (see Cohen 1995): where age represents days since metamorphosis, and size is expressed in millimetres. Under this growth model, the initial snout-vent length (SVL) of simulated metamorphs was 7.54 mm, and they took 52 days to grow into juveniles (SVL = 30 mm), and another 243 days to grow into sexually mature adults size = 124.435(1 − e −0.0041(age+15.330) ) Fig. 1 Overview of the virToad model. virToad's simulation landscape and environment are parameterized using four model inputs: (1) daily rainfall, (2) landscape structure, and optionally, (3) the location of the management area, and (4) the location of fenceline barriers. Our model input data were obtained from Beatrice Hill Farm, Northern Territory, Australia (12.621794ʹ S, 131.305016ʹ E), which we used for the simulation experiments. virToad was calibrated to reproduce species-specific empirical patterns and data; and for this reason, the model's demographic and behavioural parameters are not regarded as inputs (they should not be altered). virToad outputs include forecasts of daily cane toad population size and distribution (SVL = 90 mm). The asymptotic size of simulated adults was 124.44 mm (see Zug and Zug 1979;Cohen 1995).

Water ecophysiology
Each day, juvenile and adult cane toads lose an average of 5.76% of their body weight due to cutaneous water loss. Meanwhile, because smaller individuals lose water at higher rates, metamorphs lose 30% of their body weight daily due to cutaneous water loss (Schwarzkopf and Alford 1996). Simulated cane toads are considered dehydrated when they lose ≥ 10% of their body weight due to water loss, and they die from desiccation upon losing 40% (Jørgensen 1991). In response to dehydration, simulated terrestrial phase toads will attempt to rehydrate. Terrestrial phase cane toads can sense both the season and rainfall, and they utilise this information to adjust their response to water stress (and cue breeding, detailed below). When the landscape is sufficiently wet, typically during the wet season and days with sporadic rain events, toads will attempt to rehydrate from their immediate environment (Schwarzkopf and Alford 1996;Seebacher and Alford 2002;Brown et al. 2011). During the dry season, toads will actively disperse to rehydration sites (i.e., margins of water bodies and moist microhabitats in human habitation; Carpenter and Gillingham 1987;González-Bernal et al. 2016).

Behavioural responses
Each day, simulated cane toads are either active or, in the case of terrestrial phase toads, remain inactive in diurnal shelter sites Brown et al. 2011). Simulated toads that are active are either breeding, rehydrating, or foraging. As is the case in the wild, only simulated adults breed, and they only breed during the wet season (Yasumiba et al. 2016). Moreover, simulated adult females only breed until they spawn, whereas simulated adult males engage in breeding activities multiple times throughout a breeding season (Hearnden 1991;Cohen 1995). Simulated dehydrated cane toads decide whether to rehydrate or not, with the likelihood of rehydrating increasing with water stress (depicted by body weight, as described above; Jørgensen 1991). Active, non-breeding, and non-rehydrating simulated toads continuously forage.

Reproduction
Breeding adult cane toads disperse to spawning areas (i.e., the margins of still-moving waterbodies, such as lakes and ponds; Hagman and Shine 2006). Both male and female cane toads attempt to select spawning sites that optimise the survival of their offspring; and the prospect of cannibalism or chemical suppression of their offspring by older tadpoles will influence spawning site selection (Hearnden 1991;Clarke et al. 2015). The best spawning sites are those with no or few tadpoles, and with high densities of calling adult males (Hearnden 1991;Clarke et al. 2019;McCann et al. 2020a). Once a spawning site has been selected, simulated breeding females mate with a randomly selected breeding male, and then spawn a clutch of eggs (but see Bowcock et al. 2013;Yasumiba et al. 2015). Clutch sizes (the number of eggs produced) are computed using an equation fitted with data obtained from field experiments (see Hearnden 1991): Mother's SVLs range between 90-124.4 mm, which results in clutch sizes ranging from 7675 to 14,288.

Foraging ecology
Cane toad foraging behaviour varies based on age class. Foraging tadpoles disperse along the margins of natal sites (Hearnden 1991). Foraging metamorphs disperse to vegetated areas near their natal sites (Child et al. 2008), while post-metamorphs generally disperse to naturally or anthropogenically disturbed areas (González-Bernal et al. 2015Lettoof et al. 2018). Due to a lack of available studies examining cane toad foraging ecology in detail, virToad assumes that simulated cane toads can obtain at least their minimum required resources whenever foraging. We feel this assumption is clutch size = 1512.0642e 0.01805(mother � sSVL) highly plausible, as real-world cane toads typically have excessive energy reserves, and mortality resulting from insufficient food intake has not been observed in the wild (Brown et al. 2006;Yasumiba et al. 2016).

Mortality
Simulated cane toad mortality was evaluated daily. Survival rates capture the multiple interacting factors (except desiccation; see water ecophysiology above) that may cause a toad to die, such as density effects, resource availability, disease, and predation. Notably, density effects on cane toad survival in Australia were found to be negligible (Hearnden 1991;Greenlees et al. 2007;Lindstrom et al. 2013;González-Bernal et al. 2016) and, thus, were not separately accounted. Survival rates for each life stage were obtained from published studies (Table 1), except for metamorphs and juveniles whose rates we calibrated using a Pattern-Oriented Modelling approach (see Sect. 3.7 of Online Resource S1; Grimm and Railsback 2012) that attempted to reproduce observed adult population densities (Smart et al. 2020). Our calibrated metamorph survival rate is well within the range observed in field enclosures (1.16% to 17.64%; Cohen and Alford 1993). Meanwhile, the survival rate of juvenile toads in Australia is unknown. Bayliss (1994) estimated the survival rate of juvenile cane toads in Brazil to be about 10% of that of their adult counterparts; however, population models demonstrated that this approximation underestimates juvenile survival rates in Australia (Lampo and De Leo 1998), where the cane toad's natural predators, competitors, and diseases are absent.

Management strategies
virToad can simulate five management strategies, which can each be flexibly implemented along a spectrum of effort: Hand-capturing post-metamorphs involves the manual removal of juvenile and adult toads, and is typically performed along the margins of spawning areas (Greenlees et al. 2018;Shine et al. 2018;Greenlees et al. 2020). In virToad, hand-capturing is implemented on a scheduled basis, corresponding to a specified level of effort. On days when hand-capturing is implemented, active (breeding, foraging, or rehydrating) post-metamorphic cane toads located along the shorelines of spawning areas (0-20 m from the margin of spawning areas), within the management area, are removed from the simulation. The number of toads removed per day varies based on the detection probability of post-metamorphs on that day.
Trapping post-metamorphs involves baiting of foraging and breeding post-metamorphic toads using UV lights and recorded cane toad advertisement calls (Yeager et al. 2014;Muller and Schwarzkopf 2017). In virToad, trapping post-metamorphs is implemented on a scheduled basis, and associated with a specified level of effort. We simulated post-metamorph trapping by placing virtual traps about 240 m apart, along the margins of spawning areas within the management area. An average of six post-metamorphs are captured by each trap per day of implementation, as observed in field trials (see Muller and Schwarzkopf 2017).
Fencing waterbodies is used to exclude postmetamorphic cane toads from accessing these critical resources (Florance et al. 2011). To simulate this management strategy, we made use of HexSim barrier maps to restrict terrestrial toad movements. Hex-Sim barriers can impart mortality, deflect individuals that encounter them, or they may be crossed without incident. The probabilities of these three outcomes may vary across space and time, but must always sum to 100%. virToad mimicks waterbody fencing by imposing barriers that have a 100% probability of deflection and 0% mortality (but in simulation experiments detailed in Online Resource S2, we also tested the implications of breaches in barriers). In virToad, fencing waterbodies is implemented on a scheduled basis, which is governed by a specified level of effort.
Trapping tadpoles involves the baiting of tadpoles using chemical attractants McCann et al. 2019a, b). In virToad, tadpole trapping is simulated by placing traps at least 10 m apart on a proportion or along the entire margin of spawning areas within the management area. Trap placement frequency was associated with a specified level of effort. As is the case in the wild, trapped simulated tadpoles could not be detected by breeding adults prospecting for optimal spawning sites (McCann et al. 2020a). Moreover, trapped simulated tadpoles were only removed from the simulation after they metamorphosed, or when the virtual traps were removed from the simulation, which allowed them to continue producing pheromones that suppress the growth and development of nearby eggs.
Chemically suppressing tadpoles involves the use of chemicals that inhibit the growth, development, and survival of eggs (Crossland and Shine 2012;Clarke et al. 2015;McCann et al. 2020a, b). Because virToad already simulates the production and behavioural consequences of these compounds, it is straightforward to replicate this management strategy by adding extra virtual chemical pheromones at least 10 m apart on a proportion or along the entire margin of spawning areas within the management area. We associated the density of pheromone release sites with a specified level of effort.

Simulation experiments
We ran simulation experiments that illustrated vir-Toad's ability to predict emergent population dynamics in response to simulated cane toad management strategies. We also included a 'baseline scenario' in which population dynamics were unconstrained by management. We examined 40 competing management scenarios within which the level of effort varied across space and/or time (Online Resource S2; Table S2.1). Simulated management strategies were implemented within a management area (Fig. 1), typifying the extent to which typical on-the-ground management actions are implemented in Australia (Florance et al. 2011;Shine et al. 2018;Greenlees et al. 2020). To both remove the confounding effects of the start-up bias period, and to speed up computation, we started all scenarios using a fixed initial population obtained from a simulation that was allowed to run to steady-state in the absence of management (see Sect. 8 of Online Resource S1). All simulations were run for five years. Preliminary analysis indicated that ten replicates were sufficient to produce reliable statistical inferences (see Sect. 8 of Online Resource S1). We therefore ran ten replicates of our baseline scenario, and ten replicates of each of the 40 management scenarios. For each individual scenario, we then combined the model outputs across all replicates to analyse emergent population dynamics.

Emergent population dynamics
virToad simulations generate daily forecasts of the cane toad's population size and point distribution (see Sect. 5 of Online Resource S1). Using the model outputs from the baseline scenario, we analysed the cane toad's spatiotemporal population dynamics at stable-state, by computing (i) aerial density, (ii) linear density, (iii) detection probability in spawning areas, (iv) sex ratio in spawning areas, (v) spatial distribution, and (vi) spatial segregation. Unless otherwise specified, we stratified these measures of population dynamics into age and sex classes: metamorph, juvenile, adult male, and adult female. Likewise, using the model outputs from our management scenarios, we assessed the efficacy of competing management responses based on the cane toad's linear density on the year of implementation, and at 1-year and 5-year intervals post-implementation. These evaluations were only applied to the local toad population present within the management area. All analyses were conducted in the R statistical platform (v.4.1; R Core Team 2021), using the suite of tidyverse packages We first computed aerial density (individuals/km 2 ) of cane toads across the entire landscape. To better compare our measures of population dynamics to findings from earlier studies (see Cohen 1995;Smart et al. 2020), we also computed linear density (individuals/km of shoreline), as well as post-metamorph detection probability and adult sex ratios in spawning areas. Linear density was quantified by counting individuals within a 250 m distance from the margin of spawning areas, which comports with cane toad site fidelity and water-sensing ability (Carpenter and Gillingham 1987;Sinsch 1987;Child et al. 2008;Brown et al. 2015). We measured post-metamorph detection probability as the proportion of either foraging, breeding, or rehydrating post-metamorphs found along the shoreline, within 20 m distance from the margin of spawning areas (Smart et al. 2020). Finally, we computed adult sex ratios as the proportion of adult males vs. adult females found along the shoreline of spawning areas (Cohen 1995).
To investigate the spatiotemporal dynamics of the cane toad's distribution, we generated daily spatial density maps for each age-sex class by spatially smoothing point distribution data using a Gaussian kernel density estimator (bandwidth = 100 m). We then visually analysed density maps for spatial and temporal patterns. To investigate differences in spatial distribution during different periods, we also performed a pairwise Warren's D similarity analysis of these density maps (Warren et al. 2008). Here, we sampled the daily density maps to capture the distribution during the early-wet (November), mid-wet (December to March), late-wet (April), early-dry (May), mid-dry (June to September), and late-dry (October) periods in Beatrice Hill Farm. To investigate the spatiotemporal dynamics of the cane toad's spatial segregation, we generated pairwise bias maps of age-sex classes during these different periods. Bias maps reveal the relative differences in density between age-sex classes across space. We visually analysed the bias maps for spatial and temporal patterns. We also computed the density and relative frequency of toads in different habitat types (individuals/ km 2 of habitat) and time periods.

Sensitivity analysis
In the spirit of Marcot et al. (2015), we used our baseline scenario to conduct a sensitivity analysis of 24 model parameters. Whenever applicable, we varied parameters by ± 10% and ± 20%, relative to their default values (Table S1.17). We selected parameters involved in simulating the cane toad's behaviour, development, movement, reproduction, survival, and water ecophysiology. This resulted in 92 sensitivity analysis (SA) scenarios, each of which we initialised at steady-state conditions, ran for ten years, and replicated ten times. We quantified the model's sensitivity to each parameter by comparing the population dynamics of each SA scenario vs. that of the baseline scenario. The aspects of population dynamics we examined include: (1) adult population growth/decay rate (depicted by the rate coefficient, β, of a Gamma-Generalized Linear Regression Model); (2) the average daily adult aerial and linear population density in the 10 th simulation year; and (3) adult female and male spatial distribution in the 10 th simulation year (See Sect. 8 of Online Resource S1 for details).

Population dynamics without management
Results from our baseline scenario indicated that the cane toad's population density fluctuates seasonally (Fig. 2). The average daily aerial density of adults was 1034.24/km 2 . The lowest density was observed during the late dry period (949.87/km 2 ), but this was immediately followed by an abrupt increase to the highest observed density (1105.64/km 2 ). In contrast, the average daily linear density of adults was 398.61/ km shoreline, with the lowest density observed during the late-wet period (325.09/km shoreline) and highest during the late-dry period (516.88/km shoreline). Juvenile population density peaked during the midwet period (maximum aerial density = 394.75/km 2 ; median linear density 155.17/km shoreline), and then gradually declined. Similarly, the population density of metamorphs peaked during the early-wet season (maximum aerial density = 6102.84/km 2 ; median linear density = 3476.25/km shoreline) and abruptly declined thereafter. Peaks followed by declines in aerial densities can be explained by surges of new toads attributable to recruitment, and the constant rate of toad mortalities. In addition to these demographic rates, seasonal fluctuations in adult and juvenile linear densities can be explained by seasonal emigrationimmigration of toads from rehydration sites to foraging habitats/spawning areas, and vice versa. This does not occur in metamorphs because of their dependence on waterbodies for rehydration, and thereby, inability to emigrate.
The detection probability of post-metamorphic cane toads along shorelines of spawning areas fluctuated seasonally, reaching an average of 23.8% during the wet season, and then falling to an average of 11.34% during the dry season (Fig. 3a). Decreasing detection probability during the dry season can be explained primarily by reduced activity rate, which is a behavioural response to the risk of dying from desiccation. Meanwhile, the adult sex ratio along shorelines of spawning areas was highly male-biased Fig. 2 Daily trend in aerial (individuals/km 2 ; left panel) and linear (individuals/ km shoreline; right panel) population density of terrestrial phase cane toads. Solid black lines indicate the mean population density, taken across replicates of a baseline scenario, whose raw values are shown in grey. Red circles indicate the reference pattern (Smart et al. 2020). Broken black lines indicate one standard deviation Fig. 3 Daily post-metamorph cane toad detection probability (a) and adult cane toad sex ratio (b) along the shorelines of spawning areas. The solid black trend lines represent the mean of detection probabilities observed over all replicates of a baseline scenario. Red circles indicate the reference pattern (Hearnden 1991;Smart et al. 2020) during the wet season, with an average of 93.8% of adults being males. This ratio decreased during the dry season to 61.7% (Fig. 3b). Male-biased sex ratios during the wet season can be explained by the contrasting breeding and foraging behaviours between female and male adults. As mentioned above, simulated adult females only breed until they spawn, whereas simulated adult males engage in breeding activities multiple times throughout a breeding season (Hearnden 1991;Cohen 1995). In this sense, adult females would spend much less time along spawning areas during the wet or breeding season, relative to their male counterparts. Meanwhile, the subsequent diminishing bias in adult sex ratios can be explained by the cane toad's dependence on rehydration sites during the dry season, which causes males and females to aggregate along rehydration sites such as spawning areas.
Visual analysis of daily density maps indicated a "pulsating" pattern of distribution throughout the landscape (Fig. 4; also see animation in Online Resource S3). During the dry season, terrestrial phase toads are concentrated in or near rehydration sites, and then disperse throughout the landscape shortly after the beginning of the wet season. This seasonal fluctuation in spatial densities is also apparent in our Warren's D similarity analysis, which revealed dissimilar distributions between wet and dry periods among all age-sex classes ( Figure  S1.12). Taken together, densities, bias maps, as well as relative frequencies observed in different habitats, indicate a seasonal pattern of spatial segregation of age-sex classes, which becomes most evident during the wet season. During the wet season, adult females dominated within human habitations and woodland habitats; adult males dominated the interiors of bushland habitats where spawning areas were most abundant; juveniles were most abundant in the margins of bushlands; and metamorphs were most abundant in bushland interiors, near natal sites. Meanwhile, during the dry season, spatial segregation was less evident, with all agesex classes aggregating near and within rehydration sites (Fig. 5).

Relative effect of management strategies
Management scenarios affected linear densities across age-sex classes, with the most significant impact being observed in the year of implementation ( Fig. 6; see also Figure S2.1-2). The magnitude of these effects increased with effort.
Hand-collecting and trapping post-metamorphs had the most significant impacts on cane toad population densities. As expected, a daily effort sustained throughout a full year resulted in a near-complete eradication of toads across age-sex classes. However, we found that similar results could be obtained when hand-capturing post-metamorphs at least once per week (also over an entire year). In contrast, trapping a The mean density (individuals/km 2 habitat) and relative frequency of toad age-sex classes in different habitat types and in different time periods, measured across replicates of a baseline scenario. b The mean inter-class pairwise bias maps of cane toad age-sex classes during late-wet and mid-dry periods, measured across replicates of a baseline scenario. Red shades indicate that the area has more individuals of the age-sex class in the column facet (X-biased), and coversely, blue shades indicate more individuals of the age-sex class in the row facet (Y-biased). Polygons indicate the margins of spawning areas post-metamorphs required more effort (daily for an entire wet season) to attain a similar impact. We observed that, in either case, hand-collecting and trapping post-metamorphs must proceed for at least seven days to accrue any benefits.
The simulated trapping and chemical suppression of tadpoles impacted metamorph and juvenile linear densities, but had little consequence for adults. While these effects did increase with effort, the benefits from deploying tadpole traps in at least 25% of the shoreline were almost as great as those accrued from managing the entire shoreline. In contrast, we found that chemical suppressants had little impact unless deployed along at least 75% of the shoreline. Interestingly, we found that both management strategies slightly increased adult linear densities. The resulting reduction of detectable tadpoles by tadpole traps and chemical suppressants made the spawning area relatively attractive and thus lured prospecting breeding adults, resulting in increased linear density.
The impacts of fencing waterbodies on cane toad densities were least anticipated. To accrue any benefits, the fencing of waterbodies must be implemented throughout the entire wet season or yearlong. Notably, even the slightest breach rendered the fences completely ineffective. Both fencing scenarios reduced the recruitment rate of metamorphs to zero, but they also increased adult linear densities. Since fencing did not involve culling, and death resulting from desiccation was unlikely due to the availability of alternative rehydration sites, the sensory cues emanating from waterbodies resulted in the accumulation of adults prospecting for virgin spawning areas or rehydration sites in Beatrice Hill Farm's mosaic landscape. In this way, the addition of fencing to waterbodies produced the observed increase in adult linear densities.
We found that only the hand-capturing and trapping of post-metamorphic toads had an effect on adult linear densities that persisted beyond the year of implementation. Any long-term consequences of fencing waterbodies or trapping and chemically suppressing tadpoles were likely swamped by the high rates of post-metamorph immigration from neighbouring source habitats.

Sensitivity analysis
The results from our sensitivity analysis were fairly consistent across all of the model outputs that we examined. In summary, most of the parameters examined had little effect on cane toad population growth, density, and/or distribution. As expected, emergent measures of cane toad population dynamics were the most sensitive to variations in traits that directly impacted fitness, such as those controlling development, survival, and reproduction ( Figure  S1.18-21). The most sensitive parameters include the terrestrial phase growth rate, the adult survival rate, and the number of tadpole development days (when not exposed to chemical suppressants). Notably, these parameters are relatively well-understood. Interestingly, parameters controlling water ecophysiology, specifically the distance at which cane toads can sense waterbodies, and the rate of water loss from the skin, had a strong influence on distributions during the dry season ( Figure S1.20-21), but not on population growth rate or density.

Discussion
We developed virToad, a mechanistic and biologically nuanced individual-based model, for forecasting the cane toad's local-to landscape-scale spatiotemporal population dynamics, and its responses to management strategies implemented across a spectrum of efforts. virToad's simulated toads help identify points of vulnerability in the species' life history (e.g., water ecophysiology, mechanisms underpinning kin selection during the aquatic phase, sensory cues) that are influenced by the environment and can be exploited by management strategies. In response to mechanistic environmental constraints Fig. 6 Relative effect of management scenarios on cane toad population dynamics. The magnitude of effect of a management scenario (rows) was evaluated in the year of implementation, plus one and five years after implementation (columns). Effects were computed as the relative difference between an age-sex classes' average daily linear density under a management scenario vs. the baseline scenario. DR refers to deflectance rate. See Online Resource S2 for more details ◂ and management actions, simulated toads employ fitness-maximising behavioural responses, which collectively give rise to the emergence of population dynamics. Its capability to produce defensible forecasts of population dynamics and management outcomes makes virToad a useful decision support tool capable of meeting the established need for models that can inform on-the-ground management responses to the cane toad's invasions.
Structurally realistic framework virToad's simulated cane toad population dynamics matched key patterns observed in real-world toad populations (see Sect. 6 and 7 of Online Resource S1). The aerial and linear densities of simulated terrestrial phase toads were well within the ranges observed in the wild (Table S1. 3;Alford et al. 1995;Cohen 1995;Child et al. 2008;Greenlees et al. 2020). And seasonal fluctuations in simulated adult linear densities closely resembled findings of recent localscale surveys by Smart et al. (2020). The seasonally pulsating spatial distribution of simulated cane toads emulated a long-standing hypothesis by Phillips et al. (2007), that toads are confined near rehydration sites during the dry season, and dispersed throughout the landscape during the wet season (see also Brown et al. 2011;Tingley and Shine 2011;Greenlees et al. 2018). Simulated toads also exhibited temporally dynamic spatial segregation among age-sex classes, which was most evident in the highly skewed adult male-bias in spawning areas during the wet season, in agreement with patterns observed in wild toads (González-Bernal et al. 2015Hearnden 1991). Coincidentally, as is the case in the wild (see Smart et al. 2020), the detection probability of simulated toads in spawning areas tended to be highest during the wet season and lowest during the dry season. Under the rationale of Pattern-Oriented Modelling, in which models are considered structurally realistic if multiple patterns have been at least weakly reproduced (Grimm and Railsback 2012), we assert that virToad provides a structurally realistic platform for predicting cane toad population dynamics.
virToad's forecasts of cane toad population dynamics emerge from a suite of mechanistic fitnessmaximising behavioural responses to dynamic environmental conditions and social interactions. In the model, as is the case in the wild, water availability in the environment, and kin selection during the aquatic phase, affect the cane toad's growth, development, physiology, reproduction, and survival (Hearnden 1991;Jørgensen 1991;Schwarzkopf and Alford 2002;Tingley and Shine 2011;Clarke et al. 2015;Yasumiba et al. 2016). As a consequence of this design, our simulated cane toads are able to modify their behaviour, with the goal of maximising their fitness. The implications of these fitness-maximising behaviours were reflected in our sensitivity analysis, since parameters controlling the cane toad's water ecophysiology affected population distribution, but had a negligible impact on population growth and densities. These findings align with expectations that cane toads cope well with hydric challenges, and avoid death due to desiccation by modifying their behaviours and distributions (Carpenter and Gillingham 1987;Jørgensen 1991;Brown et al. 2011;Florance et al. 2011). Similarly, our results consistently indicated that adult linear densities increased under management scenarios that inhibited the recruitment of new toads (fencing waterbodies, and trapping and chemically suppressing tadpoles). This emergent model property matches our understanding of how kin selection during the aquatic phase motivates breeding adults to disperse and search for virgin breeding sites where they can maximising the survivability of their offspring Hudson et al. 2015). While the consequences of environmental constraints and kin selection on individual cane toad fitness and behaviour are well recognized, their population-scale implications only become evident when extrapolated across larger spatiotemporal scales within a mechanistic IBM such as virToad.
Identifying points of vulnerability in the cane toad's life history has been the holy grail of research on the species' biology, ecology, and invasions . These points of vulnerability have long inspired the development and implementation of management responses to cane toad invasions. For example, the traditional management focus on spawning areas (Florance et al. 2011;Shine et al. 2018) reflects our understanding that the cane toad's life cycle revolves around waterbodies. Managers artificially mimic acoustic and visual cues that cane toads associate with optimal breeding and foraging habitats to lure post-metamorphs into traps (Muller and Schwarzkopf 2017;Yeager et al. 2014). Chemical cues associated with cane toad eggs have been isolated and used to trap cannibalistic tadpoles, and tadpole pheromones are being examined as a tool for suppressing aquatic phase growth and development McCann et al. 2019b). Unlike previous cane toad models, virToad is able to directly simulate these and other key life history vulnerabilities, and the management strategies intended to exploit them. virToad's ability to extrapolate such a wide array of nuanced biological, ecological, and behavioural detail across time and space is precisely what makes its predictions of cane toad population dynamics plausible, defensible, and useful.
Decision support virToad's ability to flexibly simulate competing management scenarios, and produce realistic predictions at spatiotemporal scales relevant to on-the-ground management, make it an invaluable decision support tool (Matzek et al. 2015;Schuwirth et al. 2019). Our simulation experiments revealed that management strategies targeting reproductive adults, such as handcapturing and trapping post-metamorphs, had the most potential to suppress cane toad invasions, supporting hypotheses of previous studies (Lampo and De Leo 1998;González-Bernal et al. 2015;Vimercati et al. 2021). Perhaps more importantly, our simulation experiments demonstrate that well-designed management responses can generate highly impactful longterm outcomes, even at low levels of effort; meanwhile, high effort but relatively ineffective strategies can prove futile. For example, the benefit of handcapturing post-metamorphs once a week (52 days, low effort) was comparable to the results obtained from doing so every day (365 days, high effort), and it produced a better outcome than daily interventions conducted throughout the dry season (173 days, moderate effort). Meanwhile, hand-capturing postmetamorphs, when implemented for a single night in a year (1 day, extremely low effort), which typifies community-led interventions, was found to be completely ineffective (affirming Greenlees et al. 2020).
We developed virToad to inform cane toad management across different landscapes and environments. To achieve this goal, we premised virToad's design on basic principles that are not specific to any location or ecological system (i.e., water availability in the environment and kin selection during the aquatic phase ;Smith 1964;Wells 2007;Willmer et al. 2009). Additionally, virToad can simulate cane toads with variable life-history traits and behaviours characteristic of different populations in Australia. Due to its design and functionality, virToad can be used to simulate cane toad invasions and management scenarios across the species current and tentative future invaded ranges in Australia and potentially other parts of the world (see Sect. 8 of Online Resource S1 for examples). And doing so is straightforward, given virToad's minimal and accessible input data requirements (Fig. 1). For instance, users intending to control cane toad invasions in the Kimberley-Pilbara corridor of Western Australia can simulate cane toads from the invasion front population in the Kimberley region. Similarly, supposing similar traits and behaviours (Ward-Fear et al. 2016), researchers might employ virToad to simulate invasions in Hawaii using cane toads from long-established populations in Queensland.
Naturally, the effectiveness of competing management scenarios will vary in different landscapes and environments, and as a consequence of different management objectives (i.e., eradicating founding populations, suppressing established populations, or controlling spread). Moreover, simultaneously using combinations of management strategies will likely result in amplified effectiveness. For example, our simulation experiments showed that adult densities increased after fencing waterbodies. Because no direct culling was involved, and death due to desiccation was unlikely due to alternative rehydration sites in Beatrice Hill Farm's mosaic landscape, the sensory cues emanating from waterbodies would result in an accumulation of adults prospecting for virgin spawning or rehydration sites. To reduce cane toad populations in mosaic landscapes, fencing waterbodies should be implemented in tandem with other methods, particularly hand-capturing of post-metamorphs. Nevertheless, fencing waterbodies is arguably most effective in arid areas, such as the current western invasion front along the Kimberley-Pilbara corridor in Western Australia, where precluding cane toads from accessing scarce rehydration sites leaves them vulnerable to desiccation (Florance et al. 2011;Tingley et al. 2013;Southwell et al. 2017;Tingley et al. 2017). In light of these real-world concerns, virToad provides an indispensable platform enabling researchers and practitioners the ability to explore and deliberate cane toad management responses in different landscapes, environments, and management contexts.
In Australia, cane toad management is typically implemented at local scales, and tends not to be repeated at regular intervals (Florance et al. 2011;Shine et al. 2018;Greenlees et al. 2020). Our results suggest that such small-scale interventions generate negligible long-term benefits, even with a considerable level of effort. For instance, though the localised hand-capturing of post-metamorphs implemented every day for a year resulted in the near-complete eradication of the target population, this local population began re-establishing itself the following year. This exemplifies a core effect whereby reinvasions at the peripheries continually undermine local management efforts (Glen et al. 2013). virToad, given its construction within HexSim (Schumaker and Brookes 2018), provides a platform that is uniquely well-suited for studying the cane toad's demographic processes in light of landscape connectivity (e.g., immigration-emigration rates, source-sink dynamics), and for addressing the implications of these dynamics for management. For example, the identification and targeting of highly permeable and/or high traffic corridors linking source habitats to management areas will be necessary for lowering immigration rates and, thereby, having any hope of lowering reinvasion pressure (Rollins et al. 2011;Glen et al. 2013;Lurgi et al. 2016). And critically, this landscape-scale approach to cane toad management will be necessary to protect vulnerable wildlife and their habitats from future cane toad invasions (Glen et al. 2013).

Conclusion
Our virToad model is an essential new addition to the existing suite of decision support tools intended for predicting and responding to the environmental and socio-economic consequences wrought by cane toad invasions. virToad will distinguish itself from previous models due to its ability to realistically predict the cane toad's responses, across space and time, to mechanistic constraints and management actions. vir-Toad's high degree of structural realism is the source of its usefulness as a forecasting tool. Importantly, virToad is a freely-available open-source tool that researchers and practitioners can use to explore, optimise, and prioritise local-to landscape-scale management plans intended to control cane toad invasions. In the context of the mounting problem of alien species invasions globally, we hope that virToad's framework will inspire and facilitate research into other systems, especially emerging invasive alien toad species (e.g., Duttaphrynus melanosticus, Sclerophrys gutturalis, Bufo bufo; van Wilgen et al. 2018).