Paleohistological inferences of thermometabolic regimes in Notosuchia (Pseudosuchia: Crocodylomorpha) revisited

Abstract Notosuchia is a group of mostly terrestrial crocodyliforms. The presence of a prominent crest overhanging the acetabulum, slender straight-shafted long bones with muscular insertions close to the joints, and a stable knee joint suggests that they had an erect posture. This stance has been proposed to be linked to endothermy, because it is present in mammals and birds and contributes to the efficiency of their respiratory systems. However, a bone paleohistological study unexpectedly suggested that Notosuchia were ectothermic organisms. The thermophysiological status of Notosuchia deserves further analysis, because the methodology of the previous study can be improved. First, it was based on a relationship between red blood cell size and bone vascular canal diameter tested using 14 extant tetrapod species. Here we present evidence for this relationship using a more comprehensive sample of extant tetrapods (31 species). Moreover, contrary to previous results, bone cross-sectional area appears to be a significant explanatory variable (in addition to vascular canal diameter). Second, red blood cell size estimations were performed using phylogenetic eigenvector maps, and this method excludes a fraction of the phylogenetic information. This is because it generates a high number of eigenvectors requiring a selection procedure to compile a subset of them to avoid model overfitting. Here we inferred the thermophysiology of Notosuchia using phylogenetic logistic regressions, a method that overcomes this problem by including all of the phylogenetic information and a sample of 46 tetrapods. These analyses suggest that Araripesuchus wegeneri, Armadillosuchus arrudai, Baurusuchus sp., Iberosuchus macrodon, and Stratiotosuchus maxhechti were ectothermic organisms.


Introduction
Notosuchia is a group of extinct, mostly terrestrial crocodyliforms. The presence of several morphological features suggests that they had an erect posture: the prominent crest overhanging the acetabulum observed in Chimaerasuchus paradoxus (Wu and Sues 1996), Notosuchus terrestris (Pol 2005), Araripesuchus tsangatsangana (Turner 2006), Baurusuchus albertoi (Nascimento and Zaher 2010), and Stratiotosuchus maxhechti (Riff and Kellner 2011); the straight-shafted long bones described in Anatosuchus minor and Araripesuchus spp. (Sereno and Larsson 2009); the slender limb bones with muscular insertions close to the joints reported in Malawisuchus mwakasyungutiensis (Gomani 1997); and the tight/stable knee joint shown in Pissarrachampsa sera (Godoy et al. 2016). Among extant tetrapods, only endotherms (mammals and birds) show an upright stance. This last feature has been proposed to be linked to endothermy, because it contributes to the efficiency of the respiratory system (Carrier 1987). Thus, according to this morphological evidence, we can reasonably hypothesize that Notosuchia were endothermic. However,  concluded that they were primitively ectothermic using two proxies: resting metabolic rate (RMR) and red blood cell size (RBC size ).
RMR is the minimal consumption of oxygen over time per unit of body mass measured under postabsorptive conditions during the period of normal activity of the daily cycle in resting, nonreproductive specimens (Andrews and Pough 1985;Montes et al. 2007). RMRs of extant endotherms are at least one order of magnitude higher than those of extant ectotherms of similar body mass, because the mechanisms of thermogenesis operating in the former are costly in terms of energy (Clarke and Pörtner 2010;Legendre and Davesne 2020). RMRs inferred by  for Notosuchia were significantly lower than the threshold separating ectotherms from endotherms.
Within extant tetrapods, RBC size is lower in endotherms (mammals and birds) than in ectotherms (Amphibia, Squamata, Testudines, and Crocodylia) (Hartman and Lessler 1964;Snyder and Sheafor 1999;Soslau 2020). It has been suggested that the acquisition of lungs together with the subsequent evolution of the cardiovascular system was the driving force explaining the evolution of vertebrate RBC size (Snyder and Sheafor 1999). In endotherms, thermogenetic mechanisms use a huge amount of oxygen, producing high RMRs. Considering that "Smaller capillaries [and smaller RBCs] are associated with increased potential for diffusive gas exchange" (Snyder and Sheafor 1999: 189), these features may have been positively selected in endotherms. Huttenlocker and Farmer (2017) found that RBC size values are related to, and can be inferred from, bone vascular canal diameter.  inferred notosuchian RBC size values using this last relationship and concluded that they were significantly higher than the threshold separating ectotherms from endotherms.
To sum up, both proxies (RMR and RBC size ) suggest that Notosuchia were ectothermic organisms. Considering that paleohistological evidence (suggesting low RMR, large RBC size , and ectothermy) is not congruent with morphological evidence (suggesting an erect posture, cursoriality, and endothermy), the thermophysiological status (i.e., either ectothermic or endothermic) of Notosuchia deserves further analysis.
The approach used by  to perform these inferences can be improved in two ways. First, notosuchian thermophysiological status inferred using RBC size is based on the quoted relationship between RBC size and bone vascular canal diameter ). This relationship was tested by Huttenlocker and Farmer (2017) using a rather small sample size (14 extant tetrapod species). Here we tested this relationship using a more comprehensive sample of extant tetrapods (31 species) and phylogenetic generalized leastsquares regression (PGLS). Second, RBC size estimations were performed using phylogenetic eigenvector maps (PEMs), and this method excludes a fraction of phylogenetic information. This is because PEM generates a high number of eigenvectors (n − 1, with n being the number of terminal taxa analyzed), thus requiring a selection procedure to compile a subset of eigenvectors to avoid model overfitting (Guénard et al. 2013;Legendre et al. 2016).
Here we inferred the thermophysiology of Notosuchia using phylogenetic logistic regression (PLR) (Ives and Garland 2010;Tung Ho and Ané 2014), a method that overcomes this problem, because it includes all (instead of a fraction) of the phylogenetic information. Figure 1 Pol et al. (2014) for Notosuchia. Both phylogenies were dated using Time Tree of Life (http://www.timetree.org). When the ages of two successive nodes collapsed, we arbitrarily added 1 Myr in between the moreinclusive and less-inclusive nodes to facilitate the graphic visualization of the topology. For Notosuchia, nodes were dated according to the last appearance datum (LAD) of the oldest fossil included in each node taken from the Paleobiology Database (https://paleobiodb. org). The age of the node Notosuchia (113 Myr) corresponds to the LAD of Malawisuchus mwakasyungutiensis. The age of the node Armadillosuchus-Baurusuchus (100.5 Myr) corresponds to the LAD of Chimaerasuchus paradoxus. The age of the node Iberosuchus-Baurusuchus (83.6 Myr) corresponds to the LAD of Comahuesuchus brachybuccalis, Pehuenchesuchus enderi, Cynodontosuchus rothi, and Wargosuchus australis. Finally, the age of the node Stratiotosuchus-Baurusuchus (66 Myr) corresponds to the LAD of these taxa. The latter (Stratiotosuchus-Baurusuchus), as well as Armadillosuchus arrudai, come from the Adamantina Formation, the age of which is still debated. We follow the hypothesis of a Campanian-Maastrichtian age proposed by some authors (e.g., Gobbo-Rodrigues et al. 1999;Batezelli 2017).

Phylogenies in
Testing the Relationship between RBC size and Bone Vascular Canal Diameter Using PGLS.-Supplementary File 1 contains the sample (31 species of extant tetrapods) and the phylogeny (topology and divergence times) used to test the relationships between the response variables (RBC width and RBC area ) and the explanatory variables (femoral vascular canal diameter and femoral cross-sectional area including the medullary cavity). Thin sections of extant taxa are curated at the Vertebrate Hard Tissue Collection of the Museum national d'Histoire naturelle, Paris, and are available on request to the curator (D. Germain). RBC width (defined as RBC minimum diameter) and RBC area (either published values or values computed using maximum and minimum published diameters and assuming an ellipse) were taken from the literature (Supplementary File 2). Femoral vascular canal diameters (white arrowheads in Fig. 2) were computed as Can harmean and Can min , as defined by Huttenlocker and Farmer (2017). Can harmean , Can min , and femoral crosssectional area were either quantified in this study or taken from Huttenlocker and Farmer (2017) (data available in Supplementary File 2).
The method of ordinary least-squares regression makes the assumption of no covariance FIGURE 1. Phylogenetic relationships among extant taxa used to construct the thermophysiology inference model and the extinct Notosuchia for which we performed paleobiological inferences. Sources of topology and divergence times are cited in the main text. Scale on the right: geologic time in millions of years (Myr). between residuals obtained from the regression equation (i.e., the off-diagonals of the variancecovariance matrix are expected to contain zeros) (Symonds and Blomberg 2014). In analyses using interspecific data, this assumption is not verified because of the hierarchical, shared phylogenetic history among terminal taxa (i.e., closely related species are more similar than expected by chance). PGLS (Grafen 1989;Martins and Hansen 1996;Rohlf 2001;Symonds and Blomberg 2014) overcomes this problem by using a variance-covariance matrix in which off-diagonals correspond to the phylogenetic history shared by the two species under comparison. Symonds and Blomberg (2014) described PGLS as a "weighted regression" in which data points corresponding to closely related species are "downweighted." We ran PGLS using the function pgls of the package caper (Orme et al. 2013) in R (R Development Core Team 2008).
Inferring the Thermophysiology of Notosuchia Using PLR.- Figure 1 shows the phylogenetic relationships among extant taxa (46 species of tetrapods) used to construct the thermophysiology inference model (to infer the probability of being endothermic) and the extinct Notosuchia for which we performed paleobiological inferences. This model was constructed using femoral vascular canal diameter (Can harmean and Can min ) and femoral crosssectional area as explanatory variables. As noted earlier, femoral Can harmean , Can min , and femoral cross-sectional area were either quantified in this study or taken from Huttenlocker and Farmer (2017)   : Araripesuchus wegeneri, Armadillosuchus arrudai, Baurusuchus sp., Iberosuchus macrodon, and Stratiotosuchus maxhechti. Considering that the model was constructed using femora of extant species, we performed inferences only for those Notosuchia for which data for femora were available. Sample size was smaller for PGLS analyses, because data for RBC size were not available for many species analyzed in PLR analyses. PLR is a generalized linear model explaining the probability of occurrence of the state "presence" of a binary response (dependent) variable (here the "presence of endothermy") using continuous explanatory (independent) variables when residual variation of the former variable is phylogenetically structured (Ives and Garland 2010). The regression coefficients computed do account for phylogenetic correlation; when data are not phylogenetically structured, these coefficients are those of standard logistic regression (Ives and Garland 2010). PLR models contain two components. The first is controlled by parameters α (the transition rate) and μ (the asymptotic probability of being in state 1 [here the asymptotic probability of being endotherm]). Parameter α equals α 1 + α 0 ; α 1 being the probability that the response FIGURE 2. Transverse thin section, lateral side, of the femur of Araripesuchus wegeneri Buffetaut, 1981, observed in crosspolarized light with lambda wave plate. The thin section was made from a partial femur (MNHN.F.GDF660) from the Aptian of Gadoufaoua (Niger), and is curated at the Museum national d'Histoire naturelle (MNHN) (Paris, France). The cortex is made of lamellar-zonal bone. It is composed of three zones formed at moderate growth rate and containing vascular canals (white arrowheads) included in primary osteons, and three annuli formed at low growth rates and made of parallel fibered bone (black arrowheads). Periosteum is on the top and medullary cavity on the bottom. The continuous black line occurring near the medullary cavity is an artifact. Scale bar, 0.5 mm.
variable switches from 0 to 1 in each small time increment when it evolves up a phylogenetic tree, whereas the α 0 parameter is the probability that it evolves from 1 to 0 (Ives and Garland 2010). The transition rate α is a measure of phylogenetic signal, because the larger the α, the quicker the evolutionary transitions and the lower the phylogenetic structure of data (Ives and Garland 2010). In the second component, the probability of occurrence of the state "presence of endothermy" is modeled using values of the independent (explanatory) variable (here bone vascular canal diameter). Parameters α and μ are estimated using an iterative process in which μ is estimated given α, using the quasi-likelihood function, and α is estimated given μ, using least squares until convergence (Ives and Garland 2010). Analyses were performed using the package phyloglm (Tung Ho and Ané 2014) in R (R Development Core Team 2008).

Results
Testing the Relationship between RBC size and Bone Vascular Canal Diameter Using PGLS.-We used PGLS to test the relationships between RBC width and RBC area and the explanatory variables femoral vascular canal diameter (computed as Can harmean and Can min ) and femoral cross-sectional area (data available in Supplementary File 2). Shapiro-Wilk normality tests showed that residuals of PGLS regression of RBC area to Can harmean + bone cross-sectional area and the regression RBC area to Can min + bone cross-sectional area do not follow a normal distribution ( p-values of 0.0007557 and 0.001896, respectively). Thus, we performed a log transformation of all variables. After log transformation, residuals of all four PGLS regressions (RBC area and RBC area to the explanatory variables bone cross-sectional area and either Can min or Can harmean ) do follow a normal distribution. All four of these PGLS regressions were significant and, in each regression, both explanatory variables (bone crosssectional area and either Can min or Can harmean ) were significant (Table 1).
Inferring the Thermophysiology of Notosuchia Using PLR.-We used PLR to construct models aimed at computing the probability of being endothermic using paleohistological features (data available in Supplementary File 3). When Can harmean was used as the explanatory variable, we obtained a model with a transition rate α of 0.00144, an intercept estimate of 6.04 ( p-value = 0.004) and an estimate for the coefficient of Can harmean of −0.45 ( p-value = 0.001). The negative sign of the Can harmean coefficient indicates that the probability of being endothermic decreases as vascular canal diameter increases. (2) Ives and Garland (2010: p. 17) stated that "we assume that if μ i < m, then trait Y will evolve toward 0; […] Conversely, if μ i > m, then trait Y will evolve toward 1," where m is the mean probability of being endotherm in our sample. Thus, we considered that m is the cutoff probability, so that an inferred probability higher than m would be evidence for endothermy. Conversely, a probability lower than m would be evidence for ectothermy. When Can harmean was used as the explanatory variable, m = 0.59. To evaluate the predictive power of the model, we constructed a contingency table in which we inferred the thermometabolic regime of each extant species of the sample using its Can harmean . Lines contain predictions (0, inferred ectothermy; 1, inferred endothermy) and columns contain true states (0, observed ectothermy; 1, observed endothermy): state 0 state 1 prediction 0 14 2 prediction 1 3 27 The specificity (the ratio of quantity of true 1 inferred as 1 on the quantity of true 1; Sp = 27/ (27 + 2)) equals 0.931. The sensitivity (the ratio of quantity of true 0 inferred as 0 on the quantity of true 0; Se = 14/(14 + 3)) equals 0.824. The classification error (the ratio (quantity of true 0 inferred as 1 + quantity of true 1 inferred as 0)/ total; error = (3 + 2)/(14 + 2 + 3 + 27)) equals 0.109. This classification error is quite low, so we used the cutoff probability of 0.59 to perform paleobiological inferences using Can harmean as the explanatory (predictor) variable.
When Can min was used as the explanatory variable, we obtained a model with a transition rate α of 0.00052, an intercept estimate of 2.58 ( p-value = 0.032), and an estimate for the coefficient of Can min of −0.49 ( p-value = 0.018).   Again we evaluated the quality of the model by constructing a contingency table in which we inferred the thermophysiological regime of each extant species of the sample using its Can min . We used a cutoff probability of m = 0.32 (the mean probability of being endothermic in our sample), so that an inferred probability lower than 0.32 is evidence for ectothermy and a probability higher than 0.32 is evidence for endothermy: state 0 state 1 prediction 0 15 9 prediction 1 2 20 The specificity (the ratio of quantity of true 1 inferred as 1 on the quantity of true 1; Sp = 20/(20 + 9)) equals 0.690. The sensitivity (the ratio of quantity of true 0 inferred as 0 on the quantity of true 0; Se = 15/(15 + 2)) equals 0.882. The classification error (the ratio (quantity of true 0 inferred as 1 + quantity of true 1 inferred as 0)/total; error = (2 + 9)/ (15 + 9 + 2 + 20)) equals 0.239. This classification error is quite high, so we recomputed a new cutoff probability of 0.22 using the receiver operating characteristic curve. Then we constructed a new contingency table in which we inferred the thermophysiological regime of each extant species of the sample using its Can min and considering that an inferred probability higher than 0.22 is evidence for endothermy: state 0 state 1 prediction 0 15 1 prediction 1 2 28 With this contingency table, the classification error [the ratio (quantity of true 0 inferred as 1 + quantity of true 1 inferred as 0)/total; error = (2 + 1)/(15+ 1 + 2 + 28)] equals 0.0652. This classification error is quite low, so we used the cutoff probability of 0.22 to perform paleobiological inferences using Can min as the explanatory (predictor) variable.
We used these equations and cutoff probabilities and femoral Can harmean and Can min values published by  for Notosuchia to compute the probability of these taxa being endotherms (Table 2).

Discussion
Notosuchia is an extremely diversified group of crocodyliforms. This diversity is particularly striking regarding their diet, suggesting that they occupied various ecological niches (Carvalho and Bertini 1999;Iori and Carvalho 2018). Uruguaysuchidae (Pol et al. 2014), the most basal notosuchians, of which Araripesuchus wegeneri (in our sample) is a representative, range from the Aptian (Araripesuchus gomesii) to the Maastrichtian (Araripesuchus tsangatsangana) (Price 1959;Turner 2006). Several species of this group have been inferred as being omnivorous, or even insectivorous, based on their dental complexity (Sereno and Larsson 2009;Soto et al. 2011;Nieto et al. 2021) and postcranial remains suggest that they had an erect posture (see "Introduction"). Uruguaysuchids were smaller than Sphagesauridae Godoy et al. 2019). Armadillosuchus (also sampled by us) belongs to the large-bodied sphagesaurids group (Melstrom and Irmis 2019). The diagnosis of this clade is based on its peculiar dentition morphology (Price 1950). They show extremely complex manducatory systems, with evidence of "chewing" mechanisms, dental wear, and propalinal movements (e.g., see Ősi 2014; Iori and Carvalho 2018). The foraging abilities of some notosuchians, such as Armadillosuchus, Mariliasuchus, or Malawisuchus, to locate food or water have led some authors to propose the presence of burrowing habits (Gomani 1997;Nobre et al. 2008;Marinho and Carvalho 2009), a behavior that might play a role in thermoregulation (e.g., to search for a cooler shelter during dry periods, as in extant crocodilian species; Campos and Magnusson 2013). This behavior has also been proposed for the sebecosuchian Baurusuchus salgadoensis . Sebecosuchia (to which Iberosuchus and Stratiotosuchus, also sampled by us, belong) were large predators with an erect posture (see "Introduction") and cursorial abilities (Nascimento and Zaher 2010;Riff and Kellner 2011), feeding on large prey, including small sphagesaurids (Godoy et al. 2014). Indeed, their ziphodont teeth (unicuspidated, laterally compressed with serrated carinae) associated with the biomechanical performances of their skull allowed sebecosuchians to effectively handle prey after wounding it (Montefeltro et al. 2020). It is noteworthy that the inferred ectothermic sebecosuchians occupy a niche usually occupied by endothermic theropod dinosaurs (Benson et al. 2013;Zanno and Makovicky 2013).
The ectothermic condition of Notosuchia suggested by  is supported by the results of the present study using larger sample sizes of extant species and a more robust phylogenetic comparative method (PLR). First, the finding that RBC size is related to bone vascular canal diameter (Huttenlocker and Farmer 2017) is supported by our results obtained using a sample size more than twice of that used by these authors. Thus bone vascular canal diameter can be used as a proxy to infer RBC size , and then endothermy (because within tetrapods, RBC size is lower in endotherms than in ectotherms; Snyder and Sheafor 1999). Huttenlocker and Farmer (2017) included bone cross-sectional area (in addition to bone vascular canal diameter) as an explanatory variable in models aimed at explaining the variation of RBC size . However, in their study, bone cross-sectional area did not improve the explanatory power of models and was not retained (Huttenlocker and Farmer 2017; Supplemental Information, "Analysis I, Training Data Set for Extant Taxa"). Unexpectedly, our analyses (using a larger sample size) showed that bone cross-sectional area significantly improves the explanatory power of models and is retained, together with bone vascular canal diameter, in models explaining the variation of RBC size . The fact that the estimate for bone cross-sectional area is always negative (Table 1) indicates that RBC size decreases as bone cross-sectional area increases. Second, notosuchian thermometabolism is inferred here using a larger sample of extant tetrapods (more than three times the sample used by ) and logistic phylogenetic regressions (a method more reliable than those used in previous studies; see below). We are aware of the fact that Huttenlocker and Farmer (2017) showed that histological changes reflect changes in VO 2max better than changes in thermometabolism. However, we have found here that microstructural variation is linked to thermometabolism too. The models constructed to infer the probability of endothermy using vascular canal size as an explanatory (predictive) variable were highly significant, and the classification errors obtained are quite low (6.5% using Can min ). Thus we conclude that we can use these models confidently in paleobiological inference of thermometabolism.
Thermal paleophysiology is an emergent discipline (Cubo and Huttenlocker 2020). It has great potential resulting from the synergy between physiological studies aimed at deciphering the mechanisms of thermogenesis in extant amniotes (e.g., Bal and Periasamy 2020; Jastroch and Seebacher 2020; Grigg et al. 2021) and paleobiological inferences in extinct amniotes using phylogenetic comparative methods (e.g., Cubo et al. 2012Cubo et al. , 2022Legendre et al. 2016;Huttenlocker and Farmer 2017;Olivier et al. 2017;Fleischle et al. 2018;Cubo and Jalil 2019;Faure-Brac et al. 2021; Knaus et al. 2021). Logistic phylogenetic regressions are the third step in efforts performed during the last decade to carry out reliable inferences of thermometabolic status in extinct amniotes. A decade ago, Cubo et al. (2012) inferred bone growth rates using bone histological features and multiple linear regressions tested for significance using permutations in order to circumvent the nonindependence of the observations due to the phylogeny. Considering that bone growth rates are significantly related to RMR in amniotes (Montes et al. 2007), the former was used by Cubo et al. (2012) as a proxy to infer the thermometabolic status of extinct archosaurs. This method was used by Legendre et al. (2013) to infer the bone growth rate and the thermometabolic status of Euparkeria. In a second step, Legendre et al. (2016) adapted Guénard et al.'s (2013) PEMs to perform paleobiological inferences of RMR. This contribution represented significant methodological progress, because paleobiological inference models included the phylogeny (rather than circumventing its effects, as did the preceding method), assuming an evolutionary model (Molina-Venegas et al. 2018). PEMs have been widely used to infer the thermometabolic status of extinct amniotes (Legendre et al. 2016;Olivier et al. 2017;Fleischle et al. 2018;Cubo and Jalil 2019;Cubo et al. , 2022Faure-Brac and Cubo 2020;Faure-Brac et al. 2021;Knaus et al. 2021). Using logistic phylogenetic regressions is a new step in this sequence. This method improves upon the previous approach by using all of the phylogenetic information (rather than a fraction of it, as did PEMs in order to avoid model overfitting). An encouraging sign is that results are congruent in spite of the diversity of methods used to obtain them. Inferring the maximum metabolic rate of Notosuchia using the size of femoral nutrient foramina (Seymour et al. 2012) would be the next promising step to fully understand the thermophysiology of these amazing crocodylomorphs.