Statistical Estimation of Surface Heat Control and Exchange in Endotherms

In this paper, we combine polynomial functions, Generalized Estimating Equations, and bootstrap-based model selection to test for signatures of linear or nonlinear relationships between body surface temperature and ambient temperature in endotherms. Linearity or nonlinearity is associated with the absence or presence of cutaneous vasodilation and vasoconstriction, respectively. We obtained experimental data on body surface temperature variation from a mammalian model organism as a function of ambient temperature using infrared thermal imaging. The statistical framework of model estimation and selection successfully detected linear and nonlinear relationships between body surface temperature and ambient temperature for different body regions of the model organism. These results demonstrate that our statistical approach is instrumental to assess the complexity of thermoregulation in endotherms.


Introduction
Endotherms (mammals and birds) regulate body temperature primarily by balancing metabolic heat production with heat exchanged with the ambient [1].
Heat is exchanged at the interface between body surface and the ambient [1], and under standard experimental settings, the differential between surface vasoconstriction is a major physiological mechanism of vasomotor control, modulating heat exchange with the ambient [1]. Whereas there is experimental evidence for an association between nonlinearity or linearity and control of heat exchange or lack thereof [2], a theoretical justification is lacking. Recently, Boldrini et al. [5] demonstrated mathematically, based on physical first principles, that linearity or nonlinerity is expected under realistic assumptions of the thermoregulatory process and relevant thermal feedback control mechanisms. Nonlinearity or linearity thus appears to be signatures of vasomotor control. Furthermore, it is well-established that the body surface of an endotherm is not homogeneous with respect to vasomotor control [1]. Detection of signatures of cutaneous vasoreaction across the body surface is fundamental not only for understanding the dynamics of heat control, but also to quantify the relative importance of different body regions for global heat exchange in endotherms [6] [7] [8] [9].
Body surface temperature is accurately measured using infrared thermal imaging (IRT) either at the large scale of whole organisms [10] or at local, small scales such as that of the cornea [11]. Infrared thermal measurements taken at sensible gradients of ambient temperatures allow to estimate whether body surface temperature varies linearly or nonlinearly with ambient temperature, and, therefore, to determine the existence of heat exchange control [6]. The same thermographic data can be used in connection with biophysical models to estimate how different body regions contribute to total heat balance [6]. Such estimates are instrumental to assess the dynamics of heat control and exchange in endotherms [6]. Infrared thermal data are nevertheless intrinsically prone to within-individual correlations, missing values, time-varying covariates, and irregularly-timed measurements. Additionally, there is a need to discriminate between linear or nonlinear functional relationships between surface body temperature and ambient temperature. Therefore, to infer physiological process from patterns in thermographic data, we need reliable statistical methods. Here, we advance a three-pronged framework to test for signatures of linear or nonlinear relationships between body surface temperature and ambient temperature, and to quantify the relative importance of different surface body regions for the maintenance of total heat balance. First, we use polynomial functions that capture linearity and nonlinearity of the functional relationship between surface body temperature and ambient temperature. Second, we use Generalized Estimating Equations [12], a modeling formalism particularly suited to handle the idiosyncrasies of thermographic data, to estimate the relevant parameters. Third, we assess model parsimony, that is, the compromise between model complexity and residual variance, using a bootstrap-based strategy [13]. We demonstrate  [15]. Trapping and handling methods followed the guidelines of the American Society of Mammalogists (Animal Care and Use Committee 1998).

Experimental Procedures
Changes in radiative heat exchange at the surface of different body regions of G. microtarsus were monitored by video thermography as ambient temperature was varied from 8˚C to 38˚C. On the day of an experiment, each individual animal was transferred from its home-cage to an experimental container (7.5 cm × 9.0 cm × 23.5 cm). This container was framed with wire mesh (1 × 10 cm), which allowed the monitoring of body temperature surface without any considerable interference. The animal container was then placed inside a temperature-controlled chamber (FANEM Ltd., São Paulo) and the thermal camera (FLIR SC640, FLIR Systems, Inc.) was positioned under the animal cage, allowing us to track body surface temperature while manipulating ambient temperature. First, the individual was allowed one hour to habituate to the experimental conditions at 23˚C. We set the temperature change protocol for a 1˚C stepwise increment or decrement from the initial temperature (23˚C) up to 38˚C or down to 8˚C, respectively. This range of a T covers the range of a T the gracile mouse opossum is commonly exposed in its natural habitat ( Figure 2).
Upon each step change in a T , animals were kept for 12 min at the target temperature before the next step change. Previous tests indicated that this time interval assured that animals had reached a steady-state with a T . Incrementing and decrementing temperature protocols were separated by an interval of at least  4 days, in which animals were returned to their home-cages and allowed to feed and drink. Whether individual animals were first submitted to the incrementing and decrementing temperature protocols was a random decision. Thermal images were recorded at 10 frames•s −1 using an IR camera streamed to a computer where data acquisition was managed by the software Thermacam 2.9 (FLIR Systems, Inc). Thermal imaging cameras measure the amount of near infrared radiation (typically wavelengths: 8 -12 nm) emitted by a surface and then convert this measured radiation to a radiative temperature reading where,  is the emissivity of the surface and σ is the Stefan-Boltzmann constant (5.67 × 10 -8 W•m −2 •K −4 ). We assumed that the body surface of G. microtarsus has an homogenous emissivity of 0.96, which is typical of organic materials [8]. Surface temperatures were visualised as images in grey or colour scales [6].
We analyzed surface temperature ( s T ) separately for different body regions including ears (split into distal and proximal regions), feet, tail, chest, back (evaluated when the animal turned on its back) and ventral areas. Initially, we chose one frame for each a T from the last two minutes of the exposure period to any given a T . Thus, at the time used for data collection, animals had already been exposed to that particular a T for at least 10 min. For each chosen frame, the regions of interest were digitally drawn to obtain the average surface temperature of each body part. Typically, we used the same frame to analyze the s T from all body parts. In cases in which the animal positioning did not allow for the analysis of all the body regions on the same frame, we used the temporally closest frame to analyze the missing body region. T with a family of polynomial functions.

Model Estimation
A response variable that is repeatedly measured on the same subject at different time points is the key feature of longitudinal data sets. In this setting, the correlation between observations from a given individual must be accounted for.
Otherwise, downstream analyses may be affected by a number of factors, including false conclusions due to underestimated variance terms. Generalized estimating equation (GEE) models are an extension of generalized linear models devised to analyze data, which arise commonly in applied sciences [12]. Typically, in such cases, the data are collected in clusters in which observations within a cluster tend to be correlated, whereas observations in separate clusters are independent [12]. Clustered data arise from longitudinal or panel data, family studies or studies with spatial structures [12]. GEE models are used to test hypotheses about the dependent variables that are not necessarily normally distributed. The GEE framework is based on the concept of quasi-likelihood, which requires only an assumed relationship between the expected value (first moment) of the dependent variable and the covariates and between the conditional mean and variance (first and second moments) of the dependent variable. This approach uses the mean and variance of the response variable to derive the quasi-likelihood and its estimating equations [20] [21]. Because only the first two moments are needed, this approach is called semi-parametric. GEE allows for different choices of correlation structures, namely, independent, exchangeable, autoregressive, and unstructured. The method yields asymptotically unbiased and consistent estimates even if the incorrect working correlation structure is chosen. GEE accounts for within-individual correlations, allows for missing data, time-varying covariates, irregularly-timed or mistimed measurements. Fitting models using GEE requires the definition of the link function (which associates the linear predictor to the mean), an assumed distribution for the response variable, and a correlation structure of the response variable, often referred to as the working correlation matrix [12].  (Figure 3).
We estimated heat exchange (Q) in Watts for each body surface ( ears   Table 1).
where k is the thermal conductivity of the air at a particular a T (W•m −1 •˚K).
The relationship between k and a T was estimated by Tattersall et al. [2] as, , , p c Pr k µ = where a is the coefficient of thermal expansion of air, g is the acceleration of gravity, d is the characteristic dimension of the body region (Table 1), v is the kinematic viscosity of air, p c is the specific heat of air, µ is the dynamic viscosity of air, and k is the thermal conductivity of air. Experimental evidence shows that the free-convection function for heat exchange should be used when 2 16

Gr
Re > and the forced-convection function when 2 0.1 Gr Re < [24].  convection considering a wind velocity of 0.1 ms −1 according to Tattersall et al. [2], Gates [24], and Greenberg et al. [9], to obtain the values for Gr and Re. The relationship between Nu and Gr, Re, and Pr has been determined empirically for a range of geometric shapes and flow regimes (see Table 1 for the relationships for each body region). According to Monteith and Unsworth [23] and Gates [24] we assumed, Assuming that each gracile mouse opossum was in thermal balance during the experimental procedure we should expect that the relationship between total Q and a T would follow the relationship between metabolic heat production and a T . Based on the rates for oxygen uptake determined by Cooper et al. [18], we estimated the total metabolic heat production (Watts) for a gracile mouse opossum with a body mass of 40 g at 12 C a T =  , 20˚C, 28˚C, 30˚C, and 32˚C.
Thereafter, we compared these rates with the total heat exchange rate (Watts) estimated by the biophysical model ( total Q ).

Results
Our primary objective is to search for signatures of linearity or non-linearity in the relationship between the differential Within the range of a T from 8˚C to 38˚C, surface temperatures ( s T ) for all body regions of the gracile mouse opossum increased with a T at increasing rates depending on the body region ( Figure 4). For furred body regions such as the back, chest, and ventral region, we observed a linear relationship between T ∆ and a T . Conversely, for furless body regions such as the ears, feet, and tail, the relationship between T ∆ and a T was nonlinear ( Table 2)    The comparison between the rates of heat production estimated from Cooper et al. [18] and the total heat exchange (present study) for the gracile mouse opossum ( Figure 5) revealed that our estimates were consistently lower than the rates measured by Cooper et al. [18]. At  Heat exchanged by each body region varied across the range of ambient temperature as indicated by the significant interactions between i Q and a T (Table 3 and Figure 6). However, the pattern of variation was considerably different between furred and furless body regions. Furred regions such as the chest, dorsal, and ventral regions exchanged most of the heat, from 50% to 75%, dissipated by the gracile mouse opossum at lower temperatures. The furred ears

Discussion
Heat exchange in the surface of an endotherm is a process that can be modeled by the equation Therefore, in order to understand the heat exchange process as a function of ambient temperature across non-uniform surfaces in endotherms one must: 1) measure s T as function of a T across the surface of the organism; 2) assess the relationship between s T and a T across the surface of the organism; 3) understand how s T and a T interact with each body region to determine the role of each body region for total heat exchange at each a T [ compared to the dorsal region. Together, these differences indicate a better insulation of the dorsal region compared to the chest and ventral regions and, possibly, differences in the capacity to modulate body insulation by piloerection.
These interpretations also supported by the results on heat exchange.
For the range of temperatures tested, temperature equalization (i.e., for the furred body regions occurred only at the highest temperatures tested (differently from what happened with the furless body regions). At these temperatures (>35˚C), ambient temperature most likely had already surpassed internal body temperature. Therefore, we suspect that, in this case, fur insulation of the animals reduced heat gain from the environment at greater rates, shifting the point of temperature equalization to a temperature a few degrees higher than their normal body temperature. As at these temperatures energy expensive heat loss mechanisms are likely to intervene, preventing extra heat gain from the ambient might be highly relevant for the total heat balance and energetics of G. microtarsus. Finally, we should acknowledge the possibility that the atainement of equalization at higher than expected temperatures (normal b T ) could be attributable to hyperthermia. We did not monitored internal body temperature of our animals during the experiments. Nevertheless, the results from temperature equalization obtained for the furless regions, particularly the ear, suggest that the animals did not experience hyperthermia. Thus, we believe that the linear relationship between T ∆ and a T for the furred regions of G. microtarsus reveals not only the insulative properties of the fur in preventing excessive heat loss at low temperatures, but also the less commonly acknowledged role of preventing excessive heat gain at higher temperatures.
Furless body regions such as the ears, feet, and tail exhibited a nonlinear relationship between T ∆ and a T . However, there were differences among different furless regions. The relationship found for the ears was quite interesting.
For a T higher than 35˚C, T ∆ tended to zero and remained unchanged. Our interpretation of this result is that around 35˚C a T equalized with body temperature [18] [27], which induced a vasomotor response to limit of blood flow to the ears, thus preventing more heat to be gained through this region.
Interestingly, below 35˚C and even at the lowest temperature tested, we found no indication of vasoconstriction for the ears, which could have been expected as a heat conservation strategy under cooler conditions. The dynamics of heat exchange variation in response to changes in a T for the ears was quite different from that observed for the other furless regions, such as the feet and tail. Feet and tail both exhibited clear vasomotor responses (Figure 4). For feet, the threshold for vasodilation was detected at 20˚C and increased until around 30˚C, declining  [36], and ungulate horns and antlers [37] [38].
The tail, in particular, is commonly acknowledged as an important organ to modulate radiative heat dissipation via vasomotor adjustments in rats [39] [40] [41]. However, the potential for controlling heat exchange rate through the tail is probably higher in rodents than in G. microtarsus, which has semi-arboreal habits and uses its tail for climbing [14] [42]. We suspect that conflicting functional demands between thermoregulation and climbing may explain the difference in the potential to control heat exchange through the tail in rodents and the gracile mouse opossum. A similar situation seems to occur with young toco-toucans, in which the vasomotor control of heat exchange through the bill is thought to conflict with a demand associated with bill growth [2]. We are not aware of other studies focusing on the heat exchange properties of other small marsupials and/or bearers of prehensile tails. Therefore, it remains uncertain whether a functional trade-off in tail use may occur in G. microtarsus or whether it reflects a more general difference between different mammalian groups.
Our estimates of total heat exchange for the gracile mouse opossum were consistently lower than the rates of heat production estimated by Cooper et al. [18] using metabolic measurements ( Figure 6). Assuming that the gracile mouse opossum reached a thermal steady state at each a T tested, this comparison shows that our estimates of heat dissipation of gracile mouse opossum using biophysical modeling are good, given that they balance its internal heat production.
The difference observed between the two curves might be attributable to mechanisms of heat exchange, such as evaporative heat exchange, which were not considered in our model.

Conclusion
In conclusion, the process of heat exchange across an endotherm surface as a function of ambient temperature is complex and depends on as many dimensions of the organism and its surrounding ambient. Assessing this process under simplified conditions within the unifying principles of thermoregulation unveils and clarifies the rich diversity of mechanisms underlying the process. Under this approach, one is able to better describe how the dynamics of heat exchange mechanisms vary between body regions, allowing endotherms to thermoregulate across ambient temperatures. This approach is based on detailed understanding of the relative importance of each body region for the entire budget across the range of ambient temperature.