SUMMARY
Assessment of the relative importance of variation in enzyme concentration [E] and metabolic regulation in accounting for interspecific variation in metabolic rates is an unrealized area of research. Towards this end, we used metabolic flux rates during hovering and enzymatic flux capacities (V_{max} values, equal to [E]×k_{cat}, where k_{cat} is catalytic efficiency) in flight muscles measured in vitro from 14 orchid bee species ranging in body mass from 47 to 1065 mg. Previous studies revealed that, across orchid bee species, wingbeat frequencies and metabolic rates decline in parallel with increasing body mass. V_{max} values at some enzymatic steps in pathways of energy metabolism decline with increasing mass while, at most other steps, V_{max} values are massindependent. We quantified the relative importance of `hierarchical regulation' (alteration in V_{max}, indicative of alteration in [E]) and `metabolic regulation' (resulting from variation in substrate, product or modulator concentrations) in accounting for interspecific variation in flux across species. In addition, we applied the method of phylogenetically independent contrasts to remove the potentially confounding effects of phylogenetic relationships among species. In the evolution of orchid bees, hierarchical regulation completely accounts for allometric variation in flux rates at the hexokinase step while, at other reactions, variation in flux is completely accounted for by metabolic regulation. The predominant role played by metabolic regulation is examined at the phosphoglucoisomerase step using the Haldane relationship. We find that extremely small variation in the concentration ratio of [product]/[substrate] is enough to cause the observed interspecific variation in net flux at this reaction in glycolysis.
Introduction
Orchid bees (Apidea: Euglossini) belong to a diverse clade consisting of more than 190 species belonging to five genera that vary in body mass by about 20fold. Pioneering work conducted by Casey and colleagues (Casey and Ellington, 1990; Casey et al., 1992, 1985) showed that wingbeat frequencies and massspecific metabolic rates during flight decline with increasing body mass. Inspired by this work, we examined allometric scaling patterns in greater detail, and in the context of phylogenetic relationships among species (see accompanying papers, Darveau et al., 2005a,b). We found evidence supporting the idea that the allometric scaling of wingbeat frequency, one of the determinants of muscle volumespecific power output (Pennycuick and Rezende, 1984), drives the scaling of hovering metabolic rate (Darveau et al., 2005a). Orchid bees oxidize carbohydrate as the main energy source for flight (Suarez et al., 2005), so we also conducted interspecific studies of enzymatic flux capacities through carbohydrateoxidizing pathways and examined how these relate to metabolic flux rates during flight. We found that flux capacities (V_{max} values) at one step, i.e. hexokinase, scale against body mass with an allometric exponent almost identical to that of wingbeat frequency and metabolic rate. In contrast, V_{max} values at other steps either decline with lower allometric exponents than hexokinase or are massindependent, i.e. scale isometrically (Darveau et al., 2005b).
Within an individual bee, >90% of the whole body metabolic rate during hovering flight is due to O_{2} consumption by flight muscle mitochondria (Rothe and Nachtigall, 1989; Suarez, 2000). Therefore, rates of wholebody CO_{2} production (V̇_{CO2}) or O_{2} consumption (V̇_{O2}) can be used to estimate flux rates through catabolic pathways (Suarez et al., 1996) or through mitochondrial respiratory chain enzymes (Suarez et al., 1999, 2000). These fluxes represent steadystate rates whose control can be analyzed quantitatively, at least in principle, through metabolic control analysis (Fell, 1997). On the other hand, interspecific variation in flux can be considered within the framework proposed (but, originally, in an intraspecific context) by ter Kuile and Westerhoff (2001). At any step in metabolism, such variation in flux may be due to `hierarchical regulation', i.e. variation in enzyme concentration, [E], across species. An alternative is that [E] at a given step is constant across species and variation in flux at this step is due to variation in the concentrations of substrates, products or allosteric effectors, i.e. `metabolic regulation'. A third possibility is that variation in flux may result from a combination of hierarchical and metabolic regulation.
Materials and methods
Theory
We summarize the approach described by ter Kuile and Westerhoff (2001) and our application of it to metabolic scaling in orchid bees as follows:
The hierarchical regulation coefficient ρ_{h}, which is a measure of the contribution of enzyme concentration [E], and the metabolic regulation coefficient ρ_{m}, which is a measure of the contribution of mechanisms that result in the modulation of the activity of a constant concentration of enzyme, to variation in flux, are related as: (1)
At any enzymecatalyzed step i, ρ_{h} is a function of the relative change in rate, ∂v_{i}, divided by the relative change in enzyme concentration, ∂e_{i}, times the ratio of change in e_{i} to the change in flux J. Assuming that∂ lnv_{i}/∂lne_{i} equals 1 (ter Kuile and Westerhoff, 2001), (2)
Similarly, at any step i, metabolic regulation results from a change in an enzymecatalyzed rate divided by the change in concentration of its substrate, product or allosteric modulator (X), times the ratio of change in X to change in J, such that: (3)
Orthologous enzymes from organisms of similar body temperature display similar catalytic efficiencies (k_{cat} values) (Hochachka and Somero, 2002). Because (4)
it follows that interspecific variation in V_{max} measures variation in e_{i}. From Eq. 1, the slope of lnV_{max} plotted against lnJ yields the hierarchical regulation coefficient ρ_{h}. It follows that 1ρ_{h} yields the metabolic regulation coefficientρ _{m}.
We note here that this approach was not devised for interspecific comparisons; ter Kuile and Westerhoff (2001) applied it to the analysis of variation in glycolytic rates in protozoa resulting from experimental manipulation. However, the mechanistic bases for and the quantitative relationships that describe intraspecific variation in flux among protozoa apply equally well to interspecific variation in flux among orchid bees. In a given steadystate situation (in this case, high rates of carbohydrate oxidation during hovering flight), the interspecific variation in flux can be accounted for either by interspecific variation in [E] or by interspecific variation in substrate, product or modulator concentrations. However, the interspecific variation in flux may also be due to variation in other morphological, physiological and biochemical traits besides these. Thus, in the work described herein, the value of slope of lnV_{max} plotted against lnJ is potentially influenced by other traits that covary with the biochemical parameters that we measured.
We address the above problem with the use of Felsenstein's method of phylogenetically independent contrasts (PIC analysis; Felsenstein, 1985; Garland et al., 1992). Briefly, in conventional statistical analyses of interspecific data, it is assumed that species are equally related to each other. In reality, they are related to each other by their phylogenetic relationships. Because closely related species tend to display greater similarity than more distantly related ones, estimates of ρ_{h} from slopes of linear regressions may be biased by phylogenetic relationships and correlated traits other than those hypothesized to account for interspecific variation in, as in the present study, flux rates. PIC analysis provides a statistical approach that allows for phylogenetic nonindependence and makes possible the application of the approach described by ter Kuile and Westerhoff (2001) to the analysis of interspecific data. In the present work, PIC analysis was conducted using the PDAP (Midford et al., 2003) module in Mesquite (Maddison and Maddison, 2004). The analyses were performed using both gradual and speciational models of character evolution, with branchlengths obtained from genetic distances inferred from cytochrome b partial sequences (Darveau et al., 2005a) or with branch lengths set to 1.
Sources of data
Glycolytic flux rates were estimated from V̇_{CO2} during flight (Darveau et al., 2005a). Enzyme V_{max} values are from Darveau et al. (2005b). Hypothetical phylogenies used in PIC analysis were generated on the bases of partial sequences for the mitochondrial cytochrome b gene (Darveau et al., 2005a).
Calculations
Mathematical modeling of the phosphoglucoisomerase reaction was performed using Mathcad Professional, version 8 (MathSoft, Inc., Cambridge, MA, USA). Model assumptions, equations and parameters are provided in the text.
Results and discussion
Glycolytic flux rates, estimated from V̇_{CO2} during flight, ranged from 4.64 to 13.0 μmol g^{1} min^{1} across species. Fig. 1 shows lnV_{max} plotted against lnJ with regression slopes equivalent to hierarchical regulation coefficients, ρ_{h}. At this level of analysis, completely hierarchical regulation is observed at the hexokinase (HK) step, where ρ_{h}=0.98, while steadystate flux rates at cytochrome c oxidase (COX) and glycogen phosphorylase (GP) are the result of control shared by hierarchical regulation (ρ_{h}=0.46 and 0.36, respectively) and metabolic regulation (ρ_{m}= 0.54 and 0.64, respectively). The finding that hierarchical regulation accounts for interspecific variation in HK flux remains robust after PIC analysis (Table 1). However, the statistical relationships for COX and GP are no longer significant after controlling for phylogeny. All other steps examined yield ρ_{h}=0; this means ρ_{m}=1 and that metabolic regulation completely accounts for interspecific variation in flux at these steps in energy metabolism.
Our results are similar to those obtained by ter Kuile and Westerhoff (2001) in the sense that most of the variation in glycolytic flux in their intraspecific study of protozoa was due to metabolic regulation, and complete or partial hierarchical control contributed to the regulation of only a few steps in energy metabolism. It is important to point out that ρ_{h} and ρ_{m} values are not flux control coefficients; i.e. they are not measures of the contributions of individual enzymecatalyzed steps to the control of flux relative to other steps in a pathway. Rather, they allow quantitative determination of the relative contributions of two alternative mechanisms (i.e. hierarchical or metabolic regulation) in explaining pathway flux changes at individual steps.
Although allometric variation in biochemical flux capacities is widespread among animals (Suarez et al., 2004), the regulatory mechanisms underlying allometry in [E] remain largely unexplored. Two quite noteworthy exceptions to this are intraspecific studies on fish by Yang and Somero (1996) and Burness et al. (1999). In both studies are found examples wherein [E] does not correlate directly with the tissue content of the corresponding mRNA. This indicates that translation rates and/or rates of enzyme degradation may play important roles in the allometric variation in [E]. It would be worthwhile to determine whether HK mRNA content is sizedependent in orchid bees. However, given the preponderance of metabolic regulation revealed by our results, how such mechanisms lead to the allometric scaling of metabolism in orchid bees also merits further consideration.
Metabolic regulation at the PGI step
In previous work using honeybees Apis mellifera, we combined empirical and modeling approaches to try to understand the relations between flux capacities and flux rates at the PGI step, a nearequilibrium, reversible, glycolytic reaction at which glucose 6phosphate, G6P, is converted to fructose 6phosphate, F6P (Staples and Suarez, 1997). At such a nearequilibrium step, the enzyme catalyzes reactions in both directions, and net glycolytic flux equals the forward rate minus the reverse rate. Using Haldane's equations (Haldane, 1930) in a similar manner as Veech et al. (1969), we showed that the high V_{max} at this step could be rationalized on the bases of the enzyme's dual roles, i.e. to maintain a net forward flux, equal to the overall rate of glycolysis, and to maintain the reaction close to equilibrium. We now consider how substrate and product concentrations in vivo would have to change to account for the range of net flux rates observed interspecifically among orchid bees. The Haldane (1930) equation applied to PGI is: (5)
where J is the net flux, [G6P] and [F6P] are glucose 6phosphate and fructose 6phosphate concentrations, respectively, K_{f} and K_{r} are Michaelis constants for the forward (G6P is substrate) and reverse (F6P is substrate) reactions, respectively. We used K_{f}=1.07 mmol l^{1} and K_{r}=0.117 mmol l^{1}, obtained under simulated intracellular conditions using purified honeybee PGI (Staples and Suarez, 1997). It is assumed that both K_{f} and K_{r} are massindependent and conserved across orchid bee species. K_{eq} is the equilibrium constant obtained in vitro, equal to 0.3 (Staples and Suarez, 1997). V_{r} represents the V_{max} value in the direction of F6P conversion to G6P. These scale isometrically in orchid bees and average 345.6 μmol g^{1} min^{1} (Darveau et al., 2005b). V_{f} is the V_{max} in the glycolytic direction, estimated from the Haldane (1930) equation for K_{eq}, (6)
as 948.2 μmol g^{1} min^{1}. This is close to the V_{f} value measured empirically in honeybee flight muscles (Staples and Suarez, 1997). Eq. 5 is used to solve for the mass action ratios, i.e. the ratios of [F6P]/[G6P] in the flight muscles, required to account for the range of estimated flux rates. The decline in mass action ratio (Fig. 2A) is extremely small, a result that can also be demonstrated by holding [F6P] constant to 0.1 mmol l^{1}, a value within the range measured in vivo (Staples and Suarez, 1997), and solving for the range in [G6P] required to obtain the range of observed net forward flux rates (Fig. 2B). The model predicts that a change of about 0.02 mmol l^{1} in [G6P] is enough to account for the entire range of variation in flux across orchid bee species. A highly amplified response to a small change in substrate concentration, as seen here, is a characteristic of nearequilibrium reactions (Newsholme and Crabtree, 1976). Given this result, it seems unlikely that the expected range of variation in [F6P]/[G6P] in orchid bee flight muscles would be detectable in a comparative, interspecific study. At many other steps in metabolism, this problem is exacerbated by even lower steadystate concentrations of substrate and product (e.g. Kashiwaya et al., 1994). Thus, a modeling approach such as the one described here may be the only way to address the issue of metabolic regulation at such steps.
A broader view of metabolic regulation
Metabolic control analysis using skinned muscle fibers from mammals has revealed that as much as half of the control of oxidative metabolism in muscles is exerted by actomyosin ATPase (Wisniewski et al., 1995). In the context of metabolic scaling in orchid bees, a possible interpretation of our results is that, as muscle operating frequencies decline with increasing body mass (Darveau et al., 2005a), muscle ATP hydrolysis rates decline (Darveau et al., 2005b) and this drives the allometric scaling of metabolism. Such an interpretation follows from the similarity in allometric exponents of wingbeat frequency and hovering metabolic rate and is consistent with the predominant role played by metabolic regulation seen in the results reported here.
Further insights into the role played by metabolic regulation can be derived from work concerning bioenergetic scaling in mammalian hearts. Topdown control analysis has revealed that work rate (therefore ATP hydrolysis rate) dominates the control of oxidative metabolism in the heart (Diolez et al., 2002). Using ^{31}PNMR spectroscopy, Dobson and Headrick (1995) and Dobson and Himmelreich (2002) found that cytosolic free [ADP] declines as body mass increases in hearts, such that 1/[ADP] and the cytosolic phosphorylation potential [ATP]/[ADP][P_{i}] both scale allometrically with exponents close to 0.25. Dobson and colleagues propose that these result in higher `kinetic gain' such that small changes in [ADP] result in greater fractional changes in metabolic rate in the hearts of small mammals compared with larger ones. This scheme provides a plausible control mechanism by which allometry in work rate drives allometry in metabolic rates. That this applies to orchid bee flight muscles is a testable hypothesis.
Conclusion
The relative importance of the roles played by hierarchical regulation, resulting from variation in [E], and metabolic regulation in accounting for interspecific variation in metabolic rates is not well understood. Among orchid bees, the allometry in metabolic rates during hovering flight involves the hierarchical regulation of flux at the HK step and metabolic regulation at all other enzymecatalyzed reactions that we examined. Thus, although techniques for the study of gene expression can be brought to bear on the problem of metabolic scaling, the predominance of metabolic regulation in accounting for interspecific variation in flux rates also warrants the application of a systems approach to the control of flux (Fell, 1997) and the manner in which it scales (Darveau et al., 2002; Suarez et al., 2004).
ACKNOWLEDGEMENTS
This work was supported by a NSF research grant (IBN0075817) to R.K.S., a NSERC grant to P.W.H., and predoctoral fellowships from NSERC and FCAR to C.A.D. We are deeply grateful to David Roubik for collaborative efforts, Oris Acevedo for logistical support on Barro Colorado Island, and George Somero for constructive comments on our manuscript. Dedicated to the memory of Peter W. Hochachka, inspiring mentor and dear friend.
FOOTNOTES

↵† Deceased
 © The Company of Biologists Limited 2005