Spider mite populations
Two populations of the spider mite Tetranychus urticae red form (a.k.a. T. cinnabarinus; Cruz et al., 2021) that share the same genetic background, except for the presence or absence of an allele conferring resistance to the pesticide etoxazole (i.e., a single nucleotide mutation on the chitin synthase 1 gene, CHS1; Van Leeuwen et al. 2012), were used in the experiments. The procedure used to create these populations is depicted in Fig. 1 and fully described in the sections below.
Creation of the Wi.SS populations
Five Wolbachia-infected homozygous etoxazole-susceptible (“Wi.SS”) population replicates were initially created in 2015, each by merging six compatible Wolbachia-infected populations collected in the region of Lisbon, Portugal, in 2013 (Zélé et al. 2018a, 2020b). Each population replicate was initially created with a census population size of 200 mated females and subsequently maintained with discrete generations at this same population size (cf. detailed procedure for the creation and maintenance of these populations, called “Infected controls; iC”, in Rodrigues et al. 2022a). Full homozygosity for the etoxazole susceptible allele was subsequently confirmed by PCR-RFLP (as described in Van Leeuwen et al. 2012) for each of the field-derived populations and each of the Wi.SS replicates.
Creation of the Wi.RR populations
Five Wolbachia-infected homozygous etoxazole-resistant (“Wi.RR”) population replicates were created in 2018 by introgression of an etoxazole-resistance allele into each of the Wi.SS population replicates (Fig. 1 and detailed procedure in Fig. S1). The etoxazole-resistance allele originated from the homozygous etoxazole-resistant (RR) red lab strain SB9, collected in Crete, Greece, in 2006 (Van Leeuwen et al. 2012). As coinfection with Wolbachia and Rickettsia was detected in this strain (using the multiplex PCR procedure described in Zélé et al. 2018b), and because these endosymbionts may cause reproductive incompatibilities with the Wolbachia strain that is present in the Wi.SS populations, SB9 was treated with a rifampicin solution at 0.05%, w/v for one generation as described in Zélé et al. (2020a) to create the symbiont-free donor parent strain SB9.rif. We then used a backcross procedure ensuring the complete introgression of the etoxazole-resistance allele from the SB9.rif donor strain into each of the Wi.SS recurrent parent populations, with an estimated retention of 99.6% of their nuclear background, and complete retention of their mitochondrial background and Wolbachia infection.
At the onset of the backcross procedure, SB9.rif males were crossed with females from each of the Wi.SS population replicates to obtain a 100% heterozygous F1 female progeny (hence all “Wi.RS”: Wolbachia-infected with a resistant and a susceptible allele). Given haplodiploidy in spider mites and because the allele conferring resistance to etoxazole is recessive, these females were isolated as virgin to produce “Wi.S” and “Wi.R” F2 haploid males (note that this step replaces the usual self-cross step for the introgression of recessive alleles in diploid organisms). Resistant F2 males (Wi.R) were then selected via exposure to etoxazole during the development of the male offspring, and subsequently backcrossed with females from the recurrent parent population Wi.SS to produce a new generation of heterozygous females (i.e., the first generation of backcrossed females; BC1). Heterozygous BC1 females were, as their grandmothers, isolated as virgin to produce a new generation of males (i.e., the progeny of “non crossed” females; NC1), which were again selected for etoxazole resistance during their development. This procedure was repeated for six additional backcrosses (BC2NC2 to BC7NC7). At the end of the procedure, a subset of the progeny resulting from the last generation of backcross (Wi.RS females and Wi.S males; BC7) was self-crossed (SC1) to obtain Wi.RS and Wi.SS females, which were crossed with NC7 selected cousins. Their offspring (SC2) was then exposed to etoxazole to select the homozygous resistant Wi.RR females and Wi.R males that founded the replicated populations “Wi.RR”. A census population size of at least 200 females and 100 males was used at each step of the procedure (cf. detailed procedure in Fig. S1).
Creation and merge of the Wu.SS and Wu.RR populations
The five Wi.SS and the five Wi.RR population replicates (each with a population size of 200 females) were subsequently treated for one generation with rifampicin (0.05% w/v), following the procedure described in Zélé et al. (2020a). Each newly created Wu.SS and Wu.RR population replicates was subsequently maintained in discrete generations (population size of 200 females) in absence of antibiotics (i.e., under the same conditions as the untreated populations). Three generations later, the presence/absence of Wolbachia and that of the etoxazole-resistance allele were confirmed in a pool of 100 females from each Wi.SS, Wi.RR, Wu.SS and Wu.RR population replicate by multiplex PCR (see Zélé et al. 2018b) and PCR-RFLP (see Van Leeuwen et al. 2012), respectively. All replicates from either the Wu.SS or the Wu.RR populations, were then merged in September 2020 to create the susceptible and resistant base populations used in this study, respectively.
Rearing conditions
Since their creation and merging, both Wu.SS and Wu.RR were maintained in large numbers (> 2000) under continuous generations on entire plants. These populations were kept on bean leaves (Phaseolus vulgaris, variety Contender, provided by Germisem, Portugal) under controlled conditions (25°C; photoperiod of 16L:8D). Bean plants were germinated and grown in an isolated and herbivore-free room for 14 days under controlled conditions (photoperiod of 16L at 25°C: 8D at 20°C) before being infested with spider mites.
Experimental setup
Effect of high developmental temperatures on life-history traits
First, we characterized the impact of four temperatures above the control temperature in the laboratory (control temperature: 25°C, high temperatures: 33, 35, 36 and 37°C) during development on the developmental time, survival and fertility of pesticide-susceptible spider mites (from Wu.SS).
To assess the effect of temperature on developmental time and juvenile survival, mite cohorts were created by placing on a bean leaf patch (ca. 75 cm2) 100 mated females for 2 days at 25°C, so they could oviposit synchronously. Once adult, the offspring of those females was transferred to one of the tested temperatures (25, 33, 35, 36 or 37°C) for three days, after which a maximum of 32 offspring females were individually placed on bean leaf discs (2.55 cm2) to oviposit for 24 hours at the same temperature. Egg development was followed daily, and the number of days elapsed until the first adult male and female offspring emerged (i.e., developmental time), as well as the number of eggs that reached adulthood (i.e., number of surviving juveniles) were recorded.
Due to logistic constraints, three assays were performed, always including the control temperature (25°C): i) 25°C, 33°C and 37°C; ii) 25°C and 35°C; iii) 25°C and 36°C. In all cases, replicates with damaged (i.e., hurt during manipulation) mothers were excluded from the analyses. Additionally, in the analysis of all independent variables, patches with females that did not lay eggs were excluded. Overall, we analysed 67, 65 and 82 replicates (i.e., patches) of female developmental time, male developmental time and number of surviving juveniles, respectively (see Table S1 for more details). Due to the low number of adult females surviving and ovipositing at 37°C, this temperature was not included in the statistical analyses of developmental time and number of surviving juveniles (see Table S2).
To determine the impact of heat on fertility, we focused only on female fertility, as we assumed that male fertility would also be hindered at the temperature females are affected, since males generally suffer from a reduction in fertility at lower high temperatures than females in several arthropod species (David et al. 2005; Sales et al. 2018; Iossa 2019; Zwoinska et al. 2020; van Heerwaarden and Sgrò 2021). Cohorts were created with mated females as before, but they were directly maintained at one of the tested temperatures (25, 33, 35, 36 or 37°C). Before the offspring reached adulthood, deutonymph quiescent males and females were isolated separately to ensure virginity. Once adults (ca. 24 hours later), a maximum of 40 pairs of virgin females and males were isolated on bean leaf discs (2.55 cm2) to mate and oviposit. Five treatments were established: pairs where both individuals developed at 25°C (control; ♀25 x ♂25) and pairs where only the female developed at one of the high temperatures (♀33 x ♂25, ♀35 x ♂25, ♀36 x ♂25 and ♀37 x ♂25). Pairs were maintained during 24 hours, after which males were discarded. Females were left to lay eggs for two more days, after which they were also discarded, and the eggs were counted. Leaf discs were kept at the female developmental temperature from pairing until the female was removed, being placed at 25°C afterwards. The survival status of both male and female was recorded until removal. Once reaching adulthood, offspring was counted, sexed and removed 12 and 14 days after the onset of oviposition. Daily fecundity and offspring number were thus used as proxies for female fertility.
As in the previous experiment, temperature treatments were distributed in three assays: (i) 25°C, 33°C and 37°C; (ii) 25°C and 35°C; (iii) 25°C and 36°C. For the first assay, we performed two blocks, which occurred in two consecutive days. The other two assays were done in a single block. In all cases, replicates with individuals forming the mating pairs that were damaged during manipulation were excluded from the analyses. In addition, patches with females that did not lay eggs were excluded from the analysis of offspring number. Overall, we analysed 140 and 117 replicates (i.e., patches) of daily fecundity and offspring number, respectively (see Table S1 for more details).
Interaction between pesticide resistance and response to high developmental temperature
To test if the effect of high developmental temperature on developmental time, survival and fertility differed depending on the resistance status of individuals, we exposed both resistant (Wu.RR) and susceptible (Wu.SS) populations to either 25 or 36°C and repeated the experiments described above. When testing fertility, females were given four days to oviposit instead of three, and four treatments were established to separately test the response of males and females to high developmental temperature: pairs where both individuals developed at 25°C (control; ♀25 x ♂25), and pairs where each or both sexes developed at 36°C (♀36 x ♂25, ♀25 x ♂36 and ♀36 x ♂36, respectively). Pairs were maintained at 25°C when both individuals developed at that temperature. Otherwise, pairs were kept at 36°C until the sex that developed at 36°C was removed, at which time they were moved to 25°C. When both sexes developed at 36°C, the transfer to 25°C occurred when the last individual of the pair was removed. All other experimental details are as above.
Replicates were excluded according to the same criteria as described above. Overall, we analysed the developmental time of females and males from 97 and 95 patches, respectively, the number of surviving juveniles from 114 patches and the number and sex of offspring from 233 patches (see Table S1 for more details).
Statistical analyses
All statistical analyses were performed using the software R (version 4.0.3; R Core Team 2020). Depending on the data and error structure, we performed linear models (LM), linear mixed-effects models (LMM), generalized linear models (GLM) and generalized mixed-effects models (GLMM) implemented in lme4 (version 1.1.26; Bates et al. 2015) and glmmTMB (version 1.0.2.1; Brooks et al. 2017; see Table S3 and S4). Maximal models were simplified by removing non-significant terms (p > 0.05) from the highest- to the lowest-order interaction (Crawley 2012). The significance of each explanatory variable was determined by chi-squared tests for discrete distributions, and Wald F tests for continuous distributions (Bolker et al. 2009). Graphic representations of the data were produced with the software package ggplot2 (version 3.3.3; Wickham 2016).
Effect of high developmental temperatures on life-history traits
The effect of developmental high temperatures on juvenile traits was assessed by analysing two traits: the number of days from egg until adulthood (i.e., number of days the first male or female took to reach adulthood in a patch), used as a proxy for developmental time of each sex; and the number of surviving juveniles (i.e., number of juveniles surviving until the adult stage; Table S3). The effect of developmental heat stress on female fertility, was determined by analysing daily fecundity and offspring number (Table S3).
As temperatures were not tested simultaneously, we analysed the differences of each high temperature to the control temperature. For that, within each assay, the mean trait value measured at the control temperature was subtracted to the individual trait values measured at the high temperature tested (trait value at the high temperature - x̄ trait value at 25°C). When the assay was done in more than one block, the mean of the control temperature was calculated within each block, to account for that variability. The new variable was used as dependent variable in the subsequent models.
All variables were analysed using LM models with a gaussian error distribution, except for the number of surviving juveniles that was analysed using a LMM model (Table S3). The data were box cox transformed (female developmental time: λ = 1.94; male developmental time: λ = 0.73; number of surviving juveniles: λ = 0.91; MASS package, boxcox function; version 7.3.53; Venables and Ripley 2002) to improve the fit of the model in each variable, except for daily fecundity and offspring number. To analyse the impact of temperature on developmental time and number of surviving juveniles, the developmental temperature (25, 33, 35 or 36°C) was included in the models as fixed factor (Table S3). Fecundity was added as random factor in the model testing the impact of temperature on the number of surviving juveniles, to account for the variation in the number of eggs laid (Table S3). In the analyses of daily fecundity and offspring number, the type of cross between males and females from different developmental temperatures (♀25 x ♂25, ♀33 x ♂25, ♀35 x ♂25, ♀36 x ♂25 or ♀37 x ♂25) was included in the model as fixed factor (Table S3).
To determine if the difference of each temperature to the control was significant, we used the function summary and took the value that reports the difference of the first temperature listed on the dataset to an intercept of zero (which in this case corresponds to the mean trait value at control temperature). This was repeated for all temperatures by reordering this variable on the dataset. A posteriori contrasts with Bonferroni corrections of significant explanatory variables were made to test differences among the different high temperatures using the emmeans package (version 1.5.3; Length 2020; Table 2).
Interaction between pesticide resistance and response to high developmental temperature
To test whether developmental heat stress differentially affected individuals with different resistance status in a sex-specific manner, we analysed four traits: developmental time, the number of surviving juveniles, the offspring number and the offspring sex ratio (Table S4).
Developmental time and the number of surviving juveniles were analysed using GLM and GLMM models, respectively, with a Poisson error distribution (Table S4). The developmental temperature of the individuals (25°C or 36°C), their resistance status (susceptible or resistant) and the interaction between these two factors were included in all models as fixed factors. Fecundity was included as a random factor in the analysis of the number of surviving juveniles.
The offspring number was analysed using a GLMM model with a quasi-Poisson error distribution (Table S4). Offspring sex ratio, analysed as the proportion of daughters, was computed using the function cbind with the number of daughters and sons as arguments (Table S4). This variable was analysed with a GLMM model with a beta-binomial error distribution and a parameter to account for zero inflation (ziformula ~1; package glmmTMB). In the analysis of both variables, the developmental temperature of the female and of the male (25°C or 36°C, in both variables), the resistance status of the pair (susceptible or resistant), as well as all possible interactions, were included as fixed factors. We paired males and females from two different cohorts (i.e., males and females could be 9 or 10 day-old when developed at 25°C, or 7 or 8 day-old when developed at 36°C), so the age of each sex was included as a random factor (Table S4). A posteriori contrasts with Bonferroni corrections were made using the emmeans package (Length 2020; Table 4).