Skip to main content

Historical changes in thermoregulatory traits of alpine butterflies reveal complex ecological and evolutionary responses to recent climate change



Trait evolution and plasticity are expected to interactively influence responses to climate change, but rapid changes in and increased variability of temperature may limit evolutionary responses. We use historical specimens to document changes in the size and thermoregulatory traits of a montane butterfly, Colias meadii, in Colorado, USA over the past 60 years (1953–2012). We quantify forewing wing length, ventral wing melanin that increases solar absorption, and the length of thorax setae that reduces convective heat loss.


The mean of all three traits has increased during this time period despite climate warming. Phenological shifts may have extended the active season earlier at low elevations and later at high elevations, increasing exposure to cool temperatures and selecting for increases in thermoregulatory traits. Fitness increases at higher elevations due to warming could also increase thermoregulatory traits. Warmer temperatures during pupal development and later flight dates in the season are associated with decreased wing melanin, indicating a role of phenotypic plasticity in historical trait changes.


Phenotypic shifts result from a complex interplay of ecological and evolutionary responses to climate change. Environmental variability within and across seasons can limit the evolutionary responses of populations to increasing mean temperatures during climate change.


For many organisms, the rate of climate change may be sufficiently fast to preclude evolutionary adaptation [1, 2]. However, there is some evidence for evolutionary changes in body size [3], coloration relative to background matching [46], and phenology; specifically the timing of diapause [7], of nesting [8] and of flowering [9]. Evidence for evolutionary or historical changes in traits relevant to thermal physiology is limited (but see [1013]) Natural history collections offer a largely untapped resource to examine morphological responses to climate change [14]. In particular, insects are abundant in these collections and the thermal significance of morphological traits has been extensively documented for some species.

A biophysical model for Colias butterflies enables estimating how morphological thermoregulatory traits and environmental conditions interact to determine body temperatures [1518]. Specifically, the body temperatures of butterflies increase proportionally to increases in air temperature [17]. Two key thermoregulatory traits influence body temperatures: 1) Increased melanism on the posterior ventral hindwings allows the butterflies to absorb more solar radiation via closed wing basking, and 2) longer setae on the thorax reduce convective heat loss [19, 20]. Wing melanin is heritable and increases with elevation for Colias species [21], suggesting its importance in determining fitness. Moreover, the biophysical model enables linking phenotypes to function because Colias flight is restricted to a narrow range of body temperatures (28°–40 °C) [17, 20]. Thus their flight time is often limited by climate, especially at higher elevations and latitudes [2224]. Flight time is a determinant of reproductive success because flight is required for activities including mating, feeding, and oviposition [16, 19]. Conversely, butterflies may overheat, which not only limits flight time, but can also directly reduce the number of viable eggs produced by each female [25]. Past Colias research provides a mechanistic understanding of how environmental conditions interact with thermoregulatory traits, which are well preserved in collections, to determine fitness.

Application of the biophysical model suggests that recent warming in Colorado has increased both flight time and overheating risk for a montane Colias species, C. meadii [16]. Translating estimates of flight time and overheating risk into fitness based on egg-laying suggests that recent warming has driven declines in fitness at lower elevations and fitness increases at higher elevations for C. meadii with constant thermoregulatory traits [16]. An analysis of selection gradients complicated our initial expectation that climate warming would select for decreased melanism. The analysis estimated that directional selection for decreased wing melanism has become more consistent and stronger in the lower elevation subalpine portion (3.04 km site) of C. meadii’s distribution over the last 60 years in Colorado [26]. The analysis predicted consistent selection for increased wing melanism at the high elevation extent (3.51 km site) where flight time is more limited, but the prediction could be an artifact of insufficiently capturing temperature extremes that could decrease survival of individuals with high wing melanism. The model also indicates that seasonal and annual variation in climate generates substantial variation in the strength and direction of predicted selection, which would dampen the potential for phenotypic shifts in response to recent climate change [26]. To address these predictions, we measured the ventral hindwing melanism, thorax setae length and wing length (proxy for body size) in historical specimens. Available specimens were restricted to a fairly narrow elevation range (25 to 75% quantile: 3.36 km to 3.64 km) spanning the high elevation site in the model. Thus, we expect that the direction of selection may vary with elevation across our sampled populations and that environmental variability may act to dampen directional selection.

Phenotypic plasticity is also likely to play an important role in responding to climate change, particularly in environmentally variable regions [2, 27]. Plasticity may either facilitate evolution through enabling persistence or hinder evolution through buffering selection [28]. Montane Colias exhibit plasticity where increased temperature during pupal development decreases melanism [10, 29]. Thus, we expect directional selection and phenotypic plasticity to interact for montane Colias.

Climate warming is also expected to shift insect body size via selection and plasticity. Warmer temperatures accelerate development, which can shorten the duration of development and reduce adult size (temperature-size rule) [30]. Warmer environments can also select for decreased size to decrease body temperatures [17]. Consequently, predictions and observations of insect responses to climate change generally indicate a decrease in body size [3, 31]. However, spring warming can extend the time available for development and increase size [32].

