Males Resemble Females: Re-Evaluating Sexual Dimorphism in Protoceratops andrewsi (Neoceratopsia, Protoceratopsidae)

Background Protoceratops andrewsi (Neoceratopsia, Protoceratopsidae) is a well-known dinosaur from the Upper Cretaceous of Mongolia. Some previous workers hypothesized sexual dimorphism in the cranial shape of this taxon, using qualitative and quantitative observations. In particular, width and height of the frill as well as the development of a nasal horn have been hypothesized as potentially sexually dimorphic. Methodology/Principal Findings Here, we reassess potential sexual dimorphism in skulls of Protoceratops andrewsi by applying two-dimensional geometric morphometrics to 29 skulls in lateral and dorsal views. Principal Component Analyses and nonparametric MANOVAs recover no clear separation between hypothetical “males” and “females” within the overall morphospace. Males and females thus possess similar overall cranial morphologies. No differences in size between “males” and “females” are recovered using nonparametric ANOVAs. Conclusions/Significance Sexual dimorphism within Protoceratops andrewsi is not strongly supported by our results, as previously proposed by several authors. Anatomical traits such as height and width of the frill, and skull size thus may not be sexually dimorphic. Based on PCA for a data set focusing on the rostrum and associated ANOVA results, nasal horn height is the only feature with potential dimorphism. As a whole, most purported dimorphic variation is probably primarily the result of ontogenetic cranial shape changes as well as intraspecific cranial variation independent of sex.


Conclusions/Significance
Sexual dimorphism within Protoceratops andrewsiis not strongly supported by our results, as previously proposed by several authors. Anatomical traits such as height and width of the frill, and skull size thus may not be sexually dimorphic. Based on PCA for a data set focusing on the rostrum and associated ANOVA results, nasal horn height is the only feature with potential dimorphism. As a whole, most purported dimorphic variation is probably primarily the result of ontogenetic cranial shape changes as well as intraspecific cranial variation independent of sex.

