A METHOD FOR MAPPING MORPHOLOGICAL CONVERGENCE ON THREE-DIMENSIONAL DIGITAL MODELS: THE CASE OF THE MAMMALIAN SABRE-TOOTH

: Morphological convergence can be assessed using a variety of statistical methods. None of the methods proposed to date enable the visualization of convergence. All are based on the assumption that the phenotypes either converge, or do not. However, between species, morphologically similar regions of a larger structure may behave differently. Previous approaches do not identify these regions within the larger structures or quantify the degree to which they may contribute to overall convergence. Here, we introduce a new method to chart patterns of convergence on three-dimensional models using the R function conv.map . The convergence between pairs of models is mapped onto them to visualize and quantify the morphological convergence. We applied conv.map to a well-known case study, the sabre-tooth morphotype, which has evolved independently among distinct mammalian clades from placentals to metatherians. Although previous authors have concluded that sabre-tooths kill using a stabbing ‘bite’ to the neck, others have presented different interpretations for speciﬁc taxa, including the iconic Smilodon and Thylacosmilus . Our objective was to identify any shared morphological features among the sabre-tooths that may underpin similar killing behaviours. From a sample of 49 placental and metatherian carnivores, we found stronger convergence among sabre-tooths than for any other taxa. The morphological convergence is most apparent in the rostral and posterior parts of the cranium. The extent of this convergence suggests similarity in function among these phylogenetically distant species. In our view, this function is most likely to be the killing of relatively large prey using a stabbing bite.

