Does selection for behavioral and physiological performance traits alter glucocorticoid responsiveness in bank voles?

In nature, animals are frequently exposed to challenges such as encounters with predators or competitors, harsh weather conditions or food deficiency. Successfully coping with such challenges involves a spectrum of psychological and physiological processes activated or modulated as a part of a stress response to the situation (Wingfield and Ramenofsky, 1999). An adequate stress response can help an animal to overcome the challenge, and, consequently, determine its chances of survival and reproduction. Therefore, evolving a better ability to cope with a particular challenge may involve not only physiological or biomechanical performance traits but also modification of the stress response mechanisms. Here, we tested this hypothesis within the framework of a unique experimental evolution model system: lines of a common rodent, the bank vole (Myodes glareolus), from a multidirectional selection experiment (Chrząścik et al., 2014; Sadowska et al., 2015, 2008).

In vertebrates, the hypothalamic–pituitary–adrenal (HPA) axis is one of the main mediators of the stress response and is an important regulator of homeostasis (Koolhaas et al., 2011; McEwen and Wingfield, 2003). By releasing glucocorticoids, the HPA axis affects a wide spectrum of an organism's functions, such as metabolism, energy mobilization, immune system, behavior and gene expression (Bonier et al., 2011; Coppens et al., 2010; Phuc Le et al., 2005; Schmid et al., 2013; Wingfield and Ramenofsky, 1999), and hence helps the organism to cope with the stressor (Bonier et al., 2009; Koolhaas et al., 1997; Patterson et al., 2014; Sapolsky et al., 2000; Wingfield et al., 1998). Under metabolically demanding conditions, glucocorticoids stimulate energy substrate mobilization (Jimeno et al., 2017; McEwen and Wingfield, 2003). This is reflected, for instance, by increased glucocorticoid levels observed during breeding (Kenagy and Place, 2000; Romero, 2002), at low ambient temperatures (de Bruijn and Romero, 2018) or following intense physical activity (Coleman et al., 1998; Duclos and Tabarin, 2016; Hill et al., 2008). However, a strong glucocorticoid response to metabolically demanding situations may hinder rather than promote an animal's performance (Breuner et al., 1998; Lipowska et al., 2019; Munck and Náray-Fejes-Tóth, 1992; Wingfield and Ramenofsky, 1999). Moreover, prolonged elevation of glucocorticoids can compromise an organism's health (Cohen et al., 2007; Sapolsky et al., 2000). Thus, the glucocorticoid response should be finely adjusted to the duration and intensity of the challenge.

Under natural conditions, a multitude of factors affects HPA axis activity, making it nearly impossible to single out changes occurring in response to a particular selection factor. However, experimental evolution offers a powerful tool to study responses to selection for well-defined traits controlling for random processes (Henderson, 1997; Swallow et al., 2009). The activity of the HPA axis is known to be heritable (Almasi et al., 2010; Béziers et al., 2019; Odeh et al., 2003a,b). Divergent selection in Japanese quail produced lines of quail with high and low corticosterone responses to brief restraint (Cockrem et al., 2010; Satterlee and Johnson, 1988). Quails from the two selection lines also differed in fearfulness (Jones et al., 1992), which suggests an involvement of glucocorticoids in mediating fear behavior. A similar selection experiment in zebra finches generated lines with a strong corticosterone response, but selection for a weak response was not successful (Evans et al., 2006; Hodgson et al., 2007). Thus, the consequences of selection pressure on glucocorticoid levels are mixed, even if selection acts directly on this trait.

Experiments in which selection was applied to traits indirectly associated with the HPA axis also gave mixed results. Changes in the glucocorticoid response to stress were observed in great tits selected for personality traits (Baugh et al., 2012; Carere et al., 2003) and in quails selected for long or short duration of tonic immobility (Hazard et al., 2005). In rats selected for high or low running performance, the differences in glucocorticoid response were observed only when two stressors, restraint and placement in an elevated plus maze, were applied in a sequence (Waters et al., 2010). Rats from lines selected for aggressive reaction toward humans had a higher baseline corticosterone level than those from lines selected for tame reaction (Albert et al., 2008). In mice, selection for high spontaneous physical activity resulted in increased baseline corticosterone, but not in a larger corticosterone response to restraint stress (Malisch et al., 2007, 2008). Importantly, as the baseline was also elevated in the inactive phase of the circadian cycle, the increase could not be explained by immediate effects of physical activity. The elevated baseline could contribute to impaired growth and immune function, and an enhanced predisposition for depression-like behavior in the selected mice (Malisch et al., 2009a,b). The elevated baseline, not accompanied by an increased response level, could indicate a reduced scope of the corticosterone response to stress. However, because baseline and stress-induced corticosterone levels interact primarily with different receptor types (Landys et al., 2006; Reul and de Kloet, 1985), the increased baseline could be an independent modification. Moreover, to check whether the reaction scope has indeed changed, the level achieved in response to an environmental stressor should be compared with the maximum potential level, which can be measured after pharmacological stimulation with adrenocorticotropic hormone (ACTH). Such comparisons are common in studies on stress (Dickens et al., 2009; Harris et al., 2012; Jimeno et al., 2018; Romero and Wikelski, 2010), but uncommon within the framework of selection experiments. In the single study we are aware of, selection for long or short duration of tonic immobility resulted in changes of the glucocorticoid response to both restraint stress and ACTH stimulation (Hazard et al., 2005). Another important aspect, not yet tackled by any selection experiment, is the rate of recovery to the baseline level after a stressful episode. This recovery is crucial for avoiding potentially detrimental effects of prolonged elevation of glucocorticoids (Cohen et al., 2007; Sapolsky et al., 2000). To summarize, the power of selection experiments has not yet been fully employed in comprehensive studies aimed at the question whether and how selection for ecologically relevant traits affects main characteristics of the glucocorticoid stress response.

Here, we asked how HPA axis function evolved in lines of bank voles selected for three distinct traits, which have played important roles in adaptive radiation of terrestrial vertebrates: high aerobic exercise metabolism, predatory behavior and herbivorous capability (Chrząścik et al., 2014; Sadowska et al., 2015, 2008). Increased aerobic capacity has been crucial not only in the evolution of persistent locomotion but presumably also in the evolution of endothermy in birds and mammals. Predatory behavior and herbivorous capability represent opposite extremes on the axis of diet strategies, where evolutionary choice determines many aspects of the organism's anatomy, physiology and behavior. Therefore, the selection experiment provides a suitable basis to ask the broad question whether changes in stress-coping mechanisms are essential in the evolution of the major adaptive strategies.

After about 20 generations of selection (Fig. S1), voles from Aerobic (A) lines, selected for high maximum rate of oxygen consumption in an 18 min swimming trial (V̇_O2,swim), evolved a 60% higher V̇_O2,swim than that achieved by voles from unselected Control (C) lines. Importantly, voles do not have to swim vigorously to stay on the water surface. Thus, the selection favored both the eagerness to work and a high metabolic performance per se, both of which increased in the aerobic lines (Jaromin et al., 2016, 2018, 2019; Lipowska et al., 2019); the Aerobic voles spent more time actively swimming and also swam more regularly, and in only about 5% cases, compared with about 15% in control lines, the swimming trial had to be interrupted to prevent drowning. The Aerobic voles were also more active in their home cages and in open field tests (Maiti et al., 2019), and have increased basal metabolic rate (Sadowska et al., 2015), daily food consumption and thermogenic capacity (Dheyongera et al., 2016). In Predatory (P) lines, selected for hunting crickets (a protocol based on Gammie et al., 2003: time to catch crickets in 10 min tests performed after a few hours of fasting), over 85% of the voles captured the crickets, but only about 15% of Control voles showed predatory behavior, and the selection also improved the time to capture the cricket (Fig. S1). We do not know to what extent the difference is due to a generally more proactive coping style (personality) of the Predatory voles (Maiti et al., 2019) or to more specific neurobiological mechanisms underlining hunting behavior, such as those associated with perception or coordinated movement (Haller, 2018; Hoy et al., 2019; Levenets et al., 2019; Shang et al., 2019). Although analysis of the transcriptome suggests that the selection could have affected sensitivity to hunger, behavioral observations showed that this cannot be the major factor underlying the increased propensity to catch crickets by the Predatory voles (Konczal et al., 2016). Predatory aggression, assessed with a similar protocol, was increased by chronic stress in rats (Pittet et al., 2017), but not by cold-induced increase of energy demand in California mice (Andrew et al., 2019). Thus, sensitivity to stress may even have a stronger effect on this behavior than an increased foraging drive. Voles from the Herbivorous lines, selected for the ability to maintain body mass during a 4 day test with low-quality diet, lost approximately 2 g less mass (approximately 10% of the original body mass) than the Control ones (Fig. S1). The difference was associated with an increased consumption of the low-quality diet by the Herbivorous voles (M.M.L., E.T.S., B.B.-S. and P.K., unpublished), but a decreased locomotor activity (Maiti et al., 2019) indicates that a more economical energy budget could contribute to the improved mass balance, too. Molecular analyses also revealed changes in gene expression and single nucleotide polymorphism (SNP) allele frequency in Aerobic and Predatory lines (Konczal et al., 2015, 2016), and cecal microbiome composition of Herbivorous lines (Kohl et al., 2016). To summarize, several traits potentially related to HPA axis function change as correlated response to the selection.

In each of the three selection directions (Aerobic, Predatory, Herbivorous), the voles are challenged with a stressful situation, but the duration and character of the stressors differ greatly among the line types. Selection in the Aerobic lines favors proactive coping with a brief challenge: being suddenly placed in a jar with water, i.e. in an alien environment requiring specific locomotion in a confined space. In the Predatory lines, successful hunting requires maintaining alertness despite a moderate-length fast (a few hours) in an empty cage preceding the test, followed by the introduction of a novel, mobile object (a cricket), which, as we noticed, can trigger stress-related and evasion behavior in some voles. In the Herbivorous line, the challenge is mild (reduced food quality) but lasts for 4 days. Individuals that perceive the worsening of food quality as distressing and do not eat the food at the onset of the trial may have no chance to restore body mass later. Thus, selection in the Herbivorous lines is likely to promote a calm reaction to unexpected challenges.

We hypothesized that the evolution of increased performance in the respective selection trials resulted in selection-specific modulation of the HPA axis reactivity to stressors. Such a modulation should be reflected in altered levels of corticosterone, the primary glucocorticoid hormone of small rodents. We tested this general hypothesis by comparing baseline corticosterone levels across the lines, and measuring the corticosterone level (1) after exposure to a standardized restraint stressor, commonly applied in stress-related studies on a wide range of species (Baugh et al., 2013; Buynitsky and Mostofsky, 2009; Harris et al., 2012; Malisch et al., 2007; Romero and Wikelski, 2010; Torres-Medina et al., 2018), (2) after pharmacologically stimulating a maximum corticosterone response, which allows estimation of the relative scope of the response to the stressor, and (3) after a short (21 min) and a long (60 min) recovery period after the two treatments. Considering the complexity of the HPA axis and considerable discussion in the field of stress physiology concerning the actual function of corticosterone (MacDougall-Shackleton et al., 2019), specific hypotheses about the expected changes in these traits in the three selection directions would be highly speculative. Therefore, with respect to such specific hypotheses, our study has an exploratory rather than confirmatory nature.

This work was performed on laboratory-bred bank voles, Myodes glareolus (Schreber 1780), from generation 22 and 23 of a multidirectional artificial selection experiment (Sadowska et al., 2015, 2008). The selection was applied based on the following criteria: Aerobic (A) lines – the maximum 1 min rate of oxygen consumption achieved during 18 min of swimming at 38°C (V̇_O2,swim); Predatory (P) lines – ranked time to catch a live cricket in 10 min trials by a fasted vole (fasting duration was gradually reduced in subsequent generations from 14 to 4 h as the selection progressed); Herbivorous (H) lines – body mass change in a 4 day trial, during which the voles were fed a low-quality diet, ‘diluted’ with dried, powdered grass. The tests in Aerobic and Predatory lines were performed on adults (about 75–95 days old), and the test in Herbivorous lines was performed on young, still growing animals (32–36 days). All the trait values used as selection criteria were mass-adjusted (residuals from ANCOVA models, also including other covariates and cofactors). Four replicate lines were maintained for each selection direction and for unselected Control (C), with 15–20 reproducing families in each of the 16 lines (to avoid excess inbreeding). The selection was applied mostly within families, but if more than 16 families were available, families in which all individuals scored below the line mean were excluded. Details of the colony origin, the selection protocol, animal husbandry and results of selection in successive generations are presented in our earlier work (Chrząścik et al., 2014; Sadowska et al., 2015, 2008; Lipowska et al., 2019; Maiti et al., 2019).

All the breeding, selection and experimental procedures were approved by the 1st Local Ethical Committee in Krakow, Poland (decision no. 170/2014).

The experiment was performed on 192 voles from generation 22 and 199 voles from generation 23 (391 animals in total). At the age of about 30 days, approximately 12 males and 12 females were sampled from each of 16 replicate lines (up to two animals per family) and placed individually into individually ventilated cages (AERO Mouse IVC Green Line, Techniplast, Buguggiate, Italy) fitted with sawdust bedding. Unlike rats or mice, bank voles are solitary in nature (Bujalska, 1990), and therefore the social isolation was unlikely to elicit stress-related disorders. The animals were maintained at constant temperature (20±1°C) and photoperiod (16 h:8 h light:dark, light phase starting at 02:00 h), with ad libitum access to water and food (a standard rodent food: 24% protein, 3% fat, 4% fiber; Labofeed H, Kcynia, Poland). The animals used in this study were not subjected to selection tests.

The animals were divided into eight approximately balanced blocks consisting of 48–52 animals of similar age (4 blocks per generation). The animals were also divided into two balanced groups, used in two parallel tests: restraint test (193 animals) and ACTH stimulation test (198 animals; Table S1). The pool of animals used in the tests included 99 pairs of siblings. To exclude the factor of relatedness from within-test comparisons, the siblings were split to different test groups. Each block consisted of an approximately equal number of animals from the two test groups. The tests on one block were performed on two consecutive days. In each of the 2 days, one type of test was performed, but the order of tests was randomized among blocks.

All preparatory and experimental procedures on animals were performed by one person with the occasional help of technical staff. To reduce non-specific stress resulting from human activity in the vicinity of a cage, as well as to habituate animals to handling procedures, voles were handed once per day for 14 days preceding the test. The animals were removed from their cages between 08:00 h and 18:00 h and briefly immobilized by hand. The entire procedure, from taking the cage from the shelf, through capturing and handling the animal to returning it to the cage lasted about 40 s. The exact timing, degree and method of hand immobilization (animal held in a hand or on a flat surface) was diversified among the handling days to habituate the animals to a diverse daily human encounter, rather than to an event of predictable timing and scenario. During the 3 days preceding the test, the handling also included measurement of body mass. Collection of these measurements before the day of testing ensured that the results of the tests were not biased by response to the mild immobilization required for the weighing procedure (animal placed in an opaque plastic cup, approximately 18 cm high and 8 cm in diameter). During handling, eight animals were recognized as diabetic and excluded from the experiment (diabetes appears in a small percentage of bank voles both in laboratory conditions and wild populations; Bartelik et al., 2013).

At the day of testing, the animals were 74–85 days old (mean±s.d. 80.6±2.8 days). Each experimental block was tested on a separate day, between 08:00 h and 12:30 h, and during that time each individual went through a corticosterone release-stimulating procedure – either the stress of restraint or a maximal stimulation with ACTH – and a period of recovery after the event. During the test, each individual had three blood samples taken.

In the restraint test, the animals were immobilized in custom-made restrainers: transparent tubes constructed from Plexiglas structural elements and adjustable PVC walls (length: 8.0 cm, diameter: 2.5–3.5 cm). The general design of the restrainer was similar to that used in rats or mice (Malisch et al., 2007; Vahl et al., 2005), but introduction of adjustable walls allowed us to restrain the animal's movements by means other than through fixing the tail outside the restraint chamber – a method not acceptable in voles, which have slim and delicate tails. The capped end of the tube was perforated to provide ventilation, and the open end could be blocked with a fitted Plexiglas plug. Tube width was adjusted by wrapping the elastic PVC walls of the tube around the animal, not tighter than to the point that prevented the animal from turning around in the tube. The restrainers were washed, dried with paper towels and wiped with ethanol before introduction of an animal. Pilot observations indicated that in female bank voles, corticosterone level plateaus after 10 min of restraint, whereas in males, the level continuously grows for 30 min or more. However, because in the pilot trials a few animals died when restrained for 20 min, a 15 min restraint was applied, to elicit a response that is as close to peak as possible without threatening the animals' safety. The restraint test was performed on up to six animals at the same time, but subsequent individuals were spaced by at least 3 min, which granted sufficient time to perform all steps of the procedure.

In the ACTH stimulation test, the animals were given an intra-peritoneal injection of 1.00 mg ACTH kg⁻¹ body mass (Synacthen Depot 1 mg ml⁻¹, Sigma-Tau, Rome, Italy; diluted 1:3 with sterile saline), and returned to their home cages. Pilot observations suggested that administration of 0.50 mg kg⁻¹ ACTH (similar to that applied by Du et al., 2015, or Touma et al., 2004) may not be sufficient to elicit a stable peak of maximal corticosterone response. As no adverse effects of high ACTH dose were observed, we decided to double the dose to 1.00 mg kg⁻¹ (as in Kalliokoski et al., 2015). The pilot observations also indicated that in females the corticosterone level drops between the 20th and 30th minute post-injection, which indicated that the peak response occurs near the 20th minute. To allow us to perform the test on up to eight animals at the same time, with a 3 min delay between consecutive individuals necessary to perform all steps of the procedure, the period of stimulation was set to 21 min.

Both of the tests started with taking the first (baseline) blood sample from an undisturbed animal, immediately followed by either the restraint or ACTH stimulation procedure. The second (response) blood sample was taken immediately upon completion of the procedure, after which the animal was returned to its home cage and allowed to recover. Each experimental group was further divided into two balanced sub-groups, in which the third (recovery) blood sample for each individual was taken after a short or long recovery time to measure the corticosterone decrease during the recovery (Table S1). Pilot observations based on recovery rates of rats and mice (Fediuc et al., 2006; Hare et al., 2014; Vahl et al., 2005) indicated that bank voles recover faster than the standard laboratory rodents: 30 min was sufficient for corticosterone to recover to near the baseline level, whereas the corticosterone level measured after 60 min of recovery was in some cases lower than before the stressor was applied. Therefore, we decided to take the ‘recovery’ blood samples after 21 min (short) or 60 min (long) recovery periods.

Blood was taken from the retro-orbital venous sinus using 70 µl heparinized capillary tubes (Medlab Products, Raszyn, Poland). Among a number of blood sampling procedures applicable to small rodents (Joslin, 2009; Kim et al., 2018), only retro-orbital sampling allows sufficient sample volume to be obtained from bank voles. No anesthesia was applied prior to blood sampling, as it could impair the ability to obtain samples of sufficient volume through a reduction of blood flow. Moreover, application of anesthesia would prolong the blood sampling procedure, which would compromise measurement of corticosterone level (Kim et al., 2018). All procedures were performed in a separate room, and except for the restraint procedure, animals were returned to the housing room immediately after each sampling. In all cases, the duration of blood sampling (measured from the onset of disturbance, e.g. touching the animal's cage) did not exceed 3 min (range 32–166 s; Table S2). Thus, the procedure of blood sampling was presumed to have a negligible effect on the measured corticosterone level (Romero and Reed, 2005), which was later verified in statistical analyses (see below).

Eleven animals were excluded from further analysis because of methodological errors. In seven animals, the timing of recovery blood sampling after the restraint test was incorrect, and four animals did not receive a full dose of ACTH because they moved during the injection and the needle was pulled out.

The capillary tubes containing blood samples were stored on ice and centrifuged for 15 min at 14,000 g. Plasma samples were separated into Eppendorf tubes, frozen at −20°C, and transferred to the Leibniz Institute for Zoo and Wildlife Research (Berlin), where corticosterone concentration (ng ml⁻¹) was measured according to the method described by Dehnhard et al. (2003), which allows detection with a limit of 3.75 ng ml⁻¹. The samples were extracted twice with 2 ml tert-butyl methyl ether:petrol ether (30:70 v/v) for 30 min. After freezing for 20 min at −80°C, the organic phase was decanted, dried and resolved in 40% methanol. Corticosterone concentration was measured with a microtiter plate enzyme immunoassay. The samples were analyzed in two batches representing parts of the experiment performed on animals from the two generations of selection, but within each generation the order in which the samples were analyzed was randomized. In generation 23, the samples were analyzed alongside samples from another experiment (Lipowska et al., 2019); therefore, in this generation the analyses spanned more (37) assays than in generation 22 (24 assays). The intra-assay variation equaled 9.29% in generation 22 and 9.80% in generation 23, and the inter-assay coefficient of variation equaled 11.3% and 13.8%, respectively. The samples were analyzed in duplicate, and the analysis was repeated if the difference between values obtained for the duplicates exceeded 5%. The values we present derive from average values obtained from the two measurements of each sample. The analyses were not performed for 47 samples representing 43 animals, for which we failed to obtain a sufficient amount of blood or plasma was lost during the sample processing. These animals were excluded from further analyses.

The analyses were performed on 164 animals from the restraint test and 164 animals from the pharmacological stimulation test (328 animals total, 18–21 males and 20–23 females from each selection type per test; Table S1) for which a complete set of three plasma corticosterone level measurements was obtained: baseline level (C_base), response level (C_response) and the level after either short or long recovery (C_recovery). Two additional values representing corticosterone level change between consecutive blood samples were calculated: relative increase (C_increase) – the ratio of response to baseline corticosterone levels, and relative decrease (C_decrease) – the ratio of recovery to response corticosterone levels. The ratio values are mathematically equivalent to relative changes computed as percentage difference, recommended as a suitable measure for assessing the corticosterone feedback loop (Lattin and Kelly, 2020). We decided to present and analyze the relative changes as ratios (proportions), because the percentage difference changes could not be log-transformed (the values could be zero or negative), and log-transformation was needed because of severe right-skewness and heteroscedasticity of all the non-transformed corticosterone traits (strictly speaking, the distribution of residuals in models described below).

The analyses were performed with cross-nested mixed ANCOVA models, using SAS v. 9.4 (SAS Institute, Inc., Cary, NC, USA) Mixed procedure (with REML method of estimation and variance components constrained to non-negative values), separately for each test. All the models included selection direction, sex and generation as the main fixed factors, random effects of replicate line (nested within selection direction) and experimental block (nested within generation), and two covariates: time at the start of the trial and body mass (unless body mass was the subject of analysis). The hierarchical structure of the statistical model (replicate lines nested in selection direction) is required to allow a proper distinction of the effects of selection from random genetic effects, such as genetic drift (Henderson, 1997). This basic model structure was further expanded to accommodate additional factors adequate for specific analysis, as indicated in the next two paragraphs.

One type of model was applied to data containing combined results of the two tests, with additional fixed factor of test type (restraint versus ACTH) and random effect of family, correcting for the relatedness of sibling pairs (nested in replicate line and generation interaction). The analyses were performed for body mass and three corticosterone traits: baseline (including duration of baseline blood sampling as an additional covariate), response and the relative increase. To answer the question whether hypothetical differences between selection directions in the response level are independent of the baseline level, and whether the response and baseline corticosterone levels are correlated among individuals within the groups, a similar model for the response level, but with the baseline level as an additional covariate was also analyzed. Body mass used in the analyses (BM_avg) was the mean value from the daily measurements performed during the 3 days preceding the test. The value of body mass is subject to rapid changes, as it can shift by ≥1% as a result of a feeding bout, urination and defecation. Therefore, the mean mass of a few daily measurements is a better representation of the animal's size than a single record.

The other type of model was applied to the corticosterone recovery level and relative decrease, separately for data from the restraint and ACTH tests. The analyses were performed separately because the dynamics of recovery in the two tests is not comparable. The model included an additional fixed factor of recovery duration (short versus long). Because only one animal from each sibling pair was used in a given test type, the model did not include the factor of family.

The initial models included interactions among the main fixed categorical factors (selection, sex and generation) and interactions of the main factors with the additional fixed factors (test type for the response level and relative increase, and recovery duration for the recovery level and relative decrease), along with all the respective random interactions with the replicate line and experimental block. The models were then step-wise reduced by removing non-significant interactions (P>0.05). However, one-way interactions between selection direction and sex, test type or recovery duration were a priori considered meaningful predictors for either biological or technical reasons, and therefore were retained in the models (together with adequate interactions with replicate line) irrespective of their significance. In analyses performed on data containing combined results from the two tests, the models included an additional random interaction of experimental block×experiment type, corresponding to experimental procedures performed during a single day.

In all analyses, Satterthwaite's approximation was used to calculate the effective degrees of freedom (d.f.) for t-tests or the denominator d.f. for F-tests (i.e. the d.f. was computed from a combination of the d.f. of respective random grouping effects and residual term, weighted by variance contribution of the terms; SAS Institute Inc. 2011). Thus, the d.f. could take any real value between d.f. of the random factor and d.f. of the residual term. Significance of the random effects was tested with the likelihood ratio (LR) test, based on results from models with the same structure as described above, but with variance components not constrained to positive values (‘nobound’ option in SAS Mixed procedure).

Analyses of BM_avg and corticosterone response, recovery (after ACTH test) and relative decrease (after restraint test) revealed outlying individuals, one individual per model (absolute value of studentized residual ≥3.4). These individuals were excluded from analyses of the respective traits, but were retained in analyses of other traits in which their residuals did not stand out. The exclusion of these individuals from the respective analyses improved the normality of residual distribution and the model's goodness of fit (judged by the models' AIC values).

Information about group composition, complete descriptive statistics and results of LR tests for random effects in all the mixed models are presented in Tables S1–S3. Here, we provide overall descriptive statistics for body mass and corticosterone levels at each stage of response (range and mean±s.d. for pooled individuals from all groups) in the text, and the main results of statistical models: significance of the main factors of interest and least squares means with 95% confidence intervals (LSM±CI), computed for the approximate mean value of covariates, for all analyses, in tables and figures.

Body mass, averaged from 3 days preceding the tests (BM_avg), ranged from 15.2 to 38.2 g (mean±s.d. 22.8±3.8 g; Fig. 1). Males were heavier than females (P<0.0001) and body mass varied among selection directions (P=0.001), but the effect of selection was affected by sex (selection×sex interaction: P=0.065; Tables 1 and 2). Females from Herbivorous (H) lines were heavier than those from the Control (C) and Predatory (P) lines (Tukey–Kramer P≤0.024) but not heavier than those from the Aerobic (A) lines (Tukey–Kramer P=0.3), whereas Herbivorous line males were heavier than males from all other lines (Tukey–Kramer P≤0.019). Differences between Aerobic, Control and Predatory lines were not significant for either sex (Tukey–Kramer P≥0.17). Body mass varied significantly among the replicate lines within the selection directions (P<0.001) and among families within lines (P<0.0001; Table S3). Body mass of males, but not females, was higher in generation 22 than in 23 (males: P=0.031, females: P=0.8; generation×sex interaction: F_1,187=3.44, P=0.065). BM_avg did not differ between voles assigned to the two tests (P=0.97).

Before reporting results concerning the effect of main factors on corticosterone levels, we provide here a brief summary concerning secondary factors, included in the models mainly for the purpose of statistical control. Body mass did not significantly affect any of the corticosterone traits (P>0.08; Figs 1 and 4; Tables 1 and 3). The baseline (C_base), response (C_response) and recovery (C_recovery) levels were higher in generation 22 than in generation 23 (P≤0.06; Table 1, Table 3). Because the analyses of samples from the two generations were not performed at the same time, it cannot be determined whether the difference reflects real differences between conditions in the two generations or is merely a technical artefact. The corticosterone baseline and response levels tended to decrease with time of day, whereas the post-ACTH stimulation recovery level and relative decrease of corticosterone (C_decrease) decreased with time of day, and the remaining corticosterone traits were not significantly affected by time of day (Tables 1 and 3). The baseline, response, relative increase and post-restraint recovery levels varied among replicate lines (P≤0.058), but replicate lines did not contribute significantly to the variance of the other corticosterone traits. The random effect of family (P≥0.6) and the remaining random factors (P≥0.078) did not explain a significant proportion of the variance of any of the traits.

The baseline corticosterone level (C_base) varied greatly among individuals, ranging from 5 to 187 ng ml⁻¹ (43±30 ng ml⁻¹; Fig. 1A,B; Fig. S2, Table S2), but it was not affected by duration of blood sampling (F_1,288=0.02, P=0.9) and did not differ between voles assigned to the two tests (restraint versus ACTH, P=1.0). The effects of selection, sex and selection×sex interaction were not significant (P≥0.6; Fig. 2, Tables 1 and 2).

The response corticosterone level (C_resp) ranged from 45 to 614 ng ml⁻¹ (179±85 ng ml⁻¹) in the restraint test and from 56 to 561 ng ml⁻¹ (220±88 ng ml⁻¹; Fig. 1C,D; Fig. S2, Table S2) in the ACTH test. It was higher in females than in males (P=0.001; Table 2). The overall effects of selection (P=0.25) and test type (P<0.0001) were complicated by a significant interaction of the two (P=0.002; Table 1). In the ACTH test, the Predatory voles had a lower response level than the Control and Herbivorous ones (Tukey–Kramer P≤0.010; overall effect of selection: P=0.008, Fig. 2B), whereas in the restraint test the response level did not differ among selection directions (P=0.8; Fig. 2A). Animals from the Control and Herbivorous lines achieved a higher response level in the ACTH than in the restraint test (P<0.0001), but in Aerobic and Predatory lines the response levels in the two tests did not differ significantly (P≥0.16). The selection×sex interaction was not significant (P=0.2). The response levels were significantly correlated with the baseline corticosterone level (F_1,309=24.1, P<0.0001; Fig. 3). Introduction of the baseline level as additional covariate decreased significance of the effect of time of day (P=0.088) but did not affect the significance of other factors in the model (Table 1).

The relative increase of corticosterone level (C_increase: ratio of response to baseline levels) was smaller in the restraint test (6.09±4.81; range 0.79–23.45) than in the ACTH test (8.03±7.09; range 0.79–41.44; P=0.056; Fig. 4A,B, Table 1; Table S2). The effects of selection, sex, and of selection×sex and selection×test interactions were not significant (P≥0.4; Tables 1 and 2).

The recovery corticosterone level (C_recov) after a short (21 min) recovery in the restraint test (97±55 ng ml⁻¹; range 19–289 ng ml⁻¹) was higher than that after a long (60 min) recovery (36±31 ng ml⁻¹; range 6–206 ng ml⁻¹; P<0.0001; Fig. 1E, Table 3; Fig. S2A, Table S2). Animals from the Herbivorous lines tended to have higher recovery levels than those from the other lines (overall effect of selection: P=0.072; Fig. 2A), and the difference was nearly significant for comparison of the Herbivorous and Aerobic lines (Tukey–Kramer P=0.055), but not for other between-line comparisons (Tukey–Kramer P≥0.18). The effect of sex, and of selection×sex or selection×recovery duration interactions was not significant (P≥0.26).

The relative decrease of corticosterone level (C_decrease) was less profound (i.e. the ratio of recovery to response levels was larger) after a short recovery (0.62±0.35; range 0.11–2.04) than after a long recovery (0.21±0.18; range 0.02–1.09; P<0.0001; Fig. 4C, Table 3; Table S2). The decrease tended to be more pronounced in females than in males (P=0.062; Table 2). The effect of selection was complicated by a significant interaction with sex (P=0.0006): in Herbivorous females, the decrease was less pronounced (C_decrease was larger) than in females from all other lines (Tukey–Kramer P≤0.010; effect of selection P=0.0005; Fig. S3A); in males the effect of selection was not significant (P=0.16). The effect of selection×recovery duration interaction was not significant (F_1,21.5=1.27, P=0.3; Fig. 5B,C).

The recovery corticosterone level (C_recovery) after a short (21 min) recovery in the ACTH test (185±125 ng ml⁻¹; range 32–782 ng ml⁻¹) was higher than that after a long (60 min) recovery (545±37 ng ml⁻¹; range 7–210 ng ml⁻¹; P<0.0001; Fig. 1F, Table 3; Fig. S2B, Table S2). Females had higher corticosterone levels than males after the short recovery (P=0.024), but the effect of sex was not significant for long recovery (P=0.4; effect of interaction: F_1,114=7.20, P=0.008; Table 2). The effects of selection, and of selection×sex or selection×recovery duration interactions were not significant (P≥0.14; Fig. 2B).

The relative decrease of corticosterone level (C_decrease) was less profound (i.e. the ratio of recovery to response levels was larger) after a short recovery (0.87±0.70; range 0.12–5.64) than after a long recovery (0.30±0.22; range 0.03–1.33; P<0.0001; Fig. 4D, Table 3; Table S2). Despite a significant interaction between sex and recovery duration (F_1,130=5.51, P=0.020), the effect of sex was not significant for either recovery duration (P≥0.10; Table 2). The effect of selection was significant only for females (selection×sex interaction P=0.064; females: P=0.048; males: P=0.2): in Aerobic females the decrease was more profound (the relative decrease was smaller) than in Predatory females (Tukey–Kramer P=0.033, other between-line comparisons P≥0.12; Fig. S3B). The selection×recovery duration interaction was not significant (P=0.3; Fig. 5B,C).

The results of this study indicate that individuals vary up to an order of magnitude in absolute corticosterone levels at each stage of the response (Fig. 1). This finding is similar to what was observed in brown lemmings (Fauteux et al., 2017), California mice (Dlugosz et al., 2012) and our earlier study on bank voles (Lipowska et al., 2019). The variation may be associated with the pulsatile character of corticosterone release in non-challenged animals: as found in rats, plasma corticosterone levels fluctuate by one or more orders of magnitude in the course of less than an hour (Carnes et al., 1989; Sarabdjitsingh et al., 2012; Windle et al., 1998). The effects of this ultradian corticosterone level variation, circadian corticosterone level variation (Dlugosz, 2012; Harris et al., 2012; Malisch et al., 2008) and the HPA system's sensitivity to environmental conditions are likely to contribute to the low repeatability of baseline corticosterone levels in free-living animals (Taff et al., 2018). Our experimental design reduced the effect of environmental and circadian corticosterone level variation through maintenance of standardized housing conditions and sampling within a narrow timeframe. In our pilot study on bank voles maintained in conditions similar to those applied in the experiment, the baseline corticosterone level was weakly but significantly repeatable within a 6 week period (coefficient of intraclass correlation, ρ=0.15, P=0.030; M.M.L., E.T.S., B.B.-S. and P.K., unpublished data). Thus, despite the large inter-individual differences, the effect of treatments and differences between groups should also be observable.

Besides the inter- and intra-individual variation, corticosterone level also varies greatly among species (Romero et al., 2008; Taymans et al., 1997). Moreover, analytical method has a profound effect on the glucocorticoid values measured (Gatti et al., 2009). Thus, comparisons of absolute values we have obtained with those reported in the literature for this or other species may be not very informative. However, both the baseline and restraint-induced corticosterone levels of bank voles fit within the ranges reported for mice (Malisch et al., 2007; Ridder et al., 2005; Rozeboom et al., 2007) and rats (Vahl et al., 2005). We also found that female voles achieved higher corticosterone levels after stimulation than males did, a tendency that has also been observed in other rodents (Harpaz et al., 2013; Malisch et al., 2007; Rozeboom et al., 2007; Seale et al., 2004; Taymans et al., 1997).

In a previous study on our animal model, we tested the effect of selection in the Aerobic lines on the corticosterone response to the swimming trial (Lipowska et al., 2019). The baseline corticosterone levels reported in our previous study were similar, although somewhat lower than those reported here (mean±s.d. 34±20 versus 43±30 ng ml⁻¹). Similar to the restraint stress and ACTH treatment applied here, the swimming trial elicited an increase in corticosterone levels. The corticosterone level achieved after 10 min of swimming at 38°C (168±56 ng ml⁻¹) was slightly lower than that induced by restraint (179±85 ng ml⁻¹), and considerably lower than that induced by ACTH (219±88 ng ml⁻¹). Despite the correlation between the corticosterone response to swimming and the metabolic rate achieved during swimming, we did not observe an effect of selection on the corticosterone response to the trial (Lipowska et al., 2019). However, because of the metabolism-regulating role of glucocorticoids (Girard and Garland, 2002), the swim-induced increase in corticosterone level cannot be fully attributed to the stress response. The restraint stress and ACTH stimulation we applied here provide a more direct approach to test the effect of selection on the stress response.

The absolute values of baseline and restraint-induced corticosterone levels did not differ among selection directions, whereas the response to ACTH stimulation was decreased in the Predatory lines. Interestingly, the corticosterone response to restraint in the Aerobic and Predatory lines was nearly equal to the maximal response, whereas in the Control and Herbivorous lines the response to restraint was significantly lower than the maximal response. Thus, despite the lack of difference in the absolute levels of corticosterone representing response to the stressor, it can be concluded that animals from the Aerobic and Predatory lines reacted more strongly relative to the maximal scope of the response than those from the Control and Herbivorous lines. However, this conclusion is weakened by the lack of such a difference in the relative increase of corticosterone level in response to stimulation (the ratio of response to baseline levels).

As neither the baseline nor the response corticosterone levels differed markedly between the selected lines, the across-lines comparison could not answer the question whether the baseline and response levels are correlated. However, the analysis of covariance showed a positive inter-individual correlation between the baseline (C_base) and response (C_response) corticosterone (Fig. 3). Our results documented only a phenotypic correlation, but Béziers et al. (2019) showed a positive additive genetic correlation between the baseline and stress-response levels in the barn owl. The results suggest that evolution of the HPA function in the baseline and stressed state may be not independent. Another intriguing question is whether regulatory functions of corticosterone at its baseline and stress-response levels are correlated. As corticosterone at low and high concentrations interacts primarily with different receptor types (Landys et al., 2006), the downstream regulatory functions at the baseline and stressed states could be indeed independent. However, although the mean C_base value in bank voles was more than four times lower than the C_response level (Fig. 2), the huge individual variation in both of these traits resulted in considerable overlap of the ranges (Figs 1 and 3; Fig. S2), which renders their complete functional separation unlikely. To summarize, analysis of the individual variation of corticosterone in bank voles showed that neither the levels of C_base and C_response nor their regulatory effects should be considered as independent.

As expected, the pattern of corticosterone level decrease during recovery from restraint stress differed from that after pharmacological stimulation. In the restraint test, corticosterone level observed at the response phase decreased by half within 21 min of short recovery and achieved values comparable to or even lower than the baseline level within 60 min of long recovery (Fig. 2). In the ACTH test, the pattern of corticosterone level changes suggests that in this test recovery was either slower or delayed relative to that in the restraint test. Considering that in the ACTH stimulation test the animals were returned to their home cages immediately after each blood sampling event, whereas in the restraint test the animals were challenged by being placed in an unfamiliar, stressful environment, we can assume that psychological factors play a greater role in recovery from the latter form of stimulation. Conversely, stimulating the animals with a uniform dose of ACTH omits the potential inter-individual differences in ACTH released in response to stress, highlighting the role of physiological factors in contributing to the rate of corticosterone recovery after stimulation. The difference in the character of recovery in the two tests is reflected in our findings: the slower recovery of the Herbivorous line voles was observed only in the restraint test, whereas the delayed onset of recovery of female voles was observed only in the ACTH test. Thus, we can postulate that selection affected the perception of stress in the Herbivorous lines, whereas the difference between sexes is of physiological origin.

The analyses of C_decrease during the recovery period (the ratio of recovery to response levels) revealed a significant interaction between selection and sex, which was not revealed by the analyses of the raw corticosterone level. The reduced recovery capability observed in the Herbivorous lines after restraint stress was significant only for females, whereas in the ACTH test the Aerobic line females tended to recover faster than those from other lines. Therefore, it is likely that the effect of selection on the HPA axis activity can be stronger in females than in males, which can potentially be related to the sex differences in corticosterone levels found in bank voles and other rodent species (Harpaz et al., 2013; Malisch et al., 2007; Rozeboom et al., 2007; Seale et al., 2004; Taymans et al., 1997).

We supposed that the selection in the Herbivorous lines should promote a calm reaction to unexpected challenges, but our results did not confirm such a prediction. However, we observed a decrease in the rate of recovery after a stressful event, although this trait is unlikely to contribute to the animals' performance in the selection test, and therefore should not be directly affected by selection pressure. This effect might put the Herbivorous line animals in a higher risk of negative outcomes of prolonged exposure to elevated glucocorticoids (Korte et al., 2005; Romero, 2012) if they are frequently exposed to acute stressors.

Despite the corticosterone levels varying by more than an order of magnitude among individual bank voles, we were able to observe the expected differences between sexes in the parameters of the corticosterone response curve, as well as significant or near-significant effects of selection on the maximal corticosterone level (Fig. 2). Therefore, the lack of statistical significance of the effect of selection on the basal, restraint-induced and post-ACTH recovery corticosterone levels is unlikely to be explained merely by the obscuring effect of large inter-individual variation. It is then particularly surprising that this study revealed only a few statistically significant differences among the four types of lines. The selection had a striking effect on the voles' performance in the potentially stress-inducing conditions of the selection tests, including changes in metabolic rate, body mass and behavior, all of which are traits associated with HPA axis function. It was therefore reasonable to expect changes at the level of the axis' reactivity to stimulation. Thus, the results of this study partially undermined our hypothesis that alteration of the glucocorticoid stress response is essential for the evolution of increased performance in challenging situations. Taken at a face value, the results imply that the evolution of complex physiological adaptations associated with distinct mammalian lifestyles is not constrained by a conservative characteristic of the HPA axis. However, despite the lack of change in restraint stress-induced corticosterone levels, it cannot be concluded that the stress-response system was not affected by selection. Our results suggest that the maximal scope of the response decreased in two of the selected lines, probably affecting the strength with which the stress-induced corticosterone level is perceived. Moreover, perception of the corticosterone-mediated signal can be modulated through expression of glucocorticoid receptors (Veenema et al., 2004), affecting the animal's sensitivity to corticosterone and behavioral reaction to stress (Kolber et al., 2008; Reichardt et al., 2000; Ridder et al., 2005). If selection affected receptor expression it could contribute to differences in stress perception among selection directions independent of the lack of difference in absolute corticosterone levels. As the selection experiment on bank voles is continued, the characteristics of the HPA axis activity and signaling can be addressed in further studies.

In conclusion, the results demonstrate that although evolution of markedly increased performance in challenging situations does not require changes in the absolute values of baseline and challenge-induced corticosterone levels, it can affect aspects of HPA axis activity that have a less straightforward association with response to the challenge. In particular, evolution of high alertness and readiness to perform physical activity (i.e. Aerobic and Predatory lines) may be supported by lowering of the maximal glucocorticoid response, resulting in an increased relative scope of the response to a challenge. Furthermore, evolution of increased performance in a challenge involving a prolonged period of mild stress (Herbivorous lines) might increase risks associated with prolonged exposure to elevated glucocorticoid level if an animal is challenged with acutely stressful situations. However, our findings also emphasize the stability of the HPA axis activity parameters, which can remain almost unchanged despite the evolutionary changes in performance traits.

We are very grateful to the many technicians and students who helped with animal maintenance and sampling, especially to S. Kunysz. We also thank the anonymous referees for constructive comments to improve the manuscript.