DOI: 10.5061/dryad.gmsbcc2sh
title: Data from: Environmental variables influence patterns of mammal co-occurrence following introduced predator control

description: Co-occurring species often overlap in resource use and can interact in
 complex ways. However, shifts in environmental conditions or resource
 availability can lead to changes in patterns of species co-occurrence,
 which may be exacerbated by global escalation of human disturbances to
 ecosystems, including conservation directed alterations. We investigated
 the relative abundance and co-occurrence of two naturally sympatric mammal
 species following two forms of environmental disturbance: wildfire and
 introduced predator control. Using 14 years of abundance data from repeat
 surveys at long-term monitoring sites in south-eastern Australia, we
 examined the association between a marsupial, the common brushtail possum
 Trichosurus vulpecula, and a co-occurring native rodent, the bush rat
 Rattus fuscipes. We asked: Is the increase in abundance of common
 brushtail possums following control of an introduced predator associated
 with a decline in abundance of the bush rats? Using Bayesian regression
 models, we tested hypotheses that the abundance of each species would vary
 with changes in environmental and disturbance variables, and that the
 negative association between bush rats and common brushtail possums was
 stronger than the association between bush rats and disturbance. Our
 analyses revealed that bush rat abundance varied greatly in relation to
 environmental and disturbance variables, whereas common brushtail possums
 showed relatively limited variation in response to the same variables.
 There was a negative association between common brushtail possums and bush
 rats, but this association was weaker than the initial decline and
 subsequent recovery of bush rats in response to wildfires. Using
 co-occurrence analysis, we can readily infer negative relationships in
 abundance between co-occurring species, but to understand the impacts of
 such associations, and plan appropriate conservation measures, we require
 more information on interactions between the species and environmental
 variables. Co-occurrence can be a powerful and novel method to diagnose
 threats to communities and understand changes in ecosystem dynamics.
description: Study location We used data from long-term
 annual monitoring that commenced in 2003 in Booderee National Park (BNP),
 Jervis Bay Territory, south-eastern Australia. BNP is on Indigenous land
 and is jointly managed by the Wreck Bay Aboriginal Community and Parks
 Australia. The 6600 ha park has a temperate climate, with an average
 annual rainfall of 1213 mm, spread evenly across the year (Bureau of
 Meteorology, 2021). Average temperatures range from 18.6-25.1°C in summer
 (January) to 9.9-16.1°C in winter (July) (Bureau of Meteorology, 2021).
 BNP supports a range of vegetation types such as heathlands, wetlands,
 forests, and woodlands. Two major fires have occurred in BNP over the past
 20 years (in 2003 and 2017), with each burning approximately half of the
 park. A fox baiting program has been in place in BNP since 1999 and was
 intensified in 2003 to reduce the deleterious impact of this introduced
 predator on native prey species (Dexter et al., 2013).
 Data collection We surveyed small and
 medium-sized mammals annually each summer for 14 years at 109 permanent
 sites starting in 2003, with another 20 sites added in 2008. The sites
 were surveyed along 100 m transects with 2 large (30 x 30 x 60 cm) cage
 traps at the beginning and end of transects, small (20 x 20 x 50 cm) cage
 traps every 20 m between the large cage traps, and 10 Elliott traps (10 x
 10 x 30 cm; Elliott Scientific Equipment, Upwey, Victoria) every 10 m
 (Figs 2) (Lindenmayer et al., 2008, 2016). Approximately 50% of the sites
 were surveyed each year (with the other 50% being surveyed the next year),
 depending on weather conditions (Lindenmayer et al., 2016). We recorded
 the number of individuals of both species caught at each site in a given
 year.                 We collected data on
 environmental and disturbance variables at each of our 129 sites. These
 data included visual estimates of the percentage of understorey and leaf
 litter cover in four 1 x 1 m subplots within 20 x 20 m survey plots during
 semi-annual vegetation surveys (Lindenmayer et al., 2008). We selected
 understorey and leaf litter as representative variables of the primary
 bush rat habitat, which are also components of common brushtail possum
 habitat (Callander, 2018; Cruz et al., 2012). We constructed a predictive
 model to fill the data gaps for those years when sites were not surveyed
 (MacGregor et al., 2020). We used monthly rainfall averages collected at
 the nearby Point Perpendicular weather station for the trapping period at
 each site (Bureau of Meteorology, 2021). We transformed both the
 vegetation variables and rainfall were transformed into quadratic
 functions using the poly() function in R (R Core Team, 2021).
                 We used data on fire occurrence recorded
 on-ground since 2003, and fire history data collected by Booderee National
 Park over the past 50 years (Foster et al., 2017), specifically the number
 of years since the last fire at a site. To minimise possible inaccuracies
 stemming from incorrect fire dates, or the occurrence of unreported fires,
 we grouped the number of years since fire into 10-year blocks (i.e., 0-10,
 11-20, 21-30, 30+ years). Statistical
 analysis We used Bayesian regression models with a
 hurdle step to test the response of species abundances to the selected
 variables using the brms package ver. 2.16.3 (Bürkner, 2017; 2108; Feng,
 2021) implemented in R (R Core Team, 2021). These regression models
 employed Markov Chain Monte Carlo simulations, with four chains and a
 warm-up of 1000 iterations before sampling another 1000 iterations. We
 assessed model convergence by ensuring all Rhat values were &lt;1.1
 (Bürkner, 2017; 2108). The hurdle step consisted of two components: the
 first modelled the presence/absence of the response variable, and the
 second, conditional on the species being present, modelled the conditional
 abundance using a zero-truncated Poisson (Feng, 2021). We combined these
 two components in an analysis of unconditional abundance (Feng,
 2021). We created a regression model with bush rats as
 the response variable, and a regression model with common brushtail
 possums as the response variable. Both regression models included time,
 understorey cover, leaf litter cover and rainfall as covariates within the
 conditional abundance component. These variables were included to assess
 the variation in bush rat and common brushtail possum abundances with
 environmental variables (H1). Years since fire was included in the
 conditional abundance and hurdle step of both regression models as an
 explanatory variable, as it is a prominent disturbance within BNP, and
 past research has indicated that fire has a significant effect on small
 vertebrate populations (Arthur et al., 2012). The other species was also
 input into the conditional abundance and hurdle step of both regression
 models (i.e., common brushtail possums into the bush rat model, bush rats
 into the common brushtail possum model) as an explanatory variable to
 assess the co-occurrence effect between species (H2). We also included
 site as a random effect, and used the log of the number of Elliott traps
 as a control for the bush rat models, and the number of open cage traps as
 a control for the common brushtail possum models. The control variables
 account for varying trapping effort, and were selected based on the main
 trap-type that captures the relevant species (i.e., Elliott traps for bush
 rats, cage traps for common brushtail possums). We
 performed a model selection procedure for both of the regression models,
 based on the selection for explanatory variables only. We chose not to
 perform model selection on the covariates (i.e., the environmental
 variables) as we were testing variation in species abundance in relation
 to the environment, and not predicting significant changes in abundance
 that we were with the co-occurrence and disturbance variables. Using model
 selection, we assessed the relevancy of our exploratory variables to
 changes in species abundance (H3). The chosen model was the most
 parsimonious, which was based on the simplest model which was within 2
 leave-one-out cross validation (LOOIC) scores of the best fitting
 model. We created ten variations of the regression
 models for each species, and assessed the fit of each variation using
 LOOIC (Tables 1 and 2) (Vehtari et al., 2017). LOOIC estimates the
 out-of-sample predictive fit by measuring the predictive accuracy for each
 data point using a variation of the expected log pointwise predictive
 density equation (Vehtari et al., 2017). LOOIC was selected as the
 appropriate method over other model selection methods as it is informative
 and was created for Bayesian models (Burnham & Anderson, 2002;
 Vehtari et al., 2017).
