Interspecific differences in microhabitat use expose insects to contrasting thermal mortality

Ecotones linking open and forested habitats contain multiple microhabitats with varying vegetal structures and microclimatic regimes. Ecotones host many insect species whose development is intimately linked to the microclimatic conditions where they grow (e.g., the leaves of their host plants and the surrounding air). Yet microclimatic heterogeneity at these fine scales and its effects on insects remain poorly quantified for most species. Here we studied how interspecific differences in the use of microhabitats across ecotones lead to contrasting thermal exposure and survival costs between two closely-related butterflies ( Pieris napi and P. rapae ). We first assessed whether butterflies selected different microhabitats to oviposit and quantified the thermal conditions at the microhabitat and foliar scales. We also assessed concurrent changes in the quality and availability of host plants. Finally, we quantified larval time of death under different experimental temperatures (thermal death time [TDT] curves) to predict their thermal mortality considering both the intensity and the duration of the microclimatic heat challenges in the field. We identified six processes determining larval thermal exposure at fine scales associated with butterfly oviposition behavior, canopy shading, and heat and water fluxes at the soil and foliar levels. Leaves in open microhabitats could reach temperatures 3 – 10 (cid:1) C warmer than the surrounding air while more closed microhabitats presented more buffered and homogeneous temperatures. Interspecific differences in microhabitat use matched the TDT curves and the thermal mortality in the field. Open microhabitats posed acute heat challenges that were better withstood by the thermotolerant butterfly, P. rapae , where the species mainly laid their eggs. Despite being more thermosensitive, P. napi was predicted to present higher survivals than P. rapae due to the thermal buffering provided by their selected microhabitats. However, its offspring could be more vulnerable to host-plant scarcity during summer drought periods. Overall, the different interaction of the butterflies with microclimatic

butterflies (Pieris napi and P. rapae).We first assessed whether butterflies selected different microhabitats to oviposit and quantified the thermal conditions at the microhabitat and foliar scales.We also assessed concurrent changes in the quality and availability of host plants.Finally, we quantified larval time of death under different experimental temperatures (thermal death time [TDT] curves) to predict their thermal mortality considering both the intensity and the duration of the microclimatic heat challenges in the field.
We identified six processes determining larval thermal exposure at fine scales associated with butterfly oviposition behavior, canopy shading, and heat and water fluxes at the soil and foliar levels.Leaves in open microhabitats could reach temperatures 3-10 C warmer than the surrounding air while more closed microhabitats presented more buffered and homogeneous temperatures.
Interspecific differences in microhabitat use matched the TDT curves and the thermal mortality in the field.Open microhabitats posed acute heat challenges that were better withstood by the thermotolerant butterfly, P. rapae, where the species mainly laid their eggs.Despite being more thermosensitive, P. napi was predicted to present higher survivals than P. rapae due to the thermal buffering provided by their selected microhabitats.However, its offspring could be more vulnerable to host-plant scarcity during summer drought periods.Overall, the different interaction of the butterflies with microclimatic

INTRODUCTION
Vegetation cover locally modifies climatic conditions and generates a microclimatic regime that deviates from open, free-air, and standardized measurements (i.e., macroclimate ;Geiger, 1950;Stoutjesdijk & Barkman, 2014).The absorption and reflection of solar radiation and the evapotranspirative cooling of forest canopies buffer macroclimatic temperatures and reduce thermal variation in the understory (Bramer et al., 2018;De Frenne et al., 2021;Zellweger et al., 2020).In contrast, in open areas with short and sparse vegetation, temperatures near the ground can be more extreme than those recorded at 2-m and shady conditions (i.e., thermal amplification; Carnicer et al., 2021;Stoutjesdijk & Barkman, 2014;Woods et al., 2015).Narrow ecotones that link open and forested habitats (hereafter termed open-closed ecotones) generate multiple microhabitats with varying vegetal structures, potentially exposing the organisms they harbor to contrasting microclimatic conditions.Here we studied how interspecific differences in microhabitat use across ecotones shape the thermal exposure and the associated costs to survival of two closely-related butterflies.
Microclimatic variation determines the thermal experience of the organisms and has important effects on their thermal adaptations and performance (Franken et al., 2018;Kaiser et al., 2016;Kaspari et al., 2015;Pincebourde & Casas, 2019;Woods et al., 2022).The scale at which microclimatic measurements are relevant for organisms depends on their body size and mobility (Kingsolver et al., 2011;Pincebourde et al., 2021;Pincebourde & Woods, 2020).For butterflies, microclimatic variability at the landscape scale is pertinent in the adult stage, as they can easily fly and sample between alternative habitats (e.g., woodland vs. grassland; Suggitt et al., 2011Suggitt et al., , 2012)).In contrast, the area that eggs and larvae can use is much smaller (Courtney, 1986); thus, these less mobile stages will be more likely affected by the microclimatic variability inside their habitats.Larval stages may be more likely influenced by the microclimatic conditions of the air bath surrounding them (at less than ~1 m of distance, i.e., microhabitat), while the microclimate measured at the plant surfaces might be more important for eggs and small larvae (Kingsolver et al., 2011;Pincebourde et al., 2021;Potter et al., 2009;Woods, 2010Woods, , 2013)).Microclimatic heterogeneity at these fine scales can be particularly high and comparable to that recorded at macroclimatic scales, but remains poorly quantified for most species (Pincebourde et al., 2016;Pincebourde & Woods, 2020).
The selection of microhabitat and host plant during butterfly oviposition will affect the growing environment of the offspring (Figure 1a; Courtney, 1986;Doak et al., 2006;Forsberg, 1987;Gibbs & Van Dyck, 2009).On the one hand, we expect that shadier microhabitats will offer buffered microclimates, with dampened thermal variability and extreme values, while they are amplified in open areas exposed to direct radiation (Figure 1a  A-C).On the other hand, we also expect that microhabitat conditions could induce plastic shifts on many traits that define host plant quality, depending on their shade tolerance (Figure 1a D-E; Poorter et al., 2019;Scriber & Slansky, 1981).Variations in host plant quality can interact with the microclimate experienced by the feeding larvae and influence their development in complex ways (Clissold & Simpson, 2015).
Butterfly responses to microclimatic exposure will depend on both their thermal exposure and thermal sensitivity (Carnicer et al., 2017).We hypothesize that the two studied butterflies will show diverging thermal adaptive strategies according to their microhabitat preferences and the associated thermal regime (Figure 1a F-G).The probability of surviving heat stress depends on both its intensity and its duration (Rezende et al., 2014).This relationship can be synthesized by thermal death time (TDT) curves, where the critical temperature (T ko ) that an organism can tolerate linearly decreases with the logarithm of exposure time (t) following this equation: Therefore, thermal tolerance to a heat challenge not only depends on the upper critical thermal limit (CT max ) but

Extreme and more variable temperatures
Survival time (log 10 t) also on the thermal sensitivity (z), which describes the required increase of temperature to decrease time to death by one order of magnitude (Figure 1b).The analysis of the TDT curves of 56 species of insects, bivalves, and fishes has pointed out that CT max and z are positively and tightly associated (Rezende et al., 2014).Organisms with a high CT max and z are more capable of surviving extreme temperatures, but their critical temperature (T ko ) rapidly decreases with exposure time (Figure 1b, yellow).In contrast, low z values allow longer survival times at less intense but still stressful temperatures to the detriment of CT max (Figure 1b,blue).As a consequence of this trade-off, we expect that larvae growing in open microhabitats may show a thermotolerant strategy (high CT max and z) to cope with the acute stresses of their highly-variable thermal regime (extreme temperatures for short times).In the same line, we would expect the opposite, thermosensitive strategy for species selecting shadier microhabitats with less intense temperatures but longer exposures (i.e., chronic stress).
To assess whether microhabitat variability in open-closed ecotones generates interspecific differences in their thermal exposure and the associated costs to survival, we studied two model species of butterfly: Pieris napi L. 1758 and P. rapae L. 1758.These species were studied in two Mediterranean sites, which harbor populations of two host plants with contrasting tolerances to shade, Alliaria petiolata (shade-tolerant) and Lepidium draba (light-demanding).We first assessed whether the two butterfly species selected different microhabitats from the open-closed ecotones to oviposit and quantified the thermal exposure of their offspring both at the microhabitat and foliar scales.We also assessed concurrent changes in the nutritional quality and condition of host plants across ecotones.Then, we conducted ecophysiological assays of heat tolerance with the larvae of both species to estimate their TDT curves and determine their thermal sensitivity.Based on the experimental quantification of larval survival at different temperatures, we finally applied a dynamic model to predict thermal mortality in field microclimatic conditions.

