Data analyses
To study which factors affected the overall invasion extent of exotic aquatic plants in the sites, we used regression analysis to examine the relationships between the site-level mean biomass of exotic aquatic plants (MBSE ) and the site-level richness of native plants (Rnative_site ), site-level performance of native plants [Pnative_site , a metric integrating MBSN andMCSN by principal component analysis (PCA), data values were scores of PC1 (proportion variances = 0.86)], T ,Nwater , Hsite ,Lsite-town , and Ahabitat . For each regression analysis, we trialed linear regression and nonlinear regression (quadratic and cubic) models and selected the best-fitting models with the highest fitting coefficient (R2 ) and a significant P-value (P<0.05).
To study the differences in the effects of the same factors on different life-form exotic plants, we selected three exotic plants belonging to different life-forms with the highest frequency in our field surveys: emergent plant Alternanthera philoxeroides (in 483 quadrats and145 sites), free-floating plant Eichhornia crassipes (in 515 quadrats and 87 sites), submerged plant Cabomba caroliniana (in 170 quadrats and 25 sites). We used linear regression analysis to examine the relationships between the quadrat-level population biomass of the three exotic plants and the quadrat-level richness of native plants (Rnative_quadrat ), quadrat-level performance of native plants (Pnative_quadrat , a metric integrating TBQN andTCQ N by PCA, data values were scores of PC1 and proportion variances ranged from 0.74~0.87), quadrat-level performance of co-occurring exotic plants (Pexotic_quadrat , a metric integratingTBQE and TCQ E of all exotic plants except the target exotic plant by PCA, data values were scores of PC1 and proportion variances ranged from 0.82~0.98), T , Nwater ,Hquadrat , and Dleaf .
To investigate how native plants, water nutrient status, climate, habitat features, herbivory level, and anthropogenic disturbance influence the overall invasion extent of exotic aquatic plants and the invasion extent of different life-form plants, we employed structural equation modeling (SEM) to explore the relationships between the factors and MBSE and quadrat-level population biomass ofA. philoxeroides , E. crassipes , and C. carolinianaindividually. For MBSE , we constructed a priori models and hypothesized that (1) all abiotic factors (T ,Nwater , Hsite ,Lsite-town , Ahabitat ) would influence the MBSE ,Rnative_site , andPnative_site , (2)Rnative_site andPnative_site would further influence theMBSE , and Rnative_sitewould predict Pnative_site , (3)Nwater would be predicted byAhabitat , Hsite , andLsite-town (Appendix S1: Fig. S2A). For biomass of A. philoxeroides , E. crassipes , and C. caroliniana , we constructed a priori models and hypothesized that (1) T , Nwater ,Hquadrat ,Rnative_quadrat ,Pnative_quadrat , andPexotic_quadrat would influence biomass of the three exotics, (2) T , Nwater ,Hquadrat ,Rnative_quadrat , andPnative_quadrat also would influencePexotic_quadrat , (3) bothRnative_quadrat , andPnative_quadrat would be predicted by T ,Nwater , and Hquadrat , andRnative_quadrat would predict thePnative_ quadrat(Appendix S1: Fig. S2B). TheDleafwas not included in the models because it neither reduced the biomass of exotic plants nor affected other factors.
The SEM models were fitted by the PiecewiseSEM R-package (Lefcheck 2016). Before modeling of the data, we inspected the data for outliers and multi-collinearity. To facilitate interpretation and ensure model convergence, we standardized all variables (except the richness of native plants) to have a mean of zero and a standard deviation. In the sub-models of Pnative_site andMBSE , we included the random effect of the sampling area to control the effect of the sampling area on the biomass of exotic and native plants. In all sub-models of quadrat-level response variables, we included random effects of the site to account for the non-independence of quadrats within the same site. Both site-level and quadrat-level richness of native plants were modeled with generalized linear models (GLMs) or generalized mixed-effects models (GLMMs) using a Poisson distribution and log-link function. The other response variables were modeled using linear models (LMs) or linear mixed-effects models (LMEs). All LMs, LMEs, GLMs, and GLMMs were performed by the “lme4” package (Bates et al. 2014). To obtain the best model, we removed the nonsignificant paths with the highest p-value in a stepwise manner. We added any link that was initially ignored when the D-separation test revealed significant missing relationships between variables, and including the significant missing path would improve the model fit. We assessed the goodness-of-fit of the model (Fisher’s C statistic) using Shipley’s test of directional separation. Models with an adequate fit (p > 0.05) were considered as candidate models, and their AICs were computed and compared. The model with the lowest AIC value was considered as the best-fit model. Standardized coefficients for the GLM and GLMM components of the piecewise SEM (not provided by PiecewiseSEM) were calculated by the latent theoretic approach (Lefcheck 2019).
To determine the relative importance of various factors that influence the overall invasion extent of exotic aquatic plants and the invasion extent of different life-form plants, we employed hierarchical partitioning (hier.part R package, Walsh and MacNally 2015) to analyze the independent explanatory power (I) ofRnative_site ,Pnative_site , T ,Nwater , Hsite ,Lsite-town , and Ahabitatfor the MBSE , and the independent explanatory power (I) of T , Nwater ,Hquadrat ,Rnative_quadrat ,Pnative_quadrat , andPexotic_quadrat for the biomass ofA. philoxeroides , E. crassipes , and C. caroliniana . The statistical significance of each I was determined by randomizing the data matrix 200 times to produce the distribution of I values, with observed I value that exceeded the 95th percentile considered significant. The results of the significance tests were expressed as Z-scores. All analyses were conducted by R v4.2.0 (R Development Team 2022). The biomass of A. philoxeroides , E. crassipes , and C. caroliniana was transformed using Ln (x) function.