Here we examine the roles of evolution and plasticity in responding to climate change by using museum specimen collections to document shifts over the past 60 years (1953–2013) in morphological traits in C. meadii butterflies in Colorado. We evaluate several factors expected to determine trait values via both phenotypic plasticity and selection. We examine developmental temperatures, which determine phenotypic plasticity, and temperatures during the previous flight period, which determine the reproductive success of parents. We expect mean trait values to decrease with the calendar date of collection because mean temperatures increase and daily variability decreases through the season (analysis of data from accessed 4 May 2015). Given the seasonal progression of temperatures, we also consider the interactions between developmental (or pupal) temperatures and collection date. Finally, we include collection year as an indicator of evolutionary shifts in phenotypes. However, we are unable to control for other factors that vary over long time scales such as trends in temperature mean and variability omitted from our temperature measures; trends in habitat, host plant availability, wind speed, or cloud cover; and collection biases.


Trait measurements

Colias meadii occurs in subalpine and alpine meadows above 2.50 km elevation throughout the Rocky Mountains, from northern New Mexico, USA to southern Alberta, Canada [20]. They have a single generation each year (univoltine) and diapause over the winter as 3rd instar larvae; the adult flight season lasts 4 to 6 weeks in July and August (between Julian days 180–240). The subspecies C. meadii meadii is found largely in the southern part of the range at elevations between 3 and 4 km. It is genetically distinct from more northern populations of the species in Wyoming Basin [33]. In this study we focus on C. m. meadii in the southern part of the range, where summer temperatures have increased significantly over the last 60 years [16]. In this region, C. m. meadii individuals pupate in early to mid-June (between Julian days 150–170), and adult butterflies emerge in early July.

We processed specimens from museums in the United States with the largest holdings of C. m. meadii from Colorado (the Yale Peabody and the McGuire Center for Lepidoptera and Biodiversity); museums located in Colorado (the C.P. Gillette Museum of Arthropod Diversity at Colorado State University and the University of Colorado Museum); and museums with holdings from time periods absent from other museums (the Smithsonian Museum of Natural History, the Milwaukee Public Museum, and the California Academy of Science). Because museum holdings of C. m. meadii are sparse for the past 25 years, we also solicited samples from private collectors via the Yale Lepidoptera list serve and measured samples from one private collector.

Because the collections contained few females after 1990, we restricted our analyses to males (N = 379 from 20 counties). While females tend to have darker wings than males, elevation clines in melanism are similar between the sexes [21]. Model predictions for temporal trends in melanism are based on the temperature dependence of flight time and egg viability, which were translated into female reproduction [16, 28]. The flight time constraints should influence male fitness similarly, so predicted temporal trends also apply to males. The similar elevation clines between sexes suggest that they respond to climatic gradients similarly and our observations should be indicative of both sexes.

All specimens for the C.m.meadii range were included except two specimens with a Julian collection date of 170 (June 19, collected at Como, Colorado 39.31 N and 105.89 W in 1984, McGuire Center collection) and 171 (June 20, collected at Los Pinos Pass, Colorado 38.1 N and 106.97 W in 1967, Yale Peabody collection), because they strongly influenced regressions and the date was more than two standard deviations below the mean (mean ± sd = 208.4 ± 13.65). Including these specimens increased the slope of our regressions and significance of our results.

Pinned specimens were measured by removing labels with forceps, transcribing the locality data, and then placing the head of the pin in a lump of modeling clay to expose the ventral hindwing (the wing region accounting for the majority of absorption for this ventral basking species). Elevation was recorded from the locality tag in the collections if available; otherwise we used a digital elevation model to estimate elevation based on latitude and longitude. We used Google Earth to georeference locality descriptions when necessary.

A digital micrometer was used to measure the forewing of the specimen from the thoracic insertion to the apex of the wing. The specimen was then photographed through a 100-mm macro lens in RAW format with a Canon Rebel XSi mounted on a copy-stand. Each image included a black and white standard. Because the height of the animal on the pin was variable, we used auto focusing to allow for the clearest image. We measured setal length on the ventral thorax with an ocular micrometer on a Wild M5 microscope as the longest setae between the first and second leg. All measurements were taken by MacLean and specimens were prepared for measurements in groups of five to obscure each specimen’s metadata during measurements. Specimens were distributed randomly throughout collection boxes and thus were measured randomly within each collection.

We photographically assessed the degree of wing melanism on the posterior ventral hind-wing. First, we selected a triangulated region between the eyespot, hind wing insertion, and the wing margin [21, 34]. Using a MatLab program (T. Hedrick, unpublished), we then converted the RAW image to black and white to account for potential fading. Black and white is appropriate because wing scales are pigmented by either melanin (dark) or pterin (light). We digitized the region of interest and the black and white standards for each sample, and then used the standards to calculate a standardized grey-level value between 0 (white) and 1 (black). The grey-level value should be robust to butterfly size because it represents the proportion of melanic wing scales. In order to verify our measure of grey-level as a proxy for absorptivity, samples collected in the field 2012 and 2013 were photographed, and the absorption spectrum was measured from 350 to 1050 nm in a spectroreflectometer with an optical integrating sphere (FieldSpecPro FEFR 7501, ASD Inc) for the same wing region. Previous studies used absorption at 650 nm as a measure of melanism and solar absorptivity for Colias wings (Watt 1968, Kingsolver 1983). Grey level and absorption at 650 nm were highly correlated (n = 60, R2 = 0.78), confirming that grey level is an appropriate measure of solar absorptivity.