Study system
The thermal exposure of the two species, its impacts, and the concurrent variation in their host plants were studied in two protected areas of the north-eastern Iberian Peninsula, 50 km from each other (Appendix S1: Figures S1 and S2).The sites of study are part of a long-term monitoring network that provides data on butterfly abundance at a weekly resolution since 1994 (i.e., CBMS; Pollard & Yates, 1993;Stefanescu, 2000).They contain abundant populations of the green-veined white (P.napi) and the small white (P.rapae) that have been intensively studied since 2012 (Carnicer et al., 2019;Vives-Ingla et al., 2020).Site 1 is in a mid-elevation area (539 m above sea level (asl); Can Jordà, La Garrotxa Volcanic Zone Natural Park) with a heterogeneous landscape of evergreen and deciduous woodlands, meadows, and natural ponds.Site 2 is in a coastal wetland (2 m asl; El Cortalet, Aiguamolls de l'Empordà Natural Park) surrounded by riparian deciduous forests, bush and bramble thickets, reed beds, and irrigated cropland.The landscape mosaic of both sites generates spatial gradients of vegetation cover between open and closed habitats and their transition zones (open-closed ecotones).
P. napi and P. rapae are partially syntopic species (i.e., they share some of the microhabitats for breeding).They co-occur in a wide range of habitats from the sea level to alpine regions, but P. napi is usually linked to shaded, humid sites while P. rapae is more common in dry, open areas.Both butterflies lay individual eggs on Brassicaceae species, such as A. petiolata and L. draba.However, P. rapae is a more generalist species and uses a greater diversity of host plants (Ohsaki, 1979).
The dominant host plant in the mid-elevation site is Alliaria petiolata.It is a biennial herb common in damp, shaded soils at the edges of deciduous and riverine forests.It can grow in highly contrasted environmental conditions, exhibiting considerable plasticity in different habitats (Cavers et al., 1979).Seedlings emerge during spring and early summer and persist as rosettes throughout the first year, until the next growing season when inflorescences are initiated.Lepidium draba is the dominant host plant in the lowland site.It is a perennial, rhizomatous herb and can be found in open areas and field margins (de Bol os & Vigo, 1990).Its extensive, multibranched rhizomes are notably capable of producing many new shoots, which can develop into large monocultural stands (Francis & Warwick, 2008).

Oviposition behavior
We assessed whether the differences in broad habitat preferences between the two butterflies (P.napi for humid areas and P. rapae for open habitats) led to different microhabitat selection for oviposition across open-closed ecotones.We tested this hypothesis by carrying out censuses of behavior at the two study sites.Females were followed and their behaviors were recorded for periods of 45 min.The censuses fully covered the entire daily period of flight activity, between 9:00 and 19:00, and were conducted during summer.Oviposition was considered to occur when females that landed on a leaf were observed to curl their abdomen and remain in this position for at least 3 s.Egg-laying was visually confirmed in most of the cases.Species, hour, and microhabitat type (open, closed, or intermediate) were noted.Censuses were simultaneously performed in the various microhabitat types, balancing the time spent in each type.The temperature of the leaves where eggs were laid and their position (upper vs. underside of leaves) were also recorded when possible using a wire K-type thermocouple probe (SC-TT-KI-30-1M, Omega Engineering Ltd, UK) attached to a hand-held thermocouple thermometer (HH503, Omega Engineering Ltd, UK, and HI935005N, Hanna Instruments Ltd, Spain) immediately after the female left the plant.

Microclimatic and host-plant variation
We assessed microclimatic conditions at both microhabitat and host-plant scales.For the measurements at the microhabitat scale, we installed 15 standalone data loggers (EL-USB-2-LCD, Lascar Electronics, UK) in different microsites harboring host plants and where ovipositing females had been detected.Seven of the sensors were installed in the mid-elevation site, and the other eight in the lowland site (Appendix S1: Tables S1 and S2).The sensors were programmed to record temperature ( C) and relative humidity (%) at hourly resolution and were placed 25 cm above the soil protected from direct solar radiation.
To quantify thermal exposure and trait variability of host plants, we monitored cohorts of 242 individuals of A. petiolata and 362 individuals of L. draba distributed in four representative categories of microhabitats found across the ecotones (open, semi-open, semi-closed, and closed).Each microhabitat category was assigned based on detailed measurements of the dynamics of the canopy and the ground cover by herbaceous plants (see Appendix S1: Section S1 and Figures S2 and S3 for further details).Canopy closure was measured by visual inspection in the vertical and the four cardinal directions.The herbaceous layer was characterized both as herbaceous ground cover and as mean herb height using the point-intercept method.Closed and semi-closed microhabitats presented a mean canopy closure higher than 50%, while we defined semi-open and open microhabitats by their herbaceous layer, as trees and shrubs were less common there.We monitored the host plants in 2017, from March to September, and measured the same 14 microclimatic, phenological, morphological, and physiological traits every 15 days (Appendix S1: Table S2).
Every monitoring day, we selected at least 16 host plants for each cohort (4 individuals Â 4 microhabitat types) ensuring that plants were randomly chosen, without repetition to avoid pseudoreplication.Microclimatic variables included several temperatures at the soil, foliar and air level, and soil moisture.We used a penetration thermometer (HI98509, Hanna Instruments Ltd, Spain) to measure soil temperature at a depth of 10 cm.Soil surface temperature, air temperature above the host plant, and foliar surface temperature were measured using a thermocouple (see device and measurement details in the section Oviposition behavior).A minimum of three replicates were taken for each measurement.The temperatures were measured between 10:00 and 16:00, and the time, wind, and radiation conditions were recorded.We measured soil surface temperature near the host plants, replicating it in spots exposed to direct solar radiation and in the shade; air temperature, immediately above the host plant at a height of 1 m; and foliar temperature, on the upper and underside of the leaves.The volumetric water content of the soil (% by volume) was measured at three points near each plant using a DELTA-T SM150 (Delta-T Devices Ltd, UK) soil-moisture sensor.
The traits measured in host plants included phenology, stem length, and foliar dimensions, density, and chlorophyll and water contents.These traits were selected because they have been associated with oviposition behavior and offspring performance (Awmack & Leather, 2002;Gibbs & Van Dyck, 2009;Stefanescu et al., 2006;Wolfson, 1980) and are simple to measure.The foliar chlorophyll content is also an indicator of plant nutritional condition, photosynthetic capacity, and developmental stage (Curran et al., 1990;Everiit et al., 1985).We assessed plant phenological status by classifying the individuals into one of four phenological stages: spring early vegetative, reproductive, senescent, and summer late vegetative.A representative basal, medial, and apical leaf was chosen for each plant, and its state (green or senescent) was recorded.Chlorophyll content was estimated as the mean of three measurements from a SPAD-502 chlorophyll meter (Konica Minolta Sensing, Spain).Finally, leaves were severed and immediately weighed (fresh weight, FW) using a field digital scale (PJS020, PESOLA Präzisionswaagen AG, Switzerland), and then oven-dried in the laboratory at 60 C for 2 days to a stable weight (dry weight, DW).Foliar water content was defined as (FW À DW)/DW.Foliar density was calculated as the ratio between DW and foliar length.When host plants were mature, we also counted the number of fruits per plant (siliques for A. petiolata and silicules for L. draba), as a proxy of plant reproductive performance between microhabitats and shade tolerance.A minimum of seven individuals were sampled for each microhabitat type.

Ecophysiological assays of heat tolerance
We implemented a static heat tolerance experiment using larvae of P. napi and P. rapae to determine whether they differ in their thermal strategies.If P. napi oviposits in more closed and buffered microhabitats than does P. rapae, as we hypothesized, the larvae of the former would exhibit a thermosensitive strategy, with lower values of z and CT max .Females from both locations and species were captured and their offspring were reared in growing chambers at 22 C 13L:11D, with fresh and abundant host plants (L.draba and A. petiolata).The experiment was conducted on 210 larvae from 20 family lines (Appendix S1: Table S3).Before the application of the thermal treatment, larvae were acclimated for 1 h at constant 22 C and deprived of food.We recorded the larval initial weight (g) and subsequently placed the larvae in individual plastic vials (diameter: 3.5 cm, height: 7 cm) that were submerged in a water bath programmed at a constant temperature (i.e., 40, 42, or 44 C, depending on the treatment).These temperatures can be recorded in the field and are known to be stressful for both species (Kingsolver, 2000;von Schmalensee et al., 2021).The status of the larvae (alive or dead) was checked at regular time intervals (once every 30 min for the assays at 40 C; every 20 min at 42 C; and every 10 min at 44 C), trading off between accurate detection of time to death and potential thermal fluctuations associated with larval status checking.The air temperature inside the plastic vial was continuously recorded using a data logger with a 20-s resolution to have a more accurate estimate of the thermal exposure of the larvae and its temporal fluctuations during the treatment.The average temperature recorded with the data logger was used in the subsequent analyses, rather than the fixed, programmed temperature in the water bath (which was considered a less accurate proxy).

Statistical analysis
All data were analyzed using R 3.6.1 (R Core Team, 2019).We applied generalized linear mixed models (GLMMs) to determine whether the two butterfly species selected different microhabitats to oviposit.The final data set included a total of 7217 min of census, 139 ovipositions, and 43 ovipositing females (Appendix S1: Table S4).The number of ovipositions observed per female was used as the response variable (model 1).We fitted the model using the glmer function of the lme4 package (Bates et al., 2015) by maximum likelihood (Laplace approximation) and a Poisson error distribution with a log link function.The type of microhabitat, the species, and their interaction were added as fixed factors, and site, date, and period of the day of the census were treated as categorical, random factors.The logarithm of census duration was added as an offset term (i.e., a linear predictor without an estimated regression parameter to account for differences in sampling effort between censuses; Zuur et al., 2009).We repeated the same modeling procedure using the number of ovipositing females per census as the response variable (model 2).
The seasonal dynamics of microclimatic conditions and host-plant traits were assessed by regressing LOESS models against ordinal day.To examine the spatial variation of these variables across the open-closed ecotones, an analysis of variance (ANOVA) testing for microhabitat type was applied followed by a post hoc Tukey honest significant difference (HSD) test calculated using the emmeans package (Lenth, 2020).The analyses were performed for the entire sampling period and specific phenological stages and seasons.The characterization of the thermal regimes in different microhabitats and scales considered the absolute records at different levels (air, soil, leaves), thermal tendencies relative to macroclimatic conditions, and measures of dispersion (i.e., standard deviation and skewness).Daily records of macroclimatic temperatures were obtained from two meteorological stations near the study sites (Appendix S1: Figure S1).Microhabitat thermal offset was calculated as the difference in daily mean temperatures between the data logger and the standardized weather station measurements.We defined foliar thermal offset as the difference between the upper side foliar temperature and the synchronic air temperature above the host plant at 1 m height (also termed "thermal excess" in some studies; see De Frenne et al., 2021;Pincebourde & Woods, 2012).The daily standard deviation of temperatures recorded with the data logger was used as an indicator of the microclimatic temporal variability at the microhabitat scale, while we used the daily SD of foliar temperatures in the same microhabitat as a measure of spatial thermal heterogeneity at finer scales.
To assess whether the two butterflies differ in their thermal strategy, we fitted a linear model of the logarithm of the time to larval death against the mean temperature recorded with the data logger during the static treatments for each species (i.e., TDT curve).We estimated z and CT max from the regressed equations and we then tested whether the two species presented different slopes by fitting an analysis of covariance (ANCOVA) of the effect of temperature, species, and their interaction on log 10 of larval survival time.Additional general mixed linear models were also fitted for each species on log 10 of time to larval death with thermal treatment ( C), larval weight (g), and site as fixed factors and family within site as a random factor.

Predictions of thermal mortality in the field
The semilogarithmic link between knockdown times and temperature during a heat challenge (i.e., TDT curve) can be estimated experimentally for different levels of mortality (e.g., the time where 100% or 50% of larvae were dead), resulting in parallel lines with the same slope (z).From these parallel curves, we can define how mortality rate (and, hence, survival probability) changes with exposure time for any constant temperature (Rezende et al., 2020b).We obtain in this way temperature-specific survival curves, which are mathematically related and collectively define the thermal tolerance landscape of the species (Rezende et al., 2014(Rezende et al., , 2020a(Rezende et al., , 2020b).In the field, organisms are exposed to variable temperatures and their survival responses are the result of shifts between temperature-specific curves.Therefore, if we know how temperature changes during a specific period, we can predict the survival probability under field conditions by summing up the infinitesimal changes in the survival rate that occurred during that period (Rezende et al., 2020b).
Here we adapted the methods and scripts developed by Rezende et al. (2020aRezende et al. ( , 2020b) ) to numerically predict the daily thermal mortality throughout the year for the larvae of both Pieris species in the different microhabitats.Microclimatic field conditions during the whole year were obtained from the thermal records of the data loggers.For each sensor, we calculated the mean thermal series of all the recorded years (Appendix S1: Table S1) and estimated the thermal profile at 1-min resolution of each day by non-linear interpolation between consecutive hours (Rezende et al., 2020b).Daily thermal mortality was predicted from March to September (both included), to capture the period of higher thermal stress (i.e., summer).We predicted larval survival in the dynamical field conditions based on the TDT curves, by bootstrapping the time of death of 35 larvae per species in each of the three experimental treatments (i.e., 40, 42, or 44 C, total sample size = 210).Bootstrapped data were then used to estimate the tolerance landscape of each species with a resolution of 0.001 of survival probability (see Rezende et al., 2020b for further details).From this relationship between survival and time at constant temperatures, we calculated the daily mortality curve at a 1-min resolution for each species, day, and microclimatic sensor by numerical approximation (Rezende et al., 2020b).
We defined daily thermal mortality as the maximum mortality (1-survival probability) of the day and calculated, with these daily values, the cumulative mortality during development (using a rolling window of 30 days from March to September).We then averaged cumulative mortality for the periods where it exceeded 0.01 (i.e., mostly in summer) to summarize the peaks of thermal mortality for each species in each microhabitat.Finally, we assessed how thermal mortality would change if we consider that larvae can avoid acute thermal stress on the leaves by actively moving to shadier parts of the plant (thermal avoidance behavior, TAB; Carnicer et al., 2019).This was done by truncating daily thermal profiles at fixed TAB thresholds (none, 35, 37.5, 40, 42.5, and 45 C) and repeating the whole numeric procedure with the truncated thermal profile.These thresholds were obtained from previous experimental observations of thermal avoidance behaviors for the genus Pieris (Carnicer et al., 2019).We performed 100 bootstrap replicates per threshold to have an estimate of the uncertainty associated with the predictions.

Interspecific differences in microhabitat selection
Results of the GLMM analyses applied on the number of ovipositions per female (model 1) and the number of ovipositing females (model 2) were very similar (Table 1).Both P. napi and P. rapae distributed their eggs unequally across the open-closed ecotones (model 1: microhabitat χ 2 = 13.1, df = 2, p = 0.0015).However, the microhabitats selected by females differed between the two butterflies (model 2: microhabitat Â species χ 2 = 21.2, df = 2, p < 0.0001).P. napi preferentially selected host plants from the microhabitats with intermediate vegetation covering (which include the semi-open and semi-closed microhabitats, Figure 2a).Concretely, predictions from model 1 indicated that the number of eggs laid in these microhabitats was six and three times higher than in closed and open microhabitats, respectively.In sharp contrast, P. rapae mainly laid eggs in the open microhabitats (in model 1, open microhabitats received 99% of P. rapae ovipositions).Closed microhabitats were rarely selected by either species, but a few ovipositions of P. napi were observed.The oviposition pattern did not present relevant differences between sites, as the variance estimated for this random factor was much lower than the other effects in both models (model 1: sd site = 0.02, sd date = 1.47, sd period = 1.77; model 2: sd site = 0.23, sd date = 0.7, sd period = 1.26).

Interspecific differences in thermal exposure
Foliar temperatures recorded during oviposition differed between both species (with thermal differences in the range of 5-10 C in the leaf underside, where most eggs were laid, Figure 2b).The thermal regimes characterized throughout the monitoring campaign varied between microhabitats in the same line.Daily maximum temperatures recorded with the data loggers were, respectively, a mean of 6 and 2 C lower in the semi-closed and semi-open microhabitats than in the open areas, where P. rapae oviposits (Figure 2c).A similar pattern was found for temperatures recorded at the upper surface of the leaves (Figure 2d), which were measured during the period of highest insolation (10-16 h) and when most eggs are usually laid (Appendix S1: Figure S4).
Canopy cover of the closed and semi-closed microhabitats from the mid-elevation site diminished the mean daily temperatures of the understory by an average of 3 C relative to mean macroclimatic conditions (Figure 2e).In contrast, the buffering capacity of lowland vegetation was much lower, and the largest offsets from macroclimatic mean temperatures were observed in the open and semi-open microhabitats (daily mean temperatures 1-2 C higher than meteorological records).Besides the higher air temperatures at the microhabitat scale, leaves in the open microhabitats of the lowland site were additionally subject to important thermal amplification processes that could elevate foliar temperatures up to 10 C higher than the air above the host plant (Figure 2f).Temperatures in closed and semi-closed microhabitats were more constant throughout the day and were more homogenous across this type of microsite (Figure 2g,h).Open microhabitats presented instead the most variable daily thermal profiles (Figure 2g), with more positively skewed distributions (i.e., more extreme values in the upper side of the thermal distribution, Appendix S1: Figure S8a).Thermal heterogeneity at fine spatial scales was also higher in the open microhabitats, which was calculated as the SD of the temperatures recorded on the upper-and underside of the leaves in each microhabitat on a daily basis (Figure 2h).

Sources of thermal heterogeneity at fine scales
All the thermal variables measured at the host-plant level (i.e., excluding air temperatures at the microhabitat scale recorded with the data loggers) were strongly correlated in pairwise correlations for each site and microhabitat type (R 2 = 0.81 ± 0.16; and p = 0.0042 ± 0.0198, which summarize the mean ± sd of the R 2 and p values obtained in the correlation tests).Nevertheless, there were notable differences between leaves, soils, and the air in the open microhabitat of the lowland site, where thermal heterogeneity was higher (Figure 3 and Appendix S1: Figures S5 and S9).
The soil in the open microhabitat of the lowland site was warmer and drier than in the other microhabitats (Figures 3 and 4a), especially during summer (Figure 4c).Soil surface reached temperatures higher than 45 C when exposed to full radiation (i.e., soil thermal amplification, Figure 3c).Going up from the soil, air temperature rapidly decreased with height following a hyperbolic sine function (i2 in Figure 3).In this line, basal leaves of the host plants from open microhabitats were 2-7 C warmer than apical leaves (i1 in Figure 3), while basal leaves in semi-open and semi-closed microhabitats could reach inferior temperatures than apical leaves (i.e., soil cooling effect, Figure 3b).This thermal difference between leaves could be recorded in spring plants of Lepidium draba, which were >40 cm in height.In contrast, leaves of summer resprouts (<5 cm height) presented the highest T A B L E 1 Generalized linear mixed models of the number of ovipositions per female and the number of ovipositing females per census.

Model
Response Note: Results from the type III Wald χ 2 analysis of deviance of the fixed effects.Marginal R 2 for GLMMs (including only fixed effects) were calculated following Nakagawa et al. (2017) with the version 0.9.2 of the performance package (Lüdecke et al., 2021).
temperatures and the lowest thermal heterogeneity (daily SD of foliar thermal records, Appendix S1: Figure S5c,d).
Fine-scale thermal heterogeneity was also detected between the sides of the same leaves.In intermediate ranges of microhabitat air temperature (i.e., 20-35 C), foliar underside temperatures were 1-3 C cooler than upper parts of the leaves exposed to direct radiation in open microhabitats (Figure 3d).However, at higher air  temperatures (i.e., >35 C), thermal differences between the foliar upper and underside vanished (Figure 3e).In these conditions, leaves could be 10 C warmer than the air temperature above the host plant (Figure 3f).The effects of foliar height and foliar side on thermal heterogeneity were putatively associated with radiative heating and sensible heat fluxes from the soil and with processes of stomatal closure of the leaves.

Host-plant variability
Open-closed ecotones also induced changes in host-plant traits, paralleling variation in thermal exposure (Figure 4 and Appendix S1: Figures S6 and S7).
Host plants in open microhabitats had smaller leaves with lower ratios of water content (i.e., less water per mg of foliar dry weight, Figure 4d).Variations in foliar chlorophyll content and fruit production depended on the shade tolerance of the host plant.Lepidium draba (lowland site) presented higher chlorophyll contents and produced more fruits in open microhabitats, while plants in closed and semi-closed microhabitats did not reproduce sexually and had thinner leaves, with low chlorophyll contents.On the contrary, chlorophyll content and fruit production for A. petiolata in the mid-elevation site were lowest in the open microhabitat (Appendix S1: Figure S10).Foliar water and chlorophyll contents decreased in both host plants (Figure 4e,f, and Appendix S1: Figure S11) as they senesced after fructification during late spring and early summer (ordinal days 140-180).Only non-flowering first-year rosettes (A.petiolata) and summer rhizome resprouts (L.draba) remained in midsummer after senescence (Appendix S1: Figure S12).First-year rosettes of A. petiolata notably coexisted in June with the reproductive stage of second-year individuals.In contrast, there was no temporal overlap between reproductive L. draba plants and new summer resprouts, leading to a period of scarcity of fresh host plants lasting 2-3 weeks.

Interspecific differences in thermal strategies and mortality
Time to larval death and experimental temperatures were associated in semilogarithmic curves for both species (Figure 5a), presenting a steeper slope (25% more negative) for P. napi (species Â temperature p = 0.0085 in the ANCOVA model, which explained a 76% of the variance).The estimated intersection point of the two TDT curves was located at T* = 41.2 C. Above this threshold, P. rapae exhibited greater survival than P. napi, while the opposite was true for inferior temperatures (Figure 5a).As predicted, the estimates (±SE) of thermal tolerance z and CT max were higher for P. rapae (z = 5.10 ± 0.24 C; CT max = 53.48C) than for P. napi (z = 4.10 ± 0.26 C; CT max = 51.08C).In the general mixed linear models considering more factors than temperature, strong effects were also found for the larval weight (Table 2), with smaller larvae presenting lower survival times.Concretely, an increase of 0.1 g in larval weight would cancel out the effects of an increase of 1 C in temperature (i.e., the former would provoke an increase of 30%-50% in survival times, while the latter would decrease them by 30%-50%, Appendix S1: Figure S13).
Interspecific differences in the TDT curves and the predicted thermal mortality in the field were in agreement with the pattern of microhabitat selection of the two species.Daily mortality was higher for the thermosensitive P. napi (Appendix S1: Figure S14) when thermal stresses were more acute (more intense heat challenges, mainly found in open microhabitats).In contrast, more constant thermal regimes, with less extreme but longer thermal heat challenges (i.e., chronic stress), were deadlier for P. rapae.Daily thermal mortality from March to September was usually low (<0.01),although it could reach values around 0.4 in the open microhabitats during the warmest days (<15% of the days, Appendix S1: Figure S15a).
The accumulation of low daily mortalities for a period similar to larval development time (i.e., 30 days)  F I G U R E 3 A visual summary of the diverse microclimatic mosaics and processes at fine scales quantified in the lowland site.Blue areas (left) indicate microhabitats preferentially selected by P. napi, and orange areas (right), by P. rapae.Observed thermal differences between the air, host plants, and soil allowed the identification of six processes determining larval thermal exposure: Butterfly oviposition behavior (a); the influence of soil cooling effects (b) and soil thermal amplification on basal leaves (c); canopy thermal homogenization (d); leaf underside cooling by active stomatal conductance (e); and foliar thermal amplification by reduced stomatal conductance (f).Detailed thermal data are reported in the four insets providing evidence for these processes (i1-i4).(a) Air temperatures in the semi-open and semi-closed microhabitats where P. napi oviposits reach about 3-10 C inferior values and are less variable than in open microhabitats (ovipositing microsites of P. rapae).Open microhabitats are more thermally heterogeneous at fine scales, with notable thermal differences between leaves, soil, and the air.(b) Temperatures at the soil surface influence the heat balance of basal leaves, such that basal leaves in the semi-closed microhabitat can be cooler than apical leaves (i1).(c) Soil thermal amplification in the open microhabitat warms the lower air layers and creates steep thermal gradients with height (i2), which can be also detected between basal and apical leaves (i1).(d) Canopy shading in more closed microhabitats cools leaves and air, resulting in similar temperatures in the leaf upper-and undersides (i3 and i4, results in blue).(e) When temperatures in open microhabitats were moderately warm (around 25 C), we observed that stomatal conductance and evapotranspiration cooled foliar undersides 1-3 C relative to the upper sides, which are more directly exposed to radiative heating (i3, orange line).(f) However, in dry and warm conditions (>35 C), high radiative heating from the soil and reduced evapotranspiration linked to leaf stomatal closure can bring leaves 10 C warmer than the air (i4, orange).Black lines in the soil are added for illustrative purposes and represent the patterns of variation of soil temperature and humidity across the ecotone that we observed in Figure 4a and Appendix S1: Figure S11b.More detailed results of thermal heterogeneity at fine scales can be found in Appendix S1: Figures S5 and S9. could exert considerable thermal pressure on natural populations (cumulative mortality >0.8 during the summer in the open microhabitats).However, the estimates indicated that thermal mortality would be importantly reduced if larvae conducted thermal avoidance behaviors (Appendix S1: Figure S15b).Mean thermal mortality during the development was around 0-0.2 in semi-open and semi-closed microhabitats (Figure 5b), and an average of 15% lower for P. napi than for P. rapae (Figure 5c).The opposite pattern was found in the open microhabitats, where mean thermal mortality during development was 15% higher for P. napi than for P. rapae and ranged between 0.4 and 0.8.Thus, by laying its eggs on semi-open and semi-closed microhabitats, P. napi eludes a high thermal mortality.In contrast, P. rapae selects the microhabitats with the most intense heat challenges, exposing its offspring to deadly thermal stresses that it can withstand better than P. napi.

DISCUSSION
Our work quantified the microclimatic and host-plant conditions of the microhabitats where two butterfly species lay their eggs and assessed their potential impacts on the offspring.We conducted the study with two model species from an extensively-studied family of butterflies in the ecotones between forested and open habitats of two protected areas.We combined field censuses of butterfly behavior, detailed measurements of the microclimate in the selected microhabitats, and an experimental characterization of the thermal tolerance of the larvae, to computationally predict the thermal mortality in the field.Many studies stress the need to consider temperatures at the microclimatic scale to better understand the thermal ecology of the species and their responses to global warming (Bramer et al., 2018;Pincebourde & Woods, 2020;Woods et al., 2015).As far as we know, mortality associated with microclimatic variation has been very occasionally assessed (Kaiser et al., 2016;Kingsolver, 1979;Potter et al., 2009;Woods et al., 2022).
For the first time, here we predicted the mortality derived from microclimatic temperatures considering both the intensities and the duration of the thermal exposures and their cumulative effects (Rezende et al., 2020b).The two butterflies selected the same host-plant species from different microhabitats to oviposit: P. napi laid most of their eggs in semi-open and semi-closed microhabitats, and P. rapae, in open ones (Table 1 and Figure 2).These results indicate that habitat choice preceded host-plant selection, as has long been proposed (Courtney, 1986;Dennis, 2010;Porter, 1992) and as has been found for this pair of species (Friberg & Wiklund, 2019;Ohsaki, 1982;Ohsaki & Sato, 1999).Ovipositing females generally use multiple cues following a spatially-structured and hierarchical process, from coarse to finer scales.Determining which specific cues influenced oviposition decisions is beyond the scope of our study, but both microclimatic and host-plant factors could have had a role.Ohsaki (1982) associated P. napi and P. rapae oviposition decisions with the different light conditions of the microhabitat, and Forsberg (1987) suggested that P. napi in Sweden actively oviposited in small plants to favor higher microclimatic temperatures.Other studies related oviposition decisions of P. rapae with host-plant qualities.Visual stimuli (mainly the color or the greenness of the plant) have a key role when females are searching for a host plant (Myers, 1985;Tsuji & Coe, 2014), although olfactory cues likely influence pre-alightment decisions too (Renwick & Radke, 1988).After landing on host plants, the chemical and nutrient status of the plants were found to be decisive (Hern et al., 1996).For example, leaves with higher water and nitrogen content, and with high transpiration rates, are more frequently accepted to oviposit (Myers, 1985;Wolfson, 1980).In our study, all of these factors varied between the selected microhabitats (Figures 2 and 4).
Microhabitat preferences observed in the field matched interspecific differences in the TDT curves determined in the laboratory and the predicted thermal mortality in the field (Figure 5).P. napi presented a  2017) with the version 0.9.2 of the performance package (Lüdecke et al., 2021).
more thermosensitive strategy, with lower survivals under acute stresses, but higher under longer subextreme challenges.Accordingly, its mean thermal mortality during development was approximately 15% higher than P. rapae in open microhabitats, where heat challenges are more extreme, but 15% lower than P. rapae in their selected microhabitats.Petersen (1954) already proposed that, as a result of interspecific competition, P. napi would have specialized in developing in shaded environments, where thermal regimes are cooler and host plants have thinner leaves than in the dry, open habitats selected by P. rapae.But later research suggested that competition between Pierids is unlikely to be driving their habitat segregation (Courtney, 1986), and that escape from parasitism is a more likely driver in the case of P. napi and P. rapae (Ohsaki & Sato, 1999).Regardless of the ultimate cause driving this habitat differentiation, specialization in fresh and humid habitats or dry and hot habitats usually comes with different costs and benefits.For example, studies with other butterflies have shown that specializing in hot environments can select higher fecundities and adult survival at the expense of larval survival (Karlsson & Wiklund, 2005).In this line, our results also predicted that larval survival in open microhabitats is lower (Figure 5), and previous studies found that P. rapae laid more eggs but smaller ones than P. napi (Ohsaki, 1982).
In this study, we also identified six processes (Figure 3a-f) modulating the thermal exposure of the two species at very fine scales, and hence their thermal mortality.The first one is microhabitat selection by adult butterflies, which strongly determines the thermal exposure of the larvae.Thermal differences between open and intermediately-covered microhabitats were in the range of 3-10 C for diurnal temperatures at the foliar and air levels (Figures 2 and 3a).These results are in agreement with the thermal differences between more open and closed habitats of herbivorous insects reported in other studies (Ashton et al., 2009;Friberg et al., 2008;Merckx et al., 2015;Ohsaki, 1982;Suggitt et al., 2012).Daily temperatures in open microhabitats also presented higher variation, with more extreme temperatures.Then, heat and water fluxes at fine scales create very different microclimatic mosaics between the selected microhabitats (Figure 3b-f).In open microhabitats, processes like soil and foliar thermal amplification can raise foliar temperatures 10 C relative to the air (Figure 3c,f;Carnicer et al., 2021;Pincebourde et al., 2021;Woods et al., 2022), or evapotranspiration can cool the foliar underside 3 C relative to the upper side (Figure 3e).These processes create a more heterogeneous thermal mosaic than that found in semi-open and semi-closed microhabitats, where soil cooling and canopy shading homogenize temperatures (Figure 3b,d).Our predictions of thermal mortality were based on the thermal series extracted from the data loggers, which represent air temperatures at the microhabitat level.But larvae might be more dependent on temperatures at the foliar level, especially at their initial stages (Pincebourde et al., 2021;Woods, 2013).The high fine-scale thermal heterogeneity in open microhabitats could potentially expose larvae of P. rapae to more acute thermal stresses than those we predicted, but could also offer more opportunities for behavioral thermoregulation (by moving to the foliar underside [Figure 3e], or a cooler leaf or plant [Figure 3c]).Although thermal heterogeneity was lower in summer in the lowland, some leaves presented temperatures that could reduce thermal mortality.For example, if larvae avoided thermal exposures >40 C, cumulative mortality during development would be importantly reduced (Appendix S1: Figures S5  and S15).Thermal avoidance behaviors have been reported for the larvae of P. napi (Carnicer et al., 2019) and P. rapae (Kingsolver & Gomulkiewicz, 2003).
Microclimatic exposure during the development of insects will simultaneously affect other physiological and demographical rates besides larval survival (Braschler et al., 2021;Diamond et al., 2013;Kingsolver et al., 2011;von Schmalensee et al., 2021), which can influence the patterns of larval mortality described here.For example, different microclimatic regimes can lead to shorter or longer larval development times due to nonlinear changes in development rates with temperature (thermal performance curves; Greiser et al., 2022;von Schmalensee et al., 2021; Appendix S1: Figure S16).These variations in the development time would in turn modulate the time of exposure to microclimatic lethal stresses.Moreover, our results indicated that smaller larvae were more sensitive to heat challenges (Table 2 and Appendix S1: Figure S13).TDT curves in this study were built with all the larvae, so predictions of thermal mortality assumed a constant weight during the whole development.However, as larvae grow, they would become more tolerant, suggesting that plastic and evolutionary changes in growth rates and thermal performance curves could play a key role in larval thermal mortality.Further studies should consider how growth dynamically modifies thermal tolerance landscapes during development, and the parallel effects of microclimatic regimes on growth, development, and mortality.
The differences between the microclimatic regimes of the microhabitats where the two species lay their eggs led to contrasting patterns of thermal mortality of the offspring.Despite being more thermosensitive, the offspring of P. napi were predicted to present higher survivals than P. rapae (Figure 5) due to the thermal buffering provided by the microhabitats with intermediate cover (Figure 2).Our results feed the growing literature on the key role that vegetation cover may play in buffering thermal macroclimatic stresses (Carnicer et al., 2019;De Frenne et al., 2021;Zellweger et al., 2020).However, many studies on butterflies indicate that other species-specific requirements, such as host plant condition and availability, should be met in these microhabitats to successfully buffer macroclimatic impacts (Ashton et al., 2009;Bennett et al., 2015;Carnicer et al., 2019;Kaiser et al., 2016;Nieto-S anchez et al., 2015;Stefanescu et al., 2011;Suggitt et al., 2011Suggitt et al., , 2012)).For example, in summer in the lowland site, L. draba plants reduced their foliar water and chlorophyll contents until their complete senescence (Figure 4 and Appendix S1: Figure S11), which would likely interact with the microclimatic impacts on larval development (Clissold & Simpson, 2015).We suggest that the different adaptive strategies of the species and their different oviposition behaviors pose different challenges to the populations.With the current microclimatic conditions, P. napi may not be as threatened by thermal exposure as P. rapae.However, the decay and disappearance of L. draba during the summer in the lowlands might have more negative impacts on P. napi because they rely on fewer host-plant species than P. rapae, which are more likely to find other host plants in a better condition during this period (Carnicer et al., 2019;Vives-Ingla et al., 2020).
Our study assessed both the thermal exposure and the thermal sensitivity of two species to predict the associated mortality using a dynamic model calibrated with physiological information from the experiments and simulated with the microclimatic regimes recorded in the field.Previous assessments of the vulnerability of organisms to global warming usually compared their experimental upper thermal limits with maximum temperatures recorded in the field (Duffy et al., 2015;Pincebourde & Casas, 2019;Sunday et al., 2014;Woods et al., 2022).However, this approach overlooks the time-dependent effects of thermal stress and the cumulative nature of heat injury (Jørgensen et al., 2021;Rezende et al., 2020b), which may underestimate the vulnerability of organisms to global warming (Huey & Kearney, 2020).We predicted larval thermal mortalities applying the tolerance landscape framework, which allowed us to adopt a more realistic approach by considering their thermal sensitivity to both the duration and the intensity of microclimatic thermal exposures (Rezende et al., 2014(Rezende et al., , 2020b)).This new framework can offer new insights into the role of climate change in the declines reported for many insect species around the globe (Didham et al., 2020;Wagner et al., 2021, and all the references therein), also in our study area (Colom et al., 2022;Herrando et al., 2019;Melero et al., 2016;Stefanescu et al., 2011;Ubach et al., 2021).Predictive models of ecological responses to climate change should incorporate information on climatic exposure at relevant scales and capture the key processes shaping organisms' sensitivity and performance in the dynamic thermal conditions they experience in nature.
microclimatic regimes y-intercept = t max = CT max /z F I G U R E 1 (a) Narrow ecotones that link open and forested habitats generate multiple microhabitats with distinct microclimatic regimes.Insects with contrasting thermal tolerances may select different microhabitats to oviposit, leading to large interspecific differences in the thermal exposure of eggs and larvae.Moreover, host plants with varying shade tolerances may show different distributions, traits, and conditions across ecotones, synergistically affecting larval performance.Overall, concurrent gradients of both biotic and abiotic processes are produced across ecotones.(b) Thermotolerant and thermosensitive strategies are expected to evolve under different thermal regimes.The two strategies are characterized by different thermal death time (TDT) curves, reflecting an evolutionary trade-off between survival capacity at acute, extreme stresses (thermotolerant) and chronic, less intense conditions (thermosensitive).Both species exhibit equal survival times at T*, which represents the thermal threshold between these two alternative thermal strategies.Silhouettes used in this figure were obtained from rawpixel.comand phylopic.org.
Overall, the structure of the open-closed ecotones largely determined the microclimatic processes that operated in the different microhabitats, modifying the thermal exposure of the host plants.Microhabitat and foliar temperatures could exceed 40 C in the open microhabitats (Figure 2b-d), especially in the lowland site during summer (Appendix S1: Figure S5), while temperatures remained below 40 C in semi-closed and semi-open microhabitats of both sites across the seasons (Appendix S1: Figures S6 and S7).
Interspecific differences in the thermal regime that eggs and larvae are exposed to, both at microhabitat (c, e and g) and foliar (b, d, f and h) scales.(a) Relative distribution of ovipositions for each species across the open-closed ecotones.(b) Foliar temperature during oviposition at the underside.(c) Daily maximum temperatures recorded with the data loggers.(d) Foliar temperatures at the upper side during host-plant monitoring.(e) Thermal offset calculated as the difference of the mean daily temperatures between the microhabitat and the macroclimate.(f) Foliar thermal offset calculated as the instantaneous difference in foliar and air temperatures above the host plant at 1 m height.(g) Daily temporal variability (standard deviation, SD) of the temperatures recorded at the microhabitat scale.(h) Thermal heterogeneity of foliar temperatures in the same microhabitat and time.Different letters indicate the microhabitats with p < 0.05 in pairwise Tukey HSD tests of the response variable for each site (light gray for the mid-elevation site, and dark gray for the lowland site).The lower and upper hinges of the box represent the 1st and the 3rd quartiles respectively (Q 1 , Q 3 ); its inner line, the median; and the length of the box, the IQR = Q 3 À Q 1 .The lower whisker represents the smallest value ≥Q 1 À 1.5 Â IQR; and the upper whisker, the biggest value ≤Q 3 + 1.5 Â IQR.Outliers are not shown.Colored areas in panels C-H indicate the microhabitats selected by each species (i.e., blue: P. napi; orange: P. rapae) to facilitate interspecific comparisons.C, closed; SC, semi-closed; SO, semi-open; and O, open microhabitat.
Downloaded from https://esajournals.onlinelibrary.wiley.com/doi/10.1002/ecm.1561by JOSEP Penuelas -Readcube (Labtiva Inc.) , Wiley Online Library on [03/05/2023].See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions)on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License Spatial and seasonal variation in the soil humidity (a-c) and foliar water content (d-f) of the host plants measured during the monitoring campaign.(a and d) Variation across open-closed ecotones.Different letters indicate the microhabitats with p < 0.05 in pairwise HSD tests of the response variable for each site (light gray for the mid-elevation site, and dark gray for the lowland site).(b, c, e and f) Seasonal cycle of the soil humidity (b and c) and the foliar water content (e and f) measured in different microhabitats types monitored in the mid-elevation (b and e) and lowland (c and f) sites.Open microhabitats (mainly selected by P. rapae) are represented in yellow; and semi-open and semi-closed microhabitats (mainly selected by P. napi), in blue.C, closed; SC, semi-closed; SO, semi-open; and O, open microhabitats.
U R E 5 (a) Thermal death time (TDT) curves for P. napi and P. rapae.The dashed lines indicate the intersection point between the two TDT curves, representing the temperature at which both species show equal survival times.(b) Interspecific differences of mean thermal mortality during development for the periods where 30-day cumulative mortality was equal to or higher than 0.01.Colored areas represent the density functions of the mean thermal mortality in the 100 bootstrap replicates for each species and each data logger (each temporal series of microclimatic temperatures).The lower and upper limits of the boxes represent the 1st and the 3rd quartiles of mean thermal mortality and, their inner line, the median.(c) Interspecific ratio (P.napi/P.rapae) of mean thermal mortality during development observed in each microhabitat.For both species, estimated mortality was lower in their preferred microhabitats: Pieris napi showed inferior mortalities in semi-open and semi-closed microhabitats (left panels) but higher in open microhabitats (right panels).Microhabitats with cumulative mortality inferior to 0.01 during the whole period (March-September) are not shown.Me, mid-elevation; Ld, lowland sites.C, closed; SC, semi-closed; SO, semi-open; and O, open microhabitats.
Downloaded from https://esajournals.onlinelibrary.wiley.com/doi/10.1002/ecm.1561by JOSEP Penuelas -Readcube (Labtiva Inc. F I G U R E 2 Legend on next page.) , Wiley Online Library on [03/05/2023].See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions)on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License T A B L E 2 General linear mixed model of the time to larval death for each species.Results from the type III Wald χ 2 analysis of deviance of the fixed effects.Marginal R 2 for GLMMs (including only fixed effects) were calculated followingNakagawa et al. ( Note: