Phenotypic selection analysis
To quantify phenotypic selection on traits, multiple regression was used to provide estimates of selection gradients, whereas individual univariate models were used to provide estimates of selection differentials. Selection gradients (β) describe the strength of unique or direct selection acting on a trait (after controlling for inter-trait correlations), whereas selection differentials (S) describe ‘total’ selection (Lande & Arnold 1983). For the analysis of β, a multiple regression model was fitted in which relative fitness (seed number divided by its mean) was regressed on standardized traits (standardised by standard deviation) (Lande & Arnold 1983). We also added to the model the two-level factors of ‘pollination’ (control, hand pollination) and ‘herbivory’ (addition, removal), and their interaction with all traits and each other. In this way estimates of β were output for all traits in all four treatment combinations. Following Sletvold (2019), the calculated difference in trait β between treatment combinations thereby provided estimates of pollinator-mediated selection (in the presence and absence of herbivory), herbivore-mediated selection (with and without pollinator limitation), and combined pollinator- and herbivore-mediated selection. We used the ‘emtrends’ function of R package emmeans (Lenth et al. 2018) to calculate these differences in β, and to test whether the result differed significantly from zero (after adjusting p-values for multiple comparisons via Benjamini–Hochberg correction). Although the multiple regression included only 27 genotypes for which complete trait data was available, we still considered these inferences to be robust given that: 1) the full spectrum of genetic variation in defence was represented (see above), and 2) no qualitative differences were observed in univariate regressions regardless of whether some (n = 27) or all (n = 81) genotypes were used (Table S2). For the analysis of S, univariate regressions were fitted for each trait individually. Only data from the ‘control pollination/herbivore present’ treatment combination were used for the estimation of S, as this treatment combination is typically considered most representative of natural population conditions (Sletvold 2019).