Climate data

To evaluate temporal changes in climate, we selected a representative NOAA Cooperative weather station (data retrieved from the Western Regional Climate Center, on 4 May 2015) at 3.5 km in Climax, CO (39.37 N, 106.18 W). This weather station was selected for its relevant elevation, its central location among sampling sites (median distance 61 km, 10th quanitile = 23.3 km and 90th quanitile = 96.0 km), and the completeness of its climate records throughout the time period of interest. We chose to use temperature data from one representative weather station because long-term data are lacking for most stations above 3 km in Colorado and, while interpolated data sets offered more complete coverage, they can perform poorly at capturing temporal variation for high elevation sites in the Rocky Mountains [35].

Data analysis

We estimated temperatures expected to influence each focal trait based on lab observations for the closely related Colias eriphyle. Melanism (and we assume setae length) is influenced by pupal temperatures (F1,54 = 6.4, p < 0.01 [36]) and pupation lasts for 14 days at 20 °C. We thus estimated temperatures during pupal development using data from 26 days to 6 days prior to collection (with the slightly extended duration intended to account for uncertainties). Forewing length is influenced by temperature during both late larval development (determining pupal mass: F1,108 = 676.5, p < 0.0001) [36]) and pupation (determining wing length: F1,101 = 6.5, p = 0.01 [36])). We thus estimate developmental temperatures using data from 46 days to 6 days prior to collection. We consider temperatures experienced by adults during the previous flight season, which determines reproductive success via flight time and overheating. We estimated temperatures during the previous flight season as the mean across Julian days spanning 201 to 216. The flight season was determined by taking the 25th to 75th quantiles of collection date. We used a constant flight season across years because we had little basis for predicting when in the flight season a parent flew. We examined sensitivity to our selection of time windows for developmental and pupal temperatures and found the model results to be robust. We additionally considered substituting flight season temperatures from the previous year with temperatures during the current year. Temperatures from the previous year were stronger predictors, but including current year temperatures did not qualitatively alter our model.

We focus on temperature, but note that shifts in cloudiness and radiation could also influence plasticity and evolution. Increases in cloudiness on the order of 1.4% (of the sky) per decade have been estimated for the United States between 1976 and 2004 [37]. However, shifts are regionally variable. Analysis of solar radiation data proximate to our specimens suggests that solar radiation has increased over our observation period [38], consistent with other observations for the Rockies [39]. However, all these estimates have very low confidence due to a shift in techniques for cloud observations (which also influences solar radiation data derived from modeling) [37]. The radiation trends are sufficiently uncertain and of a magnitude that justifies our focus on temperature. Additionally, our focus does not include wind speed. While wind plays a role in the heat budget for these butterflies, there is little evidence (beyond a weak change in mean at high elevation) of a temporial trend at these sites [40].

We account for uneven collection intensities by restricting our analysis to 15 individuals per site per year. We then averaged the model output to perform model selection. We randomly sampled specimens each of the 50 times we repeated the analysis. Statistical model results were similar across iterations of the analysis, so we report results averaged across the 50 iterations. We performed model selection and averaging using the R package MuMIn [41] (commands: dredge and model.avg). For each iteration of the analysis, we selected the top 20 sub-models according to AICc (sample size corrected Akaike Information Criterion). All predictor terms were included among the best (delta AICc < 2) sub-models and the full model received strong support (AICc < 0.3). As a result, we focus our results on the full model (trait ~ developmental T (or Pupal T) + previous flight season T + Date of Collection + Year + developmental T (or Pupal T) x Date, where T is temperature in °C), but additionally report the percent of the top sub-models that included each predictor term. We present the full results of the model-averaging in Additional file 1: Table S1. We normalized the predictor variables using the scale function in R. We assessed collinearity using the variance inflation factor (vif function from the R car package). We did not detect problematic collinearity among our predictor variables (vif <2 for all main effects). Model selection excluded elevation as a poor predictor of trait values, likely because collections were concentrated over a fairly narrow elevation range (25 to 75% quantile: 3.36 km to 3.64 km).

To account for trait similarity due to geographic proximity, we used maximum-likelihood spatial autoregressive (SAR) models ([42], R package spdep [43]). We selected an SAR error model over a lag or mixed model following assessment using a Moran’s I test on the residuals of linear models as well as a Lagrange multiplier specification test. We used Moran’s I tests, estimates of model performance, and spatial correlograms to select a 40 km threshold distance when developing spatial neighborhoods based on latitude and longitude. Neighbors were weighted using row standardization. We used a Moran’s I test on the model residuals to confirm that the SAR models fully accounted for spatial autocorrelation in the data [42]. We report log-likelihood tests comparing the SAR to a null model consisting only of an intercept along with Nagelkerke pseudo R2 values. We compared the SAR model output to linear models and found no qualitative differences in the direction of the effects.