Introduction
Sexual dimorphism is an expected product of sexual selection and represents an important factor in breeding success within many species [1,2,3]. Morphology, behavior and size differences distinguish males and females in many extant and extinct animals [4][5][6][7][8][9][10][11][12][13][14][15][16][17][18]. However, recognizing sexual dimorphism in extinct organisms presents special problems in that it can be difficult to distinguish sexual differences from ontogenetic or intraspecific ones. This issue has plagued attempts to recognize sexual dimorphism in non-avian dinosaurs. For instance, hypothesized males and females of hadrosaurs were later shown to be stratigraphically and temporally separated species [19,20]. Purported differences between male and female tyrannosaurs in the number of chevron bones within the tail were recently attributed to interspecific variation independent of sex [21]. Recently, the finding of the medullary bone, typical of female birds, in the femur of a Tyrannosaurus rex individual has suggested that this investigated specimen was a female [22]. However, even in cases where sexual dimorphism seems plausible, limited sample sizes often prevent a rigorous test of the hypothesis.
Protoceratops andrewsi (Neoceratopsia, Protoceratopsidae) is a well-known neoceratopsian from the Campanian-aged Djadokhta Formation outcrops of the Gobi Desert, Mongolia [31]. It is a small (<2 m total body length), quadrupedal animal, characterized by a skull with a thin, bony frill projecting over the neck and lacking the prominent horns that characterized ceratopsids such as Triceratops (Fig 1). Beginning with the early Central Asiatic Expeditions of the American Museum of Natural History during the 1920s, more than one hundred well-preserved skulls and skeletons (many from various ontogenetic stages) have been unearthed, providing information on ontogeny and intraspecific variability within P. andrewsi [23,25,35,39,44].
Since the early scientific contributions on Protoceratops, several authors have hypothesized sexual dimorphism in cranial and postcranial features [23,39,40]. Dimorphism has been inferred for a variety of anatomical traits, such as the presence of a nasal horn and relative frill width, presumably used for visual display, ritualized combats and interspecific recognition [24,25,31,35]. Dodson [25] investigated allometric cranial shape changes in Protoceratops specimens (n = 24), using linear measurements of the skull analyzed with bivariate plots. Dodson hypothesized sexual dimorphism in the frill (width and height), postorbital width of the skull, and in the nasal horn (nasal height), suggesting that the frill had importance for visual display rather than simply being a structure for muscle attachments [26]. Recent descriptions of new specimens of P. andrewsi from Mongolia accepted this hypothesis and further highlighted inferred sexual dimorphism in the skull [35]. Fig 1 shows a hypothetical male and female of P. andrewsi following the work of Brown and Schlaikjer [23] and Dodson [25]. Dodson [25] specifically considered four variables as indicative of sexual dimorphism: (1) postorbital width of the skull, (2) nasal height of the skull, (3) width of the frill, and (4) height of the frill. Each variable for each specimen in the sample was assigned a score of -1 (male trait expression) or +1 (female trait expression) based on the specimen's position relative to the line of best fit when each variable was regressed on basal skull length. A score of -4 indicates a "male" and a score of +4 indicates a "female" for all four characters. Only a few specimens reflect this assumption consistently across all characters; many specimens (e.g., AMNH 6439, AMNH 6425, AMNH 6438 and AMNH 6429) only partially correspond to "male" or "female" either. Dodson assigned the sex of these remaining specimens based on their position in the principal coordinates morphospace, suggesting a posteriori distinction for some "males" and "females" within P. andrewsi. However, no multivariate test has been performed to test for significantly different morphospace occupation between groups and thus to test the original hypothesis of sexual dimorphism within P. andrewsi.
Protoceratops andrewsi is a suitable taxon for a GM study because it is represented by numerous well preserved skulls. In 1990, Chapman [38] applied a landmark-configuration technique to investigate qualitative shape changes in two skulls of P. andrewsi (a hypothetical male and hypothetical female) to re-evaluate and confirm the sexual dimorphism hypothesis. However, that was a preliminary study intended as an example of the use of geometric morphometrics, so the sample was very small, no statistical tests were undertaken, and overall results were presented in qualitative format [38]. In this work, we analyze phenotypic differences between 29 skulls, exploring cranial shape variation in P. andrewsi using a landmark configuration on skulls and statistical tools, in order to more thoroughly assess the hypothesis of cranial sexual dimorphism. We collected photos of 29 skulls in dorsal and lateral view (including juveniles, subadults and adults) of P. andrewsi (Table A in S1 File). A Canon 400D camera was used to collect most pictures, and where necessary additional images were taken from the literature. We followed the protocols of Marcus et al. [62] and Mullin and Taylor [63] to minimize parallax and measurement error in the photographs. We obtained permission from the relevant museums or institutions to access the collections for photography.

Methods
Sex Attribution. According to previous works [25], the most critical anatomical features for identifying sexual dimorphism within P. andrewsi are the postorbital width of the skull, width of the frill, nasal height of the skull, height of the frill, and length and width of external nares. We took linear measurements from photographs of skulls (dorsal and lateral view) ( Fig  2) using tpsDig2 v2.16 [64] in order to distinguish a priori "males" and "females" in the sample, according to the criteria specified in Dodson [25] (Table B in S1 File). Although direct measurements would of course be most desirable, they were not possible in many cases due to specimen access constraints (e.g., unremovable display cases). We measured basal skull length (BLS) from photos of the sex-undetermined specimens (unshared with Dodson's work) and we regressed BSL with each of the sex-discriminant traits, as indicated by Dodson [25]. Specimens which lie below the line of best fit have a score of +1 (female trait) and specimens which lie above the line have a score of -1 (male trait; Table 2 and Fig 3). We note that in all cases, Dodson's original sex determinations are confirmed here for specimens in his sample.
Geometric Morphometrics. Geometric morphometrics (GM) is suitable to quantify morphological changes and to analyze phenotypic differences among taxa [48,60]. Thirty-four landmarks and 18 semilandmarks in two dimensions were digitized on each skull in lateral view (Fig 4A), and 14 landmarks and 9 semilandmarks were digitized on each skull in dorsal view ( Fig 4B) using tpsDig2 v2.16 software [64]. Landmarks were selected to identify major sutural contacts on the skull (Table C in S1 File), as well as to capture the shapes of structures previously hypothesized to be sexually dimorphic. Scale bars were used to scale each digitized specimen.
Taphonomic distortion represents a major concern in a GM study focusing on fossil taxa, because deformation or damage obscures the original cranial shape and thus may potentially affect the results. Although no method will completely remove the effects of distortion, several strategies were used to reduce the effects within our sample. For analyses of skulls in dorsal view, we only digitized the half of the skull that was least affected by distortion.
We removed from the sample specimens such as AMNH 6443, AMNH 6428, AMNH 6477, AMNH 6417 or AMNH 6251 (holotype), because they were too strongly distorted, damaged, or had major anatomical parts restored with plaster. By contrast, some specimens presented what we consider a reasonable and minor amount of restoration. For instance, AMNH 6438, AMNH 6432, AMNH 6637 and AMNH 6430 have partially restored frills. The parietals are lost in MPC-D-100.539, but the rest of the skull is preserved in lateral view. Thus, we included these specimens in the analysis to maximize sample size. We used the function fixLMtps() from the Linear measurements measured on skull. Redrawn following the measurement standards originally in Dodson [25]. The numbers of measurements corresponds to those reported in Table B in S1 File.

Specimen #
Postorbital width Frill width Nasal height Frill height Total score "Morpho" R package [65] to estimate missing landmarks based on the three closest complete specimens to minimize errors during the digitization process that result from restoration. Even if the use of 2D representations loses some data relative to 3D geometries [66], we feel that a 2D approach is justified here due to inaccessibility of many specimens within museum exhibits. For instance, many specimens currently are obscured by mounting hardware or glass exhibit cases. Furthermore, based on the criteria presented by Dodson [25], we posit that sexual dimorphism can be assessed with a 2D approach.
Semilandmarks are useful to capture morphological information of outlines where no homologous points can be detected. Curve or contours are assumed to be homologous among specimens [67,68]. We digitized semilandmarks at equal distance along outlines drawn on the specimens. Moreover, we digitized skull photos considering two separate subsets of landmarks. We also explored the cranial shape changes focusing on the skull without the frill and the frill alone (in lateral view) to investigate if either module possesses discriminating anatomical traits [46] (Fig 4C and 4D).
Generalized Procrustes Analysis (GPA, [69]) implemented using the procSym() function in the "Morpho" R package [65] was used to analyze shape differences among specimens in the four different samples (i.e. entire skull, in lateral and dorsal view, skull without the frill and the frill alone, in lateral view). GPA scales, aligns and rotates each landmark configuration to unit centroid size (CS = the square root of sum of squared differences between landmarks from their centroid [70]). CS represents a proxy for size of the Protoceratops specimens in our sample [66].  Table B in S1 File). Black dot represents "male", white dot represents "female", gray dot represents juvenile, and red dot represents a sex-undetermined specimen. Scale of axes is logarithmic. In the sample used here, a few specimens (AMNH 6418, AMNH 6637 and AMNH 6413) exhibit ambiguous characters for a priori distinction of "males" and "females". However, after changing the sex attribution of these specimens, the overall results of the analyses were unchanged. In addition, we provide the raw TPS files (S2 File), which will allow future researchers to reanalyze the data in light of our findings.
After GPA, Principal Components Analysis (PCA) was performed on the Procrustes coordinates to identify orthogonal axes of maximal variation in the four datasets. Additionally, we performed a PCA on linear measurements calculated on the skulls, in dorsal and lateral view, to test if the anatomical traits (e.g. nasal height of skull and width of the frill) suggested by Dodson [25] as critical for sexual dimorphism in P. andrewsi, are useful to distinguish "males" and "females".
Linear Models and UPGMA. Overall nonparametric permuted MANOVA and ANOVA (npMANOVA and npANOVA; 10,000 permutations) have been performed on the four pooled datasets to highlight difference in shape and size, respectively. Pair-wise npMANOVA (using sex as a factor), and corrected after Holm correction, was performed on all shape variables as well as on the linear measurements in order to test for sexual dimorphism among specimens of P. andrewsi in the four datasets. Additionally, a pair-wise npANOVA (using sex as a factor), and corrected after Holm correction, was performed to explore differences in size (measured as CS). We performed this tests using the function adonis() of the "vegan" R package [71] that manages with unbalanced sample size.
Due to the small sample available, no linear discriminant analysis has been performed. As stated in a recent paper [72], when the number of variables is close to the number of cases, the groups appear separated after performing a linear discriminant analysis or a canonical variate analysis even if the individuals come from the same population. Therefore, a valid assessment of group separation needs many more cases than variables [72].
A npANOVA (using sex as a factor, and a Holm correction for p-values) was performed to explore differences between each PC score accounting for up to 5% of shape variation.
Additionally, a cluster analysis was performed on the shape data of the four datasets. Procrustes distances were agglomerated by means of a UPGMA (Unweighted Pair Group Method with Arithmetic mean) algorithm. The results are four dendrograms of morphological similarities among specimens included in the sample (skull in lateral and dorsal view, skull without frill and sole frill).
Lastly, in order to explore ontogenetic shape changes in the skull, we examined the relationship between shape (as dependent) and size (as independent) variables, in the cranial datasets in lateral and dorsal view, by using a Ordinary Least Squares (OLS) linear regression model. Because we do not have any histological data for the studied specimens to assess the age at death of each, we used size as proxy for age [66]. Just for sake of visualization, we used Canonical Correlation Analysis (CCA) in order to extract the vector maximally correlated with the independent variable (size).

