## Abstract

Metabolic rate underlies a wide range of phenomena from cellular dynamics to ecosystem structure and function. Models seeking to statistically explain variation in metabolic rate across vertebrates are largely based on body size and temperature. Unexpectedly, these models overlook variation in the size of gills and lungs that acquire the oxygen needed to fuel aerobic processes. Here, we assess the importance of respiratory surface area in explaining patterns of metabolic rate across the vertebrate tree of life using a novel phylogenetic Bayesian multilevel modeling framework coupled with a species-paired dataset of metabolic rate and respiratory surface area. We reveal that respiratory surface area explains twice as much variation in metabolic rate, compared to temperature, across the vertebrate tree of life. Understanding the combination of oxygen acquisition and transport provides opportunity to understand the evolutionary history of metabolic rate and improve models that quantify the impacts of climate change.

## INTRODUCTION

The power of the metabolic theory of ecology (MTE) is that it uses metabolism to explain and predict phenomena at population, community, and ecosystem scales (*1*). In this theory, organismal metabolic rate is mathematically connected to broader ecosystem attributes through its dependence on body mass and temperature (*1*, *2*). While the mechanism surrounding the body mass component of the MTE continues to be debated (i.e., the fractal distribution network), the mathematical relationship has proven useful (*3*–*5*). However, this relationship seeks to provide only a “zeroth-order” approximation; even after accounting for body mass and temperature, a considerable amount of variation in metabolic rate across species still remains to be explained statistically (*1*, *2*). Specifically, metabolic rate for organisms of the same body mass varies over five orders of magnitude, after accounting for temperature (*1*, *2*). Although the MTE acknowledges that exchange surfaces are important in metabolic scaling, the nature of these surfaces is rarely elaborated upon (*6*). One particular trait that may explain variation in the scaling of metabolic rate is the surface area of the respiratory system. Many have long recognized the importance of these respiratory surfaces to metabolism, for example, as codified in Fick’s law of diffusion (*7*–*9*).

Respiratory organs—lungs and gills—comprise the exchange surfaces that are used to acquire oxygen from the external environment, which is subsequently distributed throughout the body via the circulatory system (*10*). Two lines of inference have shown that metabolic rate and respiratory surface area are highly intertwined both within and across species—experimental manipulations and allometric comparisons [i.e., comparing body mass–scaling exponents; (*8*, *11*, *12*)]. First, experiments on rainbow trout (*Oncorhynchus mykiss*) and other organisms reveal that the physical reduction or blockage of respiratory surface area results in concomitant reductions in oxygen uptake and metabolic scope (*11*, *13*). Second, allometric inference has revealed that ontogenetic body mass scaling exponents for metabolic rate and respiratory surface area are often similar when compared within and across species (*14*, *15*). The same pattern holds when evolutionary body mass scaling exponents (i.e., estimated across different species that differ in size) are compared (*12*, *16*). A recent study found that the body mass scaling exponent of oxygen diffusion capacity (combined area and thickness) of the respiratory surfaces matches the body mass scaling exponent of metabolic rate (as measured by oxygen consumption) across differing subsets of vertebrate species (*12*). However, our understanding of the intimate relationship between metabolic rate and respiratory surface area both within and across species is largely limited to these experimental manipulations and comparisons of body mass scaling exponents. There has not yet been a robust test of whether respiratory surface area explains variation in the scaling of metabolic rate across vertebrates, beyond what can be accounted for by body mass, temperature, thermoregulatory strategy, and evolutionary history. The lack of an adequate test likely stems from the profound analytical challenges as both metabolic rate and respiratory surface area are almost never measured at the same body mass in the same species.

Here, we ask whether respiratory surface area explains additional variation in the scaling of metabolic rate across the vertebrate tree of life. To do so, we first compile a dataset with paired species’ estimates of metabolic rate and respiratory surface area that includes all major vertebrate lineages—fishes, amphibians, reptiles, birds, and mammals. Such species-paired datasets have enabled a breakthrough in our understanding of the metabolic basis of species’ responses to climate change [e.g., (*17*)]. Second, to solve the problem that traits are often measured at mismatched body sizes—an unresolved issue in many macroecological analyses, we develop a phylogenetic Bayesian multilevel modeling framework. The first level of this model estimates the residual effect of respiratory surface area when regressed against the body mass associated with respiratory surface area. The second level then examines whether residual respiratory surface area explains significant variation in the scaling of metabolic rate, while simultaneously accounting for the additional effects of temperature, thermoregulatory strategy, and evolutionary history. A strength of our quantitative framework is that it propagates uncertainty across levels of the model as each iteration happens in succession. Last, we examine the differences in the scaling relationships of metabolic rate and respiratory surface area between species that vary in thermoregulatory strategy (i.e., endotherms versus ectotherms), as well as the type of respiratory organ (i.e., lungs versus gills).

## RESULTS

We compiled a dataset of metabolic rate, respiratory surface area, body mass measurements for both metabolic rate and respiratory surface area, and the temperature associated with metabolic rate for 109 species from all major vertebrate lineages: 8 chondrichthyan and 63 teleost fishes, 10 amphibians, 4 reptiles, 6 birds, and 18 mammals. To our knowledge, this is the first extensive vertebrate-wide paired species dataset containing all species that have published estimates for both metabolic rate and respiratory surface area.

In compiling this dataset, we found that metabolic rates and respiratory surface areas have rarely been measured for individuals of the same body mass in the same species, complicating comparison of mean trait values (Fig. 1). There were only three species with both traits measured at the same body mass (Fig. 1A). The mean body masses for metabolic rate and respiratory surface area differed by more than a tolerable amount (10%) for most (85%) species (*n* = 93/109; Fig. 1A). Furthermore, for approximately one-third of species, the mean body masses for both traits differed by over an order of magnitude (*n* = 34/109; Fig. 1, A to C). Macroecology is founded on analyses of endothermic birds and mammals that grow little after fledging or weaning (“determinate growers”). However, generalizing these types of analyses to include ectotherms, resulting in fully comparative vertebrate-wide analyses, poses a problem as this group of vertebrates generally grows throughout life (i.e., “indeterminate growers”). Almost all (84 of 85) ectothermic species in our dataset had size-mismatched traits, with 34 of these species (40%) having a mean body mass mismatch greater than an order of magnitude (Fig. 1, A to C). To overcome this mismatch in body mass for metabolic rate and respiratory surface area, we developed a Bayesian multilevel analytical framework that enabled a vertebrate-wide comparison of multiple size-dependent phenomena (metabolic rate and respiratory surface area) while simultaneously accounting for additional covariates (e.g., body mass, temperature, thermoregulatory strategy, and evolutionary history).

### Does respiratory surface area statistically explain variation in metabolic rate across vertebrates?

Our results show that the surface area of lungs and gills explains substantial variation in metabolic rate across the vertebrate tree of life. First, species with greater respiratory surface areas had higher metabolic rates (Fig. 2). This was exemplified by organisms of the same body mass—species that had higher relative respiratory surface area (i.e., residual respiratory surface area) had higher metabolic rates (both observed metabolic rate and fitted metabolic rate values estimated by the model), even after differences in thermoregulatory strategy were accounted for (Fig. 2). For example, the body mass of the endothermic kowari *Dasyuroides byrnei* (a rat-like marsupial) was nearly identical to that of the ectothermic white sucker *Catostomus commersonii* (a teleost fish), yet the kowari had ~32 times greater relative respiratory surface area and ~16 times greater metabolic rate compared with the white sucker (Fig. 2, orange and purple lines). Second, the addition of respiratory surface area consistently improved our explanatory models of metabolic rate across vertebrates (compare looic and elpd_{loo} for all “MR” models and all “C” models; table S1). Third, the addition of respiratory surface area was significant in all six models that included it as a covariate [95% Bayesian credible interval (BCI) of the effect sizes for respiratory surface area did not include zero; table S2, column “residual RSA”]. Fourth, evidence ratios (i.e., the weight of evidence of one model divided by that of another) show that including respiratory surface area to explain variation in metabolic rate was, on average, 18.5 times more likely than excluding respiratory surface area, after accounting for body mass, temperature, thermoregulatory strategy, and evolutionary relatedness (this evidence ratio ranged from 12.3 to 22.3 according to model run; table S3). Fifth, the standardized effect size of residual respiratory surface area was twice as large as that of temperature, indicating that respiratory surface area is twice as important in explaining variation in metabolic rate across vertebrates compared with temperature (Fig. 3; comparing the absolute value of standardized effect sizes of residual respiratory surface area and temperature in table S4). Collectively, these results show that respiratory surface area explains substantial variation in metabolic rate even after accounting for body mass, thermoregulatory strategy, temperature, and the evolutionary relatedness among species.

### Is respiratory surface area simply a recasting of the known difference in metabolic rates between endotherms and ectotherms?

We know empirically that ectotherms have lower metabolic rates for a given size than endotherms, which retain metabolically produced heat to maintain their body temperature within a narrow thermal range. However, it is unlikely that thermoregulatory strategy alone explains the observed variation in metabolic rate that exists after body mass and temperature have been accounted for. First, the inclusion of respiratory surface area in models explaining variation in metabolic rate substantially improved the fit of the model, even after accounting for thermoregulatory strategy (see previous section). Second, the models that included respiratory organ (i.e., lungs versus gills) in place of thermoregulatory strategy (i.e., endotherm versus ectotherm) provided a poor fit to the data (table S5). Third, if respiratory surface area and thermoregulatory strategy were interchangeable in explaining the same variance in the scaling of metabolic rate across vertebrates, then we would expect to see similar body mass scaling relationships of metabolic rate and respiratory surface area across all species, regardless of thermoregulatory strategy. However, we see a mismatch in the body mass scaling of metabolic rate and respiratory surface area for endotherms (Fig. 4). For endotherms, the mean body mass scaling exponent (i.e., allometric slope) of metabolic rate was shallower than the mean body mass scaling exponent of respiratory surface area, although the 95% BCIs marginally overlapped [compare Fig. 4, C and D, and body mass scaling exponents (and their 95% BCIs) for endotherms from models “MR3” and “RSA3” in table S2]. In contrast, the body mass scaling exponent of metabolic rate and respiratory surface area was nearly identical for ectotherms [compare Fig. 4, E and F, and body mass scaling exponents (and their 95% BCIs) for ectotherms from models MR3 and RSA3 in table S2]. This mismatch in scaling for metabolic rate and respiratory surface area for endotherms persisted even when respiratory surface area was included in the model; the body mass scaling exponent for metabolic rate was still shallower than that of respiratory surface area [compare metabolic rate and respiratory surface area body mass scaling exponents (and their 95% BCIs) for endotherms and ectotherms from models “C5” and RSA3; table S2]. Together, these results suggest that respiratory surface area is not simply a recasting of thermoregulatory strategy.

### Is respiratory organ (i.e., lungs versus gills) a better characterization of the known difference in metabolic rate and respiratory surface area between endotherms and ectotherms?

The difference in the type of respiratory organ—having lungs or gills—does not explain the differences in metabolic rate and respiratory surface area between endotherms and ectotherms. Specifically, using thermoregulatory strategy (endotherm versus ectotherm) as a covariate instead of the type of respiratory organ (lungs versus gills) provided a far better fit for all models (compare the looic of models with thermoregulatory strategy to those with respiratory organ instead of thermoregulatory strategy in table S5). As such, the characterization of the differences in respiratory surface area and metabolic rate between endotherms and ectotherms is far better explained by thermoregulatory strategy than by whether an organism has lungs versus gills (tables S5 and S6). See the Supplementary Materials for further results of the respiratory organ analyses.

## DISCUSSION

We have shown here that respiratory surface area plays a critical role in understanding variation in metabolic rate across the vertebrate tree of life. This is supported by two main findings. First, respiratory surface area substantially improved our ability to explain variation in metabolic rate across 109 vertebrate species from all major lineages, while simultaneously accounting for differences in body mass, temperature, thermoregulatory strategy, and evolutionary relatedness. We found that respiratory surface area was twice as important in explaining variation in metabolic rate compared to temperature. Second, we confirmed that respiratory surface area was not simply a recasting of the differences in metabolic rate between endotherms and ectotherms. Answering these questions was only possible because of our paired dataset in which each species had estimates of both respiratory surface area and metabolic rate, as well as a novel Bayesian multilevel modeling approach that propagates uncertainty in the effect of body mass on respiratory surface area to all levels of the model. This modeling framework offers a breakthrough in dealing with multiple size-dependent phenomena while accounting for evolutionary relatedness and can be applied to many types of comparative questions. Our paired dataset and modeling framework allowed us to extend the mathematical framework of the MTE by examining whether additional size-dependent phenomena—here, respiratory surface area—explain variation in metabolic rate across species. Together, our results show that respiratory surface area, in addition to body mass, temperature, and thermoregulatory strategy, underpins the scaling of metabolic rate across vertebrates. We focus our discussion on three key issues: (i) the importance of respiratory surface area and oxygen uptake in ecological and physiological phenomena, (ii) the differences in the body mass scaling of metabolic rate and respiratory surface area between endotherms and ectotherms, and (iii) the limitations of modeling studies in uncovering mechanistic relationships. Last, we lay out a research agenda to further dissect the relationship between metabolic rate and respiratory surface area.

Respiratory surface area appears to play a central role in several ecological and physiological phenomena including symmorphosis and oxygen limitation [including the temperature-size rule (TSR) and the gill oxygen limitation theory (GOLT)]. First, symmorphosis is the hypothesis that organismal structures (e.g., respiratory surface area) are perfectly matched to their function (e.g., acquiring oxygen to meet metabolic demand) (*18*, *19*). While some work has found that the respiratory system appears to be “overdesigned” for the function of acquiring oxygen, others have found that the body mass scaling of resting metabolic rate and respiratory surface area are closely matched across a broad size range of vertebrates (*12*, *18*, *19*). In our study, we found that the body mass scaling of metabolic rate and respiratory surface area matched closely for ectotherms, but not endotherms, suggesting the potential importance of additional traits in sculpting this relationship. Direct tests of symmorphosis would ideally be conducted within and not across species and using maximum rather than resting metabolic rate. In addition, it is important to recognize that oxygen diffusion across the respiratory surface is only one step in a series of steps involved in the acquisition of oxygen for aerobic metabolism (“oxygen cascade”). Many other steps—including oxygen binding to hemoglobin, oxygen delivery to the tissues through the circulatory system, and the density of mitochondria (the final oxygen receptor)—must be considered in a direct test of symmorphosis (*18*). Second, oxygen limitation is the idea that geometric and physiological constraints on oxygen supply will affect aerobic metabolism, particularly for larger organisms or those in warmer waters (*20*, *21*). This phenomenon is one of the proposed explanations that are thought to underlie the widespread inverse relationship between rearing temperature and ectothermic body size (*20*, *22*, *23*). Specifically, the smaller maximum size and faster growth rate observed under warmer temperatures are thought to be due to the difficulty in obtaining oxygen as temperature increases due to higher metabolic demand and decreased oxygen availability, particularly for aquatic ectotherms [e.g., (*22*)]. The GOLT proposes that respiratory surface area limits metabolic rate in fishes and other water-breathing organisms because an individual’s respiratory surface area (gill surface area) cannot grow as fast as the body mass it must supply with oxygen (i.e., a hypoallometric ontogenetic scaling of respiratory surface area) (*24*, *25*). This theory—while largely empirically untested—further predicts that respiratory surface area in fishes may be related to several metabolism-related phenomena, such as the “shrinking” of fish body size with climate warming (*24*–*26*). However, the GOLT is based on an allometric relationship (the ontogenetic scaling of gill surface area with body mass for a fish) and, as such, cannot be used to determine mechanism by itself [see discussion below; (*24*, *25*)]. While the role of oxygen in the physiology, ecology, and evolution of organisms is debated, broad, cross-species studies have shown that oxygen may shape marine species’ geographic distributions and affect the relationship among metabolic rate, body mass, and temperature in fishes (*21*, *27*). However, many within-species studies show a much more complicated relationship between oxygen acquisition, distribution, and use (*28*, *29*). While our results show that respiratory surface area substantially improves our ability to explain variation in metabolic rate across species, further experimental and modeling work—especially work that is able to incorporate variation across evolutionary time scales (i.e., selection experiments and additional cross-species analyses)—is needed to assess whether the diffusion of oxygen via the respiratory structures is a valid mechanism that underlies the GOLT and TSR. A coordinated effort among organismal physiologists, macrophysiologists, and comparative evolutionary ecologists would greatly enhance our ability to understand the role that oxygen plays in ecological and physiological phenomena, both within and across species.

We found that although respiratory surface area vastly improved our understanding of metabolic rate across both endotherms and ectotherms, ectothermic organisms had a tighter coupling of the scaling of metabolic rate and respiratory surface area than that observed for endotherms. These differences in scaling of metabolic rate and respiratory surface area were not explained by the type of respiratory organ itself (i.e., lungs versus gills), as all models that included respiratory organ in place of thermoregulatory strategy fit the data less well. Instead, our results suggest that attributes related to thermoregulatory strategy likely underlie the differences in the relationship of metabolic rate and body mass between endotherms and ectotherms. For example, the body mass scaling exponent of metabolic rate found here for endothermic organisms may support the heat dissipation theory, which suggests that there is an upper limit to metabolic rate in endothermic organisms (*30*). Endothermic organisms maintain their body temperature within a target range and have to dissipate excess heat produced by metabolism across their body surface area. Thus, endothermic organisms must balance heat production and heat loss, which is constrained by the body mass–to–surface area ratio (*30*). As organisms increase in size, the ratio of body surface area to body mass decreases, resulting in a decreased heat dissipation capacity (*30*, *31*). The heat dissipation theory suggests that the body mass scaling exponent of field metabolic rate for endothermic organisms will not significantly differ from ~0.63 to 0.67, following surface area–to–volume geometry [field metabolic rate is a measure of energy expenditure in a free-living organism; (*30*)]. While the mean body mass scaling exponent for resting metabolic rate—both with and without respiratory surface area—for endotherms found in this study was higher than 0.63 to 0.67, the 95% BCIs of both models included the 0.63 value [models MR3 and C5 in table S2; these intervals also included the predicted ¾ slope of the MTE; (*1*)]. Because ectothermic organisms do not retain metabolically produced heat, dissipation is not an issue, and hence, this may explain the steeper body mass scaling exponent of metabolic rate in ectotherms. In addition, some work has shown that the evolutionary body mass scaling exponent of maximum metabolic rate, and not resting metabolic rate as used here, is more similar to the evolutionary body mass scaling exponent of respiratory surface area [e.g., (*16*)]. However, we found a match in the body mass scaling of resting metabolic rate and respiratory surface area for ectotherms and not endotherms (this is also an evolutionary allometry), and thus, our examination of resting metabolic rate versus maximum metabolic rate cannot explain the observed difference in body mass scaling of metabolic rate and respiratory surface area in endotherms.

We provide compelling evidence that—to a first approximation—respiratory surface area, in addition to body size and temperature, explains significant variation in metabolic rate across vertebrates. Yet, we have much to learn about the causal relationships between metabolic rate and respiratory surface area. Correlative or scaling studies such as ours serve to identify broad, general patterns, which can then inspire other studies that aim to understand the underlying or driving mechanisms (e.g., experimental or selection studies). While our results show that respiratory surface area (in addition to body mass, temperature, and thermoregulatory strategy) underlies patterns of metabolic rate across vertebrates, we cannot say from our results—or other scaling studies—whether organismal metabolic rate constrains or shapes organismal respiratory surface area or vice versa [e.g., (*12*)]. A major step forward in understanding the mechanistic relationship between organismal metabolic rate and organismal respiratory surface area would be to understand the relationships among ontogenetic allometries (i.e., within an individual of a single species across its lifetime or, for traits that require lethal sampling, across individuals of the same species that span the size range of the species), static allometries (i.e., across individuals of the same species of the same life stage), and evolutionary allometries (i.e., across different species that differ in size) of both metabolic rate and respiratory surface area. For example, a recent study examining the relationships among ontogenetic, static, and evolutionary brain and body size allometries suggested that developmental constraints governed scaling relationships within and across species rather than geometric/physical constraints or physiological mechanisms (*32*). What are the constraints and causal mechanisms that underlie the relationship between metabolic rate and respiratory surface area within and across species?

To this end, we outline six specific avenues of research that would help us to understand causality between organismal metabolic rate and organismal respiratory surface area. First, the advance in phylogenetic methods has opened the door to comparing evolutionary transitions of metabolic rate to transitions in respiratory mode (*33*). Second, common-garden and long-term selection experiments, particularly of aquatic organisms, offer the opportunity to understand the phenotypic and genotypic response of organismal metabolic rate and organismal respiratory surface area to food availability, temperature, and oxygen [e.g., (*23*)]. Third, a deeper understanding of allometries—including the relationships examined in this study—has been profoundly hindered by a lack of available estimates of individual (i.e., raw) data for metabolic rate, respiratory surface area, and other traits. We urge experimental scientists to publish their raw data alongside means and other summaries. This would allow the statistical propagation of uncertainty using the approach we have developed here, which can be easily modified to include data at both individual and species scales. In addition, this would also enhance datasets such as ours and facilitate the identification of patterns across broad groups of species. Fourth, activity, metabolic rate (both resting and maximum, as these two measures are correlated), respiratory surface area, and temperature are deeply intertwined (*16*, *24*, *34*). Are metabolic rate and respiratory surface area simply proxies for activity or are metabolic rate and respiratory surface area capturing total energy availability for growth and reproduction in addition to activity? In addition, which measures of metabolic rate (i.e., standard/basal, resting, routine, field, maximum) are suitable to test the interrelationships between metabolic rate, respiratory surface area, temperature, and activity? Last, there is an incredible, but widely overlooked, diversity of respiratory systems, modes, and types of ventilation. We were unable to tackle this diversity with our dataset beyond examining respiratory organ as a potential predictor of metabolic rate. Studies that explore ventilation types and the diversity of respiratory modes even within the coarse categorizations of respiratory organs (e.g., unidirectional flow of water across fish gills, unidirectional flow of air through bird lungs, and tidal air flow in mammalian and reptilian lungs) could begin to examine this question. Air breathing and cutaneous respiration in aquatic organisms and amphibians provide further contrasts to explore. Blood flow across these respiratory surfaces also differs (i.e., counter-current in fish gills, cross-current in bird lungs, etc.) providing another avenue for exploration. Combining more advanced modeling approaches such as the one presented here with detailed physiological and ecological data both within and across species will allow us to further understand the role that oxygen plays in the ecology, physiology, and evolution of organisms.

## MATERIALS AND METHODS

### Trait data

We compiled a species-paired dataset of vertebrates that had both metabolic rate and respiratory surface area data. To do so, we collated mean estimates of whole-organism aerobic metabolic rate (termed here, “metabolic rate”), as measured by oxygen consumption (mg O_{2} min^{−1}, mg O_{2} g^{−1} hour^{−1}, mg O_{2} kg^{−1} hour^{−1}, ml O_{2} hour^{−1}, ml O_{2} min^{−1}, ml O_{2} g^{−1} hour^{−1}, ml O_{2} kg^{−1} hour^{−1}, ml O_{2} kg^{−1} min^{−1}, joules hour^{−1}, and watts), body mass (grams or kilograms) associated with the metabolic rate estimates, temperature (°C) associated with the metabolic rate measurements, whole-organism respiratory surface area (cm^{2} or mm^{2}; termed here, “respiratory surface area”), and body mass (g or kg) associated with the respiratory surface area measurements for as many vertebrate species as possible. If raw data (i.e., measurements for multiple individuals of the same species) were available, then these estimates were averaged to generate a species mean. Our ability to incorporate raw data into our modeling framework was limited because most species in our dataset only had published mean estimates of metabolic rate and respiratory surface area (only nine species of the 109 vertebrates have published raw data for both metabolic rate and respiratory surface area).

Much of our data came from two existing datasets: metabolic rate data from (*35*) and respiratory surface area from (*12*). We searched the primary literature to fill gaps for species missing either metabolic rate or respiratory surface area estimates. If we found more than one estimate of either mean metabolic rate or mean respiratory surface area for a given species, then we included the value from the study with the larger sample size. Metabolic rate estimates are from individuals at rest (resting or standard for ectotherms, basal for endotherms), with the exception of four teleost species for which we could only find estimates of routine metabolic rate (oxygen consumption during volitional movement): *Anabas testudineus* climbing perch, *Brevoortia tyrannus* Atlantic menhaden, *Channichthys rhinoceratus* inicorn icefish, and *Hoplerythrinus unitaeniatus* trahira. These four species, specifically, and other species with routine metabolic rates are regularly included in metabolic allometry studies, as the variation of metabolic rate among individuals of the same species is substantially smaller than the variation across different species [e.g., (*36*)]. For the purposes of this study, the thermoregulatory strategy of five fish species that are regionally endothermic (*Carcharodon carcharias* white shark, *Euthynnus affinis* kawakawa, *Isurus oxyrinchus* shortfin mako, *Katsuwonus pelamis* skipjack tuna, and *Thunnus albacares* yellowfin tuna) was classified as ectotherms. However, rerunning the three top models (MR3, RSA2, and C5 in table S1) without the five regionally endothermic species did not significantly change any coefficient value (i.e., the effect size of any parameter in a model). For analyses, all estimates of metabolic rate were converted to watts, respiratory surface area to square centimeters, body mass for both metabolic rate and respiratory surface area to grams, and temperature to inverse temperature for model parameterization as the Boltzmann factor (see “Basic modeling framework and analysis” section). Metabolic rate, respiratory surface area, and both associated body masses were natural log transformed before analyses.

### Phylogeny

We included a phylogenetic random effect in all models that allowed for a phylogenetic signal among residuals (i.e., error). To do so, we first constructed a supertree from a database of molecular phylogenies, TimeTree (*37*), and a recently published molecular phylogeny for Chondrichthyans (*38*). As the evolutionary position of all species in our dataset has not yet been fully resolved, we opted to use a genera-level phylogeny for all species except the Chondrichthyans (as a phylogeny for this group was recently published). In the infrequent (*n* = 7) case that two species from the paired dataset were in the same genus, the branch length was split equally among those two species. This use of a genera-level tree with a few equally split branches to accommodate species from the same genus—as opposed to a tree with all species at the tips—will yield the same conclusion, as divergence times between species in the paired dataset are quite high across the phylogeny because of the number of species included in our dataset relative to all extant vertebrates.

### Modeling framework and statistical analysis

Metabolic rate and respiratory surface area are mass-dependent traits, meaning that they change as an individual grows and increases in size. However, both traits do not increase at the same rate as body mass (i.e., the body mass scaling exponent of an ontogenetic allometry for these traits does not equal one), and therefore, the body mass at which these traits were measured must be included in all models. Mass-dependent traits are typically examined in an allometric context using a power-law, or scaling, relationship such as*t* is the mass-dependent trait (in this case, either metabolic rate or respiratory surface area), β_{0} is the intercept (i.e., the value of *t* at a given body mass, often called the “normalization constant”), *M* is body mass, and β_{mass} is the body mass scaling exponent (i.e., allometric slope) (*39*). This equation is most often examined on a logarithmic scale, resulting in a linear relationship for log-transformed data*E _{i}* is the activation energy for the biochemical reactions of metabolism,

*k*is the Boltzmann constant (8.617 × 10

^{−5}eV), and

*T*is temperature in Kelvin (

*1*,

*2*). Temperature is parameterized as the Boltzmann factor (i.e., inverse temperature) for metabolic rate scaling relationships as it best approximates how temperature affects metabolic reactions (

*1*,

*2*). For endotherms, temperature data are body temperature, and for ectotherms, temperature data are the temperature at which metabolic rate was experimentally measured. This temperature dependence does not capture the fundamental differences in metabolic rates between endotherms and ectotherms as temperature has the same effect on the biochemical reactions of respiration for both groups (

*2*,

*12*).

All models were fit in a Bayesian framework in Stan with the rstan package in R v.3.5.1 and v.4.0.1 (*40*, *41*). To ensure our results were robust to model run, we ran each model a total of four times. The results of each additional model run (after the first one) are in table S7. We also ran all models without one possible outlier, but this did not significantly change any coefficient estimates. A detailed outline of all models, their parameterization, and choice of priors is included in the Supplementary Materials. The results of the models with respiratory organ (lungs versus gills) in place of thermoregulatory strategy are expanded upon in the Supplementary Materials. Below, we detail the parameterization of models specific to each research question.

### Does respiratory surface area explain variation in metabolic rate across vertebrates?

To assess whether respiratory surface area explains variation in metabolic rate across vertebrates, above and beyond that explained by the other covariates (e.g., body mass, temperature, and thermoregulatory strategy), we compared candidate models that described variation in metabolic rate with and without respiratory surface area (table S1). For models that examined variation in metabolic rate without respiratory surface area (“metabolic rate models”), we compared the classic MTE model to that with the addition of thermoregulatory strategy as a covariate (table S1). To do this, we examined models that allowed just the intercept to vary by thermoregulatory strategy (i.e., metabolic rate for a given body mass differed for endotherms and ectotherms) and models that allowed both the slope and intercept to vary by thermoregulatory strategy (i.e., metabolic rate for a given body mass differed for endotherms and ectotherms, and the effect of body mass varied between endotherms and ectotherms). In total, we parameterized three candidate models to examine the body mass scaling of metabolic rate without respiratory surface area (“MR” models; table S1). Second, we built on these three metabolic rate models above (the classic MTE model with and without thermoregulatory strategy) by adding in respiratory surface area as a covariate (“combined models”). To do this, we used a multilevel model where we first calculated the residual respiratory surface area by regressing respiratory surface area against the measurement body mass. We subsequently incorporated this residual respiratory surface area as a covariate in the next level, in addition to other covariates. This approach allows for the uncertainty in estimated residual respiratory surface area to be propagated across levels of the model, as opposed to simply including a mean estimate of residual respiratory surface area per species (which does not incorporate the uncertainty in that estimate). In total, we parameterized six candidate models to examine the body mass scaling of metabolic rate with respiratory surface area (“C” models; table S1).

We used model selection to identify a single best model that explains variation in metabolic rate without respiratory surface area (termed here, “best metabolic rate model”) and a single best model with respiratory surface area (termed here, “best combined model”). To do this, we used Pareto-smoothing leave-one-out cross validation (PSIS-LOO). This model selection framework is based on the predictive accuracy of a model, as estimated by iteratively leaving out one observation at a time and then predicting that observation based on the model fit to the remaining data (*42*). An assumption of using PSIS-LOO is that the joint likelihood of the model observed over all observations is factorizable or pairwise conditionally independent, given the model parameters (*42*, *43*, *44*). As phylogenetic models do not meet this assumption, we instead computed the pointwise log-likelihood for nonfactorizable models (*43*, *44*). We then used the loo package in R v 5.3.1 to estimate the expected log predictive density (elpd_{loo}), the LOO information criterion value (looic), the effective number of parameters (p_{loo}), the standard error of the expected log predictive density (se_{elpd_loo}), the difference in the expected log predictive density (elpd_{diff}) for a given model compared to the best model compared to the best model, and, lastly, the weight of evidence for each model as estimated by the Bayesian stacking method (*42*, *45*). The model with the lowest elpd_{loo} value is the best fit to the data. In addition, we used a *z* score standardization to standardize the predictors of the best combined model to identify and compare the relative importance of these predictors in explaining variation in metabolic rate across vertebrates [i.e., comparing standardized effect sizes; (*46*)]. We also computed evidence ratios to measure how much more likely one model is over the other(s). Evidence ratios are simply the weight of evidence of the best model divided by the weight of evidence of the other model(s) of interest (*47*).

### Is respiratory surface area simply a recasting of the difference in metabolic rate between endotherms and ectotherms?

To assess whether respiratory surface area is simply a recasting of the difference in metabolic rate between endotherms and ectotherms, we compared the scaling relationships of metabolic rate and respiratory surface area while accounting for differences in thermoregulatory strategy. First, we parameterized three candidate models (“respiratory surface area models”) to examine the body mass scaling of respiratory surface area (“RSA” models; table S1). We then selected the single best model (termed here, “best respiratory surface area model”) from these candidate models using PSIS-LOO (see above, table S1). Second, we compared the body mass scaling of the best metabolic rate model (model MR3; table S2) and the best respiratory surface area model (model “RSA2”; table S2) by comparing the 95% BCI.

### Is respiratory organ (i.e., lungs versus gills) a better characterization of the known difference in metabolic rate and respiratory surface area between endotherms and ectotherms?

To examine whether respiratory organ (i.e., lungs versus gills) was a better predictor of the differences in metabolic rates between endotherms and ectotherms instead of thermoregulatory strategy, we replaced thermoregulatory strategy in all models that included it with respiratory organ (table S5). We did this in favor of simply adding respiratory organ as covariate in addition to thermoregulatory strategy, which would not be feasible with our dataset as the ectothermic species in our dataset were largely fishes and thermoregulatory strategy is almost entirely correlated with respiratory organ (i.e., fishes have gills). Specifically, we examined the effect of lungs versus gills in the scaling of metabolic rate (models “MR2_LG” and “MR3_LG” in tables S5 and S6), the scaling of respiratory surface area (models “RSA2_LG” and “RSA3_LG” in tables S5 and S6), and how respiratory organ affected the scaling of metabolic rate with the effect of residual respiratory surface area included (models “C3_LG,” “C4_LG,” “C5_LG,”and “C6_LG”; tables S5 and S6). All model comparisons were conducted using PSIS-LOO (see above).

## SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/7/19/eabe5163/DC1

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is **not** for commercial advantage and provided the original work is properly cited.

## REFERENCES AND NOTES

**Acknowledgments:**We would like to thank D. Greenberg, C. Bedore, C. Mull, and S. Pardo for thoughtful discussions on the topic.

**Funding:**This project was funded by the Natural Sciences and Engineering Research Council of Canada and the Canada Research Chairs Program.

**Author contributions:**J.S.B. and N.K.D conceived of and designed the project and analysis. J.S.B. collected the data and performed all analyses and visualizations. L.K.M. contributed to and validated the analysis. All authors contributed to the interpretation of results. J.S.B. drafted the manuscript and supplementary information with input from all authors. N.K.D. supervised the project.

**Competing interests:**The authors declare that they have no competing interests.

**Data and materials availability:**All data and code necessary to reproduce the results in this study are archived in Figshare and are cited in the reference section (

*49*). We place no restrictions on data or code availability. All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

- Copyright © 2021 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution NonCommercial License 4.0 (CC BY-NC).