We used an ANOVA to evaluate phenological trends over time. We standardized both year and elevation to the mean. Model selection identified elevation as an important predictor of phenology, despite its being omitted from the best models for trait values.


Mean forewing length, mean wing melanism, and mean setal length have increased during the past 60 years, despite increasing mean July temperatures during this time period. Forewing length has increased by 11% over 60 years (Fig. 1, Table 1, log likelihood ratio: 150.9, P < 0.01, R2 = 0.31). Wing melanism has increased significantly (21%) over the study period (Fig. 2, Table 1, log likelihood ratio: 40.9, P < 0.01, R2 = 0.17). Similarly, thoracic setae have lengthened significantly (68%) over the study period (Fig. 3, Table 1, log likelihood ratio: 149.1, P < 0.01, R2 = 0.46).

Fig. 1
figure 1

Partial residuals (residuals of regressing the response variable on the independent variables, but omitting the independent variable of interest) for predictors of forewing length (mm). While forewing length tends increase with increasing developmental temperatures (non-significant correlation) as well as collection date, developmental temperatures and collection date interact. Forewing length increases with year

Table 1 The output of spatially autoregressive maximum likelihood models describing wing length (mm)
Fig. 2
figure 2

Partial residuals for predictors of wing melanism (grey level). Wing melanism decreases with increasing pupal and flight season temperatures as well as collection date, but pupal temperatures and collection date interact. Wing melanism increases with year

Fig. 3
figure 3

Partial residuals for predictors of setae length (mm). Setae length increases with year and decreases with the date of collection. Setae length is influenced by pupal and flight season temperatures as well as the interaction of pupal temperatures with collection date

Influences of temperature and collection date on thermoregulatory traits suggest plastic responses to environmental temperatures. Non-significant tendencies for melanism to decrease with increases in pupal and flight season temperatures (Fig. 2) have high relative importance across models (relative importance values of 0.66 and 0.99, respectively: see Additional file 1: Table S1). Melanism decreases significantly with collection date, such that lower melanism is found later in the season (relative importance of 0.92). An interaction between pupal temperatures and collection dates also receives moderate relative importance across models (relative importance of 0.38). The interaction represents a correction accounting for the tendency for butterflies collected later in the season to have experienced warmer temperatures during pupation. Developmental temperatures and temperatures during the previous flight season received support as predictors of wing length (relative importance 0.89 and 0.34, respectively) with a tendency for warmer temperatures to increase size (Fig. 1). The interaction between developmental temperatures and collection date received some support as a predictor of forewing length (relative importance of 0.27). Setae length is increasing in response to increased temperatures in the previous flight season (Table 2). Neither pupal temperatures nor their interaction with collection date are significant predictors of setae length, but they receive support for model inclusion (relative importance of 0.85 and 0.24, respectively, Fig. 3).

Table 2 The output of spatially autoregressive maximum likelihood models describing wing melanism (grey level) and thorax setae length (mm)

We also find evidence of a phenological shift, which depends on elevation. Collection date has shifted later over time, but there was a significant interaction such that low elevation individuals were collected earlier and high elevation individuals were collected later over time (Fig. 4, ANOVA for predictors of collection date: year slope(se) = 2.8(0.71), P < 0.001;elevation slope(se) = 0.92(0.67), P = 0.18; interaction: slope(se) = -1.88 (0.90), P = 0.04; F[3,375] = 4.4, P < 0.03).

Fig. 4
figure 4

Partial residual plot for normalized predictors of collection date (J). Collection date advances with year and elevation. However, the predictors interact such that collection date has shifted earlier in the season at lower elevation and later in the season at higher elevations in later years


Phenotypic plasticity is increasingly viewed as a primary response to environmental change [2, 44, 27]. Plasticity is increasingly being incorporated in predictive models of phenotypic change and warming tolerance [45] and being explicitly tested for in empirical studies [2]. The interaction between phenotypic plasticity and adaptive evolution in responding to environmental change remains unclear: does phenotypic plasticity tend to facilitate adaptive evolution by enabling persistence or does it tend to slow adaptive evolution by buffering selection [2, 45]?

One striking result from our studies is that mean forewing length, wing melanism, and setal length in C. meadii have increased, not decreased, during the past 60 years in this region. These phenotypic changes have occurred despite overall climate warming over this time period: mean July maximum temperature increased significantly from 1953 to 2013, (p < 0.001; slope estimate (se) = 0.03 (0.008) °C), resulting in a total increase of ~ 1.8 °C over 60 years (Additional file 1: Figure S1A). These historical shifts are inconsistent with simple expectations for evolutionary responses to climate warming. However, previous analyses suggested complex patterns of selection: increases in melanism are consistent with our expectations for high elevation populations that are flight limited [26]. Previous modeling of C. meadii also highlighted the importance of inter-annual variability [26]. Such environmental variation generates fluctuations in the magnitude and direction of selection that strongly reduce the rate of adaptive evolution to directional environmental change [26, 46].