Results
First, npMANOVAs performed on the four pooled samples highlight difference in shape within each dataset (p-value <<0.05), whereas npANOVAs show significant differences in size (pvalue <<0.05).
We performed a PCA on linear skull measurements (dorsal and lateral view) to test the hypothesis of sexual dimorphism within P. andrewsi. The first 3 PCs explains collectively 98.7% of total shape variance. Fig 5 shows the relationship between PC1 (93.19% of shape variance) and PC2 (3.44% of shape variance) for linear measurements measured on the skull. At positive PC1 values, juveniles cluster together. On the other hand, "males" and "females" do not differ from each other at negative PC1 values or positive and negative PC2 values.
The first 15 principal components of PCA explain 95% of total shape variance of skulls in lateral view. Fig 6 shows the relationship between PC1 (34.74% of the total shape variance) and PC2 (12.78% of the total shape variance) of the total sample. Positive PC1 values are associated with a massive skull bearing a high and dorsally expanded frill, small orbit (relative to the rest of the skull), elevated snout, pronounced nasal horn and a ventrocaudally expanded jugal. This morphology is typical of an adult. At negative PC1 values, the skull is short and narrow with a short frill, no nasal horn, large orbit (relative to the rest of the skull) and a premaxilla-maxilla upper contact rising to the same height as the lower margin of the orbit. This morphology corresponds to a juvenile. At positive PC2 values, the skull has a caudodorsally elongated squamosal, a ventrally expanded jugal, a less pronounced nasal horn and a caudally tilted external naris. At negative PC2 values the skull bears a caudally expanded jugal, a dorsally expanded parietal, a pronounced nasal horn and a nearly vertical external naris. Juveniles cluster at negative PC1 values, whereas "male" and "female" morphospace overlaps at negative and positive PC1 and PC2 values, indicating morphological similarities between groups. S1 Fig shows the relationship between PC1, PC2 and PC3 (the latter corresponding to 10.57% of the total shape variance), and Table D in S1 File reports the shape variance explained by each of first 15 Principal Components.  The first 13 principal components of PCA, performed on lateral views of skulls with the frill excluded, explain collectively 95% of total shape variance. Fig 7 shows the relationship between PC1 (33.93% of shape variance) and PC2 (16.23% of shape variance) of the total sample. A Protoceratops adult cranial morphology having a massive skull with a high snout, pronounced nasal horn and a caudoventrally expanded jugal is associated with positive PC1 values. A short skull possessing a short snout, no nasal horn, large orbit and a narrow, caudally expanded jugal characterizes Protoceratops juvenile cranial morphology, at negative PC1 values. At negative PC2 values, the skull is short with a nearly vertical external naris, a dorsally place maxilla-premaxilla contact, and a caudally expanded jugal. At positive PC2 values, the skull is relatively longer, bearing an external naris tilted backwards and a maxilla-premaxilla upper contact rising to the same height as the lower margin of the orbit. Juveniles are distinct from "males" and "females" at negative PC1 values, whereas "male" and "female" morphospaces overlap each other.
The first 9 principal components of PCA, performed on GM data of the frill in lateral view, explain collectively 95% of total shape variance. Fig 8 shows the relationship between PC1 (32.05% of shape variance) and PC2 (21.04% of shape variance) of the total sample. A frill bearing a dorsally expanded parietal and a caudally elongated squamosal is associated with positive PC1 values. A shorter frill with a dorsally expanded parietal-squamosal complex is associated with negative PC1 values. At positive PC2 values the frill possesses a dorsally expanded parietal and a large squamosal, whereas at negative PC2 values the frill bears a longer squamosal and a caudodorsally expanded parietal. Juvenile, "male" and "female" morphospaces greatly overlap each other.
The first 8 principal components of the PCA, performed on the skulls in dorsal view, explain collectively 96% of total shape variance. Fig 9 shows the relationship between PC1 (41.36% of shape variance) and PC2 (30.54% of shape variance) of the total sample. At negative PC1 values the skull possesses a broad and laterally expanded frill, small orbit and a caudolaterally expanded jugal. This cranial morphology corresponds to a Protoceratops adult. At positive PC1 values the   skull has a short and caudally expanded frill, a large orbit and a moderately laterally expanded jugal. This morphology corresponds to a Protoceratops juvenile. At positive PC2 values the skull bears a short and caudolaterally expanded frill, long snout and a jugal that is deeply expanded caudolaterally. A skull with a wide and caudally expanded frill and a shorter snout is associated with negative PC2 values. Juvenile, "male" and "female" morphospaces still greatly overlap.
S2 Fig shows the relationship between PC1, PC2 and PC3 (the latter explains 10.40% of the total shape variance), and Table E in S1 File reports the shape variation explained by each of the first 8 Principal Components.
The pair-wise npMANOVAs performed on the linear measurements and Procrustes shape variables between groups (i.e. sex; Table 3) did not highlight any significant morphological differences between "males" and "females" of P. andrewsi. After a Holm correction, only the skull shape with the frill excluded appears to differentiate males and females. Both of them are separated from juveniles.
Pair-wise npANOVAs show a clear difference in size between juveniles and adults but not between "males" and "females" (which include adults and sub-adults; Table 4) in the datasets. After a Holm correction, pair-wise ANOVAs have similar results.
Pair-wise npANOVAs were performed on each PC score that accounts for up to 5% of shape variation, highlighting a non-significant difference between "male" and "female" cranial shape associated with each PC (Table F in S1 File). The only exception is represented by PC1 of the skulls exclusive of the frill (Table 5). Thus, these results seem to identify shape differences between groups. The UPGMA dendrograms of morphological similarities of entire skull shape in lateral and dorsal view (Fig 10A and 10B) highlight no obvious morphological differences between "males" and "females". Juveniles resemble the "female" cranial shape.
Cluster analysis performed on the sample of skulls excluding the frill (Fig 10C) also shows no clear separations between "males" and "females", nor does it show other groupings that might suggest alternative groupings for sexual dimorphism. Even if the ANOVA itself seems to identify two groups, some "males" have similar rostrum shape to the "females" and vice versa. Juveniles are distinguished from "female" rostrum shape. In the UPGMA dendrogram of frill shape similarities (Fig 10D), "males", "females" and juveniles do not show any clear separation on the basis of morphological differences.

Ontogenetic shape changes
Regressing shape and size documents the ontogenetic changes in the cranial sample of P. andrewsi in lateral and dorsal view. The OLS results highlights a significant relationship between shape and size (cranial sample in lateral view, R 2 : 0.226, p-values: 0.001; cranial sample in dorsal view, R 2 : 0.126, p-values: 0.016).
At small size values of skulls, in lateral view, P. andrewsi possesses a large orbit, short rostrum, slender and small squamosal, a low frill and a nasal with no horn. This morphology corresponds to a juvenile, as noted in previous studies [23].
The skulls are characterized by several modifications with the increase of size. The frill develops dorsally, the rostrum gets longer with a developed nasal bump, the external naris tilts upward, the squamosal enlarges caudodorsally and the orbit becomes smaller. Thus, the greatest cranial modifications occur in the frill and in the rostrum during growth (Fig 11).
When regressing shape of skull in dorsal view and size, at small size values the skull is juvenile-like. It possesses a short rostrum, a narrow frill and less laterally developed jugal. At large size values, the skull is characterized by a wide frill, laterally developed jugals, and a long rostrum. This morphological variation is adult-like. The major modifications in the dorsal view of skull during growth are the extreme lateral enlargement of the frill and the elongation of the rostrum (Fig 12).

Discussion
This work is the first to use geometric morphometrics (GM) on a large sample to investigate sexual dimorphism within P. andrewsi. As reported previously [23,25,31,39,44], sexual dimorphism within P. andrewsi should hypothetically be visible in skull shape. Thus, if P. andrewsi exhibits clear sexual dimorphism, particularly in the frill and in the nasal horn, GM might detect these morphological differences between "males" and "females". Padian and Horner [3] noted that previous investigations on sexual dimorphism in P. andrewsi are insufficient to identify sexual dimorphism. The results of this work do not support the hypothesis of sexual dimorphism within P. andrewsi.
First, considering only linear measurements taken on skulls, no consistent morphological difference appears between "males" and "females" (including adults and sub-adults; Fig 5), as classified under previously used criteria.
Concerning shape, the PCAs performed on the cranial sample, in lateral and dorsal view (Figs 6 and 9, and S1 and S2 Figs), highlighted wide shape variation between juveniles, "males" Table 5. Pair-wise npANOVA performed on each PCscore variable per-group.

Male
Female and "females" and within each group. The hypothetical "males" and "females" have similar, broadly overlapping cranial morphologies. Frill width and nasal height, suggested as critical to distinguish "males" and "females" [25], exhibit a wide range of variation within P. andrewsi when looking at the whole shape of the skull. Importantly, npMANOVAs (Table 3) do not highlight any differences in shape between hypothetical "males" and "females" within P. andrewsi as reported above. Similar results are found when testing for size differences or shape differences, for the PC scores, using npANOVA ( Table 4). The sole exception to broad overlap between hypothesized morphotypes concerns rostrum shape (Table 5). Here "males" are significantly morphologically separated from "females". However this shape separation should be considered quite tentative. Cluster analysis on the rostrum shape (Fig 9C) does not support an evident difference between rostrum morphology of "males" and "females". For instance, the typical "male" AMNH 6438, as defined in previous works [23,25], resembles the morphology of AMNH 6466, a typical "female".
Moreover, the morphological changes associated to PC1, that shows significant differences among groups (Table 5), could reveal peramorphic anatomical traits developed at different growth stages.
Overall, previously hypothesized criteria are not suitable to discriminate the two sexes. Thus, we suggest that cranial shape variation previously interpreted as sexual dimorphism is primarily related to ontogenetic shape changes and perhaps also related to intraspecific shape variation (as suggested by Makovicky et al. [73] and Frederickson [74]).
Alternatively, sexual dimorphism may have been present only in the rostrum shape, as highlighted by npANOVA results performed on PC1 values (Table 5), but not evident enough to detect reliably with our current sample size. Of course, additional work and fossil material eventually may reveal evident sexual dimorphism in the skull of Protoceratops. However, this will require major efforts to constrain ontogenetic stage (e.g., Makovicky et al. [73]) and a larger, more comprehensive sample. Furthermore, although Dodson's sex discriminant traits appear unclear, exploring the rostrum and frill morphology individually in lateral view, along with the entire shape of the skull in lateral and dorsal view, shows no evidence of new sex discriminant characters in Protoceratops. Additionally, sexual dimorphism in P. andrewsi does not occur in the cranial shape but may be restricted to the postcranium as recently argued in several works [75,76].
We also must consider the possibility that recognizable sexual dimorphism does occur within our sample, but simply in criteria differing from those of Dodson and others and along sex identifications that differ from those used here. However, we were not able to discern visually any obvious alternative groupings in the sample, either through the PC plots or the cluster analysis of various elements. Some techniques to infer dimorphism statistically have had success in extant and extinct mammals (e.g., [11,77]). However, the techniques are more easily applicable in mammals, which tend to stay at relatively constant size osteologically during adulthood relative to non-avian dinosaurs and which can also be more easily separated into juveniles and adults, removing confounding effects of ontogeny. Given the nature of growth in Protoceratops, we do not feel that techniques tested in mammals work here.
Evolutionary change in morphology is another possible explanation for the variation that we see within our Protoceratops sample. Unfortunately, stratigraphic data are extremely limited for most specimens, and the hypothesis cannot be tested further. Exploring the morphological variation of skulls, when regressing shape and size, highlights the broad cranial variation of the frill and rostrum in lateral and dorsal view during ontogeny. From juveniles to adults, the frill and rostrum are characterized by major changes (enlargement of frill and elongation of rostrum) which could be previously misunderstood as distinct morphologies pertaining to two different sexes. However, we cannot exclude non-sexual intraspecific variation within this taxon at different stages of growth.
Missing data or taphonomic distortion potentially are issues with analyses such as ours. Some specimens were excluded from the performed analyses because the fossils were too distorted or damaged. Therefore, some morphological variation within Protoceratops was lost because it was not investigated in this work. Nevertheless, the morphospaces still highlight wide morphological disparity within P. andrewsi. Missing anatomical traits reconstructed using the function fixLMtps() also could affect the results. The estimation of landmarks (conducted to avoid errors during the landmark digitization) is based on the three closest complete specimens in the morphospace and represents a conservative estimation. However, only a few specimens were affected by the reconstruction process, and these only for restricted anatomical regions (e.g., MPC-D 2006.36, AMNH 6441 and MPC-D-100.539). Furthermore, specimens subject to extensive reconstruction, mainly on the parietal, were excluded from the analyses.
Lambert et al. [31] argued for sexual dimorphism within Protoceratops hellenikorhinus, following a discriminant function analysis using linear measurements from three skulls of P. hellenikorhinus and 19 P. andrewsi specimens already studied by Dodson [25]. Males and females appeared to be distinct here, although one male (AMNH6467) is close to female morphology and juveniles resemble female shape. A recent paper by Mitteroecker and Bookstein [72] argued for the need to have a large sample size to perform those type of analyses to detect morphological differences between groups. Thus, the currently known sample of P. hellenikorhinus is not sufficient for confident identification of sexual dimorphism. Additionally, if no major sexual dimorphism occurs in P. andrewsi, it is unlikely that major sexual dimorphism occurred in the skulls of P. hellenikorhinus.
In the end, even if our results do not strongly support sexual dimorphism within P. andrewsi, we note that this hypothesis is quite difficult to test using modern statistical methods, particularly for dinosaur fossils. It is probable that males and females occur in the sample, although this distinction does not appear to be reflected in the cranial morphology, and in the frill shape in particular, except for nasal horn shape.
We suggest use of a similar approach for the exploration of the same issue in distinct groups of extinct animals. Collecting a large sample and applying GM allows deep exploration of the sexual-/morphological variation of the cranium or postcranium and assessment of differences in shape between males and females with adequate statistical support. Nonetheless, previous studies [e.g., 78] have suggested that quite large samples are needed to distinguish relatively subtle forms of dimorphism. Ultimate identification of sexual dimorphism in non-avian dinosaurs may indeed be quite difficult.
Supporting Information S1 Fig. 3D plot of Principal Component Analysis of skulls in lateral view. The black hull represents the "male" morphospace. The red hull represents the "female" morphospace. The green hull represents "juvenile" morphospace. Points dimensions are proportional to specimen Centroid Size. (PDF) S2 Fig. 3D plot of Principal Component Analysis of skulls in lateral view. The black hull represents the "male" morphospace. The red hull represents the "female" morphospace. The green hull represents "juvenile" morphospace. Points dimensions are proportional to specimen Centroid Size. (PDF) S1 File. Including the following: 1) Table A. List of material directly photographed for this study and references for those species for which we used published photos or drawings. Bold specimens are shared with those of Dodson [1]. Full references list is appended below. lv = lateral view; dv = dorsal view. 3) Table B. Linear measurements (cm) calculated on skulls of each Protoceratops andrewsi specimen occurred in the sample. The numbers of measurements corresponds to those illustrated in Fig 2. Bold specimens are shared with those of Dodson [1]. NA = not available. 4) Table C. Landmark definitions for the four modules (see Fig 4).