Living on the Edge: Physiological and Kinetic Trade-Offs Shape Thermal Tolerance in Intertidal Crabs From Tropical to Sub-Antarctic South America

Temperature is an important abiotic factor that drives the evolution of ectotherms owing to its pervasive effects at all levels of organization. Although a species’ thermal tolerance is environmentally driven within a spatial cline, it may be constrained over time due to differential phylogenetic inheritance. At the limits of thermal tolerance, hemolymph oxygen is reduced and lactate formation is increased due to mismatch between oxygen supply and demand; imbalance between enzyme flexibility/stability also impairs the ability to generate energy. Here, we characterized the effects of lower (LL50) and upper (UL50) critical thermal limits on selected descriptors of aerobic and anaerobic metabolism in 12 intertidal crab species distributed from northern Brazil (≈7.8°S) to southern Patagonia (≈53.2°S), considering their phylogeny. We tested for (i) functional trade-offs regarding aerobic and anaerobic metabolism and LDH kinetics in shaping thermal tolerance; (ii) influence of shared ancestry and thermal province on metabolic evolution; and (iii) presence of evolutionary convergences and adaptive peaks in the crab phylogeny. The tropical and subtropical species showed similar systemic and kinetic responses, both differing from the sub-Antarctic crabs. The lower UL50’s of the sub-Antarctic crabs may reflect mismatch between the evolution of aerobic and anaerobic metabolism since these crabs exhibit lower oxygen consumption but higher lactate formation than tropical and subtropical species also at their respective UL50’s. LDH activity increased with temperature increase, while KmPyr remained fairly constant; catalytic coefficient correlated negatively with thermal niche. Thermal tolerance may rely on a putative evolutionary trade-off between aerobic and anaerobic metabolism regarding energy supply, while temperature compensation of kinetic performance is driven by thermal habitat as revealed by the LDH affinity/efficiency equilibrium. The overall physiological evolution revealed two homoplastic adaptive peaks in the sub-Antarctic crabs with a further shift in the tropical/subtropical clade. The physiological traits at UL50 have evolved in a phylogenetic manner while all others were more plastic. Thus, shared inheritance and thermal environment have driven the crabs’ thermal tolerance and metabolic evolution, revealing physiological transformations that have arisen in both colder and warmer climes, especially at higher levels of biological organization and phylogenetic diversity.

Temperature is an important abiotic factor that drives the evolution of ectotherms owing to its pervasive effects at all levels of organization. Although a species' thermal tolerance is environmentally driven within a spatial cline, it may be constrained over time due to differential phylogenetic inheritance. At the limits of thermal tolerance, hemolymph oxygen is reduced and lactate formation is increased due to mismatch between oxygen supply and demand; imbalance between enzyme flexibility/stability also impairs the ability to generate energy. Here, we characterized the effects of lower (LL 50 ) and upper (UL 50 ) critical thermal limits on selected descriptors of aerobic and anaerobic metabolism in 12 intertidal crab species distributed from northern Brazil (≈7.8 • S) to southern Patagonia (≈53.2 • S), considering their phylogeny. We tested for (i) functional trade-offs regarding aerobic and anaerobic metabolism and LDH kinetics in shaping thermal tolerance; (ii) influence of shared ancestry and thermal province on metabolic evolution; and (iii) presence of evolutionary convergences and adaptive peaks in the crab phylogeny. The tropical and subtropical species showed similar systemic and kinetic responses, both differing from the sub-Antarctic crabs. The lower UL 50 's of the sub-Antarctic crabs may reflect mismatch between the evolution of aerobic and anaerobic metabolism since these crabs exhibit lower oxygen consumption but higher lactate formation than tropical and subtropical species also at their respective UL 50 's. LDH activity increased with temperature increase, while K m Pyr remained fairly constant; catalytic coefficient correlated negatively with thermal niche. Thermal tolerance may rely on a putative evolutionary trade-off between aerobic and anaerobic metabolism regarding energy supply, while temperature compensation of kinetic performance is driven by thermal habitat as revealed by the LDH affinity/efficiency equilibrium. The overall physiological evolution revealed two homoplastic adaptive peaks

INTRODUCTION
Biological systems are thermodynamically open, owing to the continuous flux of energy and matter with the surrounding environment (Cannon, 1929;Recordati and Bellini, 2004). Temperature is the main abiotic factor that drives the biogeographical evolution of ectotherms (Pörtner, 2002;Johns and Somero, 2004;Somero, 2004Somero, , 2010Dong et al., 2008;Stevens et al., 2010;Tomanek and Zuzow, 2010) since it affects all levels of organization, from the subcellular to the systemic, including the intact organism. The deleterious effects of temperature often exceed the regulatory capability of an organism's homeostatic mechanisms, and thus limit the distribution of life, imposing constraints on biogeographical distributions.
In aquatic decapod crustaceans at usual habitat temperatures, the capacity for oxygen uptake by ventilation and its delivery by circulation matches the biochemical demands of the tissues for oxygen: the fully oxygenated hemolymph thus can sustain optimal aerobic metabolism in the tissues (Pörtner and Gutt, 2016;Verberk et al., 2016). However, at the critical thermal limits, hemolymph oxygen levels are reduced due to the mismatch between oxygen supply and demand, imposing limitations on temperature tolerance (Pörtner, 2002). This limitation to systemic oxygen availability creates a hypoxic state, which is unable to sustain tissue demand, leading to anaerobic metabolism at the critical thermal limits as a consequence of elevated oxygen demand or insufficient aerobic capacity (Giomi and Pörtner, 2013). The recruitment of chaperone molecules and anti-oxidant defenses may extend survival, albeit with reduced performance at the molecular level, reverberating at the whole animal level and leading to death. This unifying principle known as "oxygen and capacity limitation of thermal tolerance" (OCLTT) characterizes sensitivity to thermal stress in aquatic invertebrates (Pörtner and Gutt, 2016;Verberk et al., 2016). During aerial exposure, however, air-breathers and amphibious species can attain higher hemolymph oxygen concentrations and increase oxygen diffusion, which may enhance their thermal tolerance by extending aerobic scope (Giomi et al., 2014).
To illustrate, the spider crab Maja squinado exhibits highly oxygenated hemolymph between 8 and 17 • C, showing decreased oxygen partial pressures below 8 • C and increased hemolymph lactate above 17 • C . The shore crabs Carcinus maenas and Cancer irroratus also accumulate hemolymph lactate at temperatures above 32 and 28 • C, respectively, after direct transfer from 12 • C (Jost et al., 2012). The transition to a predominantly anaerobic metabolism in the lobster Homarus americanus takes place gradually and continuously, accompanying a progressive increase in hemolymph lactate, reaching a maximum at 30 • C (Jost et al., 2012). Clearly, the limitation of oxygen supply to the tissues is linked to the onset of anaerobiosis due to a decreased capacity for oxygen uptake and transport. Further, patent patterns of increased oxygen demand and lactate formation are evident when temperature becomes life limiting.
Considering the subcellular environment, kinetic adjustments can result from oxygen-and capacity-limited thermal tolerance (Hochachka and Somero, 2002). Enzymes exhibit an inherent structural flexibility that guarantees their molecular function. However, they must also preserve the molecular rigidity that maintains their three-dimensional structure (Zavodszky et al., 1998). The balance between flexibility and stability in enzyme kinetics is temperature dependent, and is particularly evident with regard to lactate dehydrogenase (LDH) in vertebrate species. LDH converts pyruvate to lactate during anaerobic metabolism, reflecting the ability to produce energy under critical thermal conditions as predicted by the OCLTT principle (Zakhartsev et al., 2004).
Comparative studies of LDH have explored certain lineages of ectothermal vertebrates and invertebrates from various thermal niches, although a phylogenetic approach has only been implemented in the latter. In response to increased temperature, vertebrate kinetic patterns show decreased binding capacity and increased catalytic efficiency (Fields and Somero, 1998;Hochachka and Somero, 2002;Somero, 2003Somero, , 2004. LDH activity tends to increase at higher temperatures as a consequence of the augmented energy available for the conformational changes required during catalysis (Dunn et al., 1991;Fields, 2001;Hochachka and Somero, 2002). These kinetic properties of ectothermal vertebrate enzymes reflect a highly competent or occasionally a marginally functional LDH, depending on thermal regime. In contrast, in invertebrates, particularly exemplified by species of the anomuran genus Petrolisthes, LDH thermal stability is not affected by maximum habitat temperature, suggesting likely phylogenetic constraints on interspecific variability in LDH stability (Stillman and Somero, 2001).
Considering environments that offer a wide thermal landscape, the continental shelves of the Americas are divided into sixteen zoogeographic provinces. These have been defined based on "part of the neritic zone with a relatively narrow range of temperatures where the fauna shows certain homogeneity, " and particularly on oceanographic characteristics like water temperature and ocean currents (Boschi, 2000a,b;Boschi and Gavio, 2005). The eastern coast of South America embraces a significant thermal gradient and includes three provinces (Figure 1). The Brazilian province, between the mouth of the Orinoco river, Venezuela (9 • N) and Cabo Frio, Rio de Janeiro, southeastern Brazil (23 • S) (Briggs, 1974;Boschi, 2000a,b) in which mean surface water temperatures vary between 20.7 and 30.5 • C (median = 26 • C). The Argentinian province, between Cabo Frio and Rawson, Chubut, central Argentina (43-44 • S) (Cooke, 1895;Boschi, 2000a), which exhibits temperatures between 7.0 and 27.5 • C (median = 15 • C). And the Magellanic province, between Rawson and Tierra del Fuego (55 • S), Patagonian Argentina (Carcelles and Williamson, 1951;Balech, 1954;Boschi, 2000a,b) in which temperatures vary between 4.0 and 18.0 • C (median = 9 • C). These ranges represent the minimum and maximum mean monthly temperatures registered between 2009 and 2019 at the northern and southern geographical limits of each province (National Oceanic and Atmospheric Administration [NOAA] and U.S. Department of Commerce, 2019).
We have shown (Faria et al., 2017) that the micro-habitat temperatures (MHT) of some crab species from the Brazilian and Argentinian zoogeographical provinces are similar, but are higher than those of the Magellanic species. The upper (UL 50 ) and lower (LL 50 ) critical thermal limits of these crabs are subject to different environmental pressures, and exhibit unique evolutionary histories. The UL 50 show phylogenetic signal and correlate with micro-habitat temperature, while the LL 50 are more plastic and differ markedly among all provinces, revealing an asymmetrical evolutionary history of the critical thermal limits.
In the present study, we now investigate physiological and evolutionary patterns of thermal tolerance in this crab clade (Faria et al., 2017), characterizing the effects of critical thermal limits (LL 50 and UL 50 ) on selected descriptors of aerobic and anaerobic metabolism at the systemic and kinetic levels. Using 12 intertidal crab species distributed across the three thermal provinces from northern Brazil to southern Patagonia, we tested for: (i) the influence of shared ancestry and thermal province on metabolic evolution; (ii) the existence of functional trade-offs among aerobic and anaerobic metabolism and LDH kinetics in shaping thermal tolerance; and (iii) the presence of evolutionary convergences and adaptive peaks in the crab phylogeny. This study considers that species' traits are linked to each other through space and time, embracing most dimensions of physiological evolution.

Crab Sampling and Laboratory Maintenance
Approximately 70 adult, intermolt, male or female crabs of each of the 12 species were collected during the morning, at low tide, from mangroves, salt marshes, and sandy or rocky beaches, at the ends of the summers of 2013 and 2014. Collection sites included the Brazilian, Argentinian and Magellanic zoogeographical provinces (Boschi and Gavio, 2005; Figure 1), which correspond roughly to the Tropical, Subtropical and sub-Antarctic physiographical zones.
In the laboratory, the crabs were maintained prior to experiments for 5 days in plastic boxes, each containing a thin layer of seawater (≈3 mm deep, 30 S) to allow gill wetting, in incubators or in water baths (Fanem 347 DCG; Solab SL 224; Polystat 12002, Palmer Instrument Company) at the respective field reference temperatures (FRT) of each zoogeographical province (see above). These temperatures represent the median of the mean monthly minimum and maximum temperatures registered at the northern and southern geographical limits of each province (sensu Boschi, 2000a; National Oceanic and Atmospheric Administration [NOAA] and U.S. Department of Commerce, 2019). This protocol also minimized possible intrinsic physiological variability related to the species' vertical positions within the intertidal habitat. The crabs were not fed during the experiments, and the seawater was replaced daily in the morning. The

Critical Thermal Limits
Lower (LL 50 ) and upper (UL 50 ) critical thermal limits and microhabitat temperatures for each of the species evaluated here have been published previously ( Table 1; Faria et al., 2017). Both limits define the temperatures at which "the animal loses its ability to escape from conditions that will promptly lead to its death" (Cowles and Bogert, 1944). At each limit, 50% of the crabs were unable to right themselves when placed in a supine position. This inability reflects the loss of key motor functions and/or potentially lethal, disorganized locomotory activity (Faria et al., 2017).
In brief, LL 50 and UL 50 were established (8 ≤ N ≤ 10 crabs per temperature and species) by direct transfer of the crabs to various lower or higher temperatures, respectively, from their respective FRT's (26 • C for the tropical, 15 • C for the subtropical, and 9 • C The temperatures 26, 15, and 9 • C represent the field reference temperatures (FRT) at which the crabs from each province were acclimated. The tropical species (Aratus pisonii, Cardisoma guanhumi, Goniopsis cruentata, Ocypode quadrata, Pachygrapsus transversus, Uca maracoani, and Ucides cordatus) were collected from the Brazilian province (northeastern Brazil, ≈7.8 • S, 34.8 • W, red). The subtropical species (Armases rubripes, Neohelice granulata and Leptuca uruguayensis) were collected from the Argentinian province (southern Brazil, ≈32.1 • S, 52.1 • W, green). The sub-Antarctic species (Acanthocyclus albatrossis and Halicarcinus planatus) were collected from the Magellanic province (southern Patagonia, ≈53.2 • S, 67.2 • W, blue). The phylogenetic relationships among the species investigated (right panel) were generated employing a maximum likelihood search method using the mitochondrial 16S ribosomal gene as a marker (see Faria et al., 2017).
for the sub-Antarctic species). The crabs were examined for 30 s every 30 min during a 6 h period. This duration of exposure represents "the mean emergence time of a crab at low tide" (Faria et al., 2017) and constitutes the mean duration of exposure to aerial temperatures during low tide.
This experimental design, i.e., direct exposure from the respective FRT to the LL 50 or UL 50 of each species, enabled quantitative measurement of the critical thermal limits of each species, in line with Faria et al.'s (2017) protocol. It thus allowed direct comparison and straightforward correlations between the critical thermal limits and the physiological traits evaluated.

Mass-Specific Oxygen Consumption
Individual crabs (10 ≤ N ≤ 15) were placed into acrylic respirometer chambers containing a 3 mm deep layer of seawater (30 S), constructed specifically to accommodate each species with regard to diameter and volume. The respirometers were held in incubators or water baths (Fanem 347 DCG; Solab SL 224; Polystat 12002, Palmer Instrument Company) at the respective FRT of each species for a 12 h adjustment period prior to beginning the measurements. Each respirometer was coupled by plastic Crystal tubing to a multiplexing manometer and a gas analyzer (FoxBox-C Field Oxygen & Carbon Dioxide Analysis System, Sable Systems International Inc.), and perfused with an airflow of from 150 to 600 mL min −1 , depending on respirometer volume, standardized for each species. The standardized airflow provided a baseline percentage of atmospheric oxygen (20.8% O 2 ) in the absence of ambient pressure variation owing to the automated, barometric pressure compensation system.
After the 12 h adjustment period at the FRT, the respirometers were transferred to an incubator or water bath to measure O 2 consumption at the LL 50 or UL 50 temperatures established previously for each species, using a stop-flow respirometric protocol (Steffensen, 1989). After a 4 h period of continuous airflow, the respirometers were sealed by closing their respective manometer valves for a ≈2 h period. The O 2 respired by each crab in each respirometer was then measured by recording the decrease in percentage oxygen compared to the pre-established baseline reading while renewing the airflow during a 6 min period of flushing to the gas analyzer. Since measurements were conducted at the critical thermal limits, some crabs died during the measurements. Only the lowest O 2 consumption measurements recorded from at least 5-7 crabs alive at the end of each experiment were used.
After each experiment, each crab was cryoanesthetized in a freezer and killed by immersion for 20 min in crushed ice. The crabs were then dried in an oven at 60 • C for up to 72 h, and their dry masses were measured using an analytical balance (±1 mg precision).
Data acquisition software (Expedata, Sable Systems International Inc.) was used to calculate the integral of the area corresponding to the decrease in percentage oxygen (I A ) in the flushed air. The individual mass-specific oxygen consumption rate (QO 2 ) for each crab was obtained from the formula: where I A is oxygen percentage area; f, airflow rate (150-600 mL min −1 through each respirometer during the open phase); m, MHT data are the mean ± SEM (7 ≤ N ≤ 15). Data are independent of species' vertical positions (Stillman and Somero, 2000) within their intertidal zones (pGLS, 0.77 ≤ F ≤ 3.9, 0.23 ≤ P ≤ 0.54). Bold values are the mean ± SEM for all species in each province. dry mass; t, duration of sealed phase. The factors "100" and "60" convert O 2 percentage to O 2 volume, and minutes to hours, respectively. The effect of body mass on QO 2 was removed by transforming the QO 2 measurements into log 10 values, which were then regressed against the respective dry masses. The residues of this regression were then employed, thus ensuring mass-free data.

