Research article | Open | Published:
Physiological mechanisms of adaptive developmental plasticity in Rana temporaria island populations
BMC Evolutionary Biologyvolume 17, Article number: 164 (2017)
Adaptive plasticity is essential for many species to cope with environmental heterogeneity. In particular, developmental plasticity allows organisms with complex life cycles to adaptively adjust the timing of ontogenetic switch points. Size at and time to metamorphosis are reliable fitness indicators in organisms with complex cycles. The physiological machinery of developmental plasticity commonly involves the activation of alternative neuroendocrine pathways, causing metabolic alterations. Nevertheless, we have still incomplete knowledge about how these mechanisms evolve under environments that select for differences in adaptive plasticity. In this study, we investigate the physiological mechanisms underlying divergent degrees of developmental plasticity across Rana temporaria island populations inhabiting different types of pools in northern Sweden.
In a laboratory experiment we estimated developmental plasticity of amphibian larvae from six populations coming from three different island habitats: islands with only permanent pools, islands with only ephemeral pools, and islands with a mixture of both types of pools. We exposed larvae of each population to either constant water level or simulated pool drying, and estimated their physiological responses in terms of corticosterone levels, oxidative stress, and telomere length.
We found that populations from islands with only temporary pools had a higher degree of developmental plasticity than those from the other two types of habitats. All populations increased their corticosterone levels to a similar extent when subjected to simulated pool drying, and therefore variation in secretion of this hormone does not explain the observed differences among populations. However, tadpoles from islands with temporary pools showed lower constitutive activities of catalase and glutathione reductase, and also showed overall shorter telomeres.
The observed differences are indicative of physiological costs of increased developmental plasticity, suggesting that the potential for plasticity is constrained by its costs. Thus, high levels of responsiveness in the developmental rate of tadpoles have evolved in islands with pools at high but variable risk of desiccation. Moreover, the physiological alterations observed may have important consequences for both short-term odds of survival and long term effects on lifespan.
Developmental plasticity occurs when an environmental input induces a lasting alteration in an organism’s phenotype. Developmental plasticity is generally a non-reversible process and can modulate acclimation capacity (reversible plasticity) in later developmental-stages . Evolutionary benefits of adaptive developmental plasticity have been broadly assessed by both theoretical studies [2,3,4] and empirical approaches [5, 6], reviewed in . Theoretical studies argue that developmental plasticity can increase population viability and facilitate the maintenance of genetic variation in two main ways: reducing the effect of genetic drift by moderating bottlenecks and shielding genetic variants from selection [8, 9]. This shielding effect may slow down the response to selection , but at the same time increases the odds of population persistence under environmental heterogeneity and preserves genetic variation hence granting a greater adaptive potential [9, 10]. Because plasticity confers a higher capacity for surviving new or rapidly changing environments, it can facilitate divergence among populations [11,12,13] and ultimately foster speciation [3, 8], as both empirical and theoretical studies have shown.
Environmental heterogeneity affects the degree of adaptive developmental plasticity, as predicted by theory [9, 14] and demonstrated empirically [15, 16]. Thus, heterogeneous environments are expected to select for highly plastic genotypes, whereas homogeneous environments would tend to reduce plasticity, especially if there are maintenance costs associated to such plasticity [17,18,19]. Populations diverging in the extent of their adaptive plastic responses as a consequence of selection under different environmental regimes would require qualitative or quantitative divergence in the mechanisms underlying phenotypic expression. The expression of alternative phenotypes relies on changes in gene expression [19, 20] and/or changes in the regulation of physiological pathways . Thus, understanding the physiological mechanisms regulating developmental plasticity is key to understand the origin of evolved differences among populations in contrasting environments [22,23,24].
Developmental plasticity is particularly critical for amphibians because they are typically species with low vagility and high philopatry to highly variable habitats . The majority of amphibians exhibit ancestral aquatic reproduction  and breed in temporary water bodies where the larvae develop until metamorphosis. A major larval-stage risk is pool drying and amphibian larvae can often detect fluctuations in water level and accelerate development in response to decreased water levels to achieve an early metamorphosis [27, 28]. Theory predicts that populations exposed to more fluctuating hydroperiods would exhibit greater developmental plasticity since it would allow tadpoles to adjust the timing of metamorphosis accelerating development when ponds dry up or extending the larval period if ponds persist for longer. However, the evolution of developmental plasticity might be constrained by costs of maintaining and/or producing such responses and also limited by factors such as reduced genetic variation or the impossibility of developing appropriate phenotypic responses at all points throughout development [29,30,31]. The populations of Rana temporaria on the Swedish islands are a paradigmatic example of the relationship between environmental heterogeneity and degree of adaptive plasticity . These populations were established between 70 and 800 years ago  by frogs migrating from the mainland . The islands have rocky pools where the frogs breed. These pools vary greatly and consistently in depth and size among islands, and consequently vary in average pond duration, ranging from permanent to ephemeral . Over the last decade, several studies have analyzed the adaptive divergence in developmental rate among these R. temporaria populations. These studies have found signs of developmental canalization so that populations occupying islands with only ephemeral pools show overall faster developmental rates than populations from islands with permanent pools  while also showing reduced within-population genetic variation for developmental rate . Moreover, island populations show marked differences in plasticity of their developmental rate  according to pond duration during the breeding season. In addition,  found evidence of costs of developmental plasticity, although these were only noticeable for the most plastic populations . Faster developing populations were also found to express higher levels of thyroid hormone receptors (alpha and beta) and genes associated with higher metabolic activity . This is congruent with mechanisms of developmental acceleration found in other species, namely increased thyroid hormone and corticosterone levels .
Developmental acceleration in amphibians is mediated by neuroendocrine pathways, especially by the hypothalamic pituitary adrenal axis (HPA) . The HPA-axis is activated by external environmental inputs such as decreased water height , ultimately resulting in increased corticosterone production, which together with thyroid hormone synergistically activate or repress nuclear target genes causing increased metabolic rate, and accelerated morphogenesis . Prolonged corticosterone secretion enhances metabolism and development [40, 41] but it can be at the cost of overproduction of toxic substances called reactive oxygen species (ROS)  that can inflict considerable cellular damage. Similarly,  found that bird species with faster development and smaller body size produced higher antioxidant enzymes. In our study, we expect that R. temporaria populations exposed to divergent regimes of pool drying will show adaptive differences in developmental rate. These among-population differences in developmental plasticity would be explained by altered gene expression associated with adaptive changes in metabolic activity and rate of morphogenesis, which ultimately would involve divergent corticosterone and antioxidant levels.
High ROS production also results in DNA damage, of which telomere shortening is of great importance because of its association with life-history trade-offs and lifespan [44, 45]. Telomeres are non-coding tandem repeat sequences of the terminal regions of the chromosomes with high G-C strand asymmetry . Telomeres are determinants of cell senescence and are also involved in chromosome stability by avoiding chromosome fusion . Telomere replication occurs via reverse transcriptase telomerase that adds telomeric repeats to 3′ overhang. However, telomere ends shorten over time after each cell division until critical telomere length is reached and initiates apoptosis, or programmed cell death . Enhanced telomere abrasion has been described as a consequence of compensatory growth and developmental acceleration , mainly as a direct consequence of oxidative damage. In the wild, telomere shortening studies are principally focused on ageing [50, 51] and body condition quantification [52, 53]. However, only a few evolutionary studies include telomeres as a mechanism under selection even though telomere shortening correlates with numerous variations in life-history traits , and most of these studies have been conducted on humans [55,56,57].
Here we examine whether adaptive differences in developmental rate among R. temporaria populations are associated with changes in corticosterone levels, oxidative stress, and telomere length. We expected highly plastic populations to increase their corticosterone levels to a greater extent in response to pond drying than less plastic populations. We also hypothesized that populations with high developmental plasticity may pay a cost of maintaining such plastic ability in terms of higher constitutive levels of corticosterone and oxidative stress, and likely shorter telomeres, compared with less plastic populations.
Study area and field sampling
We studied the physiological consequences of adaptive divergence in developmental rate of Rana temporaria tadpoles from six islands located in the Gulf of Bothnia (Umea, Sweden) in a 10 km section of the coastline. The size range of the islands is between 9 and 38 ha. Frogs breed in water filled pools created by rocky depression on these islands. The pools on the islands differ in their water permanence such that some have ephemeral pools, other permanent pools, and others have a mixture of permanent and ephemeral pools. There is no relationship between predator abundance and life history traits among pools on these islands , and consequently pond duration seems to be the main environmental factor determining larval development in this system.
On 5 May 2014 we collected between 2 and 4 clutches from each of six islands. We separated the islands into three type of habitats according to their pool characteristics: two had only permanent pools (Storhaddingen 63° 40′N, 20° 25′ E, four clutches; Lillhaddingen 63° 40′N, 20° 24′E, two clutches), two only ephemeral pools (Sävar-Tärnögen 63° 45′N, 20° 36′E, four clutches; Ålgrundet 63° 41′N, 20° 25′E, two clutches) and two a mixture of permanent and ephemeral pools (Petlandsskär 63° 39′N, 20° 24′E, three clutches; Bredskär 63° 39′N, 20° 18′E, three clutches). Clutches were raised until Gosner stage 25  separately in 1 L plastic containers (9.5 cm × 9.5 cm, height 10 cm) filled with 500 mL of reconstituted soft water (see  for details) in a climatic chamber at 14 °C and a light: dark cycle of 18 h: 6 h simulating natural photoperiod conditions.
One week after egg sampling tadpoles reached the free-feeding stage and started to swim actively (Gosner stage 25). At that time, 12 tadpoles per clutch (six tadpoles per treatment) for a total of 216 experimental units were individually transferred to 1 L containers filled with 750 mL of reconstituted soft water, where each tadpole was haphazardly assigned to an experimental treatment, i.e. constant water level or simulated pool drying. We renewed water every 4 days and tadpoles were fed ad libitum with lightly boiled spinach on each day of water change . Containers were placed on shelves in a random order with respect to treatment and island in a walk-in climate chamber. Temperature and light were set to 22 °C and to 18 h: 6 h of light: dark. The experiment included a water level factor with two levels: constant water and simulating pool drying conditions. In the simulating pool drying treatment we decreased the initial water volume of 750 mL (10 cm) by 33% at each water change starting on day 5 until day 25, keeping the water volume constant at 66 mL (1 cm depth) afterwards. Water temperature did not differ significantly between treatments despite the differences in water level. The experiment ended when tadpoles reached Gosner stage 42 (front legs visible). We monitored tadpoles twice a day (09.00 and 21.00) to check for metamorphs. At this stage tadpoles were photographed to determine their length using ImageJ (version 1.47 t), and were weighed on a high precision balance to the nearest 0.1 mg. Clutch size from each island population varied between 2 and 4 (see Additional file 1).
Tadpoles at 42 Gosner were euthanized by immersion in a buffered solution of MS-222. A portion of muscle from the tadpoles’ tail was removed with a surgical blade and preserved at −20 °C for telomere analyses. The rest of the tail was snap frozen in liquid nitrogen and preserved at −80 °C until corticosterone assays were conducted. The rest of the body was also snap frozen for oxidative stress assays.
Corticosterone content was determined in tails (50–60 mg) collected from tadpoles at 42 Gosner Stage. Tissue was homogenized with an Ultraturrax TP18/10 (Hanke & Kunkel; IKA-Werk) during 30 s and the hormone was extracted following an organic phase extraction with 1 mL ethyl acetate during 30 min at 4 °C and continuous shaking. Samples were then centrifuged at 5000 rpm during 15 min and a known volume of the supernatant was taken and evaporated in a speedVac. Dried elutes were re-suspended in a final volume of 120 μL of the assay buffer provided with the enzymoinmunoassay kit supplemented with ethanol (<5% of final volume) to aid the re-suspension of steroids. We assayed each sample in duplicate (50 μL per sample) for corticosterone determination through specific enzymoinmunoassays (Arbor Assays, K014-H1/H5). The corticosterone antibody has low cross-reactivities to cortisol (0.38%), 11-desoxycorticosterone (12.3%), or progesterone (0.24%). The efficiency of the extraction was checked by spiking several aliquoted samples with 100 pg of exogenous corticosterone prior to extraction and comparing them to the non-spiked aliquotes. Recovery of exogenous corticosterone was never lower than 96%. The lowest point in the corticosterone standard curve was 2.29 pg/mL. To test for assay precision and variability, we determined the coefficient of variation (CV%) for intra- and inter-assay variation. Intra-assay variation was 8.74%. Inter-assay variation was 13.23% in the highest point of the standard curve and and 6.86% in the lowest point (n = 8 assays in both cases). The intra-sample correlation coefficient was 0.97. Corticosterone concentration was calculated from the %B/B0 curve by using the 4PLC fitting routine and following the online tool from .
Oxidative stress assay
We quantified the activity of three antioxidant enzymes (catalase, glutathione peroxidase, and glutathione reductase). We also quantified malondialdehyde (MDA) concentration, a product formed during lipid peroxidation, total glutathione (GSHt) and the ratio of oxidized to reduced glutathione (GSH/GSSG ratio). After evisceration, tadpoles were individually homogenized with a Miccra homogenizer (Miccra D-1) at 35,000 rpm in a buffered solution to inhibit proteolysis (1:4; w:v; ). The homogenates were centrifuged at 14,000 rpm for 30 min at 4 °C and the resulting supernatants were aliquoted into six 0.6 mL tubes and cryopreserved at −80 °C. We determined total protein content by standard Bradford’s procedure . The coefficient of variation in total protein determinations between duplicated samples was on average 4.27%.
Catalase (CAT) catalytic activity was indirectly quantified following . According to this protocol, potassium permanganate (KMnO4) was used as an oxidizing agent that reacts with the catalase substrate hydrogen peroxide (H2O2) giving a red color that can be read at 480 nm 5 min after KMnO4 was added. We used commercial catalase (SIGMA – 60,634) for standard curves preparation. The coefficient of variation between duplicated samples was on average 3.26%. The intra-sample correlation coefficient was 0.95. Glutathione peroxidase (GPX) activity was quantified following the protocol developed by . An excess of glutathione reductase (GR) reduces continually with NADPH the oxidized glutathione (GSSG) producing a constant level of reduced glutathione (GSH). We estimated NADPH oxidation spectrophotometrically at a wavelength of 340 nm. The coefficient of variation between duplicated samples was on average 2.12%. The intra-sample correlation coefficient was 0.88. GR activity was determined assessing the decrease in absorbance at 340 nm due to NADPH oxidation , and assays showed an average coefficient of variation between duplicated samples of 2.74%. The intra-sample correlation coefficient was 0.96. Lipid peroxidation was determined according to  by measuring MDA concentration. MDA is one product of the lipid peroxidation that reacts with thiobarbituric acids, generating a red product that absorbs at 535 nm. We quantified MDA concentration in nmol / mL by measuring optical density of each sample and then subtracting the blank values and comparing with the calibration curve. The coefficient of variation between duplicated samples was on average 3.93%. The intra-sample correlation coefficient was 0.95. Total glutathione (GSHt) was determined following . Homogenates were diluted (1:10, w:v) and homogenized in a stock buffer solution (0.01 M PBS and 0.02 M EDTA). Three working solutions were prepared from the same stock buffer as follows: (a) 0.03 mM of NADPH, (b) 6 mM 5,5′-Dithiobis(2-nitrobenzoic acid) (DNTB), and (c) 50 units of GR/mL. Solution a and b were mixed using a ratio of 7:1 (A:B) and then 160 μL of this mixture was added to 40 μL of supernatant. After 15 s, 20 μL of the solution C were added and absorbance was read at 405 nm after 30 and 60 s. Total concentration of glutathione was determined by comparing changes in absorbance between the two readings, according to a standard curve generated by serial dilutions of glutathione from 1 nM to 0.031 nM. The repeatability was 90.3% (n = 10 samples).
Relative telomere length quantification
Genomic DNA for telomere measurements was isolated using a high-salt DNA extraction protocol on a portion of tail muscle. All determinations were made from muscle tissue to avoid confounding differences between tissues in telomere length [69, 70]. Relative telomere length assays were performed using quantitative PCR (qPCR) following the protocol developed by  and calculated as the ratio of telomere repeats to a single-copy gene (T/S ratio; ). As a single copy gene we used glyceraldehyde-3-phosphate dehydrogenase (GAPDH, GenBank accession no. AF255390). Forward and reverse sequences of primers used to amplify GAPDH were 5-AACCAGCCAAGTACGATGACAT-3′ (GAPDH-F) and 5′-CCATCAGCAGCAGCCTTCA-3′ (GAPDH-R), respectively. Forward and reverse primers of the target gene were 5′CGGTTTGTTTGGGTTTGGGTTTGGGTTTGGGTTTGGGTT-3′ (Tel1b) and 5′-GGCTTGCCTTACCCTTACCCTTACCCTTACCCTTACCCT-3′ (Tel2b), respectively. We performed qPCR in two separate plates for GAPDH and telomere genes by adding 20 ng of genomic DNA from each sample. The set of primers used (Tel1b/Tel2b and GAPDH-F/GAPDH-R) was at an initial concentration of 900 nM containing 10 μL of Brilliant SYBR Green QPCR Master Mix (Roche) in a final volume of 20 μL. PCR protocol consisted of 10 min at 95 °C followed by 30 cycles of 1 min at 56 °C and 1 min at 95 °C for telomere fragment amplification, and 10 min at 95 °C followed by 40 cycles of 1 min at 60 °C and 1 min at 95 °C for GADPH fragment. We conducted qPCRs on a LightCycler 480 (Roche) and we tested the efficiency of each plate by performing a control standard curve by serially diluting a pool of samples from the different treatments and islands in triplicate (160, 40, 10, 2.5 and 0.66 ng of DNA per well). We calculated the cycle threshold (Ct) value of each sample for each plate. Threshold is the basal level of fluorescence. Ct is defined as the number of cycles needed to detect a signal above the threshold. All samples were run in duplicate and relative telomere length was calculated following the formula :
where Etelomere is the qPCR efficiency of telomere fragment; EGAPDH is the qPCR efficiency of the GAPDH fragment; ΔCt telomere is the Ct deviation of control – sample of the telomere (target gene) fragment; ΔCt GAPDH is the Ct deviation of control – sample of reference of GADPH (gene of reference) fragment. The coefficient of variation for duplicated samples was 0.65% on average for GAPDH assays and 1.36% for telomere assays. The intra-sample correlation coefficient was 0.90 and 0.86 for the telomere and GAPDH assays, respectively.
Statistical analyses were conducted in R, version 3.3.1 (R Development Core Team). We checked whether residuals followed normal distributions conducting Kolgomorov-Smirnov tests (“lillie.test” in “nortest” package, version 1.0–2). We also tested for homoscedasticity using Barlett’s tests with the function “bartlett.test” (“car” package; version 2.0–22). We fitted linear and generalized mixed effect models using “lmer” function for normal data and “glmer” for non-normal data (package “lme4”; version 1.1–7). We used water level (constant or simulated pool drying) and island habitat as fixed factors, and included clutch nested within island as random factors in the models. We used likelihood ratio tests to determine the significance of each predictor. Days to metamorphosis, corticosterone, GPX, GR, GSH, and telomere data were log-transformed to meet parametric assumptions. We estimated growth rate as log (body mass at 42 Gosner stage) – log (days to reach 42 Gosner stage). Intraclass correlation coefficients were calculated using “ICCest” function (“ICC” package; version 2.3.0).
Survival was very high throughout the experiment (93.05%) and we found no significant differences between treatments, clutches, or island habitats.
Simulated pool-drying induced accelerated development in R. temporaria larvae (χ2 = 58.05, P < 0.001; Fig. 1a), but the degree of response differed among island habitats as indicated by a significant island habitat by water level interaction (χ2 = 6.73, P = 0.035; Fig. 1a). Tadpoles from islands with only ephemeral pools showed the highest developmental plasticity since they reduced the time to metamorphosis by 8.13% on average, compared to tadpoles from islands with both pool-drying regimes and with only permanent pools, which accelerated their development by 4.40% and by 4.57% on average, respectively. Pool drying decreased larval growth rate by an average of 12.03% (χ2 = 42.39, P < 0.001; Fig. 1b). However, we did not find differences in growth rate among island habitats (χ2 = 0.49, P = 0.784) or in the interaction between island habitat and water level (χ2 = 1.99, P = 0.390; Fig. 1b). Also, pool drying decreased snout-to-vent length of the emerging metamorphs (χ2 = 51.01, P < 0.001), but we did not find differences in length among island types (χ2 = 0.31, P = 0.856) or in the interaction between island type and water level (χ2 = 2.15, P = 0.341).
Simulated pool desiccation induced endocrine responses in R. temporaria tadpoles, indicated by increased corticosterone levels by an average 25.48% (χ2 = 8.71, P = 0.003; Fig. S1 in Additional file 2). However, we did not find significant differences between island habitats or their interaction with water level (all P > 0.110; Fig. S1 in Additional file 2).
We also found alterations in antioxidant responses of the different populations in relation to water level. CAT activity decreased by 9.14% on average in tadpoles responding to pool drying (χ2 = 3.84, P = 0.049; Fig. 2a). Island habitat also affected CAT activity (χ2 = 6.56, P = 0.038; Fig. 2a) so that it was 11.21% and 16.86% lower in tadpoles from islands with ephemeral pools than in tadpoles from islands with both types of pools or with permanent pools, respectively. CAT levels were not affected by the interaction between pool drying and island habitat (χ2 = 4.37, P = 0.112; Fig. 2a). Pool drying seemed to induce a slight increase in GPX levels, but the effect was not significant (χ2 = 2.84 P = 0.092; Fig. S2 in Additional file 2). Neither island habitat nor its interaction with pool drying significantly affected GPX levels (χ2 = 3.25, P = 0.197, and χ2 = 0.01, P = 0.996, respectively; Fig. S2 in Additional file 2). Pool drying did not alter GR levels (χ2 = 2.371, P = 0.124; Fig. 2b). However, we found marginally non-significant differences between island habitats (χ2 = 5.78, P = 0.056; Fig. 2b). Tadpoles from islands with only ephemeral pools and from islands with both types of pools showed lower constitutive GR levels than tadpoles from permanent pools (6.26% and 7.84% lower GR levels on average, respectively). We did not find differences in the interaction between pool drying and island habitat in GR levels (χ2 = 1.18, P = 0.555; Fig. 2b).
We found that alterations in antioxidant enzymatic activity were not necessarily associated with cellular oxidative damage, since MDA values were not affected by pool drying, and showed no differences among island habitats, or their interaction (all P > 0.109; Fig. S3 in Additional file 2). Decreased water level resulted on average in 17.31% lower GSHt values (χ2 = 9.49, P = 0.002). Levels of GSHt was unaltered by either island habitat or the island habitat by water level interaction (χ2 = 1.90 P = 0.386 and χ2 = 4.69, P = 0.096, respectively; Fig. S4 in Additional file 2). The ratio of reduced to oxidized glutathione (GSH/GSSG ratio) was not affected by pool drying (χ2 = 3.16, P = 0.076; Fig. S5 in Additional file 2), and showed no differences among island habitats (χ2 = 2.54, P = 0.281; Fig. S5 in Additional file 2), or their interaction (χ2 = 0.36, P = 0.833; Fig. S5 in Additional file 2).
Relative telomere length was not affected by water level treatment (χ2 = 0.01, P = 0.923; Fig. 3), but it varied significantly among island habitats (χ2 = 6.078, P = 0.048; Fig. 3). Tadpoles from ephemeral pools showed on average 22.31% shorter relative telomere length than tadpoles from permanent pools, and 19.89% shorter relative telomere length than tadpoles from semi-permanent pools. The interaction between relative telomere length and water level was not significant (χ2 = 1.46, P = 0.482; Fig. 3).
Relative telomere length and GSHt, did not correlate with duration of larval period or growth rate (all P > 0.110). We found, however, a slight negative correlation between GSH/GSSG ratio and larval period (R2 = 0.089; P = 0.019).
Reaction norms can evolve under selection in divergent environmental regimes if plasticity confers adaptive value to individuals developing alternative phenotypes, but also if costs of plasticity select for less plastic genotypes in less variable environments [2, 3, 17]. Here we found evidence for divergence among Rana temporaria populations in their ability to accelerate development in response to decreased water level in accordance with the predominant pool drying regimes. The studied populations have been estimated to experience a small degree of neutral genetic differentiation , although Fst estimates could have been biased downwards due to high heterozygosity of the markers used [73, 74] and hence population differentiation may be greater than previously thought. Populations from islands with only ephemeral pools showed the greatest capacity for developmental acceleration in response to pool drying. This indicates that selection on the regulation of developmental rate under different flooding regimes has resulted in increased plasticity in these populations (genetic accommodation), not in canalized fast development or genetic assimilation [2, 75]. We suspect that such a steep reaction norm is derived within the species, but phylogenetic reconstructions of the ancestral state of plasticity for these populations would be required to resolve the polarity of the changes in our Rana temporaria populations. Selection for rapid larval development can result in canalized fast development and loss of plasticity [20, 76]. In our study, populations inhabiting islands with ephemeral pools have evolved greater developmental plasticity. Our results differ somewhat from a previous study on this system , where islands with both types of pools showed higher plasticity. This discrepancy can be a consequence of broader than expected interannual variability in pool duration , among-population variation in costs of plasticity maintenance , and/or stochastic effects of sampling different genotypes over different years.
The mechanisms underlying developmental acceleration in amphibians are well known, and rely on activation of the hypothalamic-pituitary-axis [27, 38]. This neuroendocrine activation in amphibians results in higher thyroid hormone and corticosterone levels, which enhance cell replication rate  and morphogenesis . We observed that larvae from all populations increased corticosterone levels to a similar extent (around 25%) when facing decreased water levels. Such up-regulation of corticosterone explains the ability of these populations to accelerate development, but does not reflect the observed among-population differences in their degree of plasticity. In contrast, the activity of antioxidant enzymes varied significantly among populations, which would suggest that populations under divergent environments evolved metabolic differences resulting in different levels of ROS production. In this case, highly plastic Rana populations showed lower constitutive levels of both CAT and GR enzymes than less plastic populations. CAT transforms hydrogen peroxide to water and oxygen whereas GR reduces glutathione disulfide to the sulfhydryl form of glutathione, both processes being essential in protecting cells from oxidative damage. Low enzymatic levels are associated with low production of ROS, but also with exhaustion of the enzymes as a result of oxidative stress via toxic substances or ROS production [77,78,79]. In addition,  found an increase in the gene transcript of catalase in tadpoles facing pool drying. Levels of catalase transcript were higher in populations that developed faster, supporting the idea of exhausted levels of these enzymes in populations permanently exposed to pool drying conditions. Thus, lower activities of CAT and GR might indicate an enzymatic inactivation caused by an excess of ROS produced during developmental acceleration [45, 80]. However, the absence of among-population differences in MDA or GSH levels do not fully support these conclusions since it indicates a lack of intense lipid peroxidation or antioxidant responses. Therefore, an alternative explanation might be that selection favored individuals that maximized ROS production when they experienced higher metabolic rates, hence showing lower antioxidant activity [81, 82].
Telomere length varied among populations adapted to different pool-drying regimes by evolving different extents of developmental plasticity. This is a novel and intriguing finding and it may help to understand the implications of developmental plasticity to lifespan and fitness. Telomere shortening occurs after each cell replication so that telomeres usually shorten with age [83, 84]. However, telomere shortening is also linked to increased metabolism , particularly to oxidative stress originated during intense cell respiration [54, 85]. In wild populations, telomere shortening has also been described as a predictor of mortality , reproductive costs  or the impact of infections . Therefore, telomere shortening can result a reliable indicator of the biological age of individuals . In our experiment, antioxidant responses observed in Rana populations from islands with only ephemeral pools suggest a high ROS production derived from intense metabolism, which is associated with the attrition of telomeres. Thus, based on our oxidative stress and telomere length determinations, we would expect accelerated development to bear the consequence of reduced lifespan and hence possibly reduced fitness. Among-population differences in telomere length have been associated to trans-generational effects of male age at reproduction due to a progressive elongation of telomeres in sperm with age [89,90,91]. Analogous maternal age effects have not been found . Parental exposure to stressful conditions is also relevant to inheritance of telomere length, although these processes remain poorly understood . In our study, shorter telomeres in tadpoles from islands with only ephemeral pools might be related to early age of first reproduction of males, which could be a long-term consequence of accelerated development against pool drying. However, further empirical studies will elucidate underlying mechanisms in telomere inheritance. Evaluations of telomerase activity and long-term studies testing the effects of parental telomere shortening on life-history traits of offspring will help to understand telomeric dynamic across generations.
Physiological differences among populations found in this study suggest that increased developmental plasticity may be associated with metabolic alterations that compromise the health and lifespan of individuals, as populations differed in antioxidant enzymes activity and length of terminal chromosomal regions. Reduced individual lifespan could have cascading demographic effects on population viability, although this remains to be explored.
In sum, populations evolving in contrasting environments showed divergent levels of developmental plasticity and associated oxidative stress and telomere length variation, despite the slight neutral genetic differentiation previously described. These results emphasize the importance of including physiological measurements in the study of phenotypic plasticity, in order to be able to understand the underlying mechanisms of particular evolutionary and ecological processes.
- GSH/GSSG ratio:
Reduced/oxidized glutathione ratio
- GSHt :
Beaman JE, White CR, Seebacher F. Evolution of plasticity: mechanistic link between development and reversible acclimation. Trends Ecol Evol. 2016;31:237–49.
West-Eberhard MJ. Developmental plasticity and evolution. Oxford: Oxford University Press; 2003.
Pfennig DW, Wund MA, Snell-Rood EC, Cruickshank T, Schlichting CD, Moczek AP. Phenotypic plasticity’s impacts on diversification and speciation. Trends Ecol Evol. 2010;25(8):459–67.
Snell-Rood EC. An overview of the evolutionary causes and consequences of behavioural plasticity. Anim Behav. 2013;85:1004–11.
Aubret F, Bonnet X, Shine R. The role of adaptive plasticity in a major evolutionary transition: early aquatic experience affects locomotor performance of terrestrial snakes. Funct Ecol. 2007;21:1154–61.
Touchon JC, Warkentin KM. (2008). Reproductive mode plasticity: aquatic and terrestrial oviposition in a treefrog. Proc Natl Acad Sci U S A. 2008;105:7495-7499.
Moczek AP, Sultan S, Foster S, Ledón-Rettig C, Dworkin I, Nijhout HF, et al. The role of developmental plasticity in evolutionary innovation. Proc R Soc Lond B. 2011;278:2705–13.
Draghi JA, Whitlock MC. Phenotypic plasticity facilitates mutational variance, genetic variance, and evolvability along the major axis of environmental variation. Evolution. 2012;66:2891–902.
Gomez-Mestre I, Jovani R. A heuristic model on the role of plasticity in adaptive evolution: plasticity increases adaptation, population viability and genetic variation. Proc R Soc Lond B. 2013;280:20131869.
Price TD, Qvarnström A, Irwin DE. The role of phenotypic plasticity in driving genetic evolution. Proc R Soc Lond B. 2003;270:1433–40.
Levin DA. Flowering-time plasticity facilitates niche shifts in adjacent populations. New Phytol. 2009;183:661–6.
Snell-Rood EC, Van Dyken JD, Cruickshank T, Wade MJ, Moczek AP. Toward a population genetic framework of developmental evolution: the costs, limits, and consequences of phenotypic plasticity. BioEssays. 2010;32:71–81.
Svanbäck R, Schluter D. Niche specialization influences adaptive phenotypic plasticity in the threespine stickleback. Am Nat. 2012;180:50–9.
Chevin LM, Lande R. Evolution of discrete phenotypes from continuous norms of reaction. Am Nat. 2013;182:13–27.
Gianoli E, González-Teuber M. Environmental heterogeneity and population differentiation in plasticity to drought in Convolvulus chilensis (Convolvulaceae). Evol Ecol. 2005;19:603–13.
Baythavong BS. Linking the spatial scale of environmental variation and the evolution of phenotypic plasticity: selection favors adaptive plasticity in fine-grained environments. Am Nat. 2011;178:75–87.
Lind MI, Johansson F. The degree of adaptive phenotypic plasticity is correlated with the spatial environmental heterogeneity experienced by island populations of Rana temporaria. J Evol Biol. 2007;20:1288–97.
Beldade P, Mateus ARA, Keller RA. Evolution and molecular mechanisms of adaptive developmental plasticity. Mol Ecol. 2011;20:1347–63.
Schlichting CD, Wund MA. Phenotypic plasticity and epigenetic marking: an assessment of evidence for genetic accommodation. Evolution. 2014;68:656–72.
Johansson F, Veldhoen N, Lind MI, Helbing CC. Phenotypic plasticity in the hepatic transcriptome of the European common frog (Rana temporaria): the interplay between environmental induction and geographical lineage on developmental response. Mol Ecol. 2013;22:5608–23.
Ricklefs RE, Wikelski M. The physiology/life-history nexus. Trends Ecol Evol. 2002;17:462–8.
Valladares F, Balaguer L, Martinez-Ferri E, Perez-Corona E, Manrique E. Plasticity, instability and canalization: is the phenotypic variation in seedlings of sclerophyll oaks consistent with the environmental unpredictability of Mediterranean ecosystems? New Phytol. 2002;156:457–67.
West-Eberhard MJ. Developmental plasticity and the origin of species differences. Proc Natl Acad Sci U S A. 2005;102:6543–9.
Whitehead A, Roach JL, Zhang S, Galvez F. Genomic mechanisms of evolved physiological plasticity in killifish distributed along an environmental salinity gradient. Proc Natl Acad Sci U S A. 2011;108:6193–8.
Smith MA, Green DM. Dispersal and the metapopulation paradigm in amphibian ecology and conservation: are all amphibian populations metapopulations? Ecography. 2005;28:110–28.
Gomez-Mestre I, Pyron RA, Wiens JJ. Phylogenetic analyses reveal unexpected patterns in the evolution of reproductive modes in frogs. Evolution. 2012;66:3687–700.
Denver RJ. Proximate mechanisms of phenotypic plasticity in amphibian metamorphosis. Amer Zool. 1997;37:172–84.
Richter-Boix A, Llorente GA, Montori A. A comparative analysis of the adaptive developmental plasticity hypothesis in six Mediterranean anuran species along a pond permanency gradient. Evol Ecol Res. 2006;8:1139–54.
DeWitt TJ, Sih A, Wilson DS. Costs and limits of phenotypic plasticity. Trends Ecol Evol. 1998;13:77–81.
Auld JR, Agrawal AA, Relyea RA. Re-evaluating the costs and limits of adaptive phenotypic plasticity. Proc R Soc Lond B. 2010;277:503–11.
Murren CJ, Auld JR, Callahan H, Ghalambor CK, Handelsman CA, Heskel MA, et al. Constraints on the evolution of phenotypic plasticity: limits and costs of phenotype and plasticity. Heredity. 2015;115:293–301.
Johansson F, Hjelm J, Giles BE. Life history and morphology of Rana temporaria in response to pool permanence. Evol Ecol Res. 2005;7:1025–38.
Lind MI, Ingvarsson PK, Johansson H, Hall D, Johansson F. Gene flow and selection on phenotypic plasticity in an island system of Rana temporaria. Evolution. 2011;65:684–97.
Lind MI, Persbo F, Johansson F. Pool desiccation and developmental thresholds in the common frog, Rana temporaria. Proc R Soc Lond B. 2008;275:1073–80.
Lind MI, Johansson F. Testing the role of phenotypic plasticity for local adaptation: growth and development in time-constrained Rana temporaria populations. J Evol Biol. 2011;24:2696–704.
Lind MI, Johansson F. Costs and limits of phenotypic plasticity in island populations of the common frog Rana temporaria under divergent selection pressures. Evolution. 2009;63:1508–18.
Gomez-Mestre I, Kulkarni S, Buchholz DR. Mechanisms and consequences of developmental acceleration in tadpoles responding to pond drying. PLoS One. 2013;8:e84266.
Denver RJ. Stress hormones mediate environment-genotype interactions during amphibian development. Gen Comp Endocrinol. 2009;164:20–31.
Boorse GC, Denver RJ. Endocrine mechanisms underlying plasticity in metamorphic timing in spadefoot toads. Integr Comp Biol. 2003;43:646–57.
Hayes TB, Wu TH. Interdependence of corticosterone and thyroid hormones in toad larvae (Bufo boreas). II. Regulation of corticosterone and thyroid hormones. J. Exp. Zool. A. Ecol. Genet. Physiol. 1995;271:103–11.
Wack CL, DuRant SE, Hopkins WA, Lovern MB, Feldhoff RC, Woodley SK. Elevated plasma corticosterone increases metabolic rate in a terrestrial salamander. Comp Biochem Physiol A Mol Integr Physiol. 2012;161:153–8.
Costantini D, Marasco V, Møller AP. A meta-analysis of glucocorticoids as modulators of oxidative stress in vertebrates. J Comp Physiol B. 2011;181:447–56.
Cohen AA, McGraw KJ, Wiersma P, Williams JB, Robinson WD, Robinson TR, et al. Interspecific associations between circulating antioxidant levels and life-history variation in birds. Amer Nat. 2008;172:178–93.
Von Zglinicki T. Oxidative stress shortens telomeres. Trends Biochem Sci. 2002;27:339–44.
Monaghan P, Metcalfe NB, Torres R. Oxidative stress as a mediator of life history trade-offs: mechanisms, measurements and interpretation. Ecol Lett. 2009;12(1):75–92.
Blackburn EH. Structure and function of telomeres. Nature. 1991;350:569–73.
O’Sullivan RJ, Karlseder J. Telomeres: protecting chromosomes against genome instability. Nat Rev Mol Cell Biol. 2010;11:171–81.
Campisi J. Cancer and ageing: rival demons? Nat Rev Cancer. 2003;3:339–49.
Metcalfe NB, Monaghan P. Compensation for a bad start: grow now, pay later? Trend Ecol Evol. 2011;16:254–60.
Haussmann MF, Vleck CM. Telomere length provides a new technique for aging animals. Oecologia. 2002;130:325–8.
Bize P, Criscuolo F, Metcalfe NB, Nasir L, Monaghan P. Telomere dynamics rather than age predict life expectancy in the wild. Proc R Soc Lond B. 2009;276:1679–83.
Trusina A. Stress induced telomere shortening: longer life with less mutations? BMC Syst Biol. 2014;8:27.
Boonekamp JJ, Mulder GA, Salomons HM, Dijkstra C, Verhulst S. 2014. Nestling telomere shortening, but not telomere length, reflects developmental stress and predicts survival in wild birds. Proc R Soc Lond B. 2014;281:20133287.
Haussmann MF, Marchetto NM. Telomeres: linking stress and survival, ecology and evolution. Curr Zool. 2010;56:714–27.
Seluanov A, Chen Z, Hine C, Sasahara TH, Ribeiro AA, Catania KC, et al. Telomerase activity coevolves with body mass not lifespan. Aging Cell. 2007;6:45–52.
Gorbunova V, Seluanov A. Coevolution of telomerase activity and body mass in mammals: from mice to beavers. Mech Ageing Dev. 2009;130:3–9.
Eisenberg DT. An evolutionary review of human telomere biology: the thrifty telomere hypothesis and notes on potential adaptive paternal effects. Am J Hum Biol. 2011;23:49–167.
Gosner KL. A simplified table for staging anuran embryos and larvae with notes on identification. Herpetologica. 1960;16:183–90.
Räsänen KR, Laurila A, Merilä J. Geographic variation in acid stress tolerance of the moor frog, Rana arvalis. I. Local adaptation. Evolution. 2003;57:352–62.
Richter-Boix A, Orizaola G, Laurila A. Transgenerational phenotypic plasticity links breeding phenology with offspring life-history. Ecology. 2014;95:2715–22.
http://www.myassays.com/arbor-assays-corticosterone-enzyme-immunoassay-kit.assay (Accessed 19 June 2017).
Burraco P, Gomez-Mestre I. Physiological stress responses in amphibian larvae to multiple stressors reveal marked anthropogenic effects even below lethal levels. Physiol Biochem Zool. 2016;89:462–72.
Bradford MM. A rapid and sensitive method for the quantitation of microgram quantities of protein utilizing the principle of protein-dye binding. Anal Biochem. 1976;72:248–54.
Cohen G, Somerson NL. Catalase-aminotriazole method for measuring secretion of hydrogen peroxide by microorganisms. J Bacteriol. 1969;98:543–6.
Paglia DE, Valentine WN. Studies on the quantitative and qualitative characterization of erythrocyte glutathione peroxidase. J Lab Clin Med. 1967;70:158–69.
Cribb AE, Leeder JS, Spielberg SP. Use of a microplate reader in an assay of glutathione reductase using 5, 5′-dithiobis (2-nitrobenzoic acid). Anal Biochem. 1989;183:195–6.
Buege JA, Aust SD. Microsomal lipid peroxidation. Method Enzymol. 1978;52:302–10.
Galván I, Gangoso L, Grande JM, Negro JJ, Rodríguez A, Figuerola J, et al. Antioxidant machinery differs between melanic and light nestlings of two polymorphic raptors. PLoS One. 2010;5:e13369.
Friedrich U, Griese EU, Schwab M, Fritz P, Thon KP, Klotz U. Telomere length in different tissues of elderly patients. Mech Ageing Dev. 2000;119:89–99.
Dlouha D, Maluskova J, Lesna IK, Lanska V, Hubacek JA. Comparison of the relative telomere length measured in leukocytes and eleven different human tissues. Physiol Res. 2014;63:S343.
Cawthon RM. Telomere measurement by quantitative PCR. Nucleic Acids Res. 2002;30:e47.
Pfaffl MW. A new mathematical model for relative quantification in real-time RT–PCR. Nucleic Acids Res. 2001;29:e45.
Jost L. GST and its relatives do not measure differentiation. Mol Ecol. 2008;17:4015–26.
Edelaar P, Burraco P, Gomez-Mestre I. Comparisons between QST and FST—how wrong have we been? Mol Ecol. 2011;20:4830–9.
Crispo E. The Baldwin effect and genetic assimilation: revisiting two mechanisms of evolutionary change mediated by phenotypic plasticity. Evolution. 2007;61:2469–79.
Gomez-Mestre I, Buchholz DR. Developmental plasticity mirrors differences among taxa in spadefoot toads linking plasticity and diversity. Proc Natl Acad Sci U S A. 2006;103:19021–6.
Barata C, Navarro JC, Varo I, Riva MC, Arun S, Porte C. Changes in antioxidant enzyme activities, fatty acid composition and lipid peroxidation in Daphnia magna during the aging process. Comp Biochem Physiol B Biochem Mol Biol. 2005;140:81–90.
Srinivasan R, Chandrasekar MJN, Nanjan MJ, Suresh B. Antioxidant activity of Caesalpinia digyna root. J Ethnopharmacol. 2007;113:284–91.
Slos S, Stoks R. Predation risk induces stress proteins and reduces antioxidant defense. Funct Ecol. 2008;22:637–42.
D’Autréaux B, Toledano MB. ROS as signalling molecules: mechanisms that generate specificity in ROS homeostasis. Nat Rev Mol Cell Biol. 2007;8:813–24.
Salin K, Auer SK, Rudolf AM, Anderson GJ, Cairns AG, Mullen W, et al. Individuals with higher metabolic rates have lower levels of reactive oxygen species in vivo. Biol Lett. 2015;11:20150538.
Salin K, Auer SK, Rey B, Selman C, Metcalfe NB. Variation in the link between oxygen consumption and ATP production, and its relevance for animal performance. Proc R Soc B. 2015;282:20151028.
Frenck RW, Blackburn EH, Shannon KM. The rate of telomere sequence loss in human leukocytes varies with age. Proc Natl Acad Sci U S A. 1998;95:5607–10.
Asghar M, Bensch S, Tarka M, Hansson B, Hasselquist D. Maternal and genetic factors determine early life telomere length. Proc R Soc Lond B. 2015;282:20142263.
Monaghan P, Charmantier A, Nussey DH, Ricklefs RE. The evolutionary ecology of senescence. Funct Ecol. 2008;22:371–8.
Barrett EL, Burke TA, Hammers M, Komdeur J, Richardson DS. Telomere length and dynamics predict mortality in a wild longitudinal study. Mol Ecol. 2013;22:249–59.
Bauch C, Becker PH, Verhulst S. Telomere length reflects phenotypic quality and costs of reproduction in a long-lived seabird. Proc R Soc Lond B. 2013;280:20122540.
Asghar M, Palinauskas V, Zaghdoudi-Allan N, Valkiūnas G, Mukhin A, Platonova E, et al. Parallel telomere shortening in multiple body tissues owing to malaria infection. Proc R Soc Lond B. 2016;283:20161184.
Baird DM, Britt-Compton B, Rowson J, Amso NN, Gregory L, Kipling D. Telomere instability in the male germline. Hum Mol Gen. 2006;15:45–51.
Kimura M, Cherkas LF, Kato BS, Demissie S, Hjelmborg JB, Brimacombe M, et al. Offspring’s leukocyte telomere length, paternal age, and telomere elongation in sperm. PLoS Genet. 2008;4:e37.
Eisenberg DT, Hayes MG, Kuzawa CW. Delayed paternal age of reproduction in humans is associated with longer telomeres across two generations of descendants. Proc Natl Acad Sci U S A. 2012;109:10251–6.
Haussmann MF, Heidinger BJ. Telomere dynamics may link stress exposure and ageing across generations. Biol Lett. 2015;11:20150396.
We thank G. Orizaola, B. Dijk, and A. Ritcher-Boix for assistance in animal husbandry, M. Taves (University of British Columbia) for his assistance in corticosterone analysis, and F. Miranda and O. García (Ecophysiology laboratory of Doñana Biological Station) for assistance in the oxidative stress and telomere analysis.
This study was funded by Ministerio de Economía y Competitividad (grant CGL2012–40044). P.B. is supported by fellowship F.P.U.-AP2010–5373 and by the travel scholarship Est13/00128 from Ministerio de Educación. We acknowledge support of the publication fee by the CSIC Open Access Publication Support Initiative through its Unit of Information Resources for Research (URICI).
Availability of data materials
All data analysed during this study are included in the Additional file 1.
This study was approved by the Animal Bioethics Committee from Doñana Biological Station (#12-53_Gomez) and all procedures performed involving the animals were in accordance with the ethical standards of the department and were approved by the Swedish Board of Agriculture (C21/14).
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.