Abstract: Morphological convergence can be assessed using a variety of statistical methods. None of the methods proposed to date enable the visualization of convergence. All are based on the assumption that the phenotypes either converge, or do not. However, between species, morphologically similar regions of a larger structure may behave differently. Previous approaches do not identify these regions within the larger structures or quantify the degree to which they may contribute to overall convergence. Here, we introduce a new method to chart patterns of convergence on three-dimensional models using the R function conv.map. The convergence between pairs of models is mapped onto them to visualize and quantify the morphological convergence. We applied conv.map to a wellknown case study, the sabre-tooth morphotype, which has evolved independently among distinct mammalian clades from placentals to metatherians. Although previous authors have concluded that sabre-tooths kill using a stabbing 'bite' to the neck, others have presented different interpretations for specific taxa, including the iconic Smilodon and Thylacosmilus. Our objective was to identify any shared morphological features among the sabre-tooths that may underpin similar killing behaviours. From a sample of 49 placental and metatherian carnivores, we found stronger convergence among sabre-tooths than for any other taxa. The morphological convergence is most apparent in the rostral and posterior parts of the cranium. The extent of this convergence suggests similarity in function among these phylogenetically distant species. In our view, this function is most likely to be the killing of relatively large prey using a stabbing bite. With few exceptions, such as convergence on similar body plans in some fast-swimming marine vertebrates (Lingham-Soliar 2016), studies of morphological convergence have targeted specific body parts (e.g. sabres, long necks, or wings) rather than the larger structures or bodies of which they are part. Excepting a few self-evident cases, currently available methods can determine whether convergence is present, but not identify specific regions within the larger structures and quantify the degree to which they contribute to convergence (McGhee 2011). This may impose limitations when assessing whether convergence is restricted to superficial morphological resemblance, or whether it is a consequence of shared selective pressures (Wainwright 2007;Moen 2019). For instance, by combining finite element analysis and geometric morphometrics to investigate humeral shape in fossorial mammals, Sansalone et al. (2020) noted that convergence among digging moles can only be demonstrated when mechanical performance is taken into account with morphology. Almost to the contrary, shared morphologies in sabre-toothed carnivores may obscure a rich functional diversity within the group (Lautenschlager et al. 2020). These examples highlight the difference between morphological convergence, which relates to simple phenotypic similarity, and functional convergence, which may take place even without phenotypic resemblance.
We have recently developed a novel and rapid method to address morphological convergence, deployed with the R function search.conv (Castiglione et al. 2019) embedded in the package RRphylo (Castiglione et al. 2018). This approach permits the identification of the pattern between entire clades or across unrelated species sparsely occurring across a phylogeny. The search.conv function computes the angle between vectors of principal component (PC) scores retrieved from geometric morphometric (GM) data to assess whether two shapes (vectors of PC scores) are morphologically closer (i.e. have a smaller angle between them) than would be predicted by their phylogenetic distance alone. Since principal component analysis (PCA) ordination of GM data represents both affine and non-affine components of shape variation, identifying the PC axes responsible for the pattern of convergence allows us to chart it on the focal biological shapes, then map and quantify the degree to which individual regions contribute to overall convergence on the structures under study. These concepts are central to the new methodology we present here: conv.map.
To illustrate how conv.map works, we have applied the method to address a classic example of convergence: the independent evolution of sabre-tooth morphology in mammalian carnivore lineages. All sabre-tooths are defined by the possession of elongated, laterally flattened upper canines (Emerson & Radinsky 1980;Christiansen 2008), which are widely thought to have been applied in the dispatch of relatively large prey (Akersten 1985;McHenry et al. 2007).
Variability in skull and postcranial morphology, for example, relative length of the canines and robusticity of the forelimbs, have led researchers to posit differences in killing behaviour between sabre-toothed species (Duckler 1997;Christiansen 2008;Slater & Van Valkenburgh 2008;Christiansen 2011;Figueirido et al. 2018). However, most researchers, including those above, have concluded that, notwithstanding these differences, all mammalian sabretooths specialized in killing relatively large prey with slashing bites to the neck, as opposed to suffocation through a clamp-and-hold bite that typifies conical toothed cats (Wroe et al. 2013;Figueirido et al. 2018).
Mechanical modelling has demonstrated that sabretooths perform poorly relative to conical toothed predators in lateral shaking of the skull and jaw adductor driven bites, but are better adapted to resist stabbing, dorsoventral 'bites' driven by neck muscles. On the basis of finite elementbased studies that have included both conical and sabretoothed species, a continuum has been proposed wherein sabre-tooths with increasingly longer canines are characterized by an increasing capacity to resist stabbing forces, but a diminishing ability to resist laterally directed forces (McHenry et al. 2007;Wroe et al. 2013;Figueirido et al. 2018). However, some authors have proposed widely disparate killing and feeding behaviours, including killing bites to the belly as opposed to the neck, and diets comprised of internal organs, or blood, as opposed to meat. These suggestions have been proposed for both placental (Smilodon fatalis) and metatherian (Thylacosmilus atrox) taxa (Merriam & Stock 1932;Akersten 1985;Janis et al. 2020). Perhaps the most divergent hypothesis forwarded in recent times suggests that the taxon characterized by the most hypertrophied canines of any sabre-tooth, the metatherian T. atrox, was not a predator at all, but a highly specialized scavenger (Janis et al. 2020).
If killing and feeding behaviour did differ greatly between sabre-tooth taxa, then we might expect to find that similarities in cranial shape were localized and not shared across the functionally relevant regions of the cranium for all taxa. Our objectives here were therefore to determine how many and to what degree different anatomical regions of the cranium were shared across very distantly related clades, and whether these differences were significant.
To address this question, we applied our methodology to a large sample comprising two placental sabre-tooth families, the single known metatherian sabre-tooth, and a wide range of other carnivorous taxa.

Data preparation
Thirty-two homologous landmarks were sampled manually on 92 specimens. The landmark configuration defines the shape of the dorsal regions of the cranium, including the maxillary bones and the tooth sockets. We excluded the zygomatic arch from the sampling since in fossil specimens it is rarely preserved. We then placed and slid 1000 bilateral semi-landmarks (500 on each side) automatically using the 'Morpho' R package (Schlager et al. 2020). Taxa included a barbourofelid, a dirk-toothed and a scimitartoothed felid among placentals. We also sampled a range of extant conical toothed cats including Neofelis which displays the most morphological features common to extinct sabre-tooths among extant felids (Christiansen 2008). Among metatherians we included the sparassodont Thylacosmilus atrox, as well as dasyuromorphians, and the diprotodontian Thylacoleo carnifex to provide species phylogenetically close to Thylacosmilus. Altogether our data set comprised 49 extant and extinct species (see Melchionna et al. 2020, appendix S1 for details).
Taphonomic distortion was present in two fossil specimens (Barbourofelis fricki, Homotherium serum). We symmetrized these (see Melchionna et al. 2020, appendix S1 for details) using the function retroDeformMesh (Schlager et al. 2018). Procrustes superimposition was applied using generalized Procrustes analysis (GPA), implemented with the procSym function of the R package 'Morpho'. GPA rotates, translates and scales landmark configurations to the unit centroid size, that is, the square root of squared differences between landmark coordinates and Centroid coordinates. After GPA, we applied ordination to the aligned coordinates by means of PCA. The resulting PC scores were taken to represent the shape variables.

