Populations of most species exhibit temporal fluctuations in size that may or may not coincide across space. On one hand, external drivers of temporal change that are shared across spatially distant populations can operate to produce synchrony in the highs and lows of abundance (Moran , Ranta et al. ). Conversely, environmental drivers, or intrinsic differences in population traits, may be local to each population and hence result in spatially asynchronous population responses or complex dynamics (Hugueny ). Unraveling and understanding the spatial scales over which species and populations respond similarly to extrinsic factors remain key focuses of population biology. Gaining such understanding is also of practical importance for managing populations facing environmental change and for managing multiple species that share a common environment. For example, if populations behave synchronously the risk of extinction will be higher if all populations decrease to zero at the same time and rescue through colonization is prohibited (Heino et al. ).
Sub‐populations can be synchronized by dispersal, mobile predators, and the Moran effect (Liebhold et al. , Krebs ). Dispersal is limited by organisms' ability to move through the landscape and thus, if dispersal is driving synchrony, the strength of population synchrony can be expected to decline with increasing distance between sub‐populations (Engen et al. , Hugueny ). Dispersal is less likely to play a critical role in synchronizing populations over large distances (>10 km) (Krebs ), especially for sedentary taxa. It can be nonetheless an important mechanism, even if dispersal ability is limited, for driving the scale of population synchrony in populations with low levels of density dependence (Lande et al. , Hugueny ).
Ydenberg () proposed that predation may synchronize populations of microtine rodents, prompting Ims and Steen () to demonstrate—theoretically—that a nomadic predator can synchronize populations if the predation rate is high enough. Experimental evidence later confirmed the theoretical predictions of Ims and Steen () and also found that dispersal events often enhance predation risk, inducing a feedback‐loop between predation and population growth that enhances synchrony (Ims and Andreassen , Huitu et al. , Korpimäki et al. ).
Moran's theorem states that sub‐populations sharing a common structure of density dependence should be synchronized by a spatially correlated density‐independent factor, such as climate (Moran , Ranta et al. ). Thus, if the density‐dependent structure is linear, the correlation in numbers between regional sub‐populations will be equal to that between the local density‐independent factor (Koenig , Hugueny ). However, the correlation between the density‐independent factor is often not equal to that of the sub‐populations, suggesting nonlinear dynamics or that the areal scale over which the sub‐populations occur does not match the scale of operation of the density‐independent factor (Hugueny , Sæther et al. ). Nonlinear dynamics can produce complex spatial sub‐population patterns, but simple linear models still may approximate the stochastic dynamics of many populations (Hugueny ). In addition, both linear and nonlinear models show that sub‐populations occurring in regions with spatial heterogeneity in local density‐independent factors can exhibit decreases in population synchrony (Engen and Sæther ).
Spatially synchronous population dynamics have been much studied in high latitude environments where extrinsic factors such as severe winter weather and mobile predators may operate over large spatial scales (Merritt , Krebs ). In other environments, however, the spatial sub‐structuring of populations may be less predictable, with dynamics driven by different density‐dependence structures or environmental heterogeneity. In arid regions, for example, where water is a key limiting factor (Noy‐Meir ), local landscape features such as drought refuges or geomorphological structures that buffer the arid climate can ensure that sub‐populations are independent from each other and exhibit different dynamics (Dickman et al. , Ferrer et al. ). Sub‐populations near ephemeral waters thus may exhibit increased survival or greater densities of animals than those in more xeric sites (Kok and Nel ), in turn setting the local dynamics at these sites on different trajectories (Céré et al. ). In addition, intrinsic factors such as the density‐dependence structure of each sub‐population could further de‐synchronize spatial synchrony and lead to complex patterns in spatial population dynamics (Sæther et al. ). In other arid regions, by contrast, regional conditions may act as powerful population‐synchronizing factors. Sub‐populations that fall within the same rainfall zone, for example, may experience pulses of productivity that generate similar bottom‐up increases in their numbers (Haythornthwaite and Dickman , Woodman et al. ), and those that experience the same disturbance event, such as wildfire (Letnic et al. , Evans et al. ) or predation (Letnic et al. ), may also exhibit similar population trajectories. Because arid environments often exhibit strongly varied biotic and abiotic conditions in space and time, they provide ideal arenas in which to explore the factors that shape the spatial dynamics of populations.
In this study, we attempt to disentangle the suite of factors that drive the dynamics of small mammals in climatically unpredictable landscapes in arid central Australia. Many of these mammals range over large distances (>1000 km) throughout Australia's arid central regions, and sub‐populations can irrupt at similar times (Strahan , Wardle et al. ). For example, widely separated rodent populations often share similar timing in “boom” and “bust” dynamics (Dickman et al. , D'Souza et al. , Greenville et al. , Pavey and Nano ). The mechanism for this apparent synchrony is largely unknown. Dispersal alone is unlikely to explain population synchrony as the dispersal ability of many small mammals in arid Australia is limited to ~14 km (Dickman et al. , Haythornthwaite and Dickman ). Rather than dispersing large distances (>20 km), many species survive in refuges throughout the landscape during population busts (Dickman et al. , Pavey et al. ). In addition, predator populations may be limited by the generally low levels of primary production, reducing their influence on small mammal population synchrony. In contrast to many northern hemisphere small mammal populations, for example, some rodents in Australia such as the house mouse (Mus musculus) do not experience high predatory mortality during population irruptions and their numbers are driven primary by rainfall, food supply and disease (Korpimäki et al. ). Rare, yet large rainfall events (>400 mm) and subsequent rapid increases in primary productivity are suspected to drive population booms in many organisms in central Australia (Dickman et al. , Morton et al. , Greenville et al. , Pavey and Nano , Wardle et al. ). Rainfall can exhibit strong regional synchrony across distances up to 5000 km (Koenig ), and thus may drive population synchrony across large areas, suggesting the Moran effect may operate in this environment.
Using long‐term population data on five species of small mammals at nine widely separated sites, we first test the hypothesis that annual rainfall drives regional synchrony. For species that do not exhibit synchrony, we then explore whether sub‐populations of our study mammals are asynchronous, form two sub‐populations at ephemeral water sources or open desert sites (oasis hypothesis) or at burnt vs. unburnt sites (wildfire hypothesis), or form three sub‐populations that are organized by shared rainfall gradients (productivity hypothesis). We use multivariate state‐space models to identify potential drivers, such as local rainfall, food resources, plant cover and species interactions, which may influence the populations and sub‐populations of the five small mammal species. Using Moran's theorem and assuming that the density‐dependence structure for the small mammal sub‐populations is linear, we also predict that the spatial correlation between the regional sub‐populations will be equal to that between the local density‐independent conditions (annual rainfall). Finally, for species that show synchronous spatial dynamics, we use multivariate state‐space models to predict that rainfall, operating at the landscape scale (8000 km2), will be positively associated with their populations, whereas species with asynchronous sub‐populations will be influenced largely by factors operating at local scales.
Understanding the spatial dynamics of multiple species that share the common desert environment is expected to provide vital information on how best to manage populations of these species, allowing land managers to target conservation activities at the appropriate spatial scale and identify species that are most at risk of decline. This information should also assist in predicting how populations will respond to more frequent and more extreme weather events such as droughts and flood rains, which are predicted to affect many arid regions in future (Greenville et al. , IPCC ).
The study was carried out in the Simpson Desert, central Australia (Appendix S1: Fig. S1). This region occupies 170,000 km2; dune fields comprise 73% of the region, with smaller areas consisting of clay pans, rocky outcrops and gibber flats (Shephard ). The sand dunes run parallel in a north‐south direction aligned with the prevailing southerly wind, are up to 10 m high, and spaced 0.6–1 km apart (Purdie ). Vegetation in the interdune swales and on dune sides is predominantly spinifex (Triodia basedowii) grassland with small stands of gidgee trees (Acacia georginae) and other woody Acacia shrubs or mallee eucalypts; low‐lying clay pans fill with water temporarily after heavy rain (Wardle et al. ).
During summer, daily temperatures usually exceed 40°C and minima in winter fall below 5°C (Purdie ). Highest rainfall occurs in summer, but heavy rains can fall locally or regionally throughout the year. Long‐term weather stations in the study region at Glenormiston (1890–2011), Boulia (1888–2011) and Birdsville (1954–2011) have median annual rainfalls of 186 mm (n = 121 yr), 216.2 mm (n = 123 yr), and 153.1 mm (n = 57 yr), respectively (Bureau of Meteorology ). In general, higher rainfall is experienced in the north of the study region than in the south, and there is a weak rainfall gradient from east to west (Bureau of Meteorology ).
Live‐trapping was carried out at nine sites across Carlo Station, Tobermorey Station, Cravens Peak and Ethabuka Reserves, and covered an 8000 km2 area of the north‐eastern Simpson Desert in south‐western Queensland (Appendix S1: Fig. S1). Dickman et al. () found that the maximum dispersal distances for the species investigated below were: 14 km for rodents in this region over a 2.5 yr study (mean maximum movement: 6.34 km), and 12 km for dasyurid marsupials (mean maximum movement: 1.04 km). Thus, each site was set >20 km apart in spinifex grassland to eliminate the effect of dispersal between sites (Appendix S1: Fig. S1). Small mammals were live‐trapped using pitfall traps (16 cm diameter, 60 cm deep), each equipped with a 5 m drift fence of aluminum flywire to increase trap efficiency (Friend et al. ). Pitfalls were arranged on grids comprising six lines of six traps spaced 20 m apart to cover 1 ha (i.e., 36 pitfall traps per grid). The top line of traps was positioned on a dune crest and the bottom line 100 m distant in the swale so that each grid sampled the topography of the dune field. Sites contained 2–12 grids and grids within sites were set 0.5–2 km apart in randomly chosen positions.
Traps were opened 2–6 times a year from 1990–2012 at one site (Main Camp) and from 1995–2012 at eight more sites (Shitty Site, Field River South, Field River North, South Site, Kunnamuka Swamp East, Carlo, Tobermorey East and Tobermorey West). Each trap was opened for 2–6 nights per sampling occasion. To account for unequal trapping effort, live‐capture counts were standardized per 100 trap nights (TN: trap nights = traps × nights opened) and averaged for each year. Not all sites were opened every year, but were included in the time series as missing values. Live‐trapping data were log‐transformed (log+1) as the population models used below are in log‐space.
Long‐term (17–22 yr; 130 sampling trips) live‐trapping data (205 524 trap nights) yielded six species of rodents: Pseudomys hermannsburgensis (sandy inland mouse; 7878 captures), Notomys alexis (spinifex hopping mouse; 5146 captures), Pseudomys desertor (desert mouse; 1299 captures), Rattus villosissimus (long‐haired rat; 1170 captures), Mus musculus (house mouse; 1120 captures), and Leggadina forresti (short‐tailed mouse; four captures), and eight species of dasyurid marsupials: Dasycercus blythi (brush‐tailed mulgara; 853 captures), Ningaui ridei (wongai ningaui; 863 captures), Sminthopsis youngsoni (lesser hairy‐footed dunnart; 2491 captures), Sminthopsis hirtipes (hairy‐footed dunnart; 334 captures), Sminthopsis macroura (striped‐faced dunnart; 17 captures), Sminthopsis crassicaudata (fat‐tailed dunnart; 14 captures), Planigale gilesi (Giles' planigale; two captures), and Planigale tenuirostris (narrow‐nosed planigale; 12 captures). Of these species, only P. hermannsburgensis (12 g) and N. alexis (35 g), the carnivorous D. blythi (100 g), and the insectivorous S. youngsoni (10 g) and N. ridei (8 g) had sufficient time series in both length and numbers of captures across nine sites to run the population models below. Nonetheless, the nine rarer species are reported here to highlight the advantage of long‐term and intensive monitoring in detecting species richness and are the focus of other research questions (e.g., Greenville et al. ).
Large‐scale wildfires (>1000 km2) have occurred three times in the study region since 1972, and the mean wildfire return interval for the region is 26 yr (Greenville et al. , NAFI ). To investigate whether wildfire affected the trajectory of sub‐populations of the study species, the wildfire history at the site of each sub‐population was classified according to whether it had experienced a wildfire (six sites classified as burnt, see Appendix S2) or remained unburnt (three sites classified as unburnt, Appendix S2) across the 22‐yr study period.
Daily data from automatic weather stations (Environdata, Warwick, Queensland) at each site (Appendix S2: Fig. S1) were used to calculate total annual (calendar year) rainfall antecedent to capture. Weather stations were active from 1995–2012. To obtain rainfall data pre‐1995 for Main Camp, where trapping started in 1990, we averaged daily rainfall records from the closest weather stations at Glenormiston, Sandringham, Boulia, Bedourie, and Birdsville (Bureau of Meteorology ). The 1‐yr lag was used to account for the time required for small mammals to respond, via breeding, to a large rainfall event and has been used successfully to predict mammal captures in previous studies (Greenville et al. ).
To measure cover of the dominant vegetation, spinifex (T. basedowii), we scored ground cover visually as a percentage in 2.5 m radius plots around the same six traps on each small mammal trapping grid on every sampling occasion. In addition, a seed productivity index (0–5, where 0 represents no seeding and 5 represents all plants seeding profusely) was used to score the degree of spinifex seeding in each plot. We chose spinifex seed due to its dominance in the landscape and because it is a key component in the diet of the study rodents, representing up to 52% of their diet by frequency of occurrence (Murray and Dickman ). Estimates of cover and seed productivity were pooled for each site, per year, over 17–22 yr for each of the nine sites.
To measure the degree of regional synchrony in annual rainfall and mammal sub‐populations, we calculated the zero‐lag cross‐correlation for the abundance (log) of each species and the annual rainfall (mm) experienced by the nine sub‐populations, with 1000 resamples. In addition, to investigate if the cross‐correlation varied over distance between sub‐populations, we calculated the spatial cross‐correlations for the abundance (log) of each species and annual rainfall (mm). The mean cross‐correlation coefficients and spatial cross‐correlations were calculated using the ncf 1.1‐5 package (Bjornstad ), in R 3.02 (R Core Team ). If regional synchrony was found, a synchronous (1‐state) autoregressive state‐space population model was used to identify potential population drivers (see below). If there was little evidence for population synchrony, the sub‐population structures of the small mammal populations, and their potential drivers, were investigated by using multivariate autoregressive state‐space (MARSS) models as outlined below.
We used Bayesian multivariate autoregressive state‐space (MARSS) models to analyze live‐trapping data from our nine sites and identify potential population drivers (rainfall, spinifex cover, food and predator populations). We based the MARSS models on the Gompertz population growth model and assumed that sub‐population growth rate varied exponentially with sub‐population size and that metapopulations are closed to immigration and emigration within the 8000 km2 study region (Hinrichsen and Holmes ). The MARSS framework is hierarchical and allows modeling of different spatial population structures and parameters, such as density dependence, while including both process (state) and observation variability (Ward et al. ). Process variability represents temporal variability in population size due to environmental and demographic stochasticity (Ward et al. ). Observation variability includes sampling error (e.g., temporal changes in detectability or error resulting in only a sub‐sample of the population being counted) (Hinrichsen and Holmes , Ward et al. ). The MARSS framework has proven useful for modeling population data in this system, particularly for incorporating observation error in long‐term monitoring results (Nguyen et al. ).
The process component is a multivariate first‐order autoregressive process and is written in log‐space (Holmes et al. ,b) as: [Image Omitted. See PDF]
where Xt represents a vector of all m sub‐populations (up to nine sub‐populations for each species) at time t (year), and u is a vector of length m. B and Q are matrices that denote process parameters. The B diagonal elements (Bi,i) represent the coefficients of autoregression in the populations through time and represent the strength of density dependence (diagonal element Bi,i = 1 represents density independence, Bi,i < 1 = density dependence). The off‐diagonal elements in the B matrix allow for interactions between processes, such as between species populations (Holmes et al. ). The parameter u describes the trend of the sub‐population. C represents coefficients of the covariates and ct are the covariates through time t. wt denotes process errors, which we assumed to be independent and to follow a multivariate normal distribution with a mean of zero and variance‐covariance matrix Q (i.e., Q diagonals represent process variance and off‐diagonals represent how the sub‐populations are correlated with one another). We allowed Q to have different variances across sites (diagonals, q2i,i, were allowed to vary, while off‐diagonals were set to zero), as process variation was assumed to differ across space (see Appendix S3 for model specifications).
The observation component, written in log‐space (Holmes et al. ,b), is: [Image Omitted. See PDF]
where Yt represents a vector of all observations at n sites at time t, a denotes the mean bias between sites, and Z is an n × m matrix of 0s and 1s that assigns observations to a sub‐population structure. The number of sites (n) may be different from the number of sub‐populations (m) at time t (Ward et al. ). Observation errors, vt,, are assumed to be uncorrelated and follow a multivariate normal distribution, with a mean of zero and a variance‐covariance matrix R. We set the diagonals in R (r2i,i) to have equal variance across sites for each species (but differ between species), as the same trapping methods were employed throughout the study and thus we did not expect differences in trappability. In addition, the pooling of captures for each year minimized differences in trappability between seasons and age classes. The R off‐diagonals were set to zero (see Appendix S3 for detailed information on model specification). Eq. 1 and 2 comprise the MARSS model.
For each species, we tested four possible spatial population structures: asynchronous, oasis, productivity, and wildfire hypotheses (each hypothesis was modeled by assigning a matrix of 1s and 0s via the Z‐matrix; see Appendix S2 and Appendix S3 for site selection and model formulation).
Covariates that might influence population dynamics were as follows. For P. hermannsburgensis and N. alexis, spinifex cover, the seed productivity index from the prior year, annual rainfall from the prior year and capture rates of the predatory D. blythi were included in the model. For the dasyurids S. youngsoni and N. ridei, spinifex cover, annual rainfall in the prior year, and D. blythi were entered, the latter as a potential predator and competitor (Dickman ). For D. blythi, spinifex cover and annual rainfall from the prior year were entered, as well as total rodent population size in the prior year as rodents are a key prey source (Chen et al. ). The 1‐yr lag for food resources (seed or rodent prey) was used to account for the time required for small mammals to respond, via breeding, to a resource‐pulse event (Greenville et al. ). All variables were standardized by subtracting the mean and dividing by the standard deviation (z‐scored) to allow direct comparisons between covariates to be made. De‐meaning the data (i.e., via z‐score transformation) allows the parameters u and a to be set to zero, which reduces computational problems when trying to estimate both u and B simultaneously (Holmes et al. ). In addition, by de‐meaning the data and setting u = 0, we focused specifically on the effects of the covariates on sub‐population estimates (see Appendix S3 for model formulation).
Covariates (rainfall, spinifex cover and seed) measured for each sub‐population were modeled via the C matrix and sub‐populations of the predatory D. blythi and rodent prey via the B matrix. The interaction between each species sub‐population (rodent × predator (D. blythi), small dasyurid × predator/competitor (D. blythi) and D. blythi × prey [rodents]) was assessed by examining the off‐diagonals of the B matrix (Holmes et al. ). Modeling the species interactions through the off‐diagonals of the B matrix allowed for better estimates of the species populations that were used to investigate the interactions as both process and observation errors are taken into account. In addition, the off‐diagonals of the B matrix assume that predator populations can affect prey, and vice versa, providing a more accurate biological description of predator‐prey dynamics. The diagonal elements of the B matrix (Bi,i) were assessed for density dependence. Thus, four models were run for each species and the best‐fitting model was chosen by the lowest posterior predictive loss value (D). Posterior predictive loss is a measure of the predictive ability of each model, penalized for model complexity (Gelfand and Ghosh ), and can be more appropriate for spatial and temporal hierarchical models than other model selection methods (reviewed in Hooten and Hobbs ).
We calculated 95% credible intervals for each parameter and considered these significant if they did not cross zero. To allow for model convergence, each model was fitted with three Markov chains, each having 300 000 Markov Chain Monte Carlo (MCMC) iterations, thinning of 10 and the first 280 000 iterations discarded, leaving 6000 iterations saved (see Appendix S3 for description of priors, model formation and Supplement S1 for JAGS code). Inspection of diagnostic plots indicated that all models met statistical assumptions and that models converged. Analyses were conducted in R 3.03 (R Core Team ) and JAGS 3.4.0 (Plummer ), using the R2jags 0.04‐03 (Su and Yajima ) and rjags 3‐13 (Plummer ) packages.
The mean cross‐correlation coefficient for annual rainfall across the nine study sites was 0.77 (CI: 0.65–0.86) and did not change with distance between sites (Fig. a), suggesting regional synchrony in annual rainfall across our 8000 km2 study region. In general, the timing of large rainfall events (>400 mm) was similar across the nine sites (Appendix S2: Fig. S1).
Enlarge this image.Spatial cross‐correlations (solid lines) from 1995–2012 across nine sites in the Simpson Desert, central Australia (map), of (a) annual rainfall, (b) Pseudomys hermannsburgensis, (c) Notomys alexis, (d) Dasycercus blythi, (e) Sminthopsis youngsoni and (f) Ningaui ridei. Dashed lines represent 95% confidence intervals, and dotted lines represent zero correlation. Illustrations by Alison Foster.
In total, 16 617 captures of rodents and 4586 captures of dasyurids were made over the period of study, yielding trap success rates of 8.09% and 2.23%, respectively. Both P. hermannsburgensis and N. alexis exhibited similar spatial population structure. The mean cross‐correlation coefficients for these species were 0.77 (CI: 0.66–0.86) and 0.78 (CI: 0.64–0.90), respectively, suggesting regional population synchrony; this did not decrease with distance (Fig. b,c). As predicted from the Moran effect, the mean spatial cross‐correlation coefficients for the two rodent species were equal to the mean cross‐correlation for annual rainfall, suggesting that rainfall may play a role in synchronizing rodent populations over the study region. The MARSS models showed that populations of both rodent species were affected positively by spinifex seed productivity and annual rainfall in the previous year (Table , Fig. a,b). Both species also had similar process and observation errors, indicating similarity in population dynamics and trappability (Table ). The two rodent species showed some evidence of density dependence (Table ).
The results of single population MARSS models for captures (per 100 trap nights, log +1 transformed) of three species of small mammals (Pseudomys hermannsburgensis, Notomys alexis, and Dasycercus blythi) in the Simpson Desert, central Australia. Data were z‐scored so that direct comparisons can be made. Density dependence occurs if the diagonal element Bi,i is less than one. Covariates were considered significant if the 95% credible intervals (CI) did not cross zero, and are shown in bold| Model covariates | Estimate | Lower CI | Upper CI |
| P. hermannsburgensis | |||
| Density dependence (Bii) | 0.33 | 0.03 | 0.72 |
| Annual rainfall year before | 0.52 | 0.30 | 0.77 |
| Spinifex cover | −0.24 | −0.58 | 0.12 |
| Spinifex seed | 0.52 | 0.26 | 0.78 |
| Mulgara (B1,2) | −0.38 | −0.87 | 0.08 |
| Process error (q2) | 0.51 | 0.33 | 0.79 |
| Observation error (r2) | 0.48 | 0.42 | 0.54 |
| N. alexis | |||
| Density dependence (Bii) | 0.31 | 0.02 | 0.75 |
| Annual rainfall year before | 0.55 | 0.34 | 0.75 |
| Spinifex cover | −0.19 | −0.44 | 0.08 |
| Spinifex seed | 0.40 | 0.14 | 0.65 |
| Mulgara (B1,2) | 0.04 | −0.51 | 0.48 |
| Process error (q2) | 0.39 | 0.24 | 0.61 |
| Observation error (r2) | 0.47 | 0.41 | 0.53 |
| D. blythi | |||
| Density dependence (Bii) | 0.35 | 0.05 | 0.71 |
| Annual rainfall year before | 0.41 | 0.20 | 0.62 |
| Spinifex cover | 0.22 | 0.002 | 0.46 |
| Rodents (B1,2) | 1.02 | 0.64 | 1.37 |
| Process error (q2) | 0.41 | 0.120 | 0.68 |
| Observation error (r2) | 0.67 | 0.58 | 0.76 |
Enlarge this image.Predicted population sizes (lines) and captures (dots; captures/100 trap nights) predicted from MARSS models for (a) Pseudomys hermannsburgensis, (b) Notomys alexis, and (c) Dasycercus blythi synchronous population models, Simpson Desert, central Australia (map). Time‐series data were collected from nine sites (sub‐populations) monitored 2–6 times per year for 17–22 yr (1990–2012 for one site and 1995–2012 for eight sites; 130 sampling trips). Shaded areas indicate 95% credible intervals. Illustrations by Alison Foster.
The three species of dasyurids showed very different temporal and spatial dynamics. The mean cross‐correlation coefficients for D. blythi, S. youngsoni and N. ridei were 0.58 (CI: 0.44–0.74), 0.26 (CI: 0.10–0.36) and 0.08 (CI: −0.10–0.31), respectively, suggesting regional synchrony for the predatory D. blythi, but little evidence of this for the smaller S. youngsoni and N. ridei. The mean cross‐correlation coefficient did not change at larger distances for D. blythi (Fig. d), and for S. youngsoni and N. ridei the mean cross‐correlation remained low across increasing distance between sites (Fig. e,f). The mean spatial cross‐correlation 95% confidence interval for D. blythi overlapped with rainfall, suggesting that rainfall may play a role in synchronizing sub‐populations of this species across the study region (Fig. c). In contrast, there was no evidence that regional rainfall was important for local populations of S. youngsoni and N. ridei.
As predicted, the synchronous (1‐state) MARSS population model for D. blythi showed a positive relationship between D. blythi captures and annual rainfall in the prior year (Table ). For both S. youngsoni and N. ridei, the best‐fitting MARSS model suggested that sub‐populations at all nine sites were distinct (asynchronous) and exhibited density dependence (Tables and ; Figs. and ). There was no support for the oasis, wildfire or productivity models (Table ).
Posterior predictive loss (D) values for MARSS models describing four possible spatial sub‐population structures for two species of small dasyurid marsupials (Sminthopsis youngsoni and Ningaui ridei) in the Simpson Desert, central Australia. Data are based on 17–22 yr of live‐trapping from nine sites. The smallest posterior predictive loss (D), highlighted in bold, indicates the best‐fitting model| Model | Posterior predictive loss (D) |
| S. youngsoni | |
| Asynchronous | 307 |
| Oasis | 458 |
| Productivity | 445 |
| Wildfire | 419 |
| N. ridei | |
| Asynchronous | 348 |
| Oasis | 526 |
| Productivity | 497 |
| Wildfire | 480 |
| Model covariates | Rainfall | Mulgara (B1,2) | Spinifex cover | Density dependence (Bi,i) | Process error (q2) | Observation error (r2) |
| S. youngsoni | ||||||
| Population | ||||||
| Main Camp | −0.08 (−0.46–0.29) | −0.59 (−0.83–−0.37) | 0.10 (−0.32–0.52) | 0.24 (0.03–0.45) | 0.53 (0.12–0.99) | 0.56 (0.32–0.75) |
| Carlo | −0.07 (−0.47–0.34) | −0.34 (−0.72–0.04) | 0.38 (0.08–0.91) | |||
| Field River South | −0.04 (−0.42–0.33) | 0.02 (−0.38–0.46) | 0.41 (0.09–0.90) | |||
| Kunnamuka Swamp East | 0.25 (−0.16–0.66) | −0.31 (−0.73–0.23) | 0.33 (0.08–0.86) | |||
| Shitty Site | 0.28 (−0.28–0.83) | 0.40 (−0.23–0.94) | 0.59 (0.10–1.31) | |||
| South Site | −0.10 (−0.66–0.43) | 0.65 (−0.13–1.23) | 0.43 (0.08–1.08) | |||
| Field River North | 0.59 (0.16–1.02) | −0.08 (−0.48–0.35) | 0.30 (0.07–0.75) | |||
| Tobermorey East | −0.44 (−1.03–0.19) | −0.08 (−0.73–0.58) | 0.67 (0.13–1.35) | |||
| Tobermorey West | −0.28 (−0.78–0.19) | −0.53 (−1.00–−0.01) | 0.28 (0.07–0.76) | |||
| N. ridei | ||||||
| Population | ||||||
| Main Camp | −0.35 (−0.67–−0.03) | 0.33 (0.12–0.54) | 0.34 (0.03–0.63) | 0.47 (0.22–0.69) | 0.33 (0.08–0.74) | 0.64 (0.51–0.77) |
| Carlo | −0.16 (−0.63–0.29) | 0.19 (−0.21–0.54) | 0.49 (0.12–0.99) | |||
| Field River South | 0.43 (0.06–0.81) | 0.23 (−0.10–0.57) | 0.30 (0.08–0.74) | |||
| Kunnamuka Swamp East | −0.39 (−0.89–0.06) | −0.26 (−0.75–0.19) | 0.55 (0.09–1.26) | |||
| Shitty Site | 0.55 (0.07–1.01) | 0.06 (−0.37–0.54) | 0.38 (0.08–0.96) | |||
| South Site | 0.28 (−0.23–0.74) | −0.16 (−0.75–0.46) | 0.38 (0.08–0.95) | |||
| Field River North | 0.26 (−0.17–0.68) | 0.27 (−0.12–0.66) | 0.31 (0.07–0.83) | |||
| Tobermorey East | −0.11 (−0.55–0.29) | 0.29 (−0.08–0.68) | 0.29 (0.07–0.81) | |||
| Tobermorey West | −0.43 (−0.86–0.01) | 0.16 (−0.34–0.65) | 0.35 (0.08–0.89) |
Enlarge this image.Predicted population sizes (lines) and captures (dots; captures/100 trap nights) from a MARSS model for Sminthopsis youngsoni depicting asynchronous spatial population dynamics, Simpson Desert, central Australia. Time‐series data were collected from nine sites (sub‐populations) monitored 2–6 times per year for 17–22 yr (1990–2012 for one site and 1995–2012 for eight sites; 130 sampling trips). Shaded areas indicate 95% credible intervals. Illustrations by Alison Foster.
Enlarge this image.Predicted population sizes (lines) and captures (dots; captures/100 trap nights) from a MARSS model for Ningaui ridei depicting asynchronous spatial population dynamics, Simpson Desert, central Australia. Time‐series data were collected from nine sites (sub‐populations) monitored 2–6 times per year for 17–22 yr (1990–2012 for one site and 1995–2012 for eight sites; 130 sampling trips). Shaded areas indicate 95% credible intervals. Illustrations by Alison Foster.
Observation errors and process errors were similar across sub‐populations for the three species of dasyurids, suggesting similar trappability, but there were large variations in the process errors for all three species (Tables and ).
For both S. youngsoni and N. ridei the effect of spinifex cover and antecedent annual rainfall were site‐specific, suggesting that drivers for these species sub‐populations were acting on a local scale (Table ). However, populations of D. blythi had contrasting effects on the smaller dasyurids. Sminthopsis youngsoni had a negative interaction with D. blythi but populations of N. ridei were positively associated with D. blythi (Table ).
Our results provide evidence that rainfall synchronizes the population fluctuations of some, but not all, species of small mammals in our arid study system. Thus, both species of rodents exhibited highly synchronous spatial dynamics, whereas the three species of dasyurid marsupials showed very different, and quite disparate, spatial dynamics. The predatory D. blythi showed evidence of spatial population synchrony, whereas the smaller S. youngsoni and N. ridei had asynchronous dynamics. Studies in other world deserts have also described population fluctuations of rodents in response to rain‐driven pulses of productivity (e.g., Shenbrot et al. , Thibault et al. , Meserve et al. ), although most of these studies were carried out at local sites rather than at the regional scale that we surveyed here. Some work has described rain‐driven population fluctuations in insectivorous desert mammals, although the association with rainfall is often more tenuous or more distant in this group than in rodents (Chung‐MacCoubrey et al. , Thibault et al. ). Rainfall has occasionally been linked also with declines in populations of insectivorous mammals (Woolley ), and virtually all work on this group has been based at local sites.
Of the three mechanisms proposed to explain spatial synchrony in population dynamics, dispersal between sub‐populations is the least likely to apply in our study system. In the first instance, dispersal between our study sites is highly improbable as the distance between the sites was greater than the maximum dispersal distance recorded for the study species (Dickman et al. , Letnic ). In addition, the distances moved by rodents decrease markedly during and after large rainfall events as animals settle to exploit flushes in local food resources (Dickman et al. ).
Predation is another mechanism that can drive spatial synchrony (Ydenberg , Ims and Andreassen ), and potentially contributed to the pattern of population synchrony documented here for both P. hermannsburgensis and N. alexis. On one hand, there was no relationship between these two rodents and D. blythi, although this marsupial includes both species prominently in its diet (Chen et al. ), thus suggesting that predation from this small carnivore does not affect the rodent populations. On the other hand, larger and highly mobile introduced mammalian predators may have stronger synchronizing effects, especially with respect to the decline and bust phases of the rodents' population cycles when per capita predation on both species is intense. Working in the same study region, Spencer et al. (,b) showed that N. alexis does not respond to odor cues of large mammalian predators—most notably the feral cat Felis catus and European red fox Vulpes vulpes—that might otherwise provide an indication of predation risk, and also that these ubiquitous predators consume very large numbers of the study rodents during the peak through to the bust phases of their populations. In other regions mobile avian predators act as synchronizing agents (Ims and Andreassen , Pavey et al. ), but we have little evidence that predators contributed to the spatial synchrony in the rodent populations observed here.
The third mechanism producing spatial synchrony in population dynamics, the Moran effect, was almost certainly an important organizing force for the two study rodents. The spatial autocorrelation in rainfall was very similar to that for the sub‐populations of both rodent species, and the MARSS models revealed further that rainfall was a significant driver for both species. Large rainfall events in arid Australia typically cover >10,000 km2 (Greenville et al. ) and hence are large enough to generate region‐wide increases in primary productivity and produce the spatial synchrony observed here. Indeed, we found both rodent species to be associated positively with both rainfall and with spinifex seed production in the prior year. These results confirm the importance of bottom‐up effects in driving population booms in desert rodents (e.g., Dickman et al. , Greenville et al. ) and suggest that the Moran effect can operate on some species in arid systems via the effects of regional rainfall. It is likely that declining post‐rain resources and predation act to reduce populations across regional areas, supporting interpretations drawn from both arid systems (e.g., Letnic et al. ) and more mesic environments (e.g., Korpimäki et al. , ) elsewhere. Density dependence also was found in our study rodents, suggesting that both intrinsic and extrinsic factors are important in regulating their population dynamics.
Complex and varied relationships between species populations and rainfall, food resources, plant cover and species interactions were found for the three marsupial species, and provide the focus of further discussion below. In the first instance, the carnivorous D. blythi exhibited synchronous spatial dynamics, with spinifex cover and antecedent rodent populations both associated positively with D. blythi captures. Previous research has reported inconsistent effects of spinifex cover on captures of this species. Masters () found that more D. blythi were caught in areas with high spinifex cover than in recently burnt areas with low cover. However, reduced cover due to spinifex harvesting has mixed effects on this species (Masters et al. , McCarthy and Masters ), and home ranges are little affected by wildfire when spinifex cover is diminished and fragmented (Körtner et al. ). Our results suggest that spinifex cover was important for D. blythi after rainfall, perhaps as a means of reducing exposure to predatory feral cats and red foxes, with the positive response to antecedent rodent populations reflecting greater access of D. blythi to its prey.
The asynchronous hypothesis was the best‐fitting model for both S. youngsoni and N. ridei. This suggests that at the nine sites surveyed, each sub‐population exhibited its own trajectory. Antecedent annual rainfall and spinifex cover had an inconsistent effect on individual sub‐populations of S. youngsoni and N. ridei. There was a positive effect of antecedent annual rainfall on S. youngsoni at one site, but no relationship at the other eight. Rainfall had positive, negative and no effect on individual sub‐populations of N. ridei. Rainfall has been identified as an important predictor of captures in S. youngsoni (Masters ), presumably as rainfall influences the availability of invertebrate prey, but our study suggests that this relationship may be site‐ or time‐specific. Indeed, long‐term surveys of invertebrates in our study region have revealed little effect of short‐term rainfall on invertebrate activity, but disparate effects at different spatial scales between invertebrate taxa up to 18 months after rain (Kwok et al. ). It is thus possible that both N. ridei and S. youngsoni were responding to site‐specific temporal differences in their prey.
In contrast to findings in previous studies, spinifex cover had no influence, or a negative influence, on S. youngsoni (Masters , Dickman et al. ) and was associated positively with just one sub‐population of N. ridei. Although spinifex hummocks have been reported to provide shelter and foraging sites for dasyurids hunting invertebrate prey (e.g., Bos ), inconsistent relationships between invertebrates and spinifex cover (Kwok et al. ) and the availability of other microhabitat components (e.g., shrubs, leaf litter) at most sites may have reduced the dependence of small dasyurids on spinifex alone. As noted above, we also detected an association between D. blythi and spinifex, and it is plausible that the presence of this carnivorous marsupial further influenced the suitability of spinifex for shelter or foraging in the two smaller species. There is some evidence that this was so.
Captures of D. blythi had a consistently negative effect on those of S. youngsoni and a positive effect on those of N. ridei. Dasycercus blythi preys, in part, upon large (6–10 mm long) invertebrates (Chen et al. ) and thus may compete with other dasyurid species, such as S. youngsoni, for food. However, competition for food with N. ridei is likely to be minimal as this species selects small invertebrates (<5 mm long) and suffers sharply reduced rates of energy intake when hunting larger prey (Fisher and Dickman ). By virtue of its small size and cone‐shaped head, Ningaui ridei likely has exclusive access to the center of spinifex hummocks (Bos ) where small invertebrates are most abundant (Fisher and Dickman ). By contrast, S. youngsoni uses the edges of spinifex hummocks where larger prey occur but encounters with D. blythi also are more likely. The larger D. blythi is strongly dominant in encounters with S. youngsoni (Dickman ), and hence may limit the activity of S. youngsoni via encounter or exploitation competition. In addition, as S. youngsoni in turn is dominant over smaller dasyurids (Dickman ), small species such as N. ridei may benefit indirectly from the presence of D. blythi via the suppressive effects of this species on the mid‐sized dunnarts. Such effects are equivalent to the suppressive influence of apex predators in other systems, where meso‐predator activity is reduced and small prey species derive indirect benefits from the presence of the apex predator (e.g., Soulé et al. , Brashares et al. ). Our results suggest that complex interactions may occur among micro‐carnivores in the desert system, as well as between these species, their habitats and site‐specific factors, but species manipulation experiments would be useful to disentangle and explore these interactions further.
In addition to interspecific interactions, all study species showed evidence of density dependence, and this also may play a role in driving their local populations. Density dependence commonly regulates animal populations (Brook and Bradshaw , Knape and de Valpine ), including populations of rodents in desert habitats (Shenbrot et al. ). However, perhaps because of the marked fluctuations in population size and occasional tendency for populations in local sites to crash to zero, due to strong abiotic factors (Fig. a, b), it has not previously been demonstrated in any arid Australian taxa. Both species of rodents and all species of dasyurids showed moderate to strong density dependence. Several factors can affect estimates of density dependence, including the length of the time series, the stage of the population cycle when a species is surveyed (e.g., growth phase), whether the species is at carrying capacity (Brook and Bradshaw ), and observation error (Knape and de Valpine ). Estimates of density dependence may also arise as statistical artifacts in time series as growth rates often decrease when populations are large and increase when they are small, that is, they show a statistical return tendency to some point (or long‐term mean) and thus may not be tied to an underlying biological mechanism (Wolda and Dennis ). However, in this study we used multiple long‐term time series (17–22 yr) for each species, which encompassed all stages of their population cycles and, by using the MARSS statistical framework, we incorporated observation error. Thus, we suggest that density dependence occurs in our study system, although its prevalence across all taxa precludes a specific role in driving spatial synchrony.
The spatial population dynamics of small mammals in our study system differed strongly between and within families, illustrating that varied life‐history strategies are employed by mammals that share a common desert environment. Sub‐populations fluctuated synchronously if driven by large‐scale events (e.g., rainfall), or were asynchronous, most likely if driven by local events. These results suggest that the Moran effect operates for some but not all species in this arid system. Based on our findings and observations elsewhere (Kelt , Shenbrot ), desert rodents generally appear most likely to exhibit spatial synchrony owing to their direct response to rain‐driven regional increases in productivity, whereas secondary consumers may fluctuate more idiosyncratically in response to local site factors. As arid lands globally are predicted to face climatic shifts that will exacerbate rainfall‐drought cycles (IPCC ), environmental stochasticity may have profound effects on species populations that exhibit synchrony across large spatial areas. For example, the risk of extinction for populations that behave synchronously will be higher if all populations decrease to zero at the same time due to a change in an important climate variable (Heino et al. ). Future work investigating the relative contributions of both density dependence and environmental autocorrelation in driving spatial synchrony—or lack thereof—in small mammal populations will be valuable for predicting species responses to a changing climate at different spatial scales (e.g., Hugueny , Sæther et al. ).
We thank Bush Heritage Australia, H. Jukes, G. McDonald, D. Smith, and G. Woods for allowing access to the properties in the study region; B. Tamayo, D. Nelson, C.‐L. Beh and other members of the Desert Ecology Research Group, and many volunteers for valuable assistance in the field. The original manuscript was significantly improved from comments provided by Elizabeth Holmes and two anonymous referees. Funding was provided by the Australian Research Council, an Australian Postgraduate Award (to ACG) and the Australian Government's Terrestrial Ecosystems Research Network (
© 2016. This work is published under http://creativecommons.org/licenses/by/3.0/ (the “License”). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.