Data analyses
We used generalised linear mixed‐effects models to test the interactive
effects of herbivory and pollination treatments on: (a) proportion
damaged beans per plant (beans with B. rufimanus emergence holes)
(b) faba bean yield components (individual bean weight, total bean
weight per plant, number of beans per pod, number of beans per plant,
number of pods (mature, immature and unfertilized) per plant, proportion
of mature pods per plant, and yield (kg.ha-1), and (c)
plant growth components (plant height and aboveground biomass, root
length and biomass). Normal distributions were used except for number of
pods and beans per plant, where a Poisson distribution was used or a
negative binomial when data was overdispersed, and for proportion of
damaged beans and mature pods per plant where a binomial distribution
was used (see Table 1 for model structures). The explanatory variables
in all models included the H+/H- and P+/P- treatments and their
interaction term. Despite the care taken to remove B. rufimanusfrom the H- cages at the beginning of the experiment, beans with
emergence holes were also found in these cages, and there was a large
variation in damage between plants within each H+/H- cage (Fig. S3). We
therefore, in addition to the main H+/H- treatment effect (i.e.
plant-stand scale, measured by averaging individual plant scale
responses in each cage), investigated the effect of herbivory damage at
the plant scale, measured as percentage of damaged beans per plant
within cage (% Damage). Proportion of damaged beans did not vary with
pollination treatment levels (Table 1, p=0.09). We tested all plant
yield and growth variables, and used P+/P- treatment and % Damage per
plant and their interaction term as explanatory variables. The random
structure in all models included cage identity (N=28) nested within
block (N=7), except for yield (kg.ha-1) per cage where
only block was included. If significant interactions were found,
post-hoc tests using the “emmeans” package were carried out to
investigate the direction of the effect.
To test the effects of herbivory on observed pollinator behaviour
(proportion of legitimate , robbing and EFN visits) and on pollinator
visitation rate (visits per flower per time unit) we used a generalized
mixed‐effects model with a binomial and a normal distribution
respectively. The explanatory variables included H+/H- treatment and the
number of open flowers per m2 and their interaction
term. To account for addition of sugar-water to the pollinators on the
7th of July a binary factor (sugar-water: yes/no) was
included as well as its interaction with number of open flowers per
m2 and herbivory treatment. The interaction between
herbivory treatment and sugar-water was never significant and did not
improve the models as determined by the Akaike Information Criterion
(AIC), indicating that the effect of herbivory on pollinator behaviour
did not change after the addition of the sugar-water. This interaction
was therefore excluded. To investigate the effect of herbivores on
number of open flowers per m2 in cages with
pollinators, we used a generalized mixed‐effects model with H+/H-
treatment as explanatory variable. The random structure for all models
incorporated the sampling round (N=15) nested within cage identity and
block.
The residuals of all models were visually inspected to validate the
model assumptions and additionally, generalized linear models were
checked for overdispersion using “DHARMa” (Hartig & Lohse, 2020).
Multicollinearity was checked for all models (variation inflation factor
<2). All analyses were conducted in R version 3.6.3, using
packages “nlme” (Pinheiro & Bates, 2020), “lme4” (Bates et al.,
2020), “emmeans” (Lenth et al., 2021) and “ggplot2” (Wickham et al.,
2020) to plot data.