Searching for convergence
To perform convergence analysis, we implemented the felid tree embedded in the RRphylo package to add the metatherians included in the analysis (see Melchionna et al. 2020, appendix S1 for details). We used the tree and shape data to feed the RRphylo package function search.conv (Castiglione et al. 2019). This function assesses convergence by testing whether phenotypes in distant clades in a phylogenetic tree are more similar to each other than expected from their phylogenetic distance. Phenotypes are represented by vectors of PC scores derived from geometric morphometric data analysis. Since PC axes have a score equal to zero at their intersection, the phenotypic vectors are calculated in relation to the origin of PC axes (the consensus shape in geometric morphometrics) and the angle they form represents a correlation coefficient. The angle actually spans from 0°to 180°. An angle close to 0°means convergence in shapes, whereas angles around 90°means dissimilarity, and angles close to 180°indicate phenotypes evolving in an opposite direction to the consensus. As the function is also able to test for convergence within one state (or more), we applied search.conv to the sabre-tooth group to verify the convergence in the skull shape within the category.

Relative warp analysis
The use of search.conv enables us to identify species which show convergent phenotypes. For that purpose, using PCA is ideal as it decomposes the variability of the sample into orthogonal axes describing shape variation around the sample mean shape. However, convergence implies large scale, non-affine shape variation which is best inspected by means of PCA of partial warp scores (relative warp analysis, RWA) weighted by a factor a > 0 (a spans from À∞ to +∞; at a = 0 RWA corresponds to PCA so that the affine and non-affine components of shape variation are not separated, Rohlf 1993; Rohlf & Bookstein 2003). In the case study reported here, we performed RWA using the Morpho function relWarp (Schlager et al. 2020) setting the a parameter to 1. Then, we extracted the relative warp scores (RW scores) and the relative warps vectors (RWs).

