Local adaptation in thermal tolerance for a tropical butterfly across ecotone and rainforest habitats

ABSTRACT Thermal adaptation to habitat variability can determine species vulnerability to environmental change. For example, physiological tolerance to naturally low thermal variation in tropical forests species may alter their vulnerability to climate change impacts, compared with open habitat species. However, the extent to which habitat-specific differences in tolerance derive from within-generation versus across-generation ecological or evolutionary processes are not well characterized. Here we studied thermal tolerance limits of a Central African butterfly (Bicyclus dorothea) across two habitats in Cameroon: a thermally stable tropical forest and the more variable ecotone between rainforest and savanna. Second generation individuals originating from the ecotone, reared under conditions common to both populations, exhibited higher upper thermal limits (CTmax) than individuals originating from forest (∼3°C greater). Lower thermal limits (CTmin) were also slightly lower for the ecotone populations (∼1°C). Our results are suggestive of local adaptation driving habitat-specific differences in thermal tolerance (especially CTmax) that hold across generations. Such habitat-specific thermal limits may be widespread for tropical ectotherms and could affect species vulnerability to environmental change. However, microclimate and within-generation developmental processes (e.g. plasticity) will mediate these differences, and determining the fitness consequences of thermal variation for ecotone and rainforest species will require continued study of both within-generation and across-generation eco-evolutionary processes. This article has an associated First Person interview with the first author of the paper.


INTRODUCTION
It is presently widely accepted that ongoing climate warming has clear and widespread consequences for biodiversity including local extinction, population declines, shifts in community structure and composition, and changes in phenology (Scheffers et al., 2014). Physiologically based models show that tropical ectotherms may be more vulnerable to climate warming than temperate species by virtue of being adapted to lower thermal variation (Deutsch et al., 2008). However, the sensitivity of ectotherms to warming is also structured by experienced thermal variation across elevation and habitat (Sunday et al., 2014;García-Robledo et al., 2016). For species with broad distributions, populations may occur in diverse habitat types where selective pressures (e.g. temperature) during ontogeny may act differentially and have significant repercussions on physiological traits (Hoffmann et al., 2003). Local adaptation to thermal conditions can thus lead to population-specific responses (Somero, 2010;Kaspari et al., 2015;Nadeau et al., 2017).
Recent studies in tropical ecosystems have found that forest species of ectotherms tend to have lower tolerance to warming because they are restricted to more stable thermal regimes and might be especially vulnerable to climate change compared with more open habitat species (Huey et al., 2009;Frishkoff et al., 2015;Bonebrake et al., 2016;Nowakowski et al., 2017). However, ectotherms do have the ability to behaviorally avoid extreme temperatures through thermoregulatory activities and experiencing diverse microclimates at multiple spatial scales (Sunday et al., 2014;Bonebrake et al., 2014;Pincebourde and Woods, 2020). Furthermore, thermal tolerance differences in populations across habitats can sometimes be diminished by developmental conditions; Montejo-Kovacevich et al. (2020), for example, found that for tropical Heliconius erato butterflies, differences in heat tolerance across elevation were nearly erased after rearing in a common garden.
Measured thermal traits are ultimately a combination of natural selection/genetics, intergenerational plasticity (offspring traits affected by parental environmental conditions), developmental plasticity, and acclimation or reversible plasticity (Llewelyn et al., 2018). Indeed, acclimation within generations can drive CTmax variation for amphibians across habitats (Kristensen et al., 2008;Gunderson and Stillman, 2015;Simon et al., 2015). Upper thermal limits for terrestrial ectotherms appear to be fairly constrainedbut details for how such limits vary with acclimation and plasticity across habitats or microhabitats are needed (Hoffmann et al., 2013). Specifically, knowing the extent to which differences in habitatspecific thermal tolerance can be explained by processes within or across generations can help elucidate vulnerability to climate change, especially for tropical ectotherms.
Though tropical Africa is dominated by rainforest, there is important variation in ecosystem physiognomy; ecotones that consist of a mosaic of woody or herbaceous savanna and gallery forests representing a key habitat in the landscape and region thought to be important for speciation (Smith et al., 1997). Rainforest climates tend to be constant and regular due to high tree density and high annual average precipitation. Ecotones in contrast have relatively variable climates and consist of a mosaic of homogenous vegetation types leading to more open canopy spaces than rainforests. In this study we investigated the role of habitat in structuring thermal limits, critical thermal maximum (CTmax) and minimum (CTmin), for Bicyclus dorothea (Cramer, 1779), a nymphalid butterfly found in both forest and ecotone habitats in many tropical African countries (Dongmo et al., 2017). Bicyclus dorothea exhibits differential variability in wing morphology across ecotone and forest habitats ecotone in Cameroon (Dongmo et al., 2018). We hypothesized that B. dorothea populations from the ecotone would have wider thermal tolerance breadths than forest populations (reflecting the greater thermal variation in those environments) that would be preserved after lab-rearing under common conditions.
Simulated microclimatic variation demonstrated that ecotone sites exhibit lower minimum and higher maximum temperatures than forest sites (Fig. 3). Shade and height have little influence on minimum temperatures, but 1 cm maximum temperatures with 0% shade are considerably higher (by over 10°C) than both 120 cm and 1 cm with 100% shade maximum temperatures (which are very similar to one another; Fig. 3). Relating to the measured thermal tolerance limits, mean CTmin estimates were much lower than simulated minimum temperatures (about 10°C lower) while mean CTmax estimates were comparable to simulated maximum temperatures at 1 cm with 0% shade (Fig. 3).

