Statistical Analyses
All counts were first standardised to colony-forming-units (cfu) per mL. Invasion success (relative invader fitness) was calculated as proportional change, v , of the proportion of invader to resident, calculated as: v = x2. (1 - x1 )/x1. (1 - x2 ), where x1 is the initial invader proportion and x2 the final (Ross-Gillespie et al. 2007). Initial invader proportion (x1 ) was calculated as the average frequency of the introduced invader:
\(x_{1}\ =\ E\left[\frac{I_{t}}{I_{t}\ +\ R_{t}}\right]\ =\ \frac{1}{3}\sum_{t\ =\ \left\{4,\ 8,\ 12\right\}}\frac{I_{t}}{I_{t}\ +\ R_{t}}\)(1)
where It is the density of the invader introduced on day t and Rt is the density of the residents getting invaded on day t . We could not measure resident density on days 8 and 12, because it would require destructive sampling of undisturbed treatments. We therefore used the resident density on day 4 and assumed that R4, R8 andR12 were equal for 1-, 2-, and 4-days disturbance treatments.
We sampled R4 for 1-, 2-, and 4-days disturbance treatments during their transfers, but we could not sampleR4 for 8- and 16-days disturbance treatments as it is a destructive process. The disturbance history up to day 4 for 8- and 16-days treatments is identical to that for 4-days treatment. We therefore assumed the resident community dynamics are the same for these three treatments, and used R4 for 4-days treatment (before the disturbance) to calculateR4 for 8- and 16-days treatments:
\(R_{4,\ 8-days}\ =\ R_{4,\ 16-days}\ =\ \frac{R_{4,\ 4-days}}{\text{Disturbance\ mortality\ rate}}=\ \frac{R_{4,\ 4-days}}{0.01}\)(2)
where Ri,j is the density of the resident on dayi under j- days disturbance treatment. Based on this calculation, we further assumed that R8,16-days = R12,16-days = R4,16-days for 16-days disturbance treatment, where R8,16-days =R on day 8 in the 16-day disturbance treatment and so forth. For 8-days disturbance treatment, we assumedR12,8-days = R4,8-days andR8,8-days = 0.01 R4,8-days to account for the disturbance event on day 8.
In order to eliminate zero inflation, one was added to the final invader density v (post volume standardisation) and was transformed to log(v +1) to normalise the residuals. A value greater than 0.69 (log(1+1)) would indicate that the invader increased in proportion throughout the experiment, whereas a value below this would suggest that invasion was unsuccessful.
To analyse the effect of disturbance and resource abundance on invasion success, v, a linear model was used to test effects of disturbance, resource abundance, and invader morphotype, with all two-way and three-way interactions. As the different morphotypes have distinct growth strategies, we expected their invasion success to be markedly different. Given a significant three-way interaction in the most complex model, we did all further analysis on each invader morphotype (SM & WS) separately.
For each invader morph, separate linear models were used to investigate treatment (disturbance frequency and resource abundance) effects on invasion success, evolved biodiversity (calculated using Simpson’s index (Simpson 1949)) and total resident density (log10(cfu+1 mL-1). Disturbance frequency was treated as a continuous predictor, whereas resource abundance was treated as categorical due to only having three levels. Model selection was done using likelihood ratio tests.
We then tested whether treatments indirectly affected invasion success through changes in resident populations. To do this we first used a model with resident biodiversity and total resident density, plus their interaction, as predictors of invasion success. We then included treatment (disturbance, resource abundance, and their interaction), alongside resident population effects as predictors of success. The models with both treatment and resident population effects were initially tested using an ANOVA with type III sums of squares, then with type II if no significant interactions were found to account for differences in the ordering of predictors on significance testing.
Post-hoc model comparisons were used to look at significant differences between levels of resource abundances and disturbance. For pairwise comparisons of single treatments (e.g. between high, medium and low resource abundances), model estimates were averaged over other predictors in the model. Where multiple pairwise comparisons were used, p values were adjusted using Bonferroni adjustments. When comparing slopes to 0, confidence intervals overlapping zero indicated no significant effect. All statistical analyses were carried out in R version 4.0.3 (R Core Team).