Statistical analysis
All analyses were run in the R (R Core Team, 2020) with packagesMuMIn (Bartoń, 2020) and lme4 (Bates et al., 2018). We
analysed each of the response variables separately with generalized
linear mixed-effects models. We tested the effects of urbanization,
local canopy cover and their interaction on leaf damage with Gaussian
error distribution and identity link (the results were the same with a
beta-distribution and log-link), and on the incidence of leaf-miners and
gall-inducers with binomial error distribution and logit-link in
separate models. The data was not overdispersed, visual inspection of
raw-data did not call for zero-inflated models, the distribution of
residuals met model assumptions.
In each model, Urbanization (%), Local canopy cover (%), Urbanization
× Local canopy cover, Year (as a factor), Spring temperature (°C) and
Spring precipitation (mm) were included as fixed effects; and Partner ID
as a random factor to account for the fact that some partners surveyed
multiple trees and/or several years.
We analysed the data in the framework of information theory (Burnham and
Anderson, 2002). We first built three models, one for each response
variable separately (leaf damage, gall-inducer and leaf-miner
incidences). We scaled and centred all continuous predictor variables
prior to modelling to make their coefficients comparable, and verified
that uncontrolled correlations among explanatory variables were unlikely
to bias model coefficient parameter estimates (all variance inflation
factors [VIF] lower than 2) (Zuur et al., 2009). We then applied a
procedure of parsimonious model selection based on the Akaike’s
Information Criterion corrected for small sample sizes (AICc) and
considered every model in a range of 2 units of AICc to the best model
as equally likely (Arnold, 2010). We calculated the AIC weight
(wi ) – i.e., the probability that a given model
is the best model within the set of candidate models – and also the
relative variable importance (RVI) as the sum ofwi of every model including this variable. When
multiple models were competing with the best model (i.e., when several
models with ΔAICc < 2), we implemented a multi-model inference
approach, constructing a consensus model that comprised the selected
variables from the set of best models. We subsequently averaged their
effect sizes over all models in the set of best models, utilizingwi as the weighting parameter (i.e., model
averaging). A certain predictor was deemed to have a statistically
significant effect on the response variable if its 95% confidence
interval did not bracket zero.