Several lines of past and current evidence support the contribution of phenotypic plasticity to variation in wing melanism. Lab experiments show that lower temperatures experienced during pupal development increase adult wing melanism in montane Colias [29, 36, 47]. Our historical analyses indicate that greater wing melanism is associated with lower pupal temperatures, and with earlier collection dates during the flight season (when temperatures are on average lower). These patterns are consistent with the idea that phenotypic plasticity in wing melanism in response to variation in pupal temperature is contributing to the historical pattern in this trait.

Whether such plasticity is adaptive depends on whether environmental cues during pupal development are good predictors of environmental conditions experienced by adults [4850]. June and July temperatures in this region were only weakly correlated during the past 60 years (r = 0.19, p = 0.18), suggesting that pupal temperature in June is often a poor predictor of temperatures during the flight season during July (Additional file 1: Figure S2). When we consider the estimated pupal temperature for the specimens included in the analysis, there is a slight but significant increase in mean pupal temperature over the whole 60 year time period, but the regression explains only 3.2% of the variation in pupal temperature (slope = 0.015, R2 = 0.032), with tremendous variability among years (residual variance = 5.50 °C2). This interannual variability may be contributing to the observed increase in melanism and setae length by causing variation in the direction of selection and mismatching developmental and flight temperatures. Thermal extremes that cause overheating are relatively rare in our montane sites [16], their incidence may be mostly independent of mean active season temperatures, and variation in thermoregulatory traits may be relatively ineffective in avoiding overheating. Thus, selection to increase flight time in cold years may be stronger and more consistent than selection to avoid overheating in warm years.

Change in the distribution of collection dates of historic specimens can indicate shifts in phenology [51]. We detect a phenological shift toward later collection in later years that interacts with elevation: individuals were collected earlier in the season at lower elevations and later in the season at high elevation in later years (Fig. 4). Climate warming may have extended the active season of high elevation individuals. The dependence on elevation may also reflect greater opportunity for phenological advancement at lower elevations where snow melts earlier. Exposure to cool, early or late season temperatures during development could select for wing darkening and setae lengthening and extend the developmental duration, increasing body size. Phenological shifts could result in the collection of individuals with darker wings and longer setae due to plastic responses to experiencing colder conditions earlier in the season at low elevation or cool late season conditions at higher elevations. The influences of phenological shifts at low and high elevation are difficult to distinguish.

The phenological shift toward later collection in later years at high elevation may indicate individuals from high elevations increasing fitness relative to lower elevation individuals as climate warms. These shifts could result in a higher proportion of individuals with high elevation phenotypes being sampled over time. Specimens in collections are concentrated over a sufficiently narrow elevation range that we do not detect the elevational gradient in melanism that is well established [20, 52]. The only significance effect of elevation is the interaction with year as a predictor of collection date. Fitness increases at higher elevation could have resulted in dispersal of individuals with high elevation phenotypes to lower elevations.

Our results suggest that any evolutionary responses in Colias butterflies have been more complex than a steady shift toward lighter wings, shorter setae, and smaller body size in response to warming. Factors such as mate choice [53] or interspecific competition may also contribute to selection on wing coloration [54], further complicating temporal trends. Collection biases may have also contributed. Collected individuals are biased towards those in flight, but the importance of flight to fitness suggests that collections are consistently biased toward individuals with trait values associated with higher fitness.

Related studies suggest that our findings may be broadly applicable to understanding how thermoregulatory traits mediate interactions with the environment. An analysis of images in field guides for 473 European butterfly and dragonfly species indicated that cooler climates contain darker species and that distribution shifts between 1988 and 2006 resulted in assemblages becoming lighter in regions that warmed [11]. However, whether wing coloration has any thermoregulatory significance was unclear for most species in the analysis [55, 56]. Although increasing melanism in response to decreasing developmental temperatures has been widely reported in insects, this includes many cases in which melanism plays no role in thermoregulation [5760]. Our findings correspond to those of J Stamberger [34] who documented decreases in wing melanism and setal length in C. meadii with warmer temperatures, but no trend across years, for a more restricted time period and set of sites in Colorado.


Our finding that thermoregulatory traits and body size have increased over time highlights the complexity of biological responses to climate change. Environmental variability within and across seasons can limit or alter evolutionary responses. Variation in selection across spatial climatic gradients can also result in trait shifts, particularly if warming differentially impacts fitness across climatic gradients. Phenological shifts can also alter selection. Our results also highlight the importance of considering how plasticity can modulate species’ responses to climate change. Temperature variability both within and across seasons and the potential for phenological, fitness, and distributional shifts may lessen the potential for responding to climate change via adaptive evolution.



Temperature in degrees Celcius


Spatially autoregressive


Akaike Information Criterion


  1. Hoffmann AA, Sgrò CM. Climate change and evolutionary adaptation. Nature. 2011;470(7335):479–85.

    Article  CAS  PubMed  Google Scholar 

  2. Merilä J, Hendry AP. Climate change, adaptation, and phenotypic plasticity: the problem and the evidence. Evol Appl. 2014;7(1):1–14.

    Article  PubMed  PubMed Central  Google Scholar 

  3. Sheridan JA, Bickford D. Shrinking body size as an ecological response to climate change. Nat Clim Chang. 2011;1(8):401–6.

    Article  Google Scholar 

  4. Galeotti P, Rubolini D, Sacchi R, Fasola M. Global changes and animal phenotypic responses: melanin-based plumage redness of scops owls increased with temperature and rainfall during the last century. Biol Lett. 2009;5(4):532–4.

    Article  PubMed  PubMed Central  Google Scholar 

  5. Mills LS, Zimova M, Oyler J, Running S, Abatzoglou JT, Lukacs PM. Camouflage mismatch in seasonal coat color due to decreased snow duration. Proc Natl Acad Sci. 2013;110(18):7360–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  6. Karell P, Ahola K, Karstinen T, Valkama J, Brommer JE. Climate change drives microevolution in a wild bird. Nat Commun. 2011;2:208.

    Article  PubMed  PubMed Central  Google Scholar 

  7. Bradshaw WE, Fujiyama S, Holzapfel CM. Adaptation to the thermal climate of North America by the pitcher-plant mosquito, Wyeomyia smithii. Ecology. 2000;81(5):1262–72.

    Article  Google Scholar 

  8. Visser M, Van Noordwijk A, Tinbergen J, Lessells C. Warmer springs lead to mistimed reproduction in great tits (Parus major). Proc R Soc London, Ser B. 1998;265(1408):1867–70.

    Article  Google Scholar 

  9. Anderson JT, Inouye DW, McKinney AM, Colautti RI, Mitchell-Olds T. Phenotypic plasticity and adaptive evolution contribute to advancing flowering phenology in response to climate change. Proc R Soc B Biol Sci. 2012;279(1743):3843–52.

    Article  Google Scholar 

  10. Higgins JK, MacLean HJ, Buckley LB, Kingsolver JG. Geographic differences and microevolutionary changes in thermal sensitivity of butterfly larvae in response to climate. Funct Ecol. 2014;28(4):982–9.

    Article  Google Scholar 

  11. Zeuss D, Brandl R, Brändle M, Rahbek C, Brunzel S. Global warming favours light-coloured insects in Europe. Nat Commun. 2014;5:3874.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Cowie RH. Climatic selection on body colour in the land snail Theba pisana (Pulmonata: Helicidae). Heredity. 1990;65(1):123–6.

    Article  Google Scholar 

  13. Cameron RA, Cook LM. Temporal morph frequency changes in sand‐dune populations of Cepaea nemoralis (L.). Biol J Linnean Soc. 2013;108(2):315–22.

    Article  Google Scholar 

  14. Johnson KG, Brooks SJ, Fenberg PB, Glover AG, James KE, Lister AM, Michel E, Spencer M, Todd JA, Valsami-Jones E, et al. Climate Change and Biosphere Response: Unlocking the Collections Vault. Bioscience. 2011;61(2):147–53.

    Article  Google Scholar 

  15. Boggs CL, Watt WB, Ehrlich PR. Butterflies: ecology and evolution taking flight: University of Chicago Press. 2003.

    Google Scholar 

  16. Buckley LB, Kingsolver JG. The demographic impacts of shifts in climate means and extremes on alpine butterflies. Funct Ecol. 2012;26(4):969–77.

    Article  Google Scholar 

  17. Kingsolver JG. Thermoregulation and Flight in Colias Butterflies - Elevational Patterns and Mechanistic Limitations. Ecology. 1983;64(3):534–45.

    Article  Google Scholar 

  18. Kingsolver JG, Moffat RJ. Thermoregulation and the determinants of heat transfer in Colias butterflies. Oecologia. 1982;53(1):27–33.

    Article  Google Scholar 

  19. Kingsolver JG. Ecological Significance of Flight Activity in Colias Butterflies: Implications for Reproductive Strategy and Population Structure. Ecology. 1983;64(3):546–51.

    Article  Google Scholar 

  20. Watt WB. Adaptive significance of pigment polymorphisms in Colias butterflies. I. Variation of melanin pigment in relation to thermoregulation. Evolution. 1968;22(3):437–58.

    Article  Google Scholar 

  21. Ellers J, Boggs CL. Evolutionary genetics of dorsal wing colour in Colias butterflies. J Evol Biol. 2004;17(4):752–8.

    Article  CAS  PubMed  Google Scholar 

  22. Nielsen M, Watt W. Behavioural fitness component effects of the alba polymorphism of Colias (Lepidoptera, Pieridae): resource and time budget analysis. Funct Ecol. 1998;12(1):149–58.

    Article  Google Scholar 

  23. Watt WB, Wheat CW, Meyer EH, Martin JF. Adaptation at specific loci. VII. Natural selection, dispersal and the diversity of molecular-functional variation patterns among butterfly species complexes (Colias : Lepidoptera, Pieridae). Mol Ecol. 2003;12(5):1265–75.

    Article  CAS  PubMed  Google Scholar 

  24. Watt WB. Adaption, constraint, and neutrality: mechanistic case studies with butterflies and their general implications: Cambridge University Press. 2003.

    Google Scholar 

  25. Kingsolver JG, Watt WB. Thermoregulatory Strategies in Colias Butterflies: Thermal Stress and the Limits to Adaptation in Temporally Varying Environments. Am Nat. 1983;121(1):32–55.

    Article  Google Scholar 

  26. Kingsolver JG, Buckley LB. Climate variability slows evolutionary responses of Colias butterflies to recent climate change. Proc R Soc B Biol Sci. 2015;282(1802):20142470.

    Article  Google Scholar 

  27. Sgrò CM, Terblanche JS, Hoffmann AA. What Can Plasticity Contribute to Insect Responses to Climate Change?. Annual Review of Entomology 2016;61(1):433–51.

  28. Hendry AP. Key Questions on the Role of Phenotypic Plasticity in Eco-Evolutionary Dynamics. J Hered. 2016;107(1):25–41.

    Article  PubMed  Google Scholar 

  29. Hoffman RJ. Environmental uncertainty and evolution of physiological adaptation in Colias butterflies. Am Nat. 1978;112(988):999–1015.

    Article  Google Scholar 

  30. Kingsolver JG, Huey RB. Size, temperature, and fitness: three rules. Evol Ecol Res. 2008;10(2):251–68.

    Google Scholar 

  31. Bowden JJ, Eskildsen A, Hansen RR, Olsen K, Kurle CM, Høye TT. High-Arctic butterflies become smaller with rising temperatures. Biol Lett. 2015;11:10.

    Article  Google Scholar 

  32. Boggs CL, Inouye DW. A single climate driver has direct and indirect effects on insect population dynamics. Ecol Lett. 2012;15(5):502–8.

    Article  PubMed  Google Scholar 

  33. DeChaine EG, Martin AP. Historical biogeography of two alpine butterflies in the Rocky Mountains: broad‐scale concordance and local‐scale discordance. J Biogeogr. 2005;32(11):1943–56.

    Article  Google Scholar 

  34. Stamberger J. Adaptation to temporal scales of heterogeneity in the thermal environment. Dissertation Stanford University. 2006.

    Google Scholar 

  35. McGuire CR, Nufio CR, Bowers MD, Guralnick RP. Elevation-dependent temperature trends in the Rocky Mountain Front Range: changes over a 56-and 20-year record. PLoS One. 2012;7(9):e44370.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Higgins JK. Rapid Evolution and Population Diverfence in Response to Environmental Change in Colias Butterflies. Chapel Hill: Chapel Hill University of North Carolina; 2014.

    Google Scholar 

  37. Dai A, Karl TR, Sun B, Trenberth KE. Recent trends in cloudiness over the United States: A tale of monitoring inadequacies. Bull Am Meteorol Soc. 2006;87(5):597.

    Article  Google Scholar 

  38. Buckley LB, Miller EF, Kingsolver JG. Ectotherm thermal stress and specialization across altitude and latitude. Integr Comp Biol. 2013;53(4):571–81.

    Article  PubMed  Google Scholar 

  39. Wild M, Gilgen H, Roesch A, Ohmura A, Long CN, Dutton EG, Forgan B, Kallis A, Russak V, Tsvetkov A. From dimming to brightening: Decadal changes in solar radiation at Earth’s surface. Science. 2005;308(5723):847–50.

    Article  CAS  PubMed  Google Scholar 

  40. Pepin N, Losleben M, Pepin N, Losleben M. Climate change in the Colorado Rocky Mountains: free air versus surface temperature trends. International Journal of Climatology. 2002;22(3):311–29.

  41. Bartoń K. MuMIn: multi-model inference, R package, version 0122. 2009. Available at:

    Google Scholar 

  42. Kissling WD, Carl G. Spatial autocorrelation and the selection of simultaneous autoregressive models. Glob Ecol Biogeogr. 2008;17(1):59–71.

    Google Scholar 

  43. Bivand R, Anselin L, Berke O, Bernat A, Carvalho M, Chun Y, Dormann C, Dray S, Halbersma R, Lewin-Koh N. spdep: Spatial dependence: weighting schemes, statistics and models, R package version 0.5–31. 2011. URL

    Google Scholar 

  44. Charmantier A, McCleery RH, Cole LR, Perrins C, Kruuk LE, Sheldon BC. Adaptive phenotypic plasticity in response to climate change in a wild bird population. Science. 2008;320(5877):800–3.

    Article  CAS  PubMed  Google Scholar 

  45. Chevin L-M, Lande R, Mace GM. Adaptation, plasticity, and extinction in a changing environment: towards a predictive theory. PLoS Biol. 2010;8(4):e1000357.

    Article  PubMed  PubMed Central  Google Scholar 

  46. Lynch M, Lande R. Evolution and extinction in response to environmental change. In: KareivaJK P, Huey R, editors. Biotic interactions and global change. Sunderland: Sinauer Assocs., Inc; 1993. p. 234–50.

    Google Scholar 

  47. Watt WB. Adaptive Significance of Pigment Polymorphisms in Colias Butterflies. 2. Thermoregulation and Photoperically Controlled Melanin Variation in Colias eurythme. Proc Natl Acad Sci U S A. 1969;63(3):767–74.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Huey RB, Kingsolver JG. Evolution of Resistance to High-Temperature in Ectotherms. Am Nat. 1993;142:S21–46.

    Article  Google Scholar 

  49. Moran NA. The evolutionary maintenance of alternative phenotypes. Am Nat. 1992;139(5):971–89.

    Article  Google Scholar 

  50. Padilla DK, Adolph SC. Plastic inducible morphologies are not always adaptive: the importance of time delays in a stochastic environment. Evol Ecol. 1996;10(1):105–17.

    Article  Google Scholar 

  51. Kharouba HM, Paquette SR, Kerr JT, Vellend M. Predicting the sensitivity of butterfly phenology to temperature over the past century. Glob Chang Biol. 2014;20(2):504–14.

    Article  PubMed  Google Scholar 

  52. Ellers J, Boggs CL. Functional ecological implications of intraspecific differences in wing melanization in Colias butterflies. Biol J Linnean Soc. 2004;82(1):79–87.

    Article  Google Scholar 

  53. Ellers J, Boggs CL. The evolution of wing color: male mate choice opposes adaptive wing color divergence in Colias butterflies. Evolution. 2003;57(5):1100–6.

    Article  PubMed  Google Scholar 

  54. Roulin A. Melanin‐based colour polymorphism responding to climate change. Glob Chang Biol. 2014;20(11):3344–50.

    Article  PubMed  Google Scholar 

  55. Willmer PG. Microclimate and the Environmental Physiology of Insects. In: Advances in Insect Physiology. Edited by M.J. Berridge JET, Wigglesworth VB, vol. Volume 16. London: Academic Press; 1982. p 1–57.

  56. Dennis RL. Butterflies and climate change: Manchester University Press. 1993.

    Google Scholar 

  57. Gilbert SF. Ecological Developmental Biology: Developmental Biology Meets the Real World. Dev Biol. 2001;233(1):1–12.

    Article  CAS  PubMed  Google Scholar 

  58. Stoehr AM, Goux H. Seasonal phenotypic plasticity of wing melanisation in the cabbage white butterfly, Pieris rapae L. (Lepidoptera: Pieridae). Ecol Entomol. 2008;33(1):137–43.

    Article  Google Scholar 

  59. Kemp DJ, Jones RE. Phenotypic plasticity in field populations of the tropical butterfly Hypolimnas bolina (L.) (Nymphalidae). Biol J Linnean Soc. 2001;72(1):33–45.

    Article  Google Scholar 

  60. Green JP, Rose C, Field J. The Role of Climatic Factors in the Expression of an Intrasexual Signal in the Paper Wasp Polistes dominulus. Ethology. 2012;118(8):766–74.

    Article  Google Scholar 

  61. Center WRC. Cooperative Climatological Data Summaries. 2013. Retrieved from http://www.wrccdriedu/climatedata/climsum/.

    Google Scholar 

Download references


We thank A. Warren, P. Opler, V. Scott, L. Gall and R. Robins for access to and assistance with the museum collections. We also thank C. and L. Williams, J. Hensel, J. and A. Richman for assistance during museum visits. We are grateful to T. Hedrick for development of a MatLab GUI for photo analysis. We thank R. Huey, R. Dunn, J. Terblanche and anonymous reviewers for edits on the manuscript.


This research was supported in part by NSF grants DEB-1120062 to LBB and JGK.

Authors’ contributions

Design was conceived by LBB, JGK, and HJM. Data were collected by HJM, and analyzed by LBB and HJM. Manuscript was drafted by HJM, LBB, and JGK. All authors read and approved the final manuscript.

Competing interests

The authors declare that they have no competing interests.

Consent for publication

Not Applicable.

Ethics approval and consent to participate

Not Applicable.

Open Access

This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Author information

Authors and Affiliations


Corresponding author

Correspondence to Heidi J. MacLean.

Additional file

Additional file 1:

Figure S1. a) Temperature trends for Mean July Maximum Temperature (solid line) and for Mean June Mean Temperature (dashed line) for all years across the time period for which we analyze specimens. b) Temperature data and linear regressions for Mean July Maximum Temperature (solid line, open circles) and June Mean Temperature (dashed line, closed circles) only for the years for which we have specimens. Figure S2. The relationship between June mean temperatures and July maximum temperatures for each year from 1953 to 2013 at the weather station in Climax, CO (39.37 N, 106.18 W, [61]). Figure S3. The distribution of collection dates for the measured specimens. Table S1. Results of averaging the top spatially autoregressive maximum likelihood models for wing length, wing melanism, and setae length. The top (by AICc) 20 subsets of the full model are averaged across the 50 iterations of the analysis We present coefficients for predictor variables with corresponding standard errors, Z-scores, and the probability that the Z-score deviates significantly from 0. Importance indicates the percent of top models that include each predictor variable. (DOCX 435 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

MacLean, H.J., Kingsolver, J.G. & Buckley, L.B. Historical changes in thermoregulatory traits of alpine butterflies reveal complex ecological and evolutionary responses to recent climate change. Clim Chang Responses 3, 13 (2016).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: