Rostral morphology of Spinosauridae (Theropoda, Megalosauroidea): premaxilla shape variation and a new phylogenetic inference

ABSTRACT Spinosauridae presents an extensive geographical and temporal distribution, with records from Gondwana and Laurasia, and a temporal range from Barremian (~129 Ma, Lower Cretaceous) to Cenomanian (~95 Ma, Upper Cretaceous). To date, 13 species were described, besides several specimens identified at a broader taxonomic level. One of the most notable cranial features of spinosaurids is their elongated rostrum (hypertrophied premaxilla-maxilla). Leastwise five species possess preserved premaxillae: Angaturama limai, Baryonyx walkeri, Cristatusaurus lapparenti, Oxalaia quilombensis and Suchomimus tenerensis, besides materials tentatively attributed to Spinosaurus aegyptiacus. We studied the premaxillae shape of 10 specimens of the above-mentioned species and other materials through geometric morphometrics, reviewing diagnoses and morphological descriptions. Clear allometric and phylogenetic signals could be identified by ordination methods. We carried out a phylogenetic analysis to test spinosaurid relationships, by inclusion of new landmarks-characters from the premaxillae into a published tetanuran character-taxon matrix. The phylogenetic inference recovered C. lapparenti outside Baryonychinae, which was composed by B. walkeri and S. tenerensis. Spinosaurinae was recovered as (A. limai, (O. quilombensis, (MSNM V4047, Irritator challengeri))). Our results suggest that the premaxillae provide useful phylogenetic information and that the inclusion of landmarks-characters improves our knowledge of this enigmatic clade.

The highly differentiated morphology possibly relates to a habitat preference distinct from that of other theropod dinosaurs, with some species probably spending a lot of time in the water or coastal regions (Amiot et al. 2010;Sales et al. 2016), presumably to obtain a more diversified diet made up mostly of fish as proposed on direct evidence (partially digested fish scales - Charig and Milner 1986;fish vertebrae preserved in premaxillary alveolus -Dal Sasso et al. 2005) and by the neuroanatomy of some representatives of the group (Schade et al. 2020). The jaw morphology and the presence of integumentary mechanoreceptors, possibly evolutionary converges to pike conger eels, are also effective biomechanical adaptations for grabbing prey in low-light aquatic environments (Vullo et al. 2016). In addition, the tail morphology indicates that some spinosaurids probably used it as a propellant in aquatic locomotion and foraging ). However, this palaeoecological interpretation has been challenged recently (e.g. Henderson 2018;Hone and Holtz 2021). Hone and Holtz (2021), based on a morphofunctional study of spinosaurids, found that the idea of a 'highly specialized aquatic predator' is not well supported, since behavioural and palaeoecological interpretations require systematic steps and well-defined comparisons (see Hone and Faulkes 2014 for the theory on behaviour in fossil organisms). Other studies concluded that the diet of spinosaurids also included small to medium-sized ornithopod, sauropods, and pterosaurs Milner 1986, 1997;Buffetaut and Suteethorn 1999;Buffetaut et al. 2004). Therefore, they were probably a generalist or opportunistic group of theropods, and not exclusively piscivorous, or even not fully aquatic dinosaurs (Henderson 2018;Hone and Holtz 2019;Schade et al. 2020;Hone and Holtz 2021;Souza et al. In press).
Although extensively studied in several palaeobiological aspects, the spinosaurids are still considered one of the most enigmatic groups of dinosaurs (Hone and Holtz 2017) because of the poor knowledge regarding their unusual functional morphology, and the rarity of well-preserved fossils of complete or nearly complete skeletons -magnified by the loss of the fairly complete holotype of the family-type species, Spinosaurus aegyptiacus Stromer 1915, during the World War II (Smith et al. 2006). The rarity and scarcity of the spinosaurid fossil record, when compared to other theropod groups (e.g. tyrannosaurids), as well as the fragmented nature of known specimens (Bertin 2010;Hone et al. 2010;Hone and Holtz 2017), is somewhat surprising since spinosaurids were more prone to burial and fossilisation because they supposedly did forage and feed near water bodies (see discussion in Hone et al. 2010). Despite this, the relatively long-time and palaeogeographical span of spinosaurids confirm it as a successful dinosaur group (in terms of durability and distribution, Carrano et al. 2012) that occupied specialised palaeoecological niches, having clear aquatic adaptations, but probably covering coastal palaeoenvironments distinctly from other theropods and thereby avoiding competition for resources with other terrestrial large-bodied predators (Holtz 1998;Amiot et al. 2010;Hone and Holtz 2019;Sales et al. 2016;Holtz 2017, 2021).
The morphology of the premaxilla has been of great importance for understanding the spinosaurid taxonomy. Although the clade has a relatively scarce fossil record (Bertin 2010;Hone et al. 2010), there are at least 11 premaxillae described in the literature. Five spinosaurid species known so far have the anteriormost part of the rostrum preserved: Angaturama limai Campos 1996, Baryonyx walkeri Charig andMilner 1986, Cristatusaurus lapparenti Taquet and Russel 1998, Oxalaia quilombensis Kellner et al. 2011 and Suchomimus tenerensis Sereno et al. 1998. In addition, there are a premaxilla of a juvenile individual (Lakin and Longrich 2019), besides materials tentatively related to Spinosaurus cf. aegyptiacus (Milner 2003;Dal Sasso et al. 2005), and a rostrum attributed to Spinosaurus maroccanus (Taquet and Russel 1998). However, Sp. maroccanus was considered as nomem dubium by Sereno et al. (1998) and referred to Sp. aegyptiacus by Hendrickx et al. (2016). Although, these materials cannot be compared with the holotype, or even the neotype of Sp. aegyptiacus, as there is no overlapping of cranial structures (see discussion in Evers et al. (2015) and Sales and Schultz (2017), and discussion below). It is important to highlight that the holotypes of some species (i.e. Angaturama limai, Cristatusaurus lapparenti and Oxalaia quilombensis) are composed exclusively of the premaxillae, and other well-preserved taxa (e.g. Baryonyx walkeri, Suchomimus tenerensis) also preserved this structure. Notwithstanding, despite the existence of a relatively well-preserved known premaxillae fossil record, a quantitative morphology study of the rostrum shape variation in spinosaurids has not been done yet.
Here we analysed the morphological variation of the spinosaurids premaxillae, based on two-dimensional geometric morphometrics methods (GMM), through a set of anatomical landmarks, commonly used as a quantitative method to assist in the study and description of biological shape (Rohlf and Slice 1990;Slice 2007;Cooke and Terhune 2015;Palci and Lee 2019). We also provide a review and discussion of the spinosaurid rostral fossil record. Our main goal was to verify if the observed morphological variability gained by GMM can be related in a phylogenetic context. Two approaches using the same dataset (configurations of premaxillae in different perspectives and depicting different morphologies) were carried out. First, we explore the premaxillae shape variation through ordination methods and multivariate regression, decomposing the morphospace into fewer variables and reviewing some of the spinosaurid morphological descriptions of those species. Then, we used phylogenetic morphometrics (PM) to input the phylogenetic signal that the premaxillae configurations (treated as characters) add to test the phylogenetic relationships within Spinosauridae.

Anatomical nomenclature
The anatomical nomenclature and abbreviations for cranial and dental structures adopted here follow the previous terminology of Dodson (2003), Dal Sasso et al. (2005); Carrano et al. (2012), Hendrickx and Mateus (2014); Hendrickx et al. (2015b), and Sales and Schultz (2017). Table 1 summarises the ten analysed specimens (Figures S1-S3) and indicates specimens that were analysed in person. Six of the specimens represent five named spinosaurid species and one does not have a more inclusive taxonomic identification. The remaining three specimens were previously proposed to belong to the Spinosaurus genus, but will be treated and discussed in this work as 'Spinosaurinae indet' because of the controversies regarding their referral to this genus. The specimen MNHN SAM 124, previously erected by Taquet and Russel (1998) as Spinosaurus maroccanus Russel (1996) was considered nomen dubium by Sereno et al. (1998) and Carrano et al. (2012). The same specimen was also related to Sp. aegyptiacus by Hendrickx et al. (2016). The African specimens NHMUK PV R 16420 and MSNM V4047 were tentatively referred to Sp. cf. aegyptiacus respectively by Milner (2003) andDal Sasso et al. (2005), but this was challenged by Evers et al. (2015) and Sales and Schultz (2017) because they cannot be directly compared with Sp. aegyptiacus holotype.

Species and studied specimens
Seven premaxillae from the left side of the rostrum were used to digitise the LMs and sLM (Figures S1-S3). The exceptions (the mirrored right side was used) were Suchomimus tenerensis, of which only the right premaxilla is well preserved, FSAC-KK-7281, which has only the right premaxilla preserved, and MSNM V4047, whose right premaxilla was used because of the proposed closing of the sixth alveolus in the left side (see Sales and Schultz (2017), and discussion below). In our phylogenetic searches, among the three Spinosaurinae indet., only the specimen MSNM V4047 was included (further details in Table 1) since it represents a well preserved specimen, mainly the rostral portion; the codification for this specimen in the matrix is the same as provided by Sales and Schultz (2017). Since the monophyletic status of the Baryonychinae clade has been challenged (Sales and Schultz 2017), we consider the 'classic species' originally allocated in this subfamily (see Sales and Schultz 2017) in our subsequent ordination analyses as 'Baryonychinae'.

Selection of specimens and landmarks
We explored the premaxillae shape variation of both ventral and lateral views. The dorsal view was not analysed since the morphological information in this view (contour of the anterior and lateral part of the rostrum) can also be obtained by evaluating the premaxillae in the ventral view and therefore would add redundancy/ dependence in the analysed data. In the main analysis, we included seven specimens (Table 1) and all twenty-seven LMs (descriptions on Table 2). For the ventral view, besides the main analysis, an additional exploratory analysis (AEA) was performed including three more specimens. Two poorly preserved specimens, NHMUK PV R 16420 and Angaturama limai USP GP/2T-5, were not included in the main analysis because they present a high degree of diagenetic alteration on the bone surface, affecting mainly the ventral border of the premaxillae (see Kellner and Campos 1996; and Sales and Schultz 2017 for comments on A. limai). FSAC-KK-7281 was not included because, due to ontogeny or taphonomy (Lakin and Longrich 2019), it does not have all structures that allow the acquisition of the total set of LMs (described in Table 2). Thus, in our AEA, all ten specimens (Table 1) were included, but we maintained only seven LMs (detailed in Table 2) that could be digitised in the three specimens excluded from the main ventral view analysis. This AEA aimed to explore the shape variation using a small set of LMs to verify whether the poorly preserved specimens (A. limai, NHMUK PV 16420 and FSAC-KK-7281) could be putatively associated to a more inclusive taxonomic level through the structure seen in the morphospace. Table 1. Species and specimens studied in this work and details of its inclusions on landmark-based analyses. References from which some images were used to digitise landmarks are presented. Abbreviations: Additional Exploratory Analysis (AEA), Lateral View GMM analysis (LV), Phylogenetic Morphometrics (PM), Ventral View GMM analysis (VV).

Landmark
Anatomical description 1 a Anterior-most point of interpremaxillary suture 2 a Interpremaxillary contact of the secondary palate at its anterior-most edge 3 Interpremaxillary contact of the secondary palate at the anteriormost margin of the median gap that divides the two premaxillae 4 a Intersection between the sagittal suture of the anteromedial processes of the maxillae with the transverse plane related to the distal border of the Lpm 6 b (sixth left premaxillary tooth) 5 Lateral

Landmark and semi-landmark description
In both anatomical views analysed here, the LMs were defined based on the preserved structures that have topological homology (i.e. represent an anatomical locus) among the spinosaurids premaxillae. In the ventral view, 27 LMs were defined to explore the rostral geometry and to access the overall shape and position of the premaxillary alveoli (ventral view LMs descriptions are presented in Table 2 and Figure 1(a)). For the lateral view of the rostrum, six LMs were defined to capture part of the rostral morphology and general shape of the anterior-most margin of the skull (lateral view LMs descriptions are presented in Table 3 and Figure 1(b)). Additionally, we digitised twenty semi-landmarks (sLM) to capture the geometry of the anterior-most rostral morphology between LMs 5 and 6 (see digitised LMs and sLM of each specimen in Figures S1-S3).

Data analysis
We used TpsUtil v.1.76 (Rohlf 2015) software to organise the datasets and build the Thin-Plate Spline (TPS) files that were subsequently used. Then, based on specimens' photographs, the LMs and sLM were digitised, and a scale factor standardisation was carried out using TpsDig2 v.2.31 (Rohlf 2015). Twenty sLM were digitised, in the lateral view, by tracing a curve between LMs 5 and 6 and recovered as twenty equidistant points distributed along this curve. All LMs and sLMs digitisation's were performed by MBSL to standardise the data collection, avoiding individual error input.
To remove undesirable variables from the data, the Generalised Procrustes Analysis (GPA) was performed using MorphoJ v.1.7a (Klingenberg 2011) andTpsRelw v.1.7 (Rohlf 2015). The GPA is a basic procedure in GMM and minimise differences (variance) among individuals' LM configurations through three basic operations: (i) translation; (ii) rotation and (iii) scaling (Rohlf and Slice 1990;Viscosi and Cardini 2011;Zelditch et al. 2012). The set of anatomical LMs (i.e. specimens' configurations) was rotated via Generalised Least Squares (GLS) to minimise the relative distance among each anatomical LM. By doing the Procrustes fit (i.e. GPA), a new set of LMs coordinates, so-called 'Procrustes residuals', is achieved and these coordinates retain the shape information that can be properly compared (Hammer and Harper 2008;Zelditch et al. 2012).
Along with Procrustes fit, another procedure that also precedes the exploratory GMM and cladistic analyses was carried out to properly 'slide' the set of sLM (only to lateral view). We used two functionalities of TpsUtil following the steps described in Zelditch et al. (2012): (1) the 'append tps curve to landmark' option was selected to 'transform' the sLM along the curve into LM coordinates, and (2) 'make sliders file' was selected, to establish which digitised points are sLM of the digitised curve (and, therefore, sliding the curve in the rotation operation of GPA). By doing it, the entire curve is treated as homologous instead of each individual sLM, which is a proper way of fitting sLM in Procrustes residual (Perez et al. 2006;Zelditch et al. 2012;Gunz and Mitteroecker 2013). The matrixes with Procrustes residual coordinates from lateral view (using the Procrustes fit sliding sLM) were saved from TpsRelw, and all subsequent analyses (ventral and lateral) were carried out using MorphoJ.
Two exploratory analyses were performed. First, a Principal Components Analysis (PCA), calculated from the covariance matrix derived from Procrustes residuals matrixes, was performed to find new variables (i.e. Principal Components, PCs) that retain most of the variance. Then, the individual scores of each specimen in the most relevant PCs were included in a multivariate regression using premaxilla size as the independent variable. The goal of the PCA was to ordinate the data to reduce dimensionality and check whether there is any clear structure in individual projections onto morphospace decomposed by the most relevant new variables (PCs). The multivariate regression was done to explore the allometry of spinosaurids premaxillae, correlating the shape variation explained by size variation (independent variable). The choice of PCs to be considered as a proxy to the shape variation and multivariate regression was selected using the broken stick model (Jackson 1993). The multivariate regression considered the normalised size (LN Centroid Size) as an independent variable versus scores of PC as a dependent variable, separately for each configuration (ventral, ventral additional, and lateral views). During the allometric analyses, we performed a  Table 2 for LMs descriptions); (b) Premaxilla in lateral view: six landmarks (filled blue circles; see Table 3 for LMs descriptions), and twenty semi-landmarks (open circles). Photographs and drawings of MSNM V4047 modified (mirrored) and redrawn from Dal Sasso et al. (2005) Table 3. Description of the landmarks digitised on the lateral view of the premaxilla.

Landmark
Anatomical description 1 Mesial-most point of the Lpm 6 (sixth left premaxillary tooth) at the lateral surface of the premaxilla 2 Central point of Lpm 5 at the lateral surface of the premaxilla 3 Central point of Lpm 4 at the lateral surface of the premaxilla 4 Central point of Lpm 3 at the lateral surface of the premaxilla 5 Antero-ventral extremity of premaxillae edge 6 Postero-dorsal premaxilla rim at position of the mesial-most point of Lpm 6 permutation test with 10,000 rounds to quantify how the changes in shape could be explained by size. Besides that, we performed an ANOVA (Kruskal-Wallis) test and extracted a box-plot of LN centroid sizes using PAST v.4.06 (Hammer et al. 2001) to be able to visualise and test size differences between medians of 'Baryonychinae' and Spinosaurinae species. All the matrices with GMM data are detailed in the supplemental material (SM), and available in Lacerda et al. (2021). Since a strong correlation was found between PC scores (shape variable) and normalised CS (see results), the Procrustes residuals matrix was further analysed in order to explore the allometry in the data and address the following questions: are the shape differences found between both subfamilies the result of size-related differences (i.e. allometric shape differences) or are there shape differences independent of the size effect? To answer that, we used the method developed by Mitteroecker et al. (2004)to explore relative warps in a size-shape space. This method consists on creating individual scores for the log-transformed CS to perform an ordination analysis. These scores, so-called Common Allometric Component (CAC), encompass variance on shape regarding the logCS values. The regression of shape variables for calculation of CAC scores were corrected for the mean CS value disregarded of subgroups (i.e. the mean value for all spinosaurids analysed). Other components, called Residuals of Shape Coordinates (RSC), are, therefore, calculated after 'size is regressed out' (Mitteroecker et al. 2004, p. 685), following the ordination analysis of the data as a regular PCA. The RSC retain variances that describe shape variation independent of size, while the CAC illustrated most of the shape variance due to size (i.e. allometry). Thus, we compared CAC, RSC and logCS values in scatter plots similarly to PCA out plot (described above). Also, we checked for correlations between these variables in order to find (and describe) the data allometry. All these complementary analyses were performed using the software PAST.

Characters and terminal taxa
We used the theory and algorithms for the Phylogenetic Morphometrics analysis under parsimony (PM) presented by Catalano et al. (2010), Catalano (2011, 2016), and Catalano and Goloboff (2012), which provide a reliable bridge between landmark-based GMM and cladistics analysis (providing a Farris optimisation [Farris 1970], and therefore treating LM configurations as additive characters; see Catalano et al. 2010;Goloboff and Catalano 2011 for an introduction to the methods).
In our cladistic analysis, we combine two premaxillae configurations of LMs (lateral and ventral) to a PM along with standard (discrete) characters from the original matrix of Carrano et al. (2012) and Sales and Schultz (2017). Thus, each configuration was input as a new character into the character-taxon matrix. The identification of spinosaurid OTUs included in our PM analysis is presented in Table 1. Since our main question is related to the relationships and premaxillae shape variation within spinosaurids, we did not include premaxillae configurations of non-spinosaurid megalosauroids, mainly because several LMs would not apply to these species. We assume the Spinosauridae monophyly (and do not test it) based on current known evidence (e.g. Sereno et al. 1998;Allain et al. 2012;Carrano et al. 2012;Evers et al. 2015;Sales and Schultz 2017;Malafaia et al. 2020).
Both LMs configurations were inserted into a recent Tetanurae characters-taxon matrix that includes a good sampling of spinosaurid taxa provided by Carrano et al. (2012), and modifications provided by revisions of specimens and characters made by Evers et al. (2015) and Sales and Schultz (2017). Since it is the most recently published, the character matrix provided by Sales and Schultz (2017) was used to input the LM characters. As our analysis focus on Spinosauridae internal relationships, several non-spinosaurid taxa (except for Piatnitzkysaurus floresi Bonaparte 1979, Eustreptospondylus oxoniensis Walker 1964, and Torvosaurus tanneri Galton and Jensen 1979 were removed from the analysis (total of 50 OTUs removed). We chose to keep the above-mentioned three relatively well preserved non-spinosaurids megalosauroids as they represent closely related outgroups of the well-supported monophyletic spinosaurid. We did not perform an extensive review of this character-taxon matrix because we consider it was already properly done (Carrano et al. 2012;Evers et al. 2015;Sales and Schultz 2017). However, by removing 50 OTUs, 304 characters became uninformative, which were deactivated using the WinClada v.1.61 (Nixon 1999(Nixon -2002. Our searches were carried using the 50 remaining informative characters (see SM). All additional editions in the character-taxon matrix with standard characters were done using Mesquite v.3.61 (Maddison and Maddison 2015). The files with Procrustes residual coordinates of both configurations (ventral and lateral) were then imported and incorporated into the data matrix (available in SM), which were analysed in TNT v.1.5 (Goloboff and Catalano 2016).

Tree searches
The complete character-taxon matrix consists of ten OTUs (Table 1)  The 'basal' megalosauroid piatnitzkysaurid Piatnitzkysaurus floresi was used as outgroup to polarise the characters and establish the tree root, however, we also kept other two megalosauroids in our matrix (Eustreptospondylus and Torvosaurus). Before proceeding with tree searches, both LM configurations (i.e. Procrustes residual coordinates of ventral and lateral premaxilla) were realigned via RFTRA (resistant fit theta-rho analysis). This is a necessary step since the GLS alignment is more susceptive to outliers and other alignment issues, such as the Pinocchio effect (Rohlf and Slice 1990; Palci and Lee 2019) which are not recommended in a phylogenetic context (Catalano et al. 2010;Catalano and Goloboff 2012). Therefore, different OTUs were used as reference taxon in a series of RFTRA realignments, prior to tree searches to verify if by changing the OTU to perform the realignment would lead to different topologies and/or discrepant tree scores (see details in SM).
Once the matrix matching discrete characters and LM data were realigned, heuristic tree search methods were run in TNT. All discrete characters were treated as equally weighted (=1) and multistate discrete characters were assumed as non-additive (Fitch optimisation). We applied a weight cost of 5 to the landmark-based characters, since each premaxilla configuration, represents a different anatomical view and, as expected, different morphological information is achieved, encompassing several traditionally discrete characters. Hereupon the contribution of each configuration was not considered equivalent to discrete characters in our analysis. Since we use sLM in our analysis (restricted to the lateral configuration), we use TNT's lmark wts command to reduce the contribution of sLM to 1 twentieth of the LM contribution, so that the total weight of the 20 sLM curve equals the weight of a single LM.
Traditional heuristic searches were performed in TNT using a tree memory of 24,000 trees and by generating starting trees through Wagner search (Farris 1970). The random seed was defined as 10 and we performed 500 replicas followed by TBR branch-swapping algorithm, saving up to 200 trees per replica.
Zero-length trees were collapsed after the search. The individual LM optimisation was performed using a 6×6 grid (TNT default). With the MPTs obtained, we mapped the 'ancestral configurations' of premaxillae along the spinosaurids phylogeny, comparing the premaxillae shape changes along the spinosaurids evolution. Clade support was obtained through Bremer support, which was calculated including trees up to five steps worse than the best tree.

Ventral view
In the ordination plot of the premaxillae configuration in ventral view, the two most relevant components retained, combined, 69.6% of the variance (broken stick in Figure S4). The projections onto the morphospace represented by the PC1 vs PC2 plot ( Figure  2(a)) reveals no overlap among individuals of both subfamilies, with Spinosaurinae specimens retaining positive scores for PC1 and 'Baryonychinae' retaining negative scores for this axis, so that PC1 summarises evolutionary differences between the two subfamilies. Therefore, the application of Phylogenetic Morphometrics using this LM configuration is reasonable. The PC2 summarises differences among individuals within each clade, including also ontogenetic differences, but its eigenvalue is rejected by the broken stick equation ( Figure S4) and, therefore, is a less relevant source of morphological variation. More details on the shape variation related to taxonomic, ontogenetic and allometric differences are presented in the discussion section.
The PC1 (41.1% of retained variance) summarises changes in the following portions/structures of the premaxillae (Figure 2 Spinosaurinae indet. premaxilla-premaxilla medial contact, besides changes in the internal and in the distal secondary palate median gap; (3) proportion and position of the alveoli 1, 3, 4, 5 and 6, especially on the labial and distal edges and (4) maximum lateral expansion of the secondary palate. The relative position of the alveoli 4 and 5 changes mainly in the PC1 axis.
Positive PC1 scores (all spinosaurines) are associated with premaxillae that have: (1) the shape of a posteriorly tapered structure; (2) the posterior-most portion more posteriorly positioned in the rostrum; (3) the maximum premaxillary medial constriction concentrated in the portion of the fourth premaxillary alveolus; (4) smaller alveoli 1, 4, 5 and 6; and (5) pronounced diastemata between the alveoli 3 and 4 and between the alveoli 5 and 6 (warped outline drawings in Figure 2(d)). The spinosaurines MSNM V4047 and MNHN SAM 124 better reflect the above-described condition, presenting large diastema and relatively small alveoli 1, 4, 5 and 6. This condition can be observed to a lesser extent in Oxalaia quilombensis, but it projects more distant to MSNM V4047 and MNHN SAM 124 because of the different PC2 scores.
Negative PC1 scores (all 'baryonychines') are associated with premaxillae that: (1) are anteroposteriorly shorter; (2) present the posterior portion of the premaxilla contact more anteriorly positioned in the rostrum; (3) have a large mid-lateral expansion that is more evident in the portion between the alveoli 3 and 4; (4) have a more laterally expanded secondary palate; (5) have larger premaxillary alveoli 1, 4, 5 and 6; and (6) more uniformly spaced alveoli (warped outline drawings in Figure 2(d)).
The PC2 (28.5% of retained variance) summarises changes in the following portions of the premaxillae: (1) maximum lateral expansion of the premaxilla, at the outer edge of the rostrum; (2) anteroposterior position of the maximum width point of the secondary palate; (3) mediolateral position of the premaxilla outer edge of the rostrum (based on alveolus 7); (4) anterior-most point of interpremaxillary suture; and, (5) variations in the position of the edges of the alveoli 2, 3, 4 and 5, mainly in their lingual and proximal portions (alveoli 3, 4 and 5). To a lesser extent than PC1, PC2 also helps to illustrate (6) the well-developed diastema present in the spinosaurines between the alveoli 5 and 6; and (7) the position of the lingual border of the alveolus 3.
Negative score values for the PC2 axis are associated with premaxillae that have their greater lateral extension and maximum width point of the secondary palate respectively more anteriorly and more posteriorly positioned on the rostrum, besides wireframes alveoli with smaller spacing. All 'baryonychine' specimens have negative scores for PC2, which is not the case for the juvenile Cristatusaurus lapparenti MNHN GDF 366, which has the most extreme positive scores for PC2 in the entire sample, followed by the spinosaurine Oxalaia quilombensis, both of which have the greater lateral expansion of the secondary palate more posteriorly positioned. MSNM V4047 is the spinosaurine with more negatively PC2 scores, and, among 'baryonychine' specimens, the mature form of Su. tenerensis MNN GDF 501 has the more extreme negative PC2 scores. Regarding the two specimens of Cristatusaurus, the morphological changes explained by PC2 may be associated with its ontogeny, since one specimen (MNHN GDF 366) represents a juvenile individual.

Ventral view -additional exploratory analysis
In our additional GMM analysis with a subset configuration of seven LMs and inclusion of three specimens (LMs reduced to seven and specimens increased to ten, Table 1), the PCA retained 71.5% of the variance in the two major PCs (Figure 2(b)). The other components were disregarded for being well below that indicated by the broken stick equation ( Figure S4).
The PC1 (53.8%) describes the following variation in the morphology of the premaxillae (Figure 2(e)): (1) position of the anterior-most point of interpremaxillary suture; (2) changes in the anterior-most edge of the secondary palate; (3) mediolateral width of the posterior portion of the premaxilla and (4) variations in the labial position of the fourth premaxillary alveolus.
The PC2 (17.8%) summarises the following morphological variations: (1) mediolateral width of the posterior portion of the premaxilla; (2) anteroposterior length of the premaxilla; (3) middleposterior end of the secondary palate in the region of premaxillary intra-contact; and (4) position of the lingual edge of the premaxillary alveoli 4 and 5; and (5) variations in the labial edge of the premaxillary alveolus 3. The first and second transformations described above, however, should be considered with caution because they possibly result from the mediolaterally taphonomical compression/erosion observed in Angaturama limai and NHMUK PV R 16420.
The structure observed in the morphospace of the two most relevant PCs reflects what was recovered by the analysis using the complete set of LM configurations. As in the complete analysis, the PC1 summarises differences in shape due to different taxonomic groups (i.e. evolutionary pattern). However, in the PC2 no distinction is observed between spinosaurines and 'baryonychines'.
The spinosaurines Oxalaia quilombensis, MSNM V4047, and MNHN SAM 124 plotted close together in the morphospace, presenting positive values for PC1, due to their premaxilla bearing a posterior tapering, a greater anteroposterior projection, and a posterior region more anteroposteriorly elongated. The spinosaurines Angaturama and NHMUK PV R 16420 retain slightly negative scores of PC1 plotting apart from other spinosaurines and closer to the 'baryonychines', specially Baryonyx, that also plots near to the average of PC1. The premaxilla of Angaturama and NHMUK PV R 16420 shows an elevated degree of erosion, which possibly explains its negative scores both in PC1 and in PC2: Angaturama, for example, has the width of the posterior portion of the premaxilla similar to the average of the entire sample, but the premaxilla at Lpm 3 is much more compressed (due to the lack of the lateral bone surface that resulted on a probably wrong placement of LM 10). Both Angaturama and NHMUK PV R 16420 present the usual features that were associated with the spinosaurine premaxilla 'morphotype' (as observed in the analysis using the complete set of LMs): (1) the premaxilla tapers on the posterior region (more perceptible in NHMUK PV R 16420); (2) the most extreme lateral point of the premaxilla related to the premaxillary alveolus 5 is more posteriorly positioned in the rostrum; and (3) the premaxilla is more latero-medially compressed mainly towards the premaxillary alveoli 4 and 5. These features are also observed in the juvenile specimen FSAC-KK-7281, so that it projects in the morphospace closer to the spinosaurines and has the most positive score for the entire sample (Figure 2(e)).

Lateral view
In the ordination plot of the premaxillae configuration in lateral view, the two most relevant PCs retained, combined, 82.9% of the variance. As recovered in the ordination plot of the ventral view, the individual projections onto the morphospace represented by the PC1 vs PC2 plot (Figure 2(c)) reveal no overlap among individuals of both spinosaurid subfamilies. Spinosaurine specimens retain only negative scores for PC1, whereas 'baryonychines' retain only positive scores for PC1.
The PC1 (69.3% of retained variance) summarises the following morphological shape changes in the premaxillae (Figure 2 Negative PC1 scores reflect the 'morphotype' of the spinosaurine premaxillae that has a more expanded anterodorsal rim, a posteriorly positioned alveolus 6, and a gradual dorsoventral constriction towards the posterior portion of the rostrum (that results from the more ventrally positioned posterior portion of the dorsal rim in association with the more dorsally positioned ventral rim (Figure 2 (f)). In the morphospace, Oxalaia quilombensis plots almost close to the average, so that the above-mentioned 'morphotype' is more easily observed in the other spinosaurine individuals (NHMUK PV R 16420, MSNM V4047, and MNHN SAM 124). Interestingly, the unidentified spinosaurid FSAC-KK-7281 retains similar scores with spinosaurines, plotting close to MNHN SAM 124 (also in PC2).
Positive PC1 scores reflect the 'morphotype' of the 'baryonychine' premaxillae, which is anteroposteriorly shorter and dorsoventrally expanded towards the posterior portion, which markedly differs from the posteriorly constricted and anteriorly expanded premaxilla of spinosaurines (Figure 2(f)). The dorsoventral expansion is a consequence of the sagittal crest found in some 'baryonychines', which explains the almost median score values for PC1 of Baryonyx walkeri, that does not present a sagittal premaxillary crest. Suchomimus tenerensis presents the most positive PC1 score and has a sagittal crest that is more developed than that of the Cristatusaurus lapparenti specimens (MNHN GDF 365 and MNHN GDF 366). The dorsoventral expansion in the posterior portion of the 'baryonychine' premaxilla also is related to the more irregular ventral surface of the rostrum, that is strongly concave between Lpm 3 and 5, a condition that differs from the straighter to slightly concave ventral edge observed in spinosaurines posterior to Lpm 3.
The PC2 (13.6% of retained variance) summarises the following morphological changes in the premaxillae: (1) development (expansion or retraction) of the anterodorsal-most region of the rostrum (also partially explained by PC1); (2) shape (convex or concave) of the ventral rim of the rostrum anterior to the Lpm 4; (3) position (more anterior or posterior) of the posterior portion of the premaxilla (at sixth alveolus) relative to the fifth alveolus; and (4) position (more anterior or posterior) of the premaxillary alveoli 3, 4 and 5.
Positive PC2 scores are associated with (1) a more expanded anterodorsal rim of the premaxilla; (2) a more ventrally positioned anterior extremity of the premaxilla; and (3) a more anteriorly positioned alveolus 6 (which indicates a more elongated posterior portion of the rostrum), with negative PC2 scores reflecting the opposite condition. The three spinosaurine specimens often assigned to the Spinosaurus genus (MSNM V4047, MNHN SAM 124 and NHMUK PV R 16420) are widespread along the PC2 axis: MSNM V4047 has the most negative PC2 scores, followed by NHMUK PV R 16420, while MNHN SAM 124 has the most positive scores among these three specimens. The two individuals of Cristatusaurus lapparenti, MNHN GDF 365 and MNHN GDF 366, also have distinct PC2 scores, as a reflection of the MNHN GDF 366 bearing the anterodorsal rim of the premaxilla less expanded and the anterior extremity of the premaxilla more dorsally positioned. Suchomimus tenerensis and Oxalaia quilombensis have almost average PC2 scores, while the unidentified spinosaurid FSAC-KK-7281 has positive PC2 scores, plotting close to MNHN SAM 124 (as in PC1).

Premaxilla size vs. shape
The regression analysis for the ventral view configuration detected a correlation (p = 0.0125) between shape variables (PC1 scores, selected by broken stick model) and the log of centroid size (independent variable) (Figure 3(a)). The allometry was responsible for 71.0% of the total shape variation considering the different scales of the sample composed of five spinosaurid genera. In this regression, the 'baryonychines' and spinosaurines plot clearly as distinct forms, being 'baryonychines' smaller specimens, restricted to negative scores of PC1, whereas spinosaurines are large specimens, retaining positive PC1 scores (Figure 3(a)). Kruskal-Wallis test found a clear difference between the median of 'baryonychines' and spinosaurines centroid sizes (p = 0.03389). The difference in size between both is evidenced by the box-plots (Figure 3(a)).
The regression analysis for the additional ventral view configuration (AEA) does not detect correlation (p = 0.3027) between shape variables (PC1 scores, selected by broken stick model) and the log of centroid size (Figure 3(b)). The allometry was responsible for only 12.9% of the total variation of the shape when considering the different scales of the sample (Figure 3(b)). There is no clear segregation between 'baryonychines' and spinosaurines as they superimpose, largely due to the spinosaurine Angaturama limai that plots within the 'baryonychines' range of size. However, 'baryonychines' retain the most negative PC1 scores in this plot and do not reach sizes as large as those of spinosaurines (see box-plots in Figure 3(b)). Interestingly, the size of the FSAC-KK-7281 specimen, which represents a small individual, is more related to the forms of 'baryonychines', but its PC1 score is more similar to those of the spinosaurines.
The regression analysis for the lateral view configuration does not detect a correlation (p = 0.3097) between shape variables (PC1 scores, selected by broken stick model) and the log of centroid size (Figure 3(c)). The allometry was responsible for only 14.1% of the total variation of the shape when considering the different scales of the sample. As in the analysis of the ventral view, for the lateral view 'baryonychines' occupy a different part of the morphospace than spinosaurines, but Kruskal-Wallis test does not find a clear difference between the median of 'baryonychines' and spinosaurines centroid sizes (p = 0.04581). 'Baryonychines' retain positive scores for PC1 while spinosaurines retain negative scores (Figure 3(c)). As for the AEA of the ventral view, the specimen FSAC-KK-7281 retains PC1 scores similar to those of spinosaurines, but is more related to the size of 'baryonychines'.

Common allometric component
In the ordination analysis of the premaxillae configuration in ventral view, CAC was found to be correlated with size (R2 = 0.499) and structured the individual projections (Figure 4(a)) similarly to PC1 scores (see Figure 3(a)). The RSC1 account for 49.3% of the variance. The shape-size projections onto CAC vs. RSC1 morphospace (Figure 4(a)) show that there is a clear shape difference between both clades related to the allometric component (CAC), although subtle differences also can be described by residual shape differences (RSC1) as well. Using individual CAC scores as proxy for allometric shape allows to differentiate spinosaurines and 'baryonychines' (Kruskal-Wallis, p = 0.034), but RSC1 scores do not (Kruskal-Wallis, p = 0.157). Therefore, although both, allometric and non-allometric shape differences, can describe differences between spinosaurines and 'baryonychines', it is statistically clearer (according to the conjuncture described in Dushoff et al. 2019) that the allometric differences are more pronounced. Finally, we found PC1 scores (Figure 3(a)) are strongly correlated to CAC scores (R2 = 0.87; p = 0.0015) and the variation assuming a linear model is almost isometric (slope of −0.91) ( Figure S5). Therefore, most of the variance in the data is related to allometry and the major shape differences found between both subfamilies are due to size differences of these groups (i.e. because of Spinosaurinae exacerbated gigantism).
In the AEA of ventral view, CAC was found to be correlated with size (R2 = 0.499) and structured the individual projections (Figure 4 (b)) similarly to PC1 scores (see Figure 3(b)). RSC1 account for 47.3% of the shape variance. The projections onto CAC vs. RSC1 morphospace (Figure 4(b)) show that there is, in general, a clear shape difference between both subfamilies related to the allometric component, with most of the spinosaurine specimens retaining higher positive CAC scores. The exception is the A. limai USP GP/2T-5, which has negative CAC scores (score = −0.048), closer to 'baryonychines' specimens (see Figure 4(b)). Contrary to the results found using the complete LMs dataset in ventral view, using individual CAC scores as proxy for allometric shape did not find clear differences between spinosaurines and 'baryonychines' (Kruskal-Wallis, p = 0.062). As in the complete LMs dataset in ventral view, RSC1 scores do not segregate both subfamilies as well (Kruskal-Wallis, p = 0.157). We also found that PC1 scores (Figure 3(b)) are strongly correlated to CAC scores (R2 = 0.84; p = 0.00019) and the variation assuming a linear model is close to isometry (slope of −0.69) ( Figure S5).
Beyond that, the ordination analysis in lateral view found similar results to the ventral view analysis. CAC was found to be correlated to size (R2 = 0.498) and structured the individual projections (Figure 4(c)) similarly to PC1 scores (see Figure 3(c)). RSC1 account for 63.1% of the shape variance. The projections onto CAC vs. RSC1 morphospace (Figure 4(c)) show a clear shape difference between both clades related to the allometric component, with the spinosaurine specimens retaining positive CAC scores and 'baryonychine' individuals with negative scores. Using individual CAC scores as proxy for allometric shape allows to differentiate spinosaurines and 'baryonychines' (Kruskal-Wallis, p = 0.036), but RSC1 scores do not (Kruskal-Wallis, p = 0.091). We also found PC1  Figure 3. Allometric regressions (size vs. shape) and box-plot of spinosaurid premaxilla. (a) Ventral view regression, plot of log centroid size (x axis) and shape (y axis) extracted from PC1 and box-plot of ventral premaxillae centroid size; (b) Additional ventral view regression, plot of log centroid size (x axis) and shape (y axis) extracted from PC1 and box-plot of ventral premaxillae centroid size; (c) Lateral view regression, plot of log centroid size (x axis) and shape (y axis) extracted from PC1 and box-plot of lateral premaxillae centroid size. Silhouettes represent 'Baryonychinae' and Spinosaurinae subfamilies and were extracted from http://phylopic.org and available under Creative Commons licences.   (Figure 3(c)) are strongly correlated to CAC scores (R2 = 0.91; p = 0.000061) and the variation assuming a linear model is close to isometry (slope of −0.78) ( Figure S5).
In both analyses that included the indeterminate Spinosauridae FSAC-KK-7281 (AEA, Figure 4(b), and lateral view, Figure 4(c)), this specimen retained the higher negative RSC1 score and positive CAC score (similar to most Spinosaurinae individuals). However, this specimen has one of the lowest centroid sizes, closer to values seen in 'baryonychines' specimens and small spinosaurines such as A. limai USP GP/2T-5.
Finally, the general structure of the morphospace of CAC vs. RSC1 and logCS (Figure 4), as well as the correlation tests ( Figure  S5), allows to identify that most relevant component (PC1) of previously ordination analysis retain strong allometric signal. Therefore, even though subtle non-allometric shape differences are notable (RSC1 shape variation in Figure 4), most of the shape variance that account to differentiate spinosaurines and 'baryonychines' is reflect of the size effect (i.e. allometry).

Phylogenetic inference
All tree searches executed after performing RFTRA realignments of data with different terminals and different landmark-based characters (ventral and lateral, ventral only and lateral only) found similar results, most of them retaining only one MPT (see Table S1, for a detailed comparison of realignments results). The best score was obtained with the realignment with Oxalaia quilombensis as reference (score of 86.615 when landmark-based characters from both views were included) and the results were similar to other realignments realised. The same topology, but with different scores, is found using C. lapparenti, Su. tenerensis, and MSNM V4047 to perform RFTRA realignments prior to tree search. Results slightly differed only when B. walkeri was used to perform the RFTRA realignment (Table S1). Three MPTs with 75 steps were achieved when only traditional discrete characters were used ( Figure S6), resulting in a polytomy within Spinosaurinae. A single MPT (84.614) was achieved in 4 of 5 searches (Table S1) when both ventral and lateral LM configurations were included (see Figures 5 and 6 for tree topology and mapping of ventral and lateral premaxillary shape transformations of the ancestors and terminals).
In the MPT found, Cristatusaurus lapparenti was found as the first branch of divergence within Spinosauridae, with all remaining taxa nested in a well-supported clade divided into two subclades. One nested Suchomimus tenerensis and Baryonyx walkeri, consisting, therefore, in the Baryonychinae clade (sensu Sereno et al. 1998). The second one consists of Spinosaurinae and includes Angaturama limai as the first branch of divergence followed by Oxalaia quilombensis as the sister taxon of a clade composed of MSNM V4047 and Irritator challengeri (Figures 5 and 6).

First premaxillae description and taxonomic validity of Cristatusaurus lapparenti
Two articulated premaxillae from the Elrhaz Formation (Niger), illustrated and initially described as anterior portions of the mandible by Taquet (1984), are the first spinosaurid premaxillae to be published in the literature. The premaxillae were recognised as such by Kellner and Campos (1996) and used as the holotypic specimens of the new species Cristatusaurus lapparenti by Taquet and Russel (1998). The latter referred additional premaxillae and other fragmentary cranial and postcranial material to this taxon.
Shortly after, Sereno et al. (1998) considered C. lapparenti as a nomen dubium or a junior synonym of Baryonyx walkeri pointing that it does not have sufficient diagnostic features, which has been accepted (e.g. Buffetaut and Ouaja 2002;Sues et al. 2002;Dal Sasso et al. 2005;Bertin 2010;Carrano et al. 2012;Ibrahim et al. 2014;Hendrickx et al. 2016). Kellner et al. (2011) defended the validity of C. lapparenti, but speculated that it could be congeneric or even co-specific with Suchomimus tenerensis as they have similar morphology. This synonymy was defended by Hendrickx et al. (2016) arguing that these taxa have similar stratigraphic distribution and premaxilla morphology, but they pointed out the priority name of C. lapparenti, which should be maintained as senior synonym of Su. tenerensis. However, Hendrickx et al. (2016) argued that C. lapparenti does not have autapomorphies and preferred to treat the genus and species as nomina dubia. The first study to reassess the C. lapparenti premaxillae, in the light of current spinosaurid knowledge, was carried out by Sales and Schultz (2017), who pointed out autapomorphies to support the validity of C. lapparenti: a very convex secondary palate that extends below the line of the premaxillary teeth, a morphological condition that differs from both Suchomimus tenerensis and Baryonyx walkeri.
The distinctiveness of C. lapparenti is also supported by our phylogeny (Figures 5 and 6), in which it is recovered as the basalmost taxon due to the absence of crown striations (character 142 [0]) and by the presence of premaxillary teeth serrations (character 149[0]), two conditions that differ from other spinosaurids (Sales and Schultz 2017). Therefore, we agree with Sales and Schultz (2017) on the validity of C. lapparenti and point, based on our results, that it also has the following morphological differences when compared to Su. tenerensis: (1) the first premaxillary alveolus of C. lapparenti is more posteriorly positioned on the rostrum; (2) the fourth premaxillary alveolus of C. lapparenti is slightly larger, which is also (3) more spaced from the third premaxillary; (4) the edge of the third premaxillary alveolus is slightly more ventrally placed in C. lapparenti; and (5) C. lapparenti has the posterodorsal region of the premaxilla slightly less expanded, with a less prominent sagittal crest. However, if more complete materials are found that confirm that C. lapparenti is in fact a senior synonym of Su. tenerensis as previously highlighted by Hendrickx et al. (2016), such morphological differences could be related to ontogenetic or even intraspecific variations (to details see Hendrickx et al. 2016).

The 'crocodile-mimic' -Suchomimus tenerensis
Since the original description and nomenclatural act by Sereno et al. (1998), the African spinosaurid Suchomimus tenerensis has been suggested as possible congeneric to the European Baryonyx walkeri in the literature (e.g. Buffetaut and Ouaja 2002;Sues et al. 2002;Candeiro 2015). A new combination 'Baryonyx tenerensis' was assumed by some authors (e.g. Sues et al. 2002;Candeiro 2015) but is not generally accepted. Nevertheless, Hendrickx et al. (2016) note some differences between Suchomimus tenerensis and Baryonyx walkeri related, mainly, to the size of first and second premaxillary alveoli.
Indeed, our GMM exploratory analysis of the premaxillae in ventral view reveals clear shape similarities between Baryonyx walkeri, Cristatusaurus lapparenti, and Suchomimus tenerensis ( Figure  2). However, when comparing Baryonyx walkeri with Su. tenerensis, our analysis revealed the following differences: (1) the labial border of the third premaxillary alveolus is more distally positioned in the rostrum in Su. tenerensis; (2) the middle-posterior contact at secondary palate is slightly medio-posterior in B. walkeri; (3) the edge of the third premaxillary alveolus is more distal and lingually positioned in Su. tenerensis; and (4) the most extreme lateral expansion of the secondary palate is anteriorly positioned in B. walkeri. Additionally, the exploratory analysis of the lateral view shows that (5) the anterodorsal rim of the premaxilla is more elevated in B. walkeri; (6) the dorsal rim of the premaxillae in B. walkeri does not present the gradual dorsoventral expansion of the rostrum towards the posterior direction (sagittal crest) that can be observed in Su. tenerensis (also visible in C. lapparenti, although less developed); (7) in B. walkeri the ventral premaxilla border is continuously concave, with a smooth curvature, differing from Su. tenerensis and C. lapparenti (adult specimen MNHN GDF 365) that presents an abrupt curvature and a posterior constriction near the third premaxillary alveolus; (8) in B. walkeri the third and fourth premaxillary alveoli are slightly posteriorly located in the rostrum. Therefore, the distinctiveness of B. walkeri and Su. tenerensis and, hence, the validity of both species is supported by our results. Since Stromer's (1915) firstly described Spinosaurus aegyptiacus, this giant theropod dinosaur has been viewed as an enigmatic taxon due to its unique morphology and the loss of the holotype (Smith et al. 2006;Bertin 2010;Sales and Schultz 2017). Recently, after the designation of a neotype (Ibrahim et al. 2014) and the description of new materials tentatively related to this species (e.g. Ibrahim et al. 2020), more information was added about morphological aspects of spinosaurids. Two nearly complete rostra, NHMUK PV R 16420 and MSNM V4047, described respectively by Milner (2003) andDal Sasso et al. (2005), have been referred to as Spinosaurus cf. Sp. aegyptiacus. Both specimens (considered here as Spinosaurinae indet.) are based on the rostral region of the skull and have preserved a premaxillary 'rosette', a character state previously recognised by Milner (1986, 1997), Kellner and Campos (1996), and Sereno et al. (1998) in spinosaurids. However, MSNM V4047 has six teeth in the premaxilla whereas NHMUK PV R 16420 has seven (Dal Sasso et al. 2005). Kellner et al. (2011) comment that this difference could be due to different ontogenetic stages; however, the ontogenetic stage is unknown for both specimens as well as for the lost holotype and the neotype of Sp. aegyptiacus. According to Hendrickx et al. (2016), the number of teeth in spinosaurid specimens can vary not only ontogenetically, but also intra-specifically or even within the same individual, considering that Baryonyx walkeri has six premaxillary alveoli in the right portion of the rostrum, and seven in the left portion (Charig and Milner 1997). Sales and Schultz (2017) also propose that, the Lpm 6 in MSNM V4047 is actually Lpm 7, and the real Lpm 6, as well as the Rpm 7 were closed at some point during its ontogeny, which could explain the longer diastema on the left premaxilla. In our analyses, we follow this suggestion by Sales and Schultz (2017) and consider that, in the right premaxilla of MSNM V4047, the sixth alveolus is the preserved one. Lakin and Longrich (2019) comment on differences in the dorsal region of the rostrum of these individuals: the posterior portion of the MSNM V4047 rostrum is 'straighter' when compared to the dorsal concavity noticeable in NHMUK PV R 16420; the shape and the position of the external nostril are substantially different, with MSNM V4047 having a smaller and anteriorly tapered nostril whereas that of NHMUK PV R 16420 is larger and posteriorly constricted. According to Lakin and Longrich (2019), the differences observed may represent sexual dimorphism as well as intraspecific or ontogenetic variations, but, more probably, a result of taphonomy and/ or crushing.

The African spinosaurines
Other material composed of the premaxillae has already been described in the literature and referred to as Spinosaurus maroccanus (MNHN SAM 124) by Taquet and Russel (1998). Sereno et al. (1998) consider Spinosaurus maroccanus as a nomen dubium because it does not have a satisfactory diagnosis. This proposition is accepted by Buffetaut and Ouaja (2002), Milner (2003), Dal Sasso et al. (2005, Bertin (2010), Carrano et al. (2012), Ibrahim et al. (2014), and Hendrickx et al. (2016), who also considered Spinosaurus aegyptiacus as the only valid species for the genus. Hendrickx et al. (2016) further commented on the lack of a comprehensive list of materials related to Spinosaurus maroccanus and a reliable description of it, suggesting the need for a reappraisal of the specimens related to this species. A different opinion is presented by Evers et al. (2015), who consider Spinosaurus maroccanus as a synonym of Sigilmassasaurus brevicollis. Sales and Schultz (2017) comment that MNHN SAM 124 and MSNM V4047 probably refer to the same taxon, however, they do not consider linking these specimens to Spinosaurus (contra Milner 2003;Dal Sasso et al. 2005) because the type-series material did not preserve this structure for comparison. We agree with this statement, and, conservatively, prefer to not formally identify the specimens MNHN SAM 124, NHMUK PV R 16420 and MSNM V4047 as Spinosaurus because this hypothesis can only be verified after new cranial materials are described.
Our results reinforce the fact that there are morphological differences between the specimens MSNM V4047, MNHN SAM 124 and NHMUK PV R 16420. Regardless of this, we suggest caution regarding the homology of the premaxillary alveoli, as any misinterpretation can lead to an erroneous interpretation (see Sales and Schultz 2017 for comments on MSNM V4047). In general, the differences among these specimens have been suggested as variations caused by dimorphism, intraspecific variations, or other causes (Sales and Schultz 2017;Lakin and Longrich 2019). In the PCA, the individuals MNHN SAM 124 and MSNM V4047 have similar PC1 scores from both, ventral and lateral view analyses (Figures 2(a) and 2(c)), but they have markedly different PC2 scores in the lateral view analysis, being the juvenile FSAC-KK-7281 much more closely positioned to MNHN SAM 124 (Figure 2(c)). NHMUK PV R 16420 plots closer to MSNM V4047 than to MNHN SAM 124 in the lateral view analysis (Figure 2(c)) but plots far from both these specimens in the ventral additional analysis (Figure 2(b)), but this is probably influenced by the taphonomical mediolateral compression of the specimen.
The similarity between MSNM V4047 and MNHN SAM 124 as previously noted (Lakin and Longrich 2019) is evident and confirmed in our ventral view analysis (Figure 2(a) and Figure 2 (c)). However, in the lateral view, it is notable that these specimens differ as MNHN SAM 124 has (1) the posterior region of the premaxilla more dorsoventrally constricted; (2) the anterior portion of the rostrum more expanded and robust; (3) the ventral rim of the premaxilla convex (it is concave in MSNM V4047); and (4) the third premaxillary alveolus more distally positioned. The above-mentioned features are mostly explained by non-allometric variations (Figure 4(c)), and could be diagnostic or even result from sexual dimorphism; however, the taxonomic relationships of these specimens still need to be better elucidated. Interestingly, the juvenile FSAC-KK-7281, which plots close to MNHN SAM 124 in the lateral view, also has several similarities with both MSNM V4047 and MNHN SAM 124 as they plot close in the AEA of the ventral view (Figure 2(b)). It is also evident that the differences among FSAC-KK-7281 and these spinosaurines result from non-allometric components (possibly due to ontogeny, (Figure 4(b) and Figure 4(c)). Considering the comparison above, it is evident that linking those specimens to a single taxon is still very speculative.

South American spinosaurines
The spinosaurid fossil record in South America comes from several localities, but mainly from the northern portion of Brazil (e.g. Kellner and Campos 1996;Bittencourt and Kellner 2004;Medeiros et al. 2014;Sales et al. 2017;Aureliano et al. 2018). Two taxa from Brazil have preserved the premaxillae: Angaturama limai, from the Aptian of the Araripe Basin (Kellner and Campos 1996), and Oxalaia quilombensis, from the Cenomanian of the Alcântara Formation, São Luís Basin (Kellner et al. 2011).
The only formally described specimen of A. limai, USP GP/2T-5, presents a high degree of wear on the labial edges of the rostrum (Kellner and Campos 1996) and bears evidence of post-mortem compression in the left palate (Sales and Schultz 2017). Our results of the AEA in ventral view (Figure 2(b)) suggest the palatal region of A. limai suffered little or no taphonomic distortion as the width of the posterior portion is similar to the average of the sample. Therefore, we assume the width of the premaxilla of A. limai is compatible with that of other spinosaurines, which reinforce that the original diagnosis from Kellner and Campos (1996) of an anterior part of the skull strongly compressed laterally is not valid, in agreement with Sales and Schultz (2017). Results also suggest there is a notable medio-lateral compression of the A. limai premaxilla at the third premaxillary alveoli, but we interpreted it as an artefact caused by the loss of the ventral border of the premaxilla, and part of the third premaxillary alveoli, which prevented us from collecting the correct position of the labial border of the third premaxillary alveolus ( Figure S3). The ventral morphology of the premaxilla of A. limai is roughly similar to other spinosaurines, with differences being mainly size effect (Figure 4(b)). Despite not included in the lateral view analysis due to the poor preservation of the anterior and dorsal margins, the most significant difference is observed in the lateral view: the presence, in A. limai, of a prominent sagittal premaxillary crest that reaches the anterior end of the snout, which represents a reliable diagnostic character for this species (Kellner and Campos 1996;Sales and Schultz 2017).
Regarding Oxalaia quilombensis premaxilla MN 6117-V, Kellner et al. (2011) comment that it presents a diastema between the third and fourth premaxillary alveoli and another one between the fifth and sixth alveoli, which differentiates it from the non-spinosaurine B. walkeri, Su. tenerensis, and C. lapparenti in which these diastemata are absent. In the MSNM V4047 specimen, the diastema between the fifth and sixth premaxillary alveoli is larger than in O. quilombensis (Kellner et al. 2011), but similar to those of specimen MHNM SAM 124 (Sales and Schultz 2017), besides that, O. quilombensis retain other spinosaurine diagnostic features (Sales and Schultz 2017). This species has a more subtle premaxillary constriction posterior to the third premaxillary alveolus than that observed in other spinosaurines (Sales and Schultz 2017). Sales and Schultz 2017) commented that O. quilombensis has a premaxillae morphology similar to that observed in some African spinosaurines; concluding that other Brazilian spinosaurines (from Araripe Basin) somehow retains 'intermediate' morphological features with both baryonychines and spinosaurines. Our result for the ventral view, shows that Oxalaia plots as far from 'baryonychines' as any other spinosaurines, which is explained mainly by the differences in size among them (Figure 4(a) and Figure 4(c)).
O. quilombensis plots much closer to 'baryonychines' than any other spinosaurines when the premaxillae are analysed in lateral view. This results from the shape of the dorsal rim of the rostrum of O. quilombensis, which is straight and more dorsoventrally expanded (similar to Baryonyx walkeri) than that observed in MSNM V4047, MNHN SAM 124 and NHMUK PV R 16420. However, this structure is less expanded than in Su. tenerensis, and C. lapparenti, which present a sagittal premaxillary crest. Also, O. quilombensis differs from B. walkeri in the posterior region of the ventral rim of the premaxilla, which is not expanded in the former. Analysed in lateral view, O. quilombensis retains PC1 scores intermediate of those between B. walkeri and NHMUK PV R 16420.
Recently, Smyth et al. (2020) argued that O. quilombensis does not have sufficient diagnostic autapomorphies to be distinguished from 'Sp. aegyptiacus' MSNM V4047 and referred it to this species. According to Smyth et al. (2020) the morphological differences that O. quilombensis has relative to 'Spinosaurus' specimens are variations due to polymorphisms. Based on our results, we consider that there are distinct features in O. quilombensis that, in fact, differentiate this specimen from the spinosaurine specimens (i.e. MSNM V4047, MNHN SAM 124 and NHMUK PV R 16420). However, some considerations need to be made. Compared to specimen MSNM V4047, we note three differences in O. quilombensis. First, the premaxillary constriction between the third and fourth alveoli is subtly more expanded in O. quilombensis, but we do not consider this as a diagnostic feature, since it could result from phenotypic plasticity in the premaxilla shape; however, the premaxilla of O. quilombensis is also wider in the posterior portion than that of MSNM V4047 (it is also wider than that of MNHN SAM 124, and NHMUK PV R 16420 that are similar in size to O. quilombensis; Figure 3(c)). In addition, as previously demonstrated by Kellner et al. (2011), O. quilombensis has ornamentations in the most anterior part of the palate that are absent in MSNM V4047 and other spinosaurids; this feature is considered as the only putative diagnostic feature of O. quilombensis by Smyth et al. (2020). Third, the rostral most portion of the unfused laminae that form the anteromedial process of the maxilla (anteriorly to the mesial edge of premaxillary alveolus 4) does not contact medially in MSNM V4047 but do contact in O. quilombensis. This may represent a new diagnostic feature for this species.
The spinosaurine O. quilombensis was also diagnosed as having two replacement teeth in the position of the third premaxillary alveolus (medially, in contact with the edge of the palate), which is absent in other spinosaurines according to Kellner et al. (2011) and Sales and Schultz (2017); this was rejected and considered a trivial diagnosis by Smyth et al. (2020). It is interesting to note that, in others theropods, such as the tyrannosaurid Tarbosaurus bataar Maleev (1955), double replacement teeth are part of the normal tooth replacement process. As demonstrated by Hanai and Tsuihiji (2019), there is an alternate replacement pattern between odd-and even-numbered alveoli, in which the presence of the second replacement tooth indicates the alveolus is at a more advanced stage of tooth replacement. Hanai and Tsuihiji (2019) also revealed that, during the growth of the new teeth, they migrate from a more medial position (adjacent to the edge of the palate) to a more lateral position, eroding the lingual surface of the 'old tooth' and shedding it out as it grows. Based on this, we reject the diagnose based on double replacement teeth for O. quilombensis (in agreement with Smyth et al. 2020) and consider that it represents a typical stage in the process of tooth replacement in spinosaurids.
Based on the above comments, we agree with several of the interpretations of Smyth et al. (2020) regarding the similarity between O. quilombensis and specimen MSNM V4047, but, as already pointed for MSNM V4047, we disagree that it should be synonymised with Spinosaurus aegyptiacus, as there are no overlapping structures that can be directly compared. Moreover, as previously indicated, there are two potential diagnostic features for O. quilombensis: (1) roughness of the palate (Kellner et al. 2011;Sales and Schultz 2017) and (2) shape of the anteromedial maxillary process. There are also other differences, such as the less constricted posterior portion of the premaxilla, especially in lateral view. The CAC analyses also reinforce that Oxalaia, despite quite similar in size with MNHN SAM 124, differs from other spinosaurines due to non-allometric factors from both ventral (Figure 4(a)) and lateral views (Figure 4(c)). As they also have a distinct biogeographic distribution, we agree with the validity of O. quilombensis until new materials from Africa or Brazil elucidate these taxonomic issues.

Small spinosaurid premaxilla FSAC-KK-7281
Recently, Lakin and Longrich (2019) described the right premaxilla of an undetermined juvenile specimen from the Cenomanian of Kem Kem Beds (Africa). Although the authors kept the specimen with a broad taxonomic identification (Spinosauridae indet.), they have compared the morphology of this specimen with spinosaurine premaxillae (e.g. NHMUK PV R 16420 and MSNM V4047), and comment that the small spinosaurid FSAC-KK-7281 would be related to spinosaurines based on the morphotype of specimen MSNM V4047. In our AEA of the ventral view (Figure 2(c)), the indeterminate FSAC-KK-7281 is projected on the extremely positive scores of PC1, along with MNHN SAM 124 and MSNM V4047, that is, it occupies the spinosaurine morphospace. The correlation of FSAC-KK-7281 with spinosaurines could also be seen in the lateral exploratory analysis (Figure 2(b)): its PC1 score is very similar to those of MNHN SAM 124 and MSNM V4047, and its PC2 score is very close to that MNHN SAM 124. However, the allometric analyses suggests that there are differences between FSAC-KK-7281 and other spinosaurines, and indicates that those are not caused by size differences as they have similar CAC scores (Figure 4(b) and Figure 4(c)). But this result is not definitive as it may have been influenced by uncertainties on defining the position of LMs that are related to the width of the anterior and posterior portions of FSAC-KK-7281 premaxilla.
Therefore, it seems safe to affirm that the tentative allocation of FSAC-KK-7281 within spinosaurines, as pointed out by Lakin and Longrich (2019), is corroborated by the premaxillary morphology explored in our ordination analysis (Figure 2

Allometric shape changes in the spinosaurid premaxillae
Spinosaurids represent large-bodied predators, being among the largest known theropods (Dal Sasso et al. 2005;Hendrickx et al. 2015a). According to estimates from Dal Sasso et al. (2005) for MSNM V4047 rostrum, this specimen would have a total skull length of 175 cm, rivalling in size with other gigantic theropods, such as Tyrannosaurus rex Osborn (1905).
The variation of the premaxillary shape in both ventral and lateral views seems to be, in most part, related to the size of the specimens (at least for the ventral view data, 71% of the predicted variations are explained by PC1). This interpretation is reinforced by the CAC analyses in which strong correlations were recovered between the common allometric component and the specimen sizes (logCS; Figure 4). However, the allometric component scores of both Cristatusaurus lapparenti individuals are close in all CAC analysis performed (Figure 4) even though their sizes (logCS) are clearly different. Therefore, as MNHN GDF 366 is a juvenile individual whereas MNHN GDF 365 is an adult (Kellner and Campos 1996), we interpret our results as having recovered variance in the premaxillae shape associated with evolutionary allometry. This can be perceived more clearly analysing the results of the common allometric component in the ventral view configuration (Figure 4 (a)): C. lapparenti specimens (both juvenile and adult) retain negative CAC scores whereas 'baryonychine' specimens have negative CAC scores closer to zero and spinosaurine specimens retain positive CAC scores. This gradient of shape change is strongly correlated to the centroid size (Figure 4(a)) and, therefore, can be interpreted as a morphocline influenced by size (i.e. allometry). In the lateral view, the association with evolutionary allometry is less evident, because Suchomimus tenerensis has a CAC score that is closer to both C. lapparenti specimens than to Baryonyx walkeri (Figure 4(c)), possibly due to the presence of a sagittal premaxillary crest in the first two taxa. However, as observed in the analysis of the ventral dataset, it is evident that in the lateral dataset CAC also separates the 'baryonychine' specimens (negative CAC scores) from the spinosaurine specimens (positive CAC scores; Figure 4(c)).
In the ventral view, the most notable effect of size increase could be described using shape differences of CAC as the following (Figure 4(a)): (1) mesial-distal displacement and reduction of the alveolus 4; (2) labio-lingual displacement and reduction of the alveolus 5; (3) medio-lateral displacement of the interpremaxillary contact of the secondary palate at its anterior-most edge; and (4) lateral displacement of the intersection between the sagittal suture of the anteromedial processes of the maxillae. These shape changes lead to, as the size of the specimens increases, lateral retraction of the snout in the position of the alveolus 4 and 5, reduction of posterior teeth (making the premaxilla acquire the shape of a posteriorly tapered structure), and formation of diastemata between alveoli 3 and 4 and between alveoli 5 and 6. The posterior tapering of the premaxilla as specimen size increase was also recovered when the less complete set of LM configurations were analysed (AEA; Figure 4(b)).
In the lateral view, the effect of size increase described by the shape differences of CAC include (Figure 4(c)): (1) dorsal displacement of the posterior ventral margin of the premaxilla; (2) posteroventral displacement of the posterior dorsal rim of the premaxilla; and (3) distal displacement of the alveolus 6 and mesial displacement of alveolus 5. These shape changes lead to, as the size of the specimens increases, dorsoventral constriction of the posterior portion of the premaxilla and formation of a diastema between alveoli 5 and 6.
The above-mentioned ventral and lateral shape changes describe most of the main anatomical features that differentiate 'baryonychine' and spinosaurine premaxillae. Therefore, it is evident that the differences in size between 'baryonychines' and spinosaurines (Figures 3 and 4) explain much of the differences in the observed premaxillae shape. Thus, although there is a clear influence of size on premaxillae shape, it also seems to be phylogenetic informative. It is important to notice, however, that several anatomical features of the premaxilla of both 'baryonychines' and spinosaurines are also explained by shape changes not related with size. In the components of shape space (RSC1) of the ventral view we observe different transformations to those of the size shape (CAC) that, however, also lead to a premaxilla that is more constricted and elongated in the region of alveoli 4-6 and less pointed in the anterior edge (Figure 4  (a)). Reduction and lingual displacement of alveoli 4 and 5 are explained by both CAC and RSC1 positive values, although alveolus 4 displaces mesially with positive RSC1 values and distally with positive CAC values. The medio-lateral constriction and elongation of the posterior premaxilla is explained by a distal-mesial displacement of the labial edge of the premaxilla at alveolus 6 (positive RSC1) and by alterations on the sagittal plane (positive CAC).
In the lateral view, the RSC1 component also leads to a generally similar premaxilla morphology to that explained by CAC, in which the posterior portion is dorsoventrally constricted in larger specimens. However different transformations lead to these similar shape changes. Specimens that present a premaxillary sagittal crest (Cristatusaurus lapparenti and Suchomimus tenerensis) have positive RSC1 scores, explained by both a dorsal displacement of the postero-dorsal rim of the premaxilla and by a ventral displacement of the postero-ventral margin of the premaxilla. Specimens that have a premaxilla with a strong constriction (MNHN SAM 124 and FSAC-KK-7281) retain negative RSC1 scores that is explained by changes in the opposite direction to that of positive scores. Specimens that do not have a sagittal crest and have a less constricted premaxilla (Oxalaia quilombensis, NHMUK PV R 16420, and MSNM V4047) have around average RSC1 scores. The CAC also partially explains these transformations, as the most negative values are associated with crested specimens, whereas positive values are related to specimens with strong constriction in the premaxilla. The RSC1 component also explains alterations in the anterior-most rim of the premaxilla, that has a flatter contour in specimens with positive values (C. lapparenti, Su. tenerensis), and projects anteriorly in specimens that have negative RSC1 values . Kellner and Campos (1996) comment that the two Cristatusaurus lapparenti specimens (MNHN GDF 365 and MNHN GDF 366) have similar premaxillae morphology, with the perceptible shape differences related to 'softer' bone surface ornamentation in the juvenile specimen (MNHN GDF 366), as well as differences in the alveoli size. These authors related these differences to the ontogenetic stage of each individual, having suggested that the specimen MNHN GDF 366 is a premaxilla from a young specimen because the interpremaxillary suture remains open.

Premaxillary ontogenetic allometry
The differences related to the lingual and distal border of the third and fourth premaxillary alveoli are noticeable between the two specimens in the PCA as well in the CAC analysis of the ventral view. In the adult specimen MHNH GDF 365, the arrangement of the alveoli is more uniform in the premaxillae when compared to the juvenile form MNHN GDF 366, which has a small diastema between the fourth and fifth premaxillary alveoli. In addition, in the adult specimen the third and fourth alveoli are larger and the widest portion of the premaxilla is placed posteriorly to the labial border of the third alveolus, whereas it is more anterior in the juvenile individual. All those morphological changes are explained by the components of shape space (RSC1) that exclude the effect of size allometry. Therefore, we interpret our results of the ventral view as having recovered, in the RSC1, variance in the premaxillae shape associated, at least in part, with ontogenetic allometry. So, the above-mentioned differences between the two C. lapparenti specimens are not related to the size of the specimens but, instead, possibly represent ontogenetic allometry, which agree with Kellner and Campos (1996). C. lapparenti MNHN GDF 366 is the only juvenile included in the ventral view CAC analysis (Figure 4(a)). Among spinosaurines, Oxalaia quilombensis may represent a younger individual than MNHN SAM 124 and MSNM V4047 (Sales and Schultz 2017). Both C. lapparenti MNHN GDF 366 and O. quilombensis present roughly similar CAC scores to, respectively, 'older' 'baryonychines' and spinosaurines, which also reinforces our interpretation that the RSC1 recovered variance associated with ontogenetic allometry. If this interpretation is correct, the shape changes that occur as 'baryonychines' grow older include: (1) the alveolus 3 expands distal-lingually; (2) the alveolus 4 expands mesially and labially (due to ontogeny, RSC1), but also posteriorly and medially (due to size allometry, CAC); and (3) the alveolus 5 expands labially. Combined, those changes cause a mediolateral expansion of the premaxilla at the region of alveoli 4 and 5 that contrasts with the spinosaurine premaxillae, that becomes narrower at this portion. It is worthy of notice that in the adult C. lapparenti MNHN GDF 365 only the region of the alveolus 4 is enlarged, which, therefore, differs from the 'baryonychine' adult condition.
The spinosaurids C. lapparenti MNHN GDF 366 and O. quilombensis present remarkably similar positive RSC1 scores in the ventral view of the premaxilla (Figure 4(a)). More interestingly, the RSC1 scores of C. lapparenti MNHN GDF 366 in the ventral view dataset are closer to spinosaurine individual scores than to the adult C. lapparenti MNHN GDF 365 and 'baryonychine' specimens (Figure 4(a)). Based on the phylogenetic relationships recovered (Figures 5 and 6), assuming C. lapparenti could serve as a proxy to the premaxillae shape of the ancestor of Spinosaurinae + Baryonychinae clade, and taking into account that C. lapparenti MNHN GDF 366 is a juvenile, we might speculate that shape differences recovered by RSC1 in ventral view could be associated with heterochrony. This hypothesis needs to be further investigated in future studies as it would be necessary a larger sample size for determining individual and ontogenetic variations (Kellner et al. 2011). However, using shape variation retained by RSC1 of ventral view and assuming spinosaurine species are paedomorphic, we could suggest the following premaxillae morphology as juvenile characteristics retained by this subfamily: (1) posterior displacement of the lateral portion of the rostrum at the border of alveolus 3; (2) reduction of alveolus 3 by a lingual-labial displacement of its distal-lingual border; (3) relatively more anterior alveolus 4; and (4) distal expansion of the posterior most lateral portion of the premaxilla. These shape changes lead to the lateral expansion of the mid portion of the premaxilla, labial displacement of the alveolus 3, and stretching of the posterior portion of the premaxilla.
The above-mentioned premaxillae morphology is clearly observed in C. lapparenti MNHN GDF 366, O. quilombensis, and in the juvenile FSAC-KK-72815 (not included in the ventral view analysis). Some of those features change in larger spinosaurine specimens in varying degrees, possibly because they represent more advanced ontogenetic stages. Both MNHN SAM 124 and MSNM V4047 have the lateral portion of the rostrum at the border of alveolus 3 mesially (instead of distally) positioned, which is also a consequence of size, as this change is also explained by the size related allometric component (CAC, Figure 4(a)). Also, MSNM V4047 has an enlarged alveolus 3, proportionately similar in size to those of adult 'baryonychines'. Therefore, it is possible that spinosaurines develop the 'adult morphology' only when they grow very old (and very large), which agrees with the proposition of Kellner et al. (2011) and Sales and Schultz (2017) that MSNM V4047 may represents an older individual.

Spinosaurids evolution and phylogenetic hypotheses
The spinosaurids monophyly is well supported by several authors, being consistently recovered within megalosauroid theropods (e.g. Sereno et al. 1998;Holtz 2004;Allain et al. 2012;Carrano et al. 2012). However, the evolutionary interrelationships of the spinosaurids remain under debate as there is no consensus on the evolutionary relationships within the group (Malafaia et al. 2020).
According to Sereno et al. (1998), two major clades of spinosaurids were recognised, namely, Baryonychinae and Spinosaurinae. They were distinguished, among other features, by characters mapped in the premaxilla, that refer to the diastemata of the premaxillary alveoli and the size of the first premaxillary alveolus, besides the morphology and ornamentations of the teeth (see original characters in Sereno et al. 1998). In a recent analysis, performed by Sales and Schultz (2017), we noticed that five of the 354 morphological characters that were used correspond to the premaxillae morphology (premaxilla characters 7, 139, 352, 353, 354), which are particularly informative for spinosaurids. However, the shape of the premaxillae has never been analysed in a quantitative way as in the present study. Sales and Schultz (2017) reviewed several cranial and dental characters and, based on their results, suggested the subfamily Baryonychinae represents a paraphyletic group. More recently, Malafaia et al. (2020) recovered the two subclades in a polytomy in the most consensus tree. Other approaches generated results with subtle differences in the position of some taxa along the tree (e.g. Holtz 2004;Allain et al. 2012;Carrano et al. 2012;Sales and Schultz 2017).
Our phylogenetic inference recovered both Baryonychinae and Spinosaurinae, as monophyletic groups, as proposed by Sereno et al. (1998), and similar with the results obtained by Carrano et al. (2012). However, unlike previous analyses, we recovered Cristatusaurus lapparenti as a basally branching spinosaurid and the sister taxon of the clade Baryonychinae + Spinosaurinae, so that Baryonychinae is only comprised of Baryonyx walkeri and Suchomimus tenerensis. In our result, we also obtained a more resolved topology for Spinosaurinae than in previous analyses: Angaturama limai is the sister taxon of a clade that has Oxalaia quilombensis as the first clade to diverge and MSNM V4047 and Irritator challengeri as more derived sister taxa.
This result is interesting because it demonstrates a more diverse composition for the group in which Cristatusaurus lapparenti represents a different radiation, which suggests a more complex spinosaurid evolutionary history than previously imagined. This more basal position of C. lapparenti brings attention to the evolution of two characters in Spinosauridae. According to Sales and Schultz (2017), the acquisition of teeth crown striations (character 142) and the loss of premaxillary teeth serrations (character 149) occurred at the ancestor of Spinosauridae and that these characters reversed to the basal state in C. lapparenti. However, our result suggests that the derived state of these characters appeared in the ancestor of Baryonychinae + Spinosaurinae. But as there is much debate regarding C. lapparenti status (e.g. Sereno et al. 1998;Taquet and Russel 1998;Hendrickx et al. 2016;Sales and Schultz 2017), we suggest caution with these results until new materials or better descriptions become public.
In our results, Angaturama limai and Irritator challengeri are recovered as distinct taxa, which may reinforce the validity of the former taxon. Our analyses also found a clade with MSNM V4047 and I. challengeri, which indicates a relationship between the African (MSNM V4047) and the Brazilian (I. challengeri) taxa. The same result was also obtained by Sales and Schultz (2017); however, they recovered a similar topology only when Angaturama limai is excluded from the analysis.
The better-resolved topology within Spinosaurinae demonstrate that the landmark-based characters on the premaxilla configurations refine the internal resolution of Spinosauridae, especially when compared to recent phylogenies (e.g. Carrano et al. 2012;Evers et al. 2015;Sales and Schultz 2017;Malafaia et al. 2020) that generally recovered polytomies. The inclusion of landmark-based characters proved to be useful to increases the group's phylogenetic studies, since Spinosaurinae clade has better resolution in our topology than in previous studies. This is interesting and according with neontological (e.g. Catalano et al. 2014;Botero-Trujillo et al. 2017) and palaeontological (e.g. Hendrickx et al. 2016;Ezcurra et al. 2020) recently studies that applied Phylogenetic Morphometrics to their phylogenetic inferences. But it is worth highlighting the scarcity of its application with palaeontological data. Therefore, the inclusion of new configurations of landmarks seems to be promising to improve the results of standard characters datasets (see Catalano et al. 2014; and Catalano and Torres 2016 for further details).

Premaxilla 'ancestral' shape
The discrete characters analysed here are based on Carrano et al. (2012), and the revision of Sales and Schultz (2017) and references therein. These authors provide an extensive discussion on the evolution of these characters. Thus, we focus solely on the evolution of landmark-based characters. We mapped the estimated ancestral configurations of the premaxilla to each node in the phylogeny and the trajectory of individual landmarks compared to its ancestor position ( Figure 5, Figure 6, and Figure S7). These configurations consist of most parsimonious reconstructions and therefore provide useful information to understand the morphological evolution of the premaxillae in spinosaurids.
The trajectory of individual LMs is presented in our topology for both, ventral ( Figure 5) and lateral (Figure 6) configurations. In the ventral view of the spinosaurids premaxilla, the main changes in the 'ancestral' shape, although small, include the position of the anterior-most point of the interpremaxillary suture, as well as the position of the anterior portion of the secondary palate. Also, the first premaxillary tooth becomes subtly more distal and the most lingual point of the third premaxillary tooth migrates lingually. The greatest shape changes seem to be retained in the fourth tooth of the premaxilla: the labial portion expands laterally and the mesial border protrudes distally. In the lateral view, the main changes in the premaxillary 'ancestral' shape of spinosaurids are similar to some of the changes aforementioned for the ventral view: changes in the position of the fifth and sixth premaxillary teeth and changes in the shape of the anteroventral-most portion of the premaxilla.
In the most recent common 'ancestor' between Baryonychinae and Spinosaurinae, the main shape changes in the ventral view occur in the interpremaxillary contact at the secondary palate, and in the displacement of the labial and distal borders of the fifth tooth of the premaxilla. The main change in the lateral view is related to the position of the mesial border of the sixth premaxillary alveolus.

Conclusion
The morphology of the spinosaurids premaxilla presents a relevant phylogenetic signal and morphological variations that contribute to the diagnosis of some taxa. For the first time, the anatomy of the spinosaurids premaxilla was explored visually and quantitatively by GMM and two sets of landmark-based characters were included in a cladistic context. The results allowed us to quantify and describe the observed variation more clearly and provided evidence to discuss several aspects of the evolution of spinosaurids. Namely, our final remarks are: (1) The identification of Suchomimus tenerensis and Cristatusaurus lapparenti specimens could be verified and refined, allowing the following reinterpretation: (1) the first premaxillary alveolus of Su. tenerensis is more anteriorly positioned on the rostrum than that of C. lapparenti; (2) the fourth premaxillary alveolus of Su. tenerensis is slightly smaller than that of C. lapparenti; (3) the diastema between the third and fourth premaxillary alveoli is larger in Su. tenerensis; (4) the edge of the third premaxillary alveolus is slightly more dorsal, in Su. tenerensis than in C. lapparenti and (5) Su. tenerensis has a more expanded posterodorsal region of the premaxilla. Therefore, our results support that C. lapparenti is a valid species, distinct from both, Baryonyx walkeri and Su. tenerensis.
(2) The Spinosauridae indet. FSAC-KK-7281 specimen could be certainly referred to Spinosaurinae (in agreement with Lakin and Longrich 2019). However, as FSAC-KK-7281 is a juvenile individual, a more exclusive taxonomic affinity remains unclear, since some of its morphological characteristics possibly reflect ontogenetic changes, with potential differences due to allometry. (3) The premaxillae specimens referred to Spinosaurus (NHMUK PV R 16420, MSNM V4047, and MNHN SAM 124) occupy a large area in the morphospace. These specimens differ mainly in the degree of constriction of the posterior region of the premaxilla and the degree of expansion of the anterodorsal rim of the premaxilla, but also in size and placement of the third and fourth premaxillary alveoli. Most of these changes are a result of size differences, but some are possibly due to ontogeny as demonstrated.
Therefore, we reinforce that identifying these specimens as Sp. aegyptiacus is putative and should not be done before the discovery of better-preserved specimens. (4) The diagnosis of Oxalaia quilombensis is not well established as recently suggested by Smyth et al. (2020), but two possible diagnostic features of the premaxilla are: (1) the roughness of the palate (Kellner et al. 2011;Sales and Schultz 2017), and (2) the proximal portion of the anteromedial maxillary process that contact medially. Moreover, Oxalaia quilombensis is projected onto the same morphospace as spinosaurines but has the lateral morphology closer to the average of Spinosauridae, therefore, differing from other spinosaurines. Thus, we reject the synonym with Sp. aegyptiacus, as suggested by Smyth et al. (2020), because it is putative and there is clear evidence to consider O. quilombensis as a valid name-bearing taxon. (5) The premaxillae shape variation in spinosaurids seems to be correlated with size. Several features that differentiate 'Baryonychinae' and Spinosaurinae premaxillae (e.g.

MSNM V4047
Oxalaia quilombensis  diastemata, medial constriction, alveoli size, posterior lateral expansion of premaxillae) have a causal explanation related with size differences (i.e. allometry), specifically, to the increase in size observed in Spinosaurinae. (6) Shape changes unrelated to size are possibly explained by ontogenetic allometry seen in Cristatusaurus lapparenti MNHN GDF 365 and MNHN GDF 366. The following shape changes possibly represent juvenile features: (1) posterior displacement of the lateral portion of the rostrum at the labial border of the third alveolus; (2) reduction of the third alveolus by an anterolabial displacement of its posterolingual border; (3) relatively more anterior fourth alveolus; and (4) posterior expansion of the lateral portion of the premaxilla. Those features are present in the studied juvenile C. lapparenti, also in the adult Spinosaurinae specimens, so it is suggested proposed that they may be paedomorphic features of Spinosaurinae.
(7) The Spinosauridae monophyly is well supported in the literature over the years, despite the unresolved relationships within the group. Our analysis recovered Cristatusaurus lapparenti as the basal most spinosaurid. Thus, Baryonychinae is composed by Baryonyx walkeri and Suchomimus tenerensis only. Within the Spinosaurinae, Angaturama limai was recovered as the first branch of divergence, as the sister taxon of a clade recovered as (Oxalaia quilombensis, (MSNM V4047, Irritator challengeri)). (8) The topology presented here is better resolved than that previously proposed. Therefore, the inclusion of landmarkbased characters proved to be useful to increase the phylogenetic reconstruction of spinosaurids. Along with the reappraisal of standard morphological characters, as generally suggested by Catalano et al. (2014) Figure 6. Most parsimonious tree (score = 86.615) obtained by the analysis of morphological characters (standard) and landmark-based characters (configurations) of Spinosauridae species. The transformations of the ancestral shape of the premaxillary branches refer to the lateral view. Red lines: inferred configuration from ancestor; gray lines: overlapping ancestor configuration; black lines: displacement of the LMs between ancestor and descendants (comparison between pairs).
structures (other cranial and post-cranial elements) seems to be promising to improve the results and understand the evolution of this enigmatic clade of theropod dinosaurs.