Extraction of Hemolymph and Abdominal Muscle Samples
Specimens (10 ≤ N ≤ 15) of each crab species were transferred from their respective FRT's to their respective LL 50 and UL 50 temperatures in incubators or water baths. After a 6 h exposure period at each thermal limit, a 50 µL hemolymph sample was drawn through the arthrodial membrane at the base of the last pereiopod of the surviving crabs using a #27-5 gauge needle coupled to an insulin syringe. Approximately 50 mg of abdominal muscle were also dissected. All samples were frozen immediately in liquid nitrogen and transferred in dry ice to the Laboratory of Crustacean Physiology at the University of São Paulo in Ribeirão Preto where they were stored at −80 • C until analyses could be performed.

Hemolymph Lactate
Perchloric acid (6%) was added to each hemolymph sample (3:1, v/v), followed by centrifugation at 10,000 g for 5 min at 4 • C. The supernatant was neutralized with 5 M K 2 CO 3 (10:1, v/v) and recentrifuged. For microplate reading, an aliquot of each sample was added to a buffer containing 0.5 M glycine, 0.4 M hydrazine and 10 mM EDTA at pH 9.0. 1 mM NAD + (final concentration) and 20 U/mL LDH (final concentration) were added to each well, and the microplate was incubated at 25 • C for 1 h in the dark and read at 340 nm (Bergmeyer and Bernt, 1974;Santos and Keller, 1993;Henry et al., 1994). The lactate concentration (mg L −1 ) was obtained from the formula: where A, Absorbance; 0.3, final fluid volume in each well; 6.22, molar extinction coefficient of NADH at 340 nm; 0.874, well optical path; 0.015, sample volume in each well; 3.4, sample dilution factor; and 90.1, conversion factor from mmol L −1 to mg L −1 .

Temperature Coefficient
The effect of temperature on QO 2 (aerobic sensitivity) and on lactate concentration (anaerobic sensitivity) was expressed using the formula Q 10 = (k 1 /k 2 ) 10/(t1−t2) , where k 1 and k 2 are the mass specific oxygen consumption rates at temperatures t 1 and t 2 , respectively, where t 1 > t 2 . Q 10O2 and Q 10Lac represent the temperature coefficients for changes in QO 2 and hemolymph lactate, respectively.

Muscle Lactate Dehydrogenase Kinetics
Each abdominal muscle sample was homogenized (Eurostar Power-B, IKO Labortechnik) on crushed ice at 1,500 rpm for 60 s in 50 mM K 2 HPO 4 buffer (7:1 v/m, mL g −1 , unadjusted pH) (modified from Zietara et al., 1996). The homogenate was then centrifuged at 14,000 g for 40 min at 4 • C and the supernatant was separated for enzyme kinetic assays. The LDH assays (1 ≤ N ≤ 4, pools of samples from 2 to 3 individual crabs) were performed using 100 mL of supernatant diluted from 10 to 2,000 times, depending on condition and species. This assay evaluates the oxidation of NADH to NAD + assuming steady-state conditions such as the linear relationship between NAD + formation and time (Segel, 1976). Controls without added homogenate were conducted to evaluate spontaneous NADH oxidation.
LDH activity was measured continuously, accompanying the oxidation of NADH to NAD + at 340 nm in thermostatted cells (ε 340 nm,pH 6 . 9 = 5,860 mM −1 cm −1 ; Cary 60 UV-Vis Spectrophotometer, Agilent Technologies) at both 25 • C and at the respective FRT, LL 50 , and UL 50 temperatures of each crab species for up to 5 min. LDH activity was not measured at the LL 50 's of Acanthocyclus albatrossis and Halicarcinus planatus owing to the inability of the spectrophotometer to read at these low temperatures (≈-1 • C).
For the kinetic assays, up to 14 pyruvate concentrations ranging from 1 µM to 7 mM were used in 80 mM Imidazole-HCl buffer at pH 6.9, containing 150 µM NADH in a final volume of 1 mL (modified from Stillman and Somero, 2001). The pH of the Imidazole buffer was adjusted according to the temperature of each assay. The kinetic parameters, i.e., maximum LDH activity (V max , U mg −1 ), the Michaelis-Menten constant for pyruvate (K m Pyr , mM) and the catalytic coefficient (min −1 10 3 ) were calculated using the saturation curve for LDH activity against pyruvate concentration, employing the SigrafW application (Leone et al., 2005). K m Pyr describes the enzymesubstrate affinity, representing the pyruvate concentration at which LDH activity is V max /2, indicating the ability to form an enzyme-substrate complex. The catalytic coefficient was measured at the FRT and was calculated for 25 • C using a Q 10 of 2, reflecting a measure of the rates of enzyme-substrate formation and dissociation (Segel, 1976;Northrop, 1998) and, thus, the rate constant for transforming pyruvate (substrate) into lactate (product). V max also was calculated for 25 • C, using a Q 10 of 2. One enzyme unit (U) is defined as the amount of LDH that hydrolyzes 1 mmol of NADH per minute. Assays were performed in duplicate aliquots.
Total protein titers in the homogenates were measured following Read and Northcote (1981), using bovine serum albumin as the standard.

Statistical Analyses
For the intra-specific analyses, the effects of critical thermal limits on the systemic parameters QO 2 , lactate concentration and their respective Q 10 's, and on the kinetic parameters V max and K m Pyr , for each species, were evaluated using a one-way Analysis of Variance followed by the Student-Newman-Keuls (SNK) multiple means comparison procedure to detect differences among groups (P ≤ 0.05). Data are expressed as the mean ± SEM.
Inter-specific comparisons were conducted based on a molecular phylogenetic analysis of 36 species, performed by a maximum likelihood procedure using the 16S ribosomal gene (Faria et al., 2017), a marker commonly used in phylogenetic studies of decapod Crustacea (e.g., Schubart et al., 2000;Mantelatto et al., 2009a,b;Pileggi and Mantelatto, 2010). Briefly, the phylogenetic analysis was undertaken using the Maximum Likelihood search method (Felsenstein, 1973(Felsenstein, , 1981 developed in RAxML 7.2.7 (Stamatakis, 2006) and implemented on the CIPRES system, employing the GTR + ∫ + I substitution model (Tavaré, 1986). A rapid bootstrap method (1,000 replicates) (Stamatakis et al., 2008) was used to evaluate topology consistency, only confidence values greater than 50% being considered. The original tree was pruned, maintaining only those species physiologically evaluated here. Branch lengths reflect genetic divergence.
The hypothesis that zoogeographical province may have affected the evolution of the physiological traits at the FRT, LL 50 , or UL 50 was tested using Phylogenetic Generalized Least Squares (PGLS) models. PGLS assumes that the residual variation among species is correlated owing to their shared ancestry (Grafen, 1989;Garland and Ives, 2000;Lavin et al., 2008), following a hierarchical disposition through time. Since character evolution is mathematically modeled along each branch of the phylogeny, we tested two evolutionary models simultaneously for each PGLS, employing the α parameter. The α-value is a measure of the selection strength that constrains a trait toward an optimum value: unconstrained character evolution (α = 0) is best fitted by the "Brownian motion" evolutionary model, while trait evolution with weakly limited or unlimited minimum and maximum values is best represented by the "Ornstein-Uhlenbeck" (O-U) model of evolution (Revell, 2010). This procedure was performed to ensure the best-adjusted evolutionary correlation for each set of inter-specific data.
Phylogenetic signal, a measure of the tendency of closely related species to resemble each other more than those drawn randomly from a clade, was calculated using Abouheif 's test (Abouheif, 1999;Pavoine et al., 2008). An independent Monte Carlo test was conducted for each trait to provide its Moran's I index and respective P-value. Moran's I varies from −1 to +1, revealing physiological dissimilarities and similarities, respectively, between closely related species (Gittleman and Kot, 1990;Diniz-Filho, 2001).
Overall physiological variability was analyzed using a phylogenetic principal components analysis (pPCA; Jombart et al., 2010). The phylogenetic structure of the multivariate comparative data was obtained by including all measured traits simultaneously (except V max and K m Pyr for the sub-Antarctic species where data were unavailable), including microhabitat temperature and both critical thermal limits, for all species. The scores of the two main dimensions (PC 1 and PC 2 ) for each species were then plotted against the phylogeny.
We also tested for the presence of convergent evolutionary events using the SURFACE method (Ingram and Mahler, 2013) to evaluate whether the hypothesis of physiological shifts among the three zoogeographic provinces is more likely than that expected by chance. Comparisons among several simulated O-U models were performed using the stepwise Akaike Information Criterion (AIC) procedure, and all improvements in AIC values were considered.
Our findings can be considered free of bias from the species' vertical positions within their intertidal zones since this parameter had no effect on any physiological trait (pGLS ANOVA, 0.77 ≤ F ≤ 3.9, 0.23 ≤ P ≤ 0.54). Our characterization of the vertical distributions assumed the five categories distinguished by Stillman and Somero (2000).
Critical temperatures affected mean QO 2 (ANOVA 20 ≤ F ≤ 46, P ≤ 0.004) and mean hemolymph [Lac] (ANOVA, 16 ≤ F ≤ 25, P ≤ 0.05) in the tropical and Sub-Antarctic species. These traits increased at the respective mean UL 50 's of ≈37 and ≈26 • C and diminished at the mean LL 50 's of ≈14 • C and ≈0 • C, both compared to the FRT's of 26 and 9 • C, respectively (SNK, P ≤ 0.05) ( Table 2). However, for the subtropical representatives, the mean QO 2 's and mean hemolymph [Lac] were similar when comparing these rates at the FRT of 15 • C and the mean LL 50 of ≈9 • C (SNK, P = 0.23), demonstrating that metabolism was unchanged during cooling. Values were higher, however, at the mean UL 50 of ≈37 • C (SNK, P ≤ 0.05).
FIGURE 3 | Effect of direct transfer of tropical (red), subtropical (green) and sub-Antarctic (blue) crab species on hemolymph lactate concentration after 6 h exposure to their respective upper and lower critical thermal limits, from their corresponding field reference temperatures of 26, 15, and 9 • C. Data (mg L −1 ) are mean values (4 ≤ N ≤ 7), and were adjusted to a second degree polynomial curve.
FIGURE 4 | Mean mass-specific oxygen consumption (QO 2 , mL O 2 g −1 h −1 ) and mean hemolymph lactate concentrations (mg L −1 ) for tropical, subtropical and sub-Antarctic crab species at their respective UL 50 , after 6 h direct transfer from their corresponding field reference temperatures of 26, 15, and 9 • C. The horizontal line within each box indicates the median QO 2 while the box boundaries show the interquartile values. Whiskers indicate the lowest and highest QO 2 's (range). Specimens (12 ≤ N ≤ 41 per species) for all species sampled from each zoogeographical province were used to provide the plots. *Significantly different from tropical and subtropical species (phyHolm-Bonferroni, P ≤ 0.05).

Thermal Effects on Muscle LDH Kinetics: V max , K m Pyr and Catalytic Coefficient
The highest V max was 7.43 ± 0.39 U mg −1 for the tropical species Ucides cordatus at its UL 50 of 39.0 • C while the lowest activity was 0.09 ± 0.00 U mg −1 in the sub-Antarctic species Halicarcinus planatus at the FRT of 9.0 • C (Table 3 and Figure 5). At the zoogeographical province level, the mean V max 's for the tropical and subtropical species were strongly sensitive to critical temperature exposure: V max 's were higher at the mean UL 50 and lower at the mean LL 50 (ANOVA; 11.5 ≤ F ≤ 15.2; P ≤ 0.01). The mean V max for the sub-Antarctic species was thermally insensitive (ANOVA; F = 1.0; P = 0.4). At an extrapolated temperature of 25 • C, mean V max values were affected by zoogeographical province (pGLS ANOVA, F = 6.0, P = 0.05, α = 16.9) with calculated values of 1.7 ± 0.3, 1.1 ± 0.1 and 0.6 ± 0.3 U mg −1 for the tropical, subtropical and sub-Antarctic species, respectively (Table 3 and Figure 5, insert). Mean K m Pyr for the subtropical and sub-Antarctic crab species were not affected by exposure at their critical thermal limits (ANOVA, 0.62 ≤ F ≤ 0.95; P ≥ 0.43) ( Table 3). However, mean K m Pyr for the tropical species was affected by exposure at the mean UL 50 (ANOVA, F = 4.6, P ≤ 0.05; SNK, P ≤ 0.05). The mean K m Pyr for the three provinces under all thermal regimes was 0.3 mM (Figure 6). In contrast, the catalytic coefficients calculated for 25 • C (3.0-26.5 min −1 10 3 ; Figure 5, insert) did regress on microhabitat temperatures (pGLS, −0.62 ≤ slope ≤ −0.64, 5.4 ≤ F ≤ 8.8, P ≤ 0.05, α = 22.9), when the sub-Antarctic species A. albatrossis and H. planatus are excluded from the analysis (see also  Catalytic coefficient (min −1 10 3 ) is calculated for 25 • C, using a Q 10 of 2 estimated from muscle samples at the FRT of each species. Data are the mean ± SEM (1 ≤ N ≤ 3, pools of 2-3 individuals per sample). *Significantly different from the respective FRT value (ANOVA, SNK, P ≤ 0.05). n.e., not estimated. Bold values are the mean ± SEM for all species in each province.

Phylogenetic Signal and Multi-Variate/Dimensional Patterns
With regard to the upper critical limit (UL 50 ), all traits, i.e., Q O2 , [Lac], V max and K m Pyr , as well as UL 50 itself, showed significant phylogenetic signal (I > 0.28, P < 0.05), thus their inter-specific variability exhibits a phylogenetic pattern. However, all other traits regarding exposure to FRT and LL 50 were unrelated to phylogeny (−0.22 ≤ I ≤ 0.19, P ≥ 0.05).
Phylogenetic principal components analysis, grouping all traits, revealed that the first two dimensions accounted for 54.3% of the total physiological variance, specifically 38.1% for PC 1 and 16.2% for PC 2 (pPCA, Figure 7): PC 1 is more related to the sub-Antarctic species while PC 2 to the tropical/subtropical clade. This scenario is in agreement with the two adaptive peaks detected in the phylogeny: (i) at the root, and maintained in A. albatrossis and H. planatus; and (ii) a shift at the outset of the tropical and subtropical species (SURFACE analysis, Figure 7). Microhabitat temperature, LL 50 and UL 50 were the traits that most strongly correlated with the more relevant pPCA dimension (R ≥ 0.30), and also contributed most to the adaptive peaks (−14 < AICc < −8) as evaluated using SURFACE.

DISCUSSION
We found that tropical and subtropical intertidal crab species showed similar subcellular and organismic responses to acute thermal stress, differing from their sub-Antarctic relatives. Oxygen consumption (QO 2 ) and hemolymph lactate ([Lac]) increased with increasing temperature in all species from the three zoogeographical provinces, revealing an apparent tradeoff with the evolution of upper critical thermal limit (UL 50 ). While maximum LDH activity (V max ) also increased with temperature, K m Pyr tended to remain constant, whereas the catalytic coefficient is driven by micro-habitat temperature in FIGURE 5 | Effect of direct transfer of tropical (red), subtropical (green), and sub-Antarctic (blue) crab species on muscle lactate dehydrogenase (LDH) activity (V max , U mg −1 ) after 6 h exposure to their respective upper and lower critical thermal limits from their corresponding field reference temperatures (FRT) of 26, 15, and 9 • C. Data are mean values (2 ≤ N ≤ 3, pools of 3 individuals each) and were adjusted to a second degree polynomial curve. V max was measured continuously in duplicate aliquots at the same experimental temperatures. Activities could not be measured for the sub-Antarctic species at their lower critical thermal limits. Insert: Mean (± SEM) muscle LDH activity at 25 • C for all species from each province (2 ≤ N ≤ 3, pools of 3 individuals each); muscle samples taken at the FRT of each species. *Significantly different from tropical crab species (phyHolm-Bonferroni, P ≤ 0.05).

FIGURE 6 | Michaelis-Menten constants (K m
Pyr , mM) for muscle lactate dehydrogenase (LDH) in tropical, subtropical and sub-Antarctic crab species at their respective lower (LL 50 ) and upper (UL 50 ) critical thermal limits after 6 h direct transfer from the corresponding field reference temperatures (FRT) of 26, 15, and 9 • C. Data are the mean ± SEM (2 ≤ N ≤ 3, pools of 3 individuals each). LDH activity was measured continuously in duplicate aliquots at the same experimental temperatures. Measurements could not be made for the sub-Antarctic species at their lower critical thermal limits. *Significantly different from tropical and subtropical crabs (phyHolm-Bonferroni, p ≤ 0.05). Insert: Catalytic coefficient (min −1 10 3 ) calculated for 25 • C, assuming a Q 10 of 2, plotted against the microhabitat temperature.
most cases. Interestingly, all physiological traits at the UL 50 have evolved in a phylogenetic fashion, while the others linked to LL 50 and FRT were more plastic. Overall metabolic evolution at FIGURE 7 | Overall physiological evolution of thermal tolerance, grouping systemic and kinetic traits in intertidal crabs distributed from northern Brazil to southern Patagonia, mapped onto a phylogeny modified from Faria et al. (2017). Different colors represent adaptive peaks identified using the SURFACE method: the first optimum (blue) was detected for the sub-Antarctic species while the second (red) was found in the tropical and subtropical representatives. Using a phylogenetic PCA, the first two dimensions accounted for 54.3% of total physiological variance, specifically 38.1% for PC 1 and 16.2% for PC 2 . The scores of both dimensions follow the same phylogenetic pattern as the adaptive peaks. Microhabitat temperature, LL 50 and UL 50 were the most important parameters in both analyses (R ≥ 0.30 for pPCA, −14 < AIC < −8 for SURFACE).
both critical limits has taken place by two independent adaptive pathways in the sub-Antarctic lineages, however one adaptive peak has been inherited by the northernmost species. These physiological patterns and evolutionary history appear to reflect the species' ability to function at higher temperatures, together with a lesser capacity to confront lower temperatures, since the crabs were likely acclimatized to the southern summer. The discussion that follows should be considered in this light.

Systemic Metabolism
The systemic metabolism complied with the paradigm that higher temperatures lead to greater energy demands through increased QO 2 and [Lac]. The highest values for these parameters were seen at the UL 50 while the lowest values occurred at the LL 50 , both compared to the FRT. At the FRT, the tropical crab species exhibited a mean QO 2 6-fold greater than the mean value for the subtropical crabs, and 3-fold greater than the sub-Antarctic representatives. The same phenomenon has occurred in the Echinodermata, Annelida, other Crustacea, and Mollusca: warmer water species show higher rates of oxygen consumption than do colder water species when measured at their respective usual environmental temperatures (Fox, 1936), reaching up to 10-fold greater in crustaceans from tropical waters than in polar species (Scholander et al., 1953). Further, the sub-Antarctic crab species showed a 2-fold higher QO 2 at the FRT of 9 • C than did the subtropical crabs at the FRT of 15 • C, which reveals temperature compensation of aerobic metabolism in cold regions, an adaptive peak probably maintained by strong stabilizing selection as revealed by the higher α-value. The evolution of temperature compensation of aerobic capacity as observed here may result from the differential regulation of hemolymph Mg 2+ as seen in some decapod crustaceans Wittmann et al., 2010Wittmann et al., , 2012 or may be due to selective pressure on cellular energy demand, as demonstrated in fish (Peck et al., 2002;Pörtner, 2002;Morley et al., 2009). Hemolymph Mg 2+ hypo-regulation is associated with higher activity levels in decapods from low temperatures since higher Mg 2+ concentrations diminish activity by inhibiting Ca 2+ influx across neuronal membranes, slowing down neuromuscular transmission and decreasing muscle contraction (Parnas et al., 1994;Dunn and Mercier, 2003). Elevated activity levels in decapods from naturally or experimentally cold temperatures are due to fairly high heart rates, ventilation frequencies and oxygen consumption, as a consequence of lower Mg 2+ levels . Increased mitochondrial densities also contribute to elevated aerobic metabolism in fish acclimated to cold temperatures revealing ultrastructural adjustments (Johnston and Maitland, 1980;Tyler and Sidell, 1984). Both hypotheses, in combination or independently, may explain the higher QO 2 values seen in the sub-Antarctic species.
With regard to phylogenetic signal for oxygen consumption, closely related species do not show similar QO 2 's at the FRT and LL 50 , which is also supported by the effect of zoogeographical province on QO 2 under both thermal regimes together with a significant effect of selection strength (17 ≤ α ≤ 20). Given the physical dependence of QO 2 on the kinetic energy of the molecules, this scenario suggests adaptive evolution strongly driven by thermal regime, regardless of the phylogenetic position of each species. To illustrate, the fiddler crabs sampled, i.e., the tropical species U. maracoani and the subtropical species L. uruguayensis, although closely related, exhibit distinct QO 2 's under their natural temperature conditions and on cold challenge. Conversely, aerobic metabolism at the UL 50 is similar among closely related species, also reinforced by the absence of a province effect, signifying that inter-specific variability in QO 2 at higher temperatures is constrained by shared inheritance, regardless of a significant effect of thermal regime. The maximum aerobic capacities at higher temperature challenge vary among species in a phylogenetic fashion.
With respect to anaerobic metabolism, mean [Lac] in the subtropical crab species did not vary as a function of their respective FRT and LL 50 , while the tropical and sub-Antarctic crabs showed a reduction of ≈50 and 70%, respectively, in mean [Lac] at the mean LL 50 compared to FRT. Since the subtropical species encounter a eurythermic environment owing to the mixing of warm water from the Brazilian Current with colder water from the Malvinas Current (Boschi, 2000a), their stable anaerobic metabolism suggests the capability for cold acclimatization and energy production even during rapid temperature changes such as might occur during exposure at low tide. In contrast, the tropical and sub-Antarctic crabs hail from more stenothermic provinces, dominated respectively by the South Equatorial Current and by sub-Antarctic waters (Boschi, 2000a;Thurman et al., 2013), suggesting that rapid cooling effects are more pervasive and more able to reduce anaerobic energy production. However, species from all three zoogeographical provinces showed an increase in mean [Lac] during exposure at their UL 50 compared to the FRT. Mean [Lac] was 1.4-and 1.8-fold higher at the mean UL 50 of ≈39 • C for the tropical and subtropical species, respectively, but 17-fold higher at the mean UL 50 of 26 • C for the sub-Antarctic species, revealing a stronger effect of higher temperatures in species acclimatized to a colder province. The lack of mobilization of lactate at the lower critical limits is possibly due to a lower cellular energy demand, which is compatible with the lower kinetic energy of the molecules.
The high [Lac] at the upper critical temperatures indicates recruitment of anaerobiosis owing to high tissue energy demand, while the high QO 2 's suggest some cardiorespiratory capability in supplying tissues with oxygen. For example, the spider crab Maja squinado shows changes in [Lac] in some tissues in tandem with a reduction in aerobic capacity. High (33 • C) and low (0 • C) temperatures provoke decoupling between O 2 supply and demand, since mitochondrial O 2 consumption is consistent with oxygen delivery by the cardiorespiratory system . Further, in decapod crustaceans subjected to hypoxia or exercise, lactate affects hemocyanin affinity for oxygen, increasing oxygen extraction by the gills and transport by the hemolymph to the tissues as seen in the brachyuran crabs Cancer, Carcinus, and Hyas (Truchot, 1980;Morris and Bridges, 1989), Homarus lobsters (Nies et al., 1992) and Macrobrachium shrimps (Ern et al., 2013).
Since aerobic and anaerobic metabolism are both linked to the energy supply of a biological system (see Pörtner et al., 2017), the lower mean UL 50 for the sub-Antarctic crabs (≈26 • C) compared to the mean UL 50 for the tropical and subtropical species (≈37 • C) may reflect a mismatch between aerobic and anaerobic energy production at the upper critical limits. In fact, under high thermal challenge, the sub-Antarctic species have evolved a lower aerobic metabolism (0.4-fold), albeit circumvented by an augmented anaerobic capacity (9.1-fold) when compared to the mean values for the tropical and subtropical species (Figure 4). Since aerobic capacity ensures higher ATP production (potentially up to 38 ATP molecules per glucose molecule) than does anaerobic metabolism (only 2 ATP molecules), the lower energy budget available to the sub-Antarctic species could impair overall cellular activities such as anti-oxidant defenses and heat-shock responses, reducing thermal tolerance. The acute stress encountered at the upper critical temperatures constrains aerobic capacity and shifts the energy generating system to a predominantly anaerobic metabolism in crabs from colder climates, suggesting an evolutionary metabolic trade-off during thermal niche diversification.

LDH Kinetics
For each zoogeographical province, mean V max tends to become reduced at the LL 50 and to increase at the UL 50 , both compared to the respective FRT. Interestingly, at the fixed temperature of 25 • C, there was an effect of zoogeographical province: mean V max for the sub-Antarctic species was 3-fold lower than that for tropical species. This tendency for lower activities at the same temperature for species from colder regions contrasts with previous data for tropical/subtropical and Antarctic fish (Kawall et al., 2002) that show thermal compensation of LDH specific activity. In this case, brain LDH activity was up to 2.5-fold higher in the Antarctic species, while muscle LDH activities were similar between the two thermal groups, both measured at the same temperature of 10 • C (Kawall et al., 2002).
The absence of thermal compensation in LDH activity in the intertidal, neotropical crabs evaluated here corroborates findings for another intertidal decapod clade distributed from Northern to Central America, anomuran crabs Petrolisthes. No effect of maximum habitat temperature was seen on the evolution of LDH thermal stability in 22 congeneric species sampled over a ≈35 • range of latitude, suggesting that the diversity of LDH stabilities also is not adaptive (Stillman and Somero, 2001). Further, similar stabilities at maximum temperatures for some sibling Petrolisthes groups were suggested, which is reinforced here for neotropical crabs, owing to the phylogenetic signal seen for V max at UL 50 . Although maximum and critical temperatures, respectively, do not constitute a selective pressure on the evolution of LDH stability (anomurans) and activity (brachyurans), the interspecific variability in these mechanistic bases demonstrates the importance of ancestral inertia in the evolutionary history of different intertidal decapod lineages.
While V max does not suggest LDH thermal compensation when comparing activities at the same temperature, the catalytic coefficient (also at 25 • C) correlated negatively with micro-habitat temperatures (MHT). Species sampled from higher MHT tend to have lower catalytic coefficients and vice-versa (except for the sub-Antarctic species A. albatrossis and H. planatus). In fact, the lower MHT values (21-26 • C) show a mean catalytic coefficient 2.5-fold higher than the mean value for the species from warmer habitats (28-36 • C), which means that species from colder MHT's show more flexible LDH molecules, while conversely, those from warmer habitats show more rigid enzymes. In contrast, mean K m Pyr tends to be constant (≈0.27 mM), reflecting retention of LDH affinity for pyruvate and, consequently, its binding efficiency.
The temperature compensation of LDH activity has evolved through the transformation of catalytic efficiency, which is in agreement with maintenance of binding capacity for pyruvate (see Hochachka and Somero, 2002). Such a biochemical tradeoff has been documented for several fish species from various thermal niches, such as Sphyraena barracudas (Holland et al., 1997), Chromis Pacific fish (Johns and Somero, 2004) and four families of Antarctic notothenioids (Fields and Somero, 1998). The maintenance of the LDH binding property demonstrates conservation of the catalytic site. Variability in enzyme efficiency results from alterations to the amino acid sequence in regions that undergo conformational changes during formation of the enzyme-substrate complex (see Hochachka and Somero, 2002). The reasons why LDH in tropical species has lost its affinity (higher K m Pyr ) during exposure at UL 50 , and why sub-Antarctic crabs have not retained a higher catalytic efficiency at their natural MHT, remain conjectural and require investigation at the molecular level.
Overall physiological evolution, when considered simultaneously for all systemic and kinetic traits, including MHT and both critical limits, consists of two steps. The first took place at the outset of the present crab phylogeny, and the second at the outset of the tropical + subtropical lineage, owing to the shift between the two main dimensions, PC 1 for the sub-Antarctic species and PC 2 for the tropical and subtropical crabs (Figure 7). Interestingly, these findings are in agreement with the two adaptive peaks detected in this crab phylogeny: one at the root, and maintained in A. albatrossis and H. planatus, and another at the outset of the tropical and subtropical species. However, since the sub-Antarctic A. albatrossis and H. planatus are a paraphyletic grouping, also belonging to different eubrachyuran families (Bellidae and Hymenosomatidae, respectively, see Faria et al., 2017), the adaptive peak revealed here should be interpreted as the result of adaptive convergence. The adaptive peaks detected using SURFACE are also in concert with most evolutionary associations detected using PGLS since they were better adjusted by the O-U model, with higher selection strengths that limit minimum and maximum physiological capabilities. The reason why the most important traits that contribute to this macroevolutionary landscape are the thermal niche, and the upper and lower critical thermal limits, may depend on their underlying biological organization, since such traits integrate lower hierarchical levels, including lesser physiological differences, following the synergistic principle of the "part-whole" relationship.
In summary, the metabolic physiology of crabs living at the edge of their thermal windows accompanies the asymmetry of their critical limits. The extension of life to higher temperatures depends on a putative trade-off between the evolution of aerobic and anaerobic metabolism with respect to energy supply, while temperature compensation of kinetic performance is driven by thermal habitat, as revealed by the LDH affinity/efficiency equilibrium. In fact, shared inheritance and the thermal environment both affect the overall physiological history of these intertidal crabs, revealing that certain evolutionary transformations have arisen both in warmer and in colder climes, particularly at higher levels of biological organization and phylogenetic diversity.

DATA AVAILABILITY STATEMENT
The datasets generated for this study are available on request to the corresponding author.

ETHICS STATEMENT
Crab sampling in Brazil was authorized by the Ministry for the Environment (MMA/ICMBio permit #29594-8 to JM), and in Argentina by the Ministry for the Environment, Sustainable Development and Climate Change of Tierra del Fuego (permit #0116/2014 to MR).