description: # Title of Dataset: Trapping data for co-occurrence analysis *** The
 trapping information per site per year, including the number of common
 brushtail possums and bush rats captured at that site for that year, and
 the site information for that year. ## Description of the data and file
 structure The information is organised by site, with each line
 representing the year that information was collected. The column named
 'T.vulpecula' represents the number of common brushtail possums
 captured, and the column named 'R.fuscipes' represents the
 number of bush rats captured. These are occumpanied by the columns
 'T.vulpeculaPA' and 'R.fuscipesPA' which are the
 presence/absence columns. The columns 'Elliotts.Open' and
 'Cages.Open' were the total number of each trap type used for
 each site that year. 'Elliott' is the logged number of Elliott
 traps used (Cage was an unused variable). The 'Start.Date',
 'Month' and 'Survey.Season' columns include
 information about when the trapping took place each year. The
 'veg' column shows the broad veg class of the site. The
 'Rainfall', 'USpred' and 'LLpred' are the
 covariates used in the co-occurrence model, with 'USpred'
 representing the predicted understorey cover, and 'LLpred'
 representing the predicted leaf litter cover. The
 'year.of.last.fire' and 'no..since.fire' were used to
 create the 'years.since.fire' variable that was used as an
 explanatory variable in the model. The 'na' within the
 'year.of.last.fire' and 'no..since.fire' columns
 represent fires that were so old, they do not have a accurate date, these
 have been translated to greater than 30 years in the
 'years.since.fire' variable. The 'Year0' column is the
 number of years since trapping commenced (in 2002). ## Sharing/Access
 information Links to other publicly accessible locations of the data: Was
 data derived from another source? No ## Code/Software All code was from R.
 The package used to create and run the models was the brms package version
 2.16.3 #Bush rat models There were 11 variations of R.fuscipes models
 created. All model variations included R.fuscipes as the response
 variable. All model variations had a hurdle step, characterised by
 'hu = '. The family was hurdle\_poisson(), and the number of
 cores was 4. A null model was created using 1 \+ \(1\|site\.number\) in
 both the population level and hurdle sections of the model\. Model
 variations were defined by the use of the T.vulpecula and years.since.fire
 variables. The Rainfall, USpred, and LLpred variables were inputted into
 models as quadratic variables using the 'poly(x,2)' function.
 The Elliott variable was inputted as a control for effort\, and the
 site\.number was inputted as a repeated function;
 '\(1\|site\.number\)'\. Model 4 (where common brushtail possums
 were only inputted into the conditional abundance step) was the most
 parsimonious model, and therefore used in this study. #Common brushtail
 possum models There werer 11 variations of T.vulpecula models created. All
 model variations included T.vulpecula as the response variable. All model
 variations had a hurdle step, characterised by 'hu = '. The
 family was hurdle\_poisson(), and the number of cores was 4. A null model
 was created using 1 \+ \(1\|site\.number\) in both the population level
 and hurdle sections of the model\. Model variations were defined by the
 use of the R.fuscipes and years.since.fire variables. The Rainfall,
 USpred, and LLpred variables were inputted into models as quadratic
 variables using the 'poly(x,2)' function. The Elliott variable
 was inputted as a control for effort\, and the site\.number was inputted
 as a repeated function; '\(1\|site\.number\)'\. Model 10 (which
 modelled bush rats and common brushtail possums in both the conditional
 abundance and hurdle step, but did not model years since fire) was
 selected the most parsimonious model, and therefore used in this study.
 #model selection We ran a model selection over the 11 variations of each
 species models using a LOOIC calculation. This was run using the
 'sapply(list())' function from the base package, and looping in
 the 'loo()' function. #Visualising the results We firstly
 visualised the results using forest plots, created using the
 'plot\_model()' function from the sjPlot package. We then model
 the conditional effects of the chosen plots for relevant variables using
 the 'plot(conditional\_effects())' from the brms package. All
 visualisation was done using the ggplot package.