Mapping convergence
We developed the conv.map function to visualize the relative intensity of convergence on 3D surfaces. Although we recommend using RWA with a parameter set to 1, we note that RW scores or PC scores could be used. Given two phenotypic vectors q 1 and q 2 (i.e. vectors of RW scores or PC scores for any two species or group of species found to converge) of length n, the angle a between them is geometrically equivalent to the correlation coefficient between the vectors (Zelditch et al. 2012; Castiglione et al. 2019). Removing a pair i of corresponding RW scores from both vectors produces the angle a i between the remaining scores q 1[Ài] and q 2[Ài] of length nÀ1. If the removed pair of scores is important to phenotypic similarity a i < a, and vice versa. In conv.map, pairs of corresponding RW scores are removed one pair at a time, and the angle between the vectors of remaining RWs computed each time. After the entire procedure is accomplished, the resulting angles a [1. . .n] are collated into a vector, from the largest to the smallest. This vector would be flat if all RWs are equally responsible for the phenotypic distance between the two original shape vectors. However, RWA decomposes shape variance in orthogonal axes of maximum variation of decreasing importance, so that the vector of ordered a [1. . .n] typically presents one or two inflection points. The first inflection point sets apart RW axes which contribute the most to the pattern of convergence (so that removing any one of them provides an angle a i ( a). The second inflection point, if present, identifies the RW axes responsible for the most important shape differences between the two phenotypic vectors. To locate the first inflection point and therefore select the PC axes responsible for convergence, we applied the function ede in the R package 'inflection' (Christopoulos 2019). ede performs an extreme distance estimator (Christopoulos 2012; Christopoulos 2016) to efficiently locate the inflection points along a curve. By finding the first inflection point the RW conv1,2 matrix of k 9 2 RW set of scores (one set of k corresponding RW axes for each species) is extracted from the q 1 and q 2 vectors. This procedure is analogous to the scree plot criterion commonly used for the selection of 'relevant' axes in PCA. It differs from the scree plot in that rather than selecting the PC axes explaining (cumulatively) most of the variance, it selects the RW conv1,2 set of scores in q 1 and q 2 that maximizes their similarity.
To evaluate the statistical significance of the RW axes selection procedure, conv.map computes the angle a conv between the two vectors in RW conv1,2 . Then, 10 000 angles are computed by randomly selecting from q 1 and q 2 pairs of corresponding RW scores of length k and testing whether a conv is smaller than 5% of the random angles, which is expected to occur if q 1 and q 2 represent convergent shapes and RW conv1,2 effectively represents the subset of RW axes that best account for convergence.
We used RW conv1,2 to retrieve a new landmark configuration using the showPC function in 'Morpho'. The new configurations (one for each compared group or species) are weighted on the variance responsible for convergence. The function automatically reconstructs a 3D mesh by using the ball-pivoting algorithm (Bernardini et al. 1999) as embedded in the vcgBallPivoting function in the R package Rvcg (Schlager & Girinon 2017). Starting from the new surfaces, conv.map estimates the area differences between corresponding triangles of each 3D mesh and the consensus shape mesh of the original RWA (including all species). If the selected groups (or species) are convergent, they should present the same pattern of variation in the same regions of the 3D surface, as compared to the consensus shape. Convergent areas will therefore show small variation (plotted in colour shades), whereas non convergent regions of the 3D surface will be plotted in white. The same procedure could be generalized to >2 shape vectors at once. In this case, the user has to supplement conv.map with a q shape vector for each species, and indicate which species were found to converge. Given j species, the function will calculate all the RW conv1..j matrices (one for each pairwise comparison), and selects the RW axes that appear more than once in the RW conv1..j matrices. By default, if j > 2 shape vectors are provided, comparisons of convergence mapping are plotted against the consensus shape, alongside pairwise comparisons.
To summarize, the conv.map function works as follows: 1. The RW scores responsible for the morphological convergence are selected from q 1 and q 2 shape vectors. 2. The landmark configurations of the selected species (or means of species groups) are reconstructed using only the RW scores and RWs (the RW conv1,2 matrix) responsible for convergence. 3. Triangle meshes of both the landmark configurations referring to q 1 and q 2 ; and consensus shape are interpolated using the ball-pivoting algorithm. 4. Each surface referring to q 1 and q 2 is compared to the consensus shape, and 3D mesh triangles areas differences are computed. 5. The mesh triangle areas referring to RW conv1,2 are calculated and plotted on the 3D surfaces. The function also provides the opportunity to exclude some RW axes from the analysis. That is because, for example, in most cases RW1 mostly captures high-order morphological differences driven by phylogeny and size variation in the sample.
As input data, conv.map needs: the data frame with the RW (or PC) scores of each group or species to be compared (q 1 and q 2 dataset); the matrix of RW (or PC) vectors; the consensus configuration (mshape); the number of the RW (or PC) that will be excluded from the comparison, if needed (exclude = NULL as default setting) ( Table 1).
The function returns the index of the RW axes selected in RW conv1,2, the angle a between q 1 and q 2 , the angle a conv between the two vectors in RW conv1,2, the angle difference a À a conv , the p value for a conv; and plots the 3D surfaces coloured according to the convergence pattern represented by RW conv1,2 ( Table 2).

Applications of conv.map to sabre-tooths
To chart convergence on sabre-tooth mammals, we first computed mean RW scores for all species in the tree. We treated Barbourofelis, Homotherium, Smilodon and Thylacosmilus as sabre-tooths and ran search.conv using all PCs as the multivariate dataset representing shape. Alternative classifications of sabre-tooths failed to find convergence for Neofelis and Yoshi, despite the fact that they are frequently cited as either showing traits shared with sabretooths (the former) or belong to the machairodontinae family (the latter). Then, we used conv.map starting from a RWA with a = 1 to perform a pairwise comparison between sabre-tooths and the consensus shape. We further compared Barbourofelis fricki against machairodont cats (averaging the shapes of Smilodon and Homotherium) and against Thylacosmilus atrox.

RESULTS
The search.conv analysis revealed that among carnivores, only species within the sabre-tooth category displayed significant convergence (mean angle = 43.88°; p = 0.017). The pairwise angle comparison between Barbourofelis, Homotherium, Smilodon and Thylacosmilus is reported in Table 3A. The angle differences (angles computed between scores vectors of the selected RWs against angles between vectors of the non-selected RWs) of all the convergent groups comparisons are significant, while the comparison with the consensus shape is marginally or not significant for each of the sabre-tooths (Table 3A). Barbourofelis fricki and Homotherium serum have the lowest distances among all pairwise comparisons, which is also evident in the convergence plots (Fig. 1). All sabretoothed carnivores are equally distant from the consensus shape (Table 3B). The 3D surfaces comparison reveals a marked similarity in the premaxillary and carnassial areas, and in the occipital region on and around the nuchal crest. The difference pattern against the consensus is similar for all sabre-tooth taxa (Fig. 1).
When Barbourofelis, the average machairodont cat skull, Thylacosmilus, and the consensus shapes are compared to each other, the angle differences for the convergent species are the only statistically significant (or marginally significant) example detected (Table 4A). The average area difference between Barbourofelis fricki and the machairodont cats is the smallest, with surfaces showing great affinity in overall shape. Thylacosmilus atrox is closer to Barbourofelis fricki than to the machairodonts, with a marked similarity in the muzzle area. All groups are distant from the consensus shape (Table 4B, Fig. 2). and convergence may fail to apply when the objects' functional performance is ignored (Sansalone et al. 2020) meaning that the structural properties of the study objects must be acknowledged (Janis et al. 2020). Revealing such particularities is possible by inspecting how different areas of the study objects contribute to the pattern of convergence. This is the aim of the current study. By rendering a visual representation of the relative contribution to The value used to set the colour gradient, computed as the maximum of all differences between each surface and the mean shape convergence of different areas of the biological objects under consideration, conv.map helps to identify the morphological regions which may have important implications for functional convergence. Because of its complexity and multipurpose functioning, the vertebrate skull represents an ideal study object to identify and quantify convergence. In the present study the only example of significant morphological convergence we detected within our sample was the sabre-tooth ecomorph. These four taxa share more anatomical features than any of other carnivore groups, suggesting that the influence of a strong selective pressure drove the iterated evolution of shared morphological features. The conv.map method revealed a range of shared anatomical features of particular importance. These were in the premaxillary area, the carnassial region, and in the occipital region around the nuchal crest, which were common to all sabretooth carnivores despite considerable phylogenetic distances, particularly with respect to Thylacosmilus. These areas are likely to enable the peculiar sabre-tooth killing behaviour, requiring masticatory muscle reorientation, strong neck musculature, low condyles and large gape to effectively use the infamous sabres. In our view this strongly supports the consensus view that despite some anatomical differences and possible functional diversification within sabre-tooths (Lautenschlager et al. 2020), the sabre-tooth morphotype universally confers a broadly comparable capacity to hunt and rapidly kill relatively large prey by applying a stabbing bite to the throat assisted by powerful neck muscles (Emerson & Radinsky 1980;Wroe et al. 2013), although this specialization may have led to their extinction at different times and locations when large prey became less abundant (Piras et al. 2018). Similarities in the rostral and dental areas are likely to represent functional adaptation related to relatively high loadings to which the muzzle was exposed when delivering these stabbing bites. In keeping with this, we found that convergence in sabre-tooths involves the facial region of the skull (particularly in the premaxillary and carnassial area). We also found close similarities in the morphology of occipital area, which was involved in extensive neck muscle depression of the cranium and pull in all sabretooths relative to other morphotypes (Duckler 1997;Christiansen 2008). The neurocranium, nasals, and maxillary regions show no bearing on convergence among sabre-tooths and set Thylacosmilus apart (Table 4, Fig. 2). Intriguingly, the barbourofelid Barbourofelis sits closer to the felid Homotherium than it does to Smilodon, probably because of similarity in their incisor arcades (Biknevicius et al. 1996). We detected greater similarity between Thylacosmilus and Barbourofelis, suggesting that the metatherian sabre-tooth converged more completely on the highly specialized barbourofelid. Previous work has suggested that these taxa converge on the presence of very long canines, a postorbital bar, mandibular flanges and a number of postcranial characters (Prevosti et al. 2012). Real angle, the real angle between the given score vectors; Selected, the angles computed between scores vectors of the selected relative warp vectors (RWs); Others, the angles computed between scores vectors of the non-selected RWs; Angle difference, the difference between Selected and Other angle values; p-value, the statistical significance of the difference in angle between the selected and nonselected RWs.

DISCUSSION
Our findings support the hypothesis that the patterns of convergence on the sabre-tooth cranium provided by the conv.map method indicate first-order functional significance, although secondary functional diversity may take place among sabre-tooths (Lautenschlager et al. 2020). The absence of convergence in the neurocranium and nasal area is consistent with previously identified features that are not indicative of prey choice or killing method (e.g. the smaller brain of the metatherian Thylacosmilus).
The most obvious example here is that the most phylogenetically distant of sabre-tooth taxa, Thylacosmilus, is strongly convergent on placental sabre-tooths across regions of the cranium widely considered to be functionally significant. However, a recent study addressing the mechanical performance of Thylacosmilus (Janis et al. 2020) proposes a very different ecology. Their FEA-based analysis suggested a greater capacity to exert pulling forces in Thylacosmilus than in Smilodon. It was inferred F I G . 1 . Visualization of the pairwise comparisons of Barbourofelis, Homotherium, Smilodon, Thylacosmilus, and the consensus configuration. The colour gradient indicates area differences between the two surfaces. Blue indicates no difference (scale bar rescaled into the range 0-1). In each case, differences between the two taxa are displayed on a reconstruction of the taxon named on the left. on this basis that Thylacosmilus did not deploy killing behaviour common to other sabre-tooths that distinguished them from conical-toothed predators and that it may have been a scavenger. We do not contest the possibility that Thylacosmilus was better adapted to exert pulling behaviour, because the only other FEA-based study to include these two taxa did not apply a neck-pulling load case (Wroe et al. 2013). On the other hand, we do note that the study of Janis et al. (2020) did not include any conical-toothed taxa by way of comparison, whereas the earlier study (Wroe et al. 2013) included the conicaltoothed Panthera pardus. Wroe et al. (2013) found that Thylacosmilus was better adapted to deliver a stabbing bite that was reliant on head depressing neck musculature than was Smilodon, but that both were superior in this respect to Panthera pardus. We suggest that Thylacosmilus was better adapted to deliver both neck-driven head pulling and depressing functions, but without applying a head-pull load case to a conical-toothed cat as well as to other sabre-tooths there is no mechanics-based evidence to indicate whether a capacity to deliver a head pull is a further defining feature of sabre-tooths or not. As well as FEA-based analysis, an intriguing feature of Thylacosmilus identified by Janis et al. (2020), determined on the basis of dental microwear patterns, was that Thylacosmilus did not habitually consume meat or bone. This raises the possibility that the metatherian sabre-tooth's diet was largely restricted to soft internal organs such as the heart, lungs and liver. These characteristics would make Thylacosmilus truly unique among mammalian scavengers, which typically do consume both bone and meat, as well as internal organs. However, large extant mammalian predators commonly consume these internal organs preferentially over meat and bone (Schaller 1972). We consider it more likely that Thylacosmilus may have concentrated on the internal organs of prey it had killed, rather than the carcasses of animals killed by other large South American carnivores (e.g. phorusrhacid birds), which were unlikely to have eaten only meat and/or bone and left the internal organs untouched. Sabre-tooths are characterized by relatively weak jaw closing muscles but large head depressors (Wroe et al. 2005;Christiansen 2011). Thylacosmilus is extreme with regard to both features (Wroe et al. 2013), and we contend that together with the evidence for a diet comprised largely of soft internal organs, these characteristics point to an even more extreme example of sabretooth killing and feeding behaviour, as opposed to an entirely divergent ecomorph and diet.

CONCLUSION
The use of three-dimensional models has revolutionized the study of fossils in both palaeontology and palaeoanthropology (Cunningham et al. 2014). In morphometric analyses, the implementation of geometric morphometric techniques based on 3D configurations has driven the development of new protocols and software suitable for 3D surfaces. These new technologies allow researchers to virtually restore and realign disarticulated elements (Gunz