DISCUSSION
The higher CTmax (by ∼3°C) in ecotone B. dorothea populations relative to forest populations, assessed using second generation individuals reared in a common environment, indicates a level of habitat-specific adaptation in thermal tolerance for this species. We also found that CTmin was lower for ecotone populations, but to a smaller extent (by ∼1°C). At least for large spatial gradients then (hundreds of kilometers, Fig. 1), our results show that local adaptation to habitat climatic conditions may result in differential thermal tolerances for species which may have consequences for vulnerability to environmental and climatic change.
Our results support the evidence of local adaptation playing a role in variable thermal tolerance limits across habitat, as has been shown in a variety of systems including tropical beetles by García-Robledo et al. (2016) and tropical lizards by Moritz et al. (2012). However, variation in thermal limits can also be a consequence of developmental effects, ontogenetic changes, and phenotypic plasticity (Bowler and Terblanche, 2008;Kellermann and Sgrò, 2018). As an example, and in contrast with our own findings, Montejo-Kovacevich et al. (2020) found little difference in thermal tolerance (knockdown time) for a tropical butterfly species across elevation following rearing in a common garden environment, despite finding significant differences in wild-caught individuals from high versus low elevation populations.
Following the terminology of Hoffmann and Sgro (2018), we have, through our experimental set up, been able to infer a genetic (across-generation) basis of thermal traits in B. dorothea across habitats, which may be an indication of local adaptation; however, we have not been able to determine what the field fitness consequences of these differences might be. Through simulations of microclimatic variation across habitat (Fig. 3), there is evidence that the CTmax differences between ecotone and forest are very similar to the maximum temperature differences at 1 cm with no shade. The observed CTmax differences therefore may have some ecological relevance to the experienced thermal variation. However, the microclimatic data also demonstrated how temperatures higher off the ground (120 cm) and in the shade lower maximum temperatures considerably (Fig. 3). These are conditions under which B. dorothea are commonly active. This result again emphasizes the importance of microclimates, and solar radiation in particular, in determining vulnerability or resilience of small ectotherms to extreme temperatures (Bonebrake et al., 2014;Pincebourde and Suppo, 2016). CTmin was much lower than modeled minimum temperatures, which may explain why the observed differences across habitats were minimal. This result also suggests that CTmin may not have great consequences for B. dorothea for our study sites.
We know from previous study of B. dorothea in Cameroon in the field, that habitat-specific differences in other traits (wing phenotypes) are clear between ecotone and forest (Dongmo et al., 2018). This then further supports the evidence for local adaptation within habitats. However, to date, no study has explored dispersal capacity, though the species does not appear to fly long distances regularly (Dongmo et al., 2017). In-depth genetic study would be ideal for better revealing the underlying processes structuring thermal tolerance and other trait differences across habitats. Research across the same ecotone-forest gradient in Cameroon, has shown evidence of morphological and genetic differentiation in the lizard Trachylepis affinis (Freedman et al., 2010) in addition to thermal performance variation similar to our findings, i.e. ∼3°C higher thermal optimum and ∼1°C higher CTmax for wild-caught ecotone lizards (Landry Yuan et al., 2018).
Other than the lack of field fitness study in our system, other limitations of this study are important for interpretation of our results. While we found different thermal limits across habitats for populations reared under a shared environment, all colonies were reared at a single temperature (26°C). The magnitude of responses of all populations to heat and cold stresses may have been less differentiated if populations were acclimated at different temperatures. We also used a fast-ramping rate (0.5°C per minute), which may increase the relevance of our results to warming (Rezende et al., 2011), but the precise ecological relevance of this ramping rate in our system remains unknown. Furthermore, heat hardening (short-term acclimation to warm temperatures) has been demonstrated in Bicyclus anynana (Fischer et al., 2010) such that short-term plastic responses to thermal variation may be key in B. dorothea and other tropical insects. Ultimately, insect responses to climate change will be the result of a complex interplay between behavior, phenology, evolution and plasticity in response to thermal variation across multiple temporal scales (Bonebrake et al., 2014;Kingsolver and Buckley, 2017).

Sampling sites
Bicyclus dorothea sampling took place in four localities in Cameroon (Fig. 1). Two sites represented tropical rainforest: one in a forest near the locality of Mbalmayo (N 3.388, E 11.47, alt. 768 m above sea level) characterized by a degraded secondary forest, and another site in Somalomo (N 3.37405, E 12.7332, alt. 638 m above sea level) with a mature forest located in the Dja faunal reserve. Two other sites represented ecotone habitats: Ako (N 6.68783, E 10.70687, alt. 706 m) and Ndikiniméki (N 4.76986, E 12.7332, alt. 812 m), both characterized by rainforestsavanna mosaic (Fig. 1). For each site we modelled local microclimatic variation using microclim simulations  at variable heights and shade conditions. At each locality (from a ∼15 km resolution), estimates of hourly minimum and maximum temperatures on a long-term climatology were downloaded (Hijmans et al., 2005;Kearney et al., 2014). Twenty-four layers (one for Base map represents forest cover as estimated from land cover GLC2000 (http:// www.diva-gis.org/gdata). each hour of the day from 0:00 h to 23:00 h) of hourly temperature were downloaded from which maximum and minimum temperature values for each latitude/longitude of the study sites were extracted with the packages 'raster' (Hijmans, 2020) and 'ncdf4' (Pierce, 2019). We then took the mean monthly values for multiple microclimatic conditions: (1) 1 cm from the ground with 0% shade, (2) 1 cm from the ground with 100% shade, and (3) 120 cm from the ground . As B. dorothea is an understory species typically flying low to the ground under variable levels of canopy cover, then these conditions are likely to be well representative of the microclimatic variation experienced by individuals in the sites.

Establishment of B. dorothea laboratory colonies
Laboratory colonies from each site were established with adult female butterflies collected with overripe-banana baited traps and hand-net captures. Sampling of wild adult butterflies was conducted during the wet season, which runs from April to October while the dry season generally starts in early November and ends in March (for all sites). The wet season was chosen because B. dorothea adults are scarce and generally cease reproduction during the dry season (Dongmo et al., 2017), while during the wet season, resources are abundant, matings are common, and females can lay large quantities of fertile eggs. For logistical reasons, sampling was then carried out in July 2016 for the localities of Ako and Ndikiniméki and in May 2017 for Mbalmayo and Somalomo.
Adult butterflies collected at each locality were brought to the laboratory at the International Institute of Tropical Agriculture (IITA) in Yaoundé, Cameroon. For each population, at least 50 gravid females were captured from the field and kept in cages (ten individuals per cage) with mashed banana, distilled water and potted millet [Pennisetum glaucum (L.) R. Br., 1810] was used as the egg-laying host plant. Eggs (F1) laid by the fieldcollected females were collected and kept in Petri dishes lined with a moistened filter paper until hatching. Hatchlings (from all populations) were reared in population cages (24×24×24 cm, made of polyester screen to facilitate air circulation) and fed on potted lawn [Axonopus compressus (Sw.) P. Beauv., 1812] in a room where temperature was maintained at a constant 26°C and high relative humidity at 80%, and a 12D:12L photoperiod. Pupae were collected daily and transferred to individual cages where they developed to eclosion. Under these room conditions the mean development time (from egg to adult) was about 45 days.
To avoid maternal effects on thermal tolerance, all adult butterflies used in this experiment were from the second generation (F2) of the laboratory colony. Second generation butterflies of each population were obtained by allowing first generation adults to mate randomly with other adults who also originated from the same initial population cage exposed in a room where the conditions were the same as in the rearing of the first-generation individuals, i.e. 26°C and 80% relative humidity and a photoperiod of 12D:12L. Each of these cages (24×24×24 cm) contained ten adults (five males and five females) from the first generation, mashed ripe banana, cotton soaked with distilled water and potted millet P. glaucum as the egglaying host plant, which was changed as needed. Multiple cohorts of butterflies (from all four sites) were reared between August 2016 and November 2017.

Thermal tolerance
To measure critical thermal minimum (CTmin), 1-day-old second-generation adult butterflies belonging to each of the population cages (four in total for each locality) were placed individually in small plastic 250 ml cups with about 25 holes (about 5 mm diameter to allow good air circulation with the environmental chamber) in their wall and were placed in an environmental chamber (I-36VL Percival Scientific Inc., Perry, IA, USA) initially set at 26°C (Fischer et al., 2010). After 1 h, the climate cabinet was set to a ramping mode with temperature decreasing at rate of 0.5°C per minute (Piyaphongkul et al., 2012). We observed butterflies' reaction to temperature variation through a window incorporated into the main door of the climate cabinet (which remained closed during the procedure). To monitor the real-time temperature inside the climate cabinet, a water proof thermometer probe (DE:30W, DER EE, New Taipei City, Taiwan) was displayed inside the climate cabinet in such a way that values on the thermometer could be read easily through the window of the climate cabinet. The critical thermal minimum was the coldest temperature at which adult butterflies were not able to flap their wings or make any movement with their appendages.
For CTmax, 1-day-old butterflies were also used. The same method was used but the climate cabinet was set in a ramping mode with increasing temperature at the same rate (0.5°C per minute). The critical thermal maximum was the high temperature at which each butterfly was not able stand on their legs (Huey and Stevenson, 1979;Lutterschmidt and Hutchison, 1997). For each individual butterfly, we first measured the critical thermal minimum and after that, it was taken back from the climate cabinet to the room maintained at 26°C for recovery. The critical thermal maximum was assessed on each recovered individual 24 h later. Individuals that died after the assessment of the critical thermal minimum were excluded from the analysis. From the CTmin and CTmax values of each individual, we calculated the thermal tolerance breath of each individual as the difference between CTmax and CTmin.

Data analysis
The effects of habitat, sex and sampling site on the critical thermal minimum, maximum, and tolerance breadth were analyzed using a nested ANOVA model, with sampling sites nested in habitat (forest and ecotone). In order to meet ANOVA requirements, data were log-transformed. Pairwise comparisons were performed with Tukey's HSD. All statistical analyses were done in R version 3.5.1 (R core team, 2019). Thermal tolerance data is available through the Data Dryad repository.

Ethical statement
All research was carried out in accordance with permission by the Ministry of Scientific Research and Innovation of Cameroon (MINRESI). and minimum (B) temperatures under variable conditions (1 cm 0% shade, 1 cm 100% shade, and 120 cm) across habitats using the microclim dataset . Mean CTmax (A) and CTmin (B) for forest versus ecotone populations represented by horizontal lines, assessed using second-generation individuals reared under common garden conditions in the laboratory.