Body temperature is a repeatable trait in a free-ranging passerine bird

Reaction norms can be used to simultaneously quantify both the level of phenotypic flexibility across an environmental gradient (i.e. reaction norm slope) and the phenotypic expression of a trait in the population-average environment (i.e. reaction norm intercept). Reaction norm components of several flexible physiological traits have been quantified, many of which are related to energy metabolism (e.g. basal metabolic rate, glucocorticoids: reviewed in Biro and Stamps, 2010) and immunity (e.g. Schreier and Grindstaff, 2020; Tieleman et al., 2010). However, in endotherms, individual variation in the phenotypic flexibility of body temperature (T_b) has been understudied by evolutionary physiologists, despite its importance in shaping organismal performance (Angilletta et al., 2010). T_b exhibits substantial plasticity in endotherms (Boyles et al., 2013), as evidenced by adjustments in T_b in response to different environmental challenges (e.g. food availability, thermal stress), activity level, daily fluctuations in rhythmicity (e.g. circadian and ultradian rhythms) and regional heterothermy within individuals as a consequence of heat flow across the body (Angilletta et al., 2010). Yet, studies investigating variation in T_b (and other thermoregulatory traits) have primarily focused on differences among species and populations (e.g. Noakes and McKechnie, 2019; Noakes and McKechnie, 2020; Smit et al., 2013), or among groups within populations (e.g. age and sex: Cai et al., 2016; Gagnon and Kenny, 2012; Andreasson et al., 2020a,b) rather than among individuals within a single population (e.g. Pessato et al., 2020; Bozinovic, 2007; Møller, 2010).

Evidence from laboratory studies, both on wild-caught and lab-raised individuals, indicates that among-individual variation in T_b exists at the individual level, and may be related to differences in metabolic traits that vary across environments (Briga and Verhulst, 2017). For instance, in birds, individuals with higher basal metabolic rate (a trait with reported genetic basis: Rønning, et al., 2007) also maintain higher T_b during exposure to temperatures below thermoneutrality (Briga and Verhulst, 2017; Stager et al., 2020), and in mice, individual consistency in T_b variation predicts torpor use (Nespolo et al., 2003; Boratyński et al., 2019). Further, variation among individuals in resting and basal metabolic rate has been widely reported in laboratory and free-ranging studies (e.g. Boratyński et al., 2017; Broggi et al., 2009; Hayes et al., 1998; Nespolo and Franco, 2007; Rønning et al., 2005; White et al., 2013). Finally, genetic variation in T_b has been reported in both domestic and laboratory species, including chickens (h² range, 0.10–0.20: Li et al., 2020; van Goor et al., 2015), dairy cattle (h²=0.17: Dikmen et al., 2012), sows (h²=0.35: Gourdine et al., 2017) and mice (h²=0.36: Connolly and Lynch, 1981).

While there is evidence demonstrating consistent individual differences in thermoregulation in laboratory and domestic animals, far fewer studies have formally examined among-individual variation in T_b in free-ranging species (e.g. Bozinovic, 2007; Dammhahn et al., 2017; Møller, 2010; Nespolo et al., 2003, 2010; Pessato et al., 2020). Free-living animals often experience multiple environmental challenges simultaneously, and the degree of individual consistency and plasticity in T_b might be adjusted differently depending on the environmental context, physiological condition and life-history strategy of the individual (e.g. Cornelius et al., 2017). For example, metabolic rate in free-living animals has been reported to be less repeatable than in laboratory animals, which may be a consequence of how animals differ in their response to the same conditions over time (Auer et al., 2016). Further, exposure to certain environmental conditions during development, such as high environmental temperatures, may induce epigenetic changes in aspects of thermoregulatory physiology (see Ruuskanen et al., 2021, for a review), such as heat shock protein regulation (Kisliouk et al., 2017) and morphology (Burness et al., 2013; see Nord and Giroud, 2020, for a review). Given, that wild animals experience fluctuating environmental conditions both intra- and inter-annually, there may be greater variation in the magnitude and type of epigenetic and parental effects experienced by wild versus laboratory animals.

To understand the degree to which T_b consistently varies among free-living birds, we examined the repeatability of T_b in tree swallows, Tachycineta bicolor, across differing ambient temperatures (T_a). We asked: (1) do birds consistently differ from each other in their average T_b (i.e. is T_b repeatable at the reaction norm intercept)?; and (2) do birds consistently differ from each other in their T_b response across ambient temperature T_a (i.e. is the reaction norm slope repeatable)? Here, we used previously published data from our research group (Tapper et al., 2020a,b) to examine repeatability in swallows.

Because inter-specific variation exists among thermoregulatory traits (Pessato et al., 2020; Stager et al., 2020), we predicted that in both female and male swallows, T_b would be repeatable, both in terms of the average trait expression (i.e. intercept) and how individuals vary T_b across T_a (i.e. slope). However, because female tree swallows have a brood patch, whereas males do not, there may be differences in heat dissipation capacity between the sexes (Nord and Nilsson, 2019). Therefore, we may expect males to forage closer to their maximally allowable T_b relative to females, which could lead to differences in repeatability between the sexes. This could occur if, for example, the upper threshold on T_b does not vary between individuals for biochemical reasons. Alternatively, if both sexes reduce their activity level at high temperatures, meaning that hyperthermia may be a limiting factor, there may be similar repeatability in T_b between male and female swallows.

In our previous studies (Tapper et al., 2020a,b), we experimentally removed the ventral feathers from a sample of breeding females to increase heat dissipation capacity (i.e. trimmed treatment). Although testing for a treatment effect on T_b repeatability was not an explicit goal of this study, previously we found that control birds had higher T_b than trimmed birds (by ∼0.5°C) across temperature (6–31°C) (Tapper et al., 2020b), suggesting that differences in heat dissipation capacity could lead to differences in the repeatability of T_b. Following the same reasoning as above, control and trimmed birds may differ in T_b repeatability if control birds forage relatively closer to their maximally allowable T_b, and if the upper threshold on T_b does not vary between individuals. Alternatively, T_b repeatability may be similar between experimental treatments if individuals tend to adjust their activity level to avoid overheating (Tapper et al., 2020b). In this case, the effects of heat dissipation capacity on T_b should be outweighed by behavioural (or physiological) differences among individuals that directly relate to individual control of T_b.

Our dataset is derived from Tapper et al. (2020b), in which a detailed description of our study sites, field methods and reproductive parameter definitions can be found. However, given that environmental variables such as wind speed, relative humidity and precipitation could affect individual variation in thermoregulation, we present these summary statistics for our study sites here (Table 1). We performed this study in May–July 2018, on two nest-box breeding populations of tree swallows, Tachycineta bicolor (Vieillot 1808), in Peterborough, ON, Canada. We captured female tree swallows during incubation (day 7–10, n=16) and, when possible, their respective male partners during early nestling provisioning (day 4–6, n=10). Upon capture, we implanted individuals with a thermal-sensitive passive integrated transponder (PIT) tag (accuracy ±0.5°C, Bio-Therm13, Biomark, Boise, IA, USA) subcutaneously into the nape of the neck to remotely monitor T_b. PIT tags can shift position post-implantation and if such shifts occurred in our study, this could change among-individual differences in T_b (i.e. increase or decrease repeatability estimates). This will depend on whether the tag was measuring T_b from different body regions within the same individual at different times. However, several studies have implanted PIT tags in the nape and reported minimal issues with respect to tag movement (e.g. Nord et al., 2013, 2016; Oswald et al., 2018), and thus variation in subcutaneous T_b from PIT tag placement seems unlikely. Similarly, in the unlikely event that a tag was placed too close to the skin surface, estimates of T_b could be underestimated as a result of the influence of environmental temperature. Finally, repeatability could be affected by differences in individual tag accuracy, but these differences are likely to be small, as reported in other studies (e.g. 0.2°C: O'Connor et al., 2021). We read PIT tags with Biomark HPR Plus readers (with the minimum time between consecutive reads of the same tag set to 10 s: Tapper et al., 2020b), which we connected to a loop antenna and positioned so that they encircled the nest box entrance.

We recaptured individual female swallows during early nestling provisioning (range: day 1–2 post-hatch) and performed our experimental manipulation (experimental trimming versus handling) to allow for increased heat loss. We assigned females to each treatment by flipping a coin between ‘pairs’ of birds, such that the first female captured received one treatment, while the second female captured received the alternative treatment. We repeated this process for each new dyad. In the case of individuals that had been experimentally manipulated in the previous year, they were typically assigned to the opposite treatment from that which they had experienced previously. Therefore, while we assigned pairs of females to each treatment based on a coin flip, this was not always possible, and on occasion we had to selectively alternate between treatment types. For each female, we either trimmed the contour and downy feathers covering the brood patch to expose the bare skin underneath (∼7% of the surface area of the bird; trimmed treatment) or ascribed them a control condition where they were handled but released with their feathers intact. All research was approved by the Trent University Animal Care Committee, in accordance with the Canadian Council on Animal Care (AUP no. 24747).

For data compilation and statistical analyses, we used R (version 4.0, R Core Team 2020; http://www.R-project.org/). To maximize our sample size (i.e. assess the repeatability of T_b under a broad range of T_a), we included all T_b measurements made between 2 and 16 days post-hatch, and between 05:00 h and 21:00 h (when swallows are diurnally active). We attempted to record individuals for 24 h at least 3 times throughout the nestling period, once between the periods of day 2–5, day 6–10 and day 11–16. While most individuals had balanced data across the nestling period (Fig. S1A), the T_b data of some individuals covered a wider range of T_a (6–31°C) than others (Fig. S1B). To control for differences in the number of T_b measurements among birds, and to minimize bias in our estimates of individual T_b variance, we calculated the weighted average hourly T_b (°C) for each bird. To do this, we used maximum likelihood to approximate the distribution of the number of reads within a single nest visit (gamma distribution, shape=0.56, scale=37.58) and the number of visits within an hour (normalized between 0 and 1, beta distribution, α=1.02, β=4.63). Approximating the distributions helped to smooth the probabilities and therefore avoid giving a certain number of reads (within a visit) or visits (within an hour) too much or too little weight. The mean number of reads per bird across females and males was 4984 (range 1559–15,917) and 1247 (range 393–3427), respectively. The mean number of visits per bird across females and males was 626 (range 198–1855) and 544 (range 193–1460), respectively. As in Tapper et al. (2020b), we excluded data from one bird that had abnormally low T_b (∼2°C) between 05:00 h and 13:00 h on one day relative to its T_b between 05:00 h and 13:00 h on all other days (i.e. 13% of all observations for that bird).

To test whether T_b was repeatable in our tree swallows, we used Bayesian hierarchical mixed effects models, constructed using the package ‘brms’ (version 2.14, Bürkner, 2017). First, we quantified the repeatability and phenotypic flexibility of T_b of each sex, based on control bird data only. While removing the trimmed birds reduced our sample size from 29 (n_female=16, n_male=13) to 16 individuals (n_female=8, n_male=8), any differences in variance between treatments may have been confounded with sex differences on T_b if trimmed birds were included in the analysis. Thus, removing trimmed individuals increases the interpretability of our results.

In our model, we used a ‘heterogeneous error structure’, partitioning the different variance components (i.e. intercept, slope, intercept–slope co-variance, and residual error terms), into female (n=8) and male (n=8) groups, respectively. Because estimates of variance in reaction norm components (i.e. intercept and slope) differ depending on whether potential confounding factors are controlled for in the analysis, we included several covariates (see below) in our model to calculate the adjusted repeatability (Nakagawa and Schielzeth, 2010). We also mean-centred and standardized all covariates to two standard deviation units (Gelman, 2008). Our response variable was therefore mean hourly T_b, and our fixed predictors were treatment, T_a, lay date, brood size, nestling provisioning rate, and a ‘sex×T_a’ interaction. We note that in our preceding study on the effects of treatment on T_b (Tapper et al., 2020b), we modelled T_b as a polynomial function of T_a. Here, we chose to model T_b as a linear function of T_a to increase interpretability of the reaction norm components (both intercepts and slopes) and to avoid overfitting our models (as we had insufficient data to properly model non-linear relationships per bird). We also ran our models with wind speed and relative humidity included as fixed effects (we did not have hourly precipitation data), but the results, and consequently our conclusions, remained similar to the model with these factors excluded. Therefore, we excluded these variables from the final analysis to reduce the number of model terms in our analysis. To determine the degree of phenotypic flexibility in T_b in response to changes in T_a, we also included random slopes for individual identity by T_a. To compare the degree of phenotypic flexibility in T_b between treatments, we ran a Levene's test on the slope estimates from the model. Lastly, we interpreted differences between sexes in the mean repeatability intercept (R_intercept) and slope (R_slope) using Bayes factors (K), as derived from one-way hypothesis tests (calculated with the Savage–Dickey density ratio method) in ‘brms’. To keep our tests conservative, we asked whether the sexes differed from each other in R_intercept and R_slope by a repeatability of >0.1, or 10% (e.g. R_{intercept trimmed}−R_{intercept control}>0.1), meaning that we considered a difference in R≤0.1 as practically equivalent to 0 (see Kruschke, 2018, for defining regions of practical equivalence). Support for a difference in repeatability would be given by K>1, with larger K values indicating stronger support and smaller K values weaker support (see Jeffreys, 1961, and Lee and Wagenmakers, 2013, for a guideline on interpretation of Bayes factors). While the choice of 0.1 is arbitrary, repeatability estimates from other avian physiological traits (e.g. basal metabolic rate, glucocorticoid levels) have been reported to cluster between 0.1 and 0.4 (Auer et al., 2016; Schoenemann and Bonier, 2018), and we therefore based our boundary on the low end of the repeatability estimates reported in the literature. Our priors for the intercept repeatabilities were beta distributed (α=12, β=3) for both sexes, and for female and male slope repeatabilities they were beta distributed (female, α=1, β=4; male, α=4, β=3).

Next, we quantified the repeatability and phenotypic flexibility of T_b in females (n_control=8, n_trimmed=8) to determine whether heat dissipation capacity (i.e. experimental trimming treatment) affected among-individual variation in T_b. Our second model was structured similarly to our first, except that population-level effects of sex and sex×T_a were replaced with the population-level effects of treatment and treatment×T_a. Additionally, instead of partitioning variance into female and male groups, we partitioned the variance components into control and trimmed groups, respectively. As above, we compared the phenotypic flexibility in T_b between treatments by running a Levene's test on the slope estimates from the model. We also tested for differences in the mean repeatability intercept and slope, with a statistical difference between groups in repeatability as >0.1. For both treatments, priors were beta distributed for the intercept (α=12, β=3) and slope (α=1, β=4).

For both the ‘sex’ and ‘treatment’ model, we used informative prior distributions for the population intercept (i.e. T_b, °C), T_a and feeding rate coefficients, which were based on model estimates in Tapper et al. (2020b). We used gamma distributions for the population intercepts (shape=84, scale=2; producing a peak density of ∼42) and T_a coefficients (shape=1.2, scale=0.75; producing a peak density of ∼0.3) to fix T_b, T_a and feeding rate estimates above 0. For feeding rate and treatment, we used normal distributions (mean=0 and s.d.=1) and for the random effect terms, we used gamma distributions (shape=1.75, scale=0.75; yielding a peak density of ∼0.5). For the remaining model parameters (i.e. population-level fixed effects and residual error terms), we used the ‘brms’ default non-informative priors, which follow Student's t distributions (with degrees of freedom=3, location=0, scale factor=2.5). We corrected for autocorrelation between adjacent T_b measurements by adding an AR1 correlation structure (ρ=0.64). We ran models for 10,000 iterations, discarding the first 1000 iterations as burn-in, and sampling every 10 iterations. for all model parameters fell between 0.99 and 1.01, and the effective sample size was greater than 2000 for all parameters. Autocorrelation was low among consecutive thinned observations (r<0.10 in all models). We viewed trace plots for fixed and random effects to ensure appropriate sampling of the posterior distribution, and we ran three chains of each model to check for model convergence using the Gelman–Rubin diagnostic (Gelman and Rubin, 1992). We report our model estimates as the mean of the posterior distributions ±95% highest density credible intervals (CI).

The grand mean T_b averaged across all control female and male swallows (n=16) was 42.1°C at the intercept (T_a=0, or ∼21°C), with T_b ranging from 41.07 to 43.8°C (Fig. 1A). On average, the mean T_b of female and male swallows increased with T_a (β=0.38, 95% CI [0.24–0.51]) and feeding rate (β=0.11 [0.07–0.15]), and there was strong evidence that females were ∼0.8°C warmer than males (i.e. sex; Table 2; 97% posterior probability distribution<0), as indicated by the strength of the coefficients, associated credible intervals, and the percentage of posterior probability distribution that fell below 0. However, there was no evidence that the sexes differed in mean T_b as a function of T_a, at the population level (i.e. sex×T_a, β=−0.03 [−0.30–0.24]).

Repeatability of average T_b was high, but similar, for the two sexes (i.e. R_intercept; Table 3), and hypothesis tests confirmed no evidence for a difference between female and male T_b repeatability (β=−0.07 [−0.18, 0.05], K=0.18). However, there was substantial evidence that females and males differed in their mean slope repeatability (Table 3, R_slope,female=0.14, R_slope,male=0.58; Fig. 2), as indicated by the Bayes factor and posterior probability (hypothesis test, β=0.34 [−0.07, 0.68], K=11.8, posterior probability=92%). Additionally, the degree of phenotypic flexibility (i.e. variance in reaction norm slopes) varied between sexes (Levene's test, F_1,14=7.63, P=0.02). Taken together, these results provide strong evidence that males displayed greater among-individual variation in their T_b across T_a compared with females. Relatedly, there was a negative relationship between intercepts and slopes in males (Fig. 2), but no evidence for a relationship in females (i.e. V_cov; Table 3). For males, this means that individuals with lower average T_b tended to be more responsive to changes in T_a compared with individuals with a higher average T_b (i.e. negative covariance).

The grand mean T_b averaged across all female swallows (n=16) was 42.1°C at T_a≈21°C, but individuals differed from each other by as much as 2.49°C (Fig. 1B; range 40.80–43.29°C). The mean T_b of female swallows increased with T_a (β=0.36 [0.23–0.47]) and feeding rate (β=0.14 [0.10–0.18]), with some evidence that T_b was lower in trimmed compared with control birds (i.e. treatment, β=−0.54 [−1.43, 0.33]; 88% posterior probability distribution<0) (Table 4). Further, while control and trimmed birds maintained relatively similar T_b across T_a (i.e. treatment×T_a, β=−0.10 [−0.29, 0.11]; Fig. 3), 85% of the posterior probability distribution was <0, suggesting that trimmed birds had, on average, shallower slopes across T_a than control birds.

Repeatability of mean T_b in control and trimmed females was high (R_{intercept control}=0.89 [0.77, 0.99]; R_{intercept trimmed}=0.79 [0.59, 0.96]; Table 5), but there was no evidence for differences in repeatability between the treatments (hypothesis test, β=−0.01 [−0.19, 0.20], K=0.88). There was low repeatability in individual slopes in both treatments (R_{slope control}=0.12 [0.00, 0.41]; R_{slope trimmed}=0.22 [0.00, 0.55]; Table 5), indicating that individuals did not consistently maintain differences in T_b across T_a, irrespective of treatment. Additionally, we found no evidence for a difference in treatment slope repeatability (hypothesis test, β=−0.20 [−0.56, 0.15], K=0.16). However, the degree of phenotypic flexibility (i.e. variance in reaction norm slope) varied between treatments (Levene's test, F_1,14=12.38, P=0.003).

We found no support for a relationship between individual identity intercepts and slopes (V_{cov control}=−0.05 [−0.97, 0.88]; V_{cov trimmed}=−0.24 [−1.00, 0.59], Table 5), in both treatments, meaning that an individual's average T_b and the plasticity of T_b (across changes in T_a) were not related.

We found strong evidence that average T_b is a repeatable trait in female and male tree swallows, meaning that in both sexes, individuals consistently differed from each other in their mean daytime T_b. We also found evidence that T_b was flexible across T_a, and this flexibility was more repeatable in males than in females. The reason for the greater repeatability in the male reaction norm slope is unclear but may be related to sex-specific differences in life-history strategies. While the feeding rate, and variation in feeding rate, was similar between males and females in our study (mean±s.d. 10.4±2.4 females, 10.9±2.7 males), males may be more variable in their overall contribution to parental care than females, which may be reflected in other parental care metrics not measured here (e.g. food load: du Plessis et al., 2012; nest cleaning). Male tree swallows are also known to exhibit polygyny (Dunn and Robertson, 1992), which may affect the time spent attending to their primary nest. Finally, males may be more variable in their body condition, or other aspects of morphology that affect insulation (and consequently T_b), such as feather density or length.

We also found that males with lower average T_b tended to have sharper increases in T_b with increasing T_a, compared with males with higher average T_b. One possible explanation is that males with high T_b may be unable (or unwilling) to raise their T_b further because of costs associated heat stress (i.e. heat dissipation limitation hypothesis: Speakman and Król, 2010). This would be consistent with our previous findings showing that limits to sustained energy expenditure in tree swallows may be a consequence of heat dissipation ability at high temperatures (Tapper et al., 2020b).

That we did not see this same effect in females, despite their higher average T_b than males at the population level, may be due to experimental error; for example, if the PIT tag of the female with the highest T_b (Fig. 1A, ID 8F) was implanted in deeper tissue (muscle) than in the other birds, it could raise the population-average T_b in females. While excluding this individual from the analysis reduced the magnitude of the difference between the sexes' average T_b (sex coefficient: 0.85 included, 0.60 excluded), it does not change our conclusions, given that most of the posterior distribution still excluded 0 (96% <0 included, 90% <0 excluded). A biological reason could be that because males lack a brood patch, they are more limited in their ability to dissipate heat relative to females. Therefore, as work rate increases, T_b may rise quicker in males than in females. Alternatively, females and males may differ in their tolerance of thermally induced cellular damage (reviewed in Hansen, 2009; Iossa, 2019; Pérez et al., 2008; Romero-Haro et al., 2016; but see Walsh et al., 2019) and, consequently, females may be less restricted in their ability to raise T_b. This, however, remains speculative in light of research showing no sex-specific differences in oxidative stress in small passerine birds (Costantini et al., 2006, 2007). A final explanation may be that male and female T_b are differentially constrained by other physiological traits. Regardless of the reason, differences in trait correlation between the sexes may imply historical differences in the strength of selection acting on T_b, because selective pressures tend to constrain trait variation within populations across time (Baker et al., 2018; Bell and Sih, 2007; DiRienzo et al., 2016).

Experimentally manipulating heat loss did not influence the repeatability of average T_b (i.e. reaction norm intercept, 21°C) or T_b across T_a (i.e. reaction norm slope), although it did influence the degree of phenotypic flexibility in T_b, as indicated by different variance in slopes between treatments. That there was greater variance in the T_b slope estimates among trimmed individuals relative to control individuals, but no difference in slope repeatability, suggests that the greater variance in slope steepness observed among trimmed birds is likely due to differences in how individuals respond to environmental factors. For example, sources of variation likely to cause within-individual differences, such as variable environmental conditions (e.g. solar radiation, wind speed, insect abundance), daily fluctuations in internal state (e.g. hormonal levels), short-term behavioural adjustments (e.g. trailing legs during flight, fluffing feathers) and sources of experimental bias (e.g. measurement error, missing data), may have influenced the T_b slope response more in trimmed compared with control birds. Alternatively, it is possible that we did not detect a difference in slope repeatability between treatments because of small sample sizes and a lack of statistical power.

While there were no differences in the repeatability of T_b between treatments, control females were warmer than trimmed females in average T_b, and their average increase in T_b across T_a (i.e. black line, Fig. 3) was also steeper than in trimmed females. Because of similar repeatability in T_b between treatments, the population-level differences seen in this study are likely to be a consequence of individuals within each group responding consistently to the effect of treatment on T_b. However, it is possible that we may have seen greater inter-individual variation in response to the treatment had we trimmed a larger surface area (i.e. increased the area of heat transfer), because of the difficulty in regulating T_b with increased surface area exposure.

Despite observing low slope repeatability in both groups, there were consistent differences in the way individuals responded to changes in T_a. For instance, some individuals increased their T_b more steeply in response to increasing T_a than others (Fig. 3), suggesting there are different strategies of regulating T_b (Tapper et al., 2020b), or differences in the tolerance of higher T_a (Noakes and McKechnie, 2019). Additionally, the flexibility in the T_b by T_a response could have been influenced by variation in individual energy expenditure, which could arise because of differences in food availability, or the interaction between food availability and T_a (Briga and Verhulst, 2017). Regardless of the reason, such phenotypic flexibility within populations opens the possibility that inter-individual differences in thermoregulatory traits are heritable, and that phenotypic flexibility in T_b may itself be heritable. However, we recognize that even heritability does not necessarily imply trait responsiveness to natural selection (Vatka et al., 2020), and that a variety of phenotypically plastic traits have been reported to be non-responsive to selection by T_a (Arnold et al., 2019).

We note that the repeatability estimates reported here may be not representative of the among-individual variation in swallow core T_b, given that our metric of T_b was peripheral to the core, and consequently is likely to be influenced by both changes in environmental conditions and heat flow across the body, at least to some degree (Nord et al., 2016; Andreasson et al., 2020b). Future studies should consider simultaneously measuring individual variation in both deep and peripheral T_b, remotely, to determine the degree of variability that arises from recording T_b closer to the skin surface. Additionally, it remains to be seen whether among-individual differences in the flexibility of T_b that we report are heritable, and ultimately adaptive.

We would like to thank A. Schubert, A. Dain and J. Baici for help with data collection; A. Gerson for loaning RFID equipment; J. Robertson for advice on code and statistics; and S. Morin for edits. We also thank S. Verhulst and an anonymous referee for many comments that improved the manuscript.