The systematic relationships and biogeographic history of ornithischian dinosaurs

Introduction

In recent years, considerable controversy arose surrounding the systematic relationships of taxa traditionally considered to be basal members of Ornithopoda (i.e., heterodontosaurids and hypsilophodontids). Once considered to be a stable, well supported portion of the ornithischian evolutionary tree (e.g., Sereno, 1999), Ornithopoda (sensu Sereno, 1998) is now rarely recovered as a monophyletic group in phylogenetic analyses of ornithischian relationships owing to the recovery of heterodontosaurids near the base of Ornithischia (e.g., Spencer, 2007; Butler, Upchurch & Norman, 2008). That situation prompted Butler, Upchurch & Norman (2008) to redefine Ornithopoda (see Table 1), removing heterodontosaurids as internal specifiers, and restricting the contents of Ornithopoda to only those taxa more closely related to Iguanodontia than to Marginocephalia. Despite this attempt to provide stability to use of the taxon name Ornithopoda, the exact contents of the clade remain poorly understood (Liu, 2004; Butler, 2005; Butler, Upchurch & Norman, 2008; Boyd et al., 2009).

The recognition of Hypsilophodontidae as a paraphyletic set of taxa (Scheetz, 1999; Butler, Upchurch & Norman, 2008; Boyd et al., 2009; Brown, Boyd & Russell, 2011) raised the question of whether all of these taxa belong within Ornithopoda (sensu Butler, Upchurch & Norman, 2008), or if some represent non-cerapodan, basal neornithischian taxa. Efforts to address this question via phylogenetic analyses have proven extremely difficult, with the position of former ‘hypsilophodontid’ taxa remaining fluid between analyses, with little consensus reached (e.g., Scheetz, 1999; Weishampel et al., 2003; Butler, 2005). Given the high level of confusion regarding their systematic position, some authors (i.e., Boyd et al., 2009) chose conservatively to refer to all non-marginocephalian, non-iguanodontian neornithischian taxa as ‘basal neornithischians’ until this question is adequately addressed. However, the majority of researchers continue to refer to these taxa as basal ornithopods or basal cerapodans, despite the fact that the use of those names implies resolved placement of taxa relative to Marginocephalia that is lacking in most recent phylogenetic analyses of ornithischian relationships (e.g., Butler, Upchurch & Norman, 2008; Butler et al., 2011; Makovicky et al., 2011).

Given the difficulties outlined above, phylogenetic analyses of basal ornithischian and/or neornithischian relationships tended to include few basal neornithischian taxa (e.g., Spencer, 2007), focusing on the most complete, well known taxa (e.g., Hypsilophodon) and ignoring less complete, but morphologically informative taxa (e.g., Zephyrosaurus). Moreover, smaller scale analyses of basal neornithischian relationships tended to include a noticeable level of geographic bias among the included taxa. For example, the dataset published by Scheetz (1999) and its subsequent modifications (e.g., Varricchio, Martin & Katsura, 2007; Boyd et al., 2009; Brown, Boyd & Russell, 2011) largely sample North American neornithischian taxa, with a few Asian taxa included. Along those same lines, analyses of South American taxa tend to heavily sample endemic taxa, while largely ignoring taxa from outside the continent (e.g., Coria, 1999; Novas, Cambiaso & Ambrosio, 2004; Calvo, Porfiri & Novas, 2007). Although these analyses may individually give the impression that the relationships of basal neornithischian taxa are well resolved, in truth the broader interrelationships of these taxa relative to each other and to the major ornithischian subclades (e.g., Marginocephalia) remain ambiguous.

The most extensive analysis of ornithischian relationships yet conducted sought to address, among other issues, the interrelationships of fifteen basal neornithischian taxa (Butler, Upchurch & Norman, 2008). That analysis met with limited success, ultimately requiring the incorporation of a combination of less-than-strict consensus methods and the removal of six basal neornithischian taxa to resolve the relationships of the remaining taxa. Butler, Upchurch & Norman (2008) conclude their discussion of these ‘hypsilophodontid’ taxa (their usage) by commenting on the need for further work on the relationships of these important but enigmatic taxa.

During the past decade there was a sharp increase in the number of new basal neornithischian taxa described from across the globe, including new taxa from Asia (e.g., Zan et al., 2005; Huh et al., 2011; Makovicky et al., 2011; Zheng et al., 2012), North America (e.g., Varricchio, Martin & Katsura, 2007; Brown, Boyd & Russell, 2011; Brown et al., 2013), South America (Novas, Cambiaso & Ambrosio, 2004; Calvo, Porfiri & Novas, 2007), and Africa (Butler, 2005). Those new taxa provide a wealth of information regarding basal ornithischian evolutionary trends and patterns, though most have yet to be included in a large-scale analysis of basal ornithischian relationships. The aim of this study is to robustly assess basal neornithischian dinosaur relationships using a newly constructed species-level dataset that is the largest yet assembled for this purpose both in the number of terminal taxa and characters. The goals of this study include assessment of the systematic relationships of Australian basal neornithischians, which were never before included in a broad analysis of basal ornithischian relationships, determination of the position of Marginocephalia within Neornithischia to clarify the contents of the clade Ornithopoda, clarification of the interrelationships of those taxa generally referred to as ‘hypsilophodontids’ and their placement relative to the major ornithischian subclades, and comparison of the results of this analysis to those of other recent phylogenetic analyses of basal ornithischian relationships (i.e., Buchholz, 2002; Spencer, 2007; Butler, Upchurch & Norman, 2008: see Fig. 1). The results of this phylogenetic analysis provide new insight into the evolutionary and biogeographic history of basal ornithischian dinosaurs and broader relationships within the clade. See Table 1 for a list of phylogenetic definitions used in this study.

Materials and Methods

Evaluation of stratigraphic congruence

The strict consensus phylogenetic hypothesis generated by this analysis was compared to the phylogenetic hypotheses of ornithischian relationships of Buchholz (2002), Spencer (2007) and Butler, Upchurch & Norman (2008) using stratigraphic consistency metrics. These metrics assume that as our understanding of the fossil record increases, phylogenetic hypotheses should become increasingly congruent with the stratigraphic record (Pol, Norell & Siddall, 2004). Under that assumption, the phylogenetic hypothesis that exhibits the closest fit to the fossil record best estimates the topology of the true tree. For this investigation the stratigraphic consistency measures minimum implied gap (MIG: Benton & Storrs, 1994; Wills, 1999), modified manhattan stratigraphic measure (MSM*: Pol & Norell, 2001) and the gap excess ratio (GER: Wills, 1999) were selected because those metrics are least affected by variations in tree size and shape (Pol, Norell & Siddall, 2004). The metric modified gap excess ratio (GER*; Wills, Barrett & Heathcote, 2008) was not calculated because the software used in this analysis (see below) does not provide those values. Additionally, accurately comparing stratigraphic congruence values calculated from different tree topologies requires that each tree includes an identical set of terminal taxa (Gauthier, Kluge & Rowe, 1988; Wills, Barrett & Heathcote, 2008). Therefore, when conducting these comparisons each tree topology was trimmed to include only those taxa that are present in both trees.

Calculations were conducted using the program Assistance with Stratigraphic Consistency Calculations v.4.0.0a (ASCC: Boyd et al., 2011a). That program provides the user with an interactive framework for designing an analysis and entering the required data (e.g., tree topology, taxon ages) and then calculates the final values. In situations where the tree topology being analyzed was incompletely resolved (i.e., polytomies were present), that systematic uncertainty was incorporated into the calculations using the ComPoly approach (Boyd et al., 2011a), which allows the full range of variation this uncertainty imparts in stratigraphic consistency values to be described. The presence of uncertainty in the age of the oldest known record for each taxon was addressed using the methods outlined by Pol & Norell (2006), which allow the full range of possible dates to be defined rather than having to select a single date for each terminal taxon. Incorporating all of these methods into this analysis ensured that the conclusions drawn from comparing the resulting stratigraphic consistency values are as accurate as possible.

Six comparisons were conducted during this study. The strict consensus tree topology generated by this analysis was compared to the tree reported by Buchholz (2002), the strict consensus tree by Spencer (2007), and the strict consensus, majority rule consensus, maximum agreement, and derivative strict reduced consensus trees by Butler, Upchurch & Norman (2008). The topology of three of these trees can be seen in Fig. 1. Values were also calculated for the unaltered strict consensus tree topology generated by this analysis (i.e., prior to being trimmed for comparison with other tree topologies). All of the resulting values are shown in Table 3.

Table 3:

Results of the stratigraphic consistency calculations.

Values of MIG are reported in millions of years.

	# of taxa	MIG	GER	MSM*	Result
Full ornithischian dataset	65	1908-1388	0.82-0.74	0.13-0.09	–
This analysis	19	480-356	0.80-0.67	0.37-0.26	More congruent
Buchholz (2002)	19	573-400	0.76-0.59	0.32-0.23
This analysis	16	289-212	0.92-0.80	0.76-0.57	Equally congruent
Spencer (2007)	16	260-208	0.94-0.84	0.80-0.63
This analysis	35	852-611	0.86-0.76	0.28-0.19	More congruent
Butler, Upchurch & Norman (2008) SCC	35	1582-844	0.77-0.49	0.20-0.10
This analysis	35	852-611	0.86-0.76	0.28-0.19	More congruent
Butler, Upchurch & Norman (2008) MR	35	1170-877	0.77-0.65	0.19-0.14
This analysis	28	723-489	0.86-0.75	0.34-0.23	More congruent
Butler, Upchurch & Norman (2008) MAS	28	816-620	0.81-0.70	0.27-0.20
This analysis	30	772-531	0.86-0.74	0.32-0.21	More congruent
Butler, Upchurch & Norman (2008) DSRC	30	967-661	0.81-0.67	0.25-0.17

DOI: 10.7717/peerj.1523/table-3

Notes:

Abbreviations

DSRC

derivative strict reduced consensus tree

GER

gap excess ratio

MAS

maximum agreement subtree

MIG

minimum implied gap

MSM*

modified manhattan stratigraphic measure

MR

majority-rule consensus tree

SCC

strict component consensus tree

Results of Phylogenetic Analysis

Analysis of the dataset as outlined above resulted in the recovery of thirty-six most parsimonious trees of length 868 (consistency index (CI) = 0.37; retention index (RI) = 0.65; rescaled consistency index (RCI) = 0.24). The strict consensus of these thirty-six trees and the resulting bootstrap support values are shown in Fig. 2. The details of the strict consensus tree topology are discussed in detail below. It should be noted that in the following discussion existing phylogenetic nomenclature was utilized whenever possible (see Table 1 for a list of phylogenetic definitions). Numbers given below in parentheses refer to the character number:character state being discussed. Descriptions of the characters cited below can be found in Table S1, and a list of unambiguous character state changes within the tree is given in Table S4. All characters discussed below are unambiguously optimized synapomorphies.

Elasmaria

This analysis recovers a slightly more inclusive Elasmaria clade than was originally proposed (see Calvo, Porfiri & Novas, 2007 and Table 1), though not as inclusive as in other recent analyses (Rozadilla et al., 2016). In addition to Talenkauen santacrucensis and Macrogryphosaurus gondwanicus, the Patagonian taxon Notohypsilophodon comodorensis is placed within this clade. All three taxa possess a highly reduced deltopectoral crest on the humerus (168:2; convergently present in basal iguanodontian Anabisetia saldiviai). T. santacrucensis and M. gondwanicus both retain the primitive condition of an epipophysis on cervical vertebra three (145:0; present in all ornithischians positioned below Agilisaurus louderbacki), which Calvo, Porfiri & Novas (2007) used to diagnose this clade, and an ovoid, or subcylindrical, ischial shaft (205:1; convergently present in Zephyrosaurus schaffi and many iguanodontians). Because the presence or absence of both of these characters cannot be assessed in Notohypsilophodon due to preservational issues, it is unclear if they unite all elasmarians, or if they are diagnostic for a more restricted clade composed of T. santacrucensis and M. gondwanicus. The presence of thin mineralized plates on the anterior portion of the thoracic ribcage (=intercostal plates: Butler & Galton, 2008; Boyd, Cleland & Novas, 2011) was also proposed to diagnose this clade (e.g., Calvo, Porfiri & Novas, 2007; Rozadilla et al., 2016). However, these structures are more widely distributed among basal neornithischian and basal ornithopod taxa (i.e., Hypsilophodon foxii, Othnielosaurus consors, Parksosaurus warreni, and Thescelosaurus neglectus: (Butler & Galton, 2008; Boyd, Cleland & Novas, 2011b) and do not diagnose this clade.

Historical Biogeography of Ornithischia

The results of the three analyses of basal ornithischian historical biogeography are presented in Fig. 3 (parsimony) and Fig. 4 (likelihood). To simplify comparisons between these sets of results, they will be referred to as follows: parsimony-based analysis (PB); likelihood-based analysis with equal branch lengths (LEB); and, likelihood-based analysis with time calibrated branch lengths set equal to implied missing fossil records (LFR). Additionally, a time-calibrated version of the strict consensus tree presented from Fig. 2 is shown in Fig. 5. By combining the temporal and geographic information presented in these figures, the biogeographic history of basal Ornithischia can be reconstructed and discussed in detail.

The LFR analysis reconstructs the ancestral area of the common ancestor of Dinosauria + Silesauridae (sensu Nesbitt et al., 2010) as Africa, while the other two analyses place the origin in South America. This difference reflects the fact that the African taxon Asilisaurus kongwe Nesbitt et al., 2010 is the oldest taxon included in this analysis, giving it more weight in the LFR analysis. All three analyses reconstruct the origin of Dinosauria and Ornithischia in South America, which is consistent with some prior proposals (e.g., Sereno, 1997). Considering two of the three basal theropods included in this study and the basal-most ornithischian taxon, Pisanosaurus mertii, are all from South America, this result is unsurprising. Ornithischia diverged from its sister taxon Saurischia by the early Late Triassic at the very latest (Fig. 5).

The ancestral area of the most recent common ancestor of the clade consisting of Heterodontosauridae + Genasauria is optimized as Africa by all three analyses, with this split likely taking place during the Late Triassic. Likewise, all three analyses reconstruct a period of rapid diversification of Heterodontosauridae to have occurred in Africa during either the Late Triassic or Early Jurassic (contra the results of Pol, Rauhut & Becerra (2011)), with the lineages leading to Echinodon becklesii, Fruitadens haagarorum and Tianyulong confuciusi later dispersing into Europe, North America, and Asia, respectively. These dispersals could have occurred anytime during the Jurassic.

The origin of Genasauria is hypothesized by all three analyses to have occurred in Africa during the Early Jurassic at the latest, and possibly during the Late Triassic, and the early diversification of Thyreophora also transpired in Africa, assuming the placement of Lesothosaurus diagnosticus at the base of this clade is accurate. There exists disagreement regarding the pattern of dispersal within Thyreophora. The LFR and LEB analyses slightly favor a scenario where basal thyreophorans dispersed from Africa into North America (giving rise to the Scutellosaurus lawleri lineage) and then migrating into Europe. The PB analysis is equivocal as to whether this is the case or if basal thyreophorans dispersed from Africa directly into Europe, with Scutellosaurus lawleri dispersing separately into North America. Either way, the diversification of these basal members of Thyreophora was completed before the late Early Jurassic. However, if L. diagnosticus is not a basal thyreophoran, then these reconstructed patterns could change dramatically.

The species Stormbergia dangershoeki constrains the origin of the Neornithischia to the Early Jurassic at the latest, and all three analyses agree that this clade arose in Africa. Sometime before the late Middle Jurassic there is an extensive radiation of neornithischian taxa, though the poor fossil record of neornithischian taxa during the Early and early Middle Jurassic make it impossible to determine precisely how rapidly this radiation occurred. However, all three biogeographic analyses agree that this radiation occurred in Asia (Figs. 3 and 4). All three analyses remain in agreement regarding the diversification of Neornithischia occurring within Asia until the most recent common ancestor of the clade consisting of Othnielosaurus consors + (Parksosauridae + Cerapoda). At this node, the LEB analysis slightly favors North America as the ancestral area (50.7% versus 42.9% for Asia), while the LFR analysis strongly favors Asia (81.1% versus 8.8% for North America). The PB analysis is equivocal. The situation is similar for the most recent common ancestor of the clade consisting of Parksosauridae + Cerapoda.

The earliest known parksosaurid taxa, Changchunsaurus parvus and Zephyrosaurus schaffi, are present in the middle Early Cretaceous (Fig. 5). However, a long ghost lineage is present for Parksosauridae, stretching from at least the Bathonian until the Aptian, a time span of at least 40 myr. The LEB and LFR analyses agree that the basal split within Parksosauridae that gave rise to the clades Orodrominae and Thescelosaurinae occurred in North America by the Aptian (the PB analysis is undecided between North America and Asia). Both the LEB and LFR analyses also agree that the diversification of orodromine taxa occurred in North America during the Cretaceous, with the lineage leading to Koreanosaurus boseongensis splitting from Oryctodromeus cubicularis either during or prior to the Cenomanian, with the former taxon eventually dispersing into Asia by the Santonian (Figs. 4 and 5). The PB analysis largely agrees with this interpretation, though it is equivocal as to whether at least some of the diversification of orodromine taxa occurred in Asia (Fig. 3).

Substantial disagreement exists between all three analyses regarding the pattern of geographic dispersals present within Thescelosaurinae. In the PB analysis (Fig. 3), thescelosaurines originated in either North America or Asia. The most recent common ancestor of the clade composed of Changchunsaurus parvus and Haya griva was located in Asia, and this clade arose by the Aptian. The ancestral area for most of the remaining thescelosaurines was North America, though a single lineage dispersed to South America by the Cenomanian, giving rise to Elasmaria. The LEB analysis largely agrees with this interpretation (Fig. 4A), though it sets the origin of Thescelosaurinae within North America, with the clade consisting of C. parvus and H. griva dispersing into Asia by the Aptian. The results of the LRF analysis strongly contrast with both of the other analyses. The LRF analysis (Fig. 4B) places the basal split within Thescelosaurinae in Asia prior to the Aptian. The sister taxon to the C. parvus + H. griva clade then migrates into South America (possibly by way of North America). Prior to the Cenomanian, two thescelosaurine lineages disperse into North America from South America. The first gives rise to Parksosaurus warreni, while the second gives rise to the Thescelosaurus clade.

The LFR and LEB analyses place the origin of Cerapoda within Asia prior to the late Middle Jurassic (Fig. 4), though the PB analysis finds North America, Europe, and Asia equally likely areas (Fig. 3). The diversification of Marginocephalia occurred most likely in Asia by the Late Jurassic (the PB analysis is uncertain if this occurs in Asia or Europe owing to the basal position of Stenopelix valdensis).

Extensive disagreement exists between each of the three analyses concerning the biogeographic history of basal iguanodontians, and each set of results will be discussed separately. The PB analysis (Fig. 3) places the origin of Iguanodontia in Europe. The ancestral location of the newly recognized Gondwanan clade of iguanodontians was either in Europe, Australia, or South America. The pattern of geographic dispersals within the Gondwanan iguanodontian clade is not sufficiently resolved in this analysis to permit further comment. The majority of the remaining iguanodontian taxa were endemic to Europe, with a single lineage dispersing to North America by the Aptian that gave rise to Tenontosaurus. The pattern of geographic dispersals involving Dryosaurus altus and Dysalotosaurus lettowvorbecki are unresolved. The LEB analysis (Fig. 4A) also places the origin of Iguanodontia in Europe; however, the basal split within the Gondwanan iguanodontian clade is placed in Australia. After this split, the ancestral area of the remaining members of the clade moves to South America, with Qantassaurus intrepidus dispersing back to Australia before the Valanginian. The remaining diversification follows that recovered by the PB analysis, though it recovers Africa as the ancestral area for the most recent common ancestor of D. altus and D. lettowvorbecki, with the former migrating into North America by the Late Jurassic.

The results of the LRF analysis contrast sharply with those of both the PB and LEB analyses (Fig. 4B). The LRF analysis places the origin of the Iguanodontia in Asia, requiring the lineage leading to Hypsilophodon foxii to disperse into Europe by the Aptian. The most recent common ancestor of the Gondwanan iguanodontian clade and Dryomorpha was also situated in Asia. The basal divergence within the Gondwanan iguanodontian clade occurred within Australia, as did all subsequent diversification within the clade, requiring two separate dispersals from Australia into South America. Above this clade, the ancestral area changes to Europe, with the lineage leading to the Tenontosaurus clade later migrating to North America and diversifying. The LFR analysis contradicts the LEB analysis in that the most recent common ancestor of the clade containing D. altus and D. lettowvorbecki migrates from Europe into North America, with the latter taxon then migrating to Africa (Fig. 4B versus Fig. 4A). One additional difference between the LFR analysis and the others is that the origin of Ankylopollexia is hypothesized to have occurred in North America and not Europe owing to the basal placement of the Jurassic taxon Camptosaurus dispar.

Results of Stratigraphic Congruence Analysis

The results of the stratigraphic congruence analysis are shown in Table 3. Comparisons were limited in some cases by the necessity of trimming each tree topology to only include congruent sets of taxa. In the most extreme case, the strict consensus topology generated by this analysis was trimmed from sixty-five terminal taxa to sixteen to facilitate comparison with the strict consensus topology from Spencer (2007), limiting the amount of data available to compare these tree topologies (see Fig. 1B versus Fig. 2). An additional complicating factor was the high number of taxa placed within polytomies in each tree topology (e.g., 14 out of 35 taxa are placed in unresolved positions in the strict consensus tree of Butler, Upchurch & Norman, 2008). As a result, the minimum and maximum recovered values for each metric tend to be are highly disparate, lowering the chances of being able to select one tree topology as more congruent with the stratigraphic record of first appearances than another.

Despite these methodological difficulties, most of these comparisons resulted in the selection of one topology as more congruent with the stratigraphic record (Table 3). In five of the six comparisons made, the strict consensus tree produced by this analysis was found to be more stratigraphically congruent than the alternative topology, and in the sixth case the two trees were found to be equally congruent (Table 3). This latter result may be at least in part due to the small number of taxa shared between these two analyses (sixteen shared taxa); however, the topology from Buchholz (2002) only shares nineteen taxa in common with the strict consensus topology of this analysis, and in that case the topology from this analysis is clearly more congruent with the stratigraphic record of first appearances (Table 3). Most importantly, the strict consensus tree topology recovered in this study is found to be more congruent with the stratigraphic record of first appearances than any of the tree topologies put forth by Butler, Upchurch & Norman (2008), which was the most comprehensive analysis of basal ornithischian relationships prior to this study.

Discussion

The strict consensus topology produced by this analysis (Fig. 2) is the most inclusive and well-resolved phylogenetic hypothesis of ornithischian relationships to date. This tree topology is equally congruent or more congruent with the stratigraphic record of first appearances than any other ornithischian phylogeny published in the last decade and a half (Fig. 5; Table 3). Comparing the results of this study to those of other recently published ornithischian phylogenetic hypotheses (i.e., Buchholz, 2002; Spencer, 2007; and Butler, Upchurch & Norman, 2008) provides important insights into those areas of the ornithischian evolutionary tree where our understanding is improving, where a consensus is beginning to be reached on contentious relationships, and where further improvement is needed.

Both this analysis and that of Butler, Upchurch & Norman (2008) recover a monophyletic Heterodontosauridae positioned outside of Genasauria at the base of Ornithischia, though more derived than Pisanosaurus mertii. This is in strong contrast to the traditional placement of Heterodontosauridae within Ornithopoda, a placement that has not been recovered since Sereno (1999). Thus, support is building for the removal of Heterodontosauridae from Genasauria. However, heterodontosaurids do convergently share some features with basal neornithischian taxa more closely related to Cerapoda than to Agilisaurus louderbacki (e.g., loss of a ventral acetabular flange on the ilium). This may account for the recovery of Heterodontosauridae at the base of Neornithischia outside of Cerapoda by the more restricted analysis conducted by Spencer (2007). Buchholz (2002) included only one heterodontosaurid, Heterodontosaurus tucki, in his analysis of ornithischian relationships, recovering it as the sister taxon to a supraspecific terminal taxon representing Marginocephalia. The only other phylogenetic analysis to recover this set of relationships also included H. tucki as the only representative of Heterodontosauridae (Xu et al., 2006). These unconventional results are likely a result of the fact that H. tucki is a relatively derived member of Heterodontosauridae (Pol, Rauhut & Becerra, 2011; this analysis), and is not an ideal exemplar species for representing Heterodontosauridae in phylogenetic analyses, at least not by itself.

This analysis recovers Eocursor parvus as a non-genasaurian ornithischian, as did Butler, Smith & Norman (2007) and Pol, Rauhut & Becerra (2011), contrasting with its placement as a basal neornithischian by Spencer (2007). Unlike Butler, Smith & Norman (2007) and Pol, Rauhut & Becerra (2011), this analysis identifies E. parvus as the basal-most heterodontosaurid. Butler (2010) provides a detailed list of features that separate E. parvus from heterodontosaurids, many of which are included as characters in this analysis (e.g., distribution of denticles on the tooth crowns and development of the coronoid process). Despite the inclusion of this evidence, Eocursor is positioned at the base of Heterodontosauridae based in part on the presence of some of the same features that are convergently shared between other heterodontosaurids and basal neornithischians more closely related to Cerapoda than to Agilisaurus louderbacki (e.g., presence of a ventral acetabular flange on the ilium). As such, it seems more plausible that Eocursor was a basal heterodontosaurid, which requires only two losses of these features within Ornithischia, as opposed to interpreting three independent losses near the base of Ornithischia.

This analysis recovers a very restricted Ornithopoda, which contains only Hypsilophodon foxii and Iguanodontia. Such a restricted Ornithopoda has never been recovered before in a published analysis, though the largely unpublished analysis of ornithischian relationships summarized in Liu (2004) recovered an even more restricted Ornithopoda that included an identical set of taxa as Iguanodontia. The reduced size of Ornithopoda in the study presented here is a result of the relatively high placement of Marginocephalia on the tree relative to other analyses (e.g., Sereno, 1999; Butler, Upchurch & Norman, 2008). As a result, most taxa previously referred to the Hypsilophodontidae are now non-cerapodan basal neornithischians, with the exception of H. foxii. Despite not being strongly supported in the bootstrap analysis (Fig. 2), the placement of Marginocephalia on the tree is by far the most parsimonious placement given the character data analyzed. Moving Marginocephalia down the tree a single node to a position below Parksosauridae adds seven steps to the total tree length. Positioning Marginocephalia further down below Jeholosaurus shangyuanensis (the location recovered by Butler, Upchurch & Norman (2008)) increases the tree length by nine steps. Thus, the recovered position of Marginocephalia is relatively well supported by the character data used in this study.

A clade composed solely of North American basal neornithischians was first recovered by Boyd et al. (2009) in their analysis of specimens previously referred to Thescelosaurus. This analysis recovers a similar clade, here termed Parksosauridae, though it now also contains Asian and South American taxa that were not included as terminal taxa by Boyd et al. (2009). The recovery of a monophyletic Parksosauridae significantly reduces the length of the inferred ghost lineages of many of its constituent members. For example, Thescelosaurus neglectus was once inferred to possess one of the longest ghost lineages in all of Dinosauria (∼105 myr: Weishampel & Heinrich, 1992). Based on its position in the strict consensus tree, the inferred ghost lineage for this taxon is reduced by more than two-thirds (Fig. 5). Thus, not only is Parksosauridae well-supported by the character evidence, it also greatly improves the stratigraphic congruence of that subsection of the tree topology. However, a sizeable ghost lineage still exists at the base of Parksosauridae, extending ∼40 myr from the Early Cretaceous back into the Middle Jurassic (Fig. 5). This implies that there is still much to learn regarding the early evolution and diversification of parksosaurids.

Recent analyses of the relationships of basal neornithischian taxa from Asia (e.g., Changchunsaurus parvus, Haya griva) have shown some support for a clade of basal neornithischian taxa endemic to Asia (Butler et al., 2011; Makovicky et al., 2011; Han et al., 2012). In its most inclusive form (Makovicky et al., 2011; Han et al., 2012) this clade consists of Haya griva as the sister taxon to a subclade composed of Jeholosaurus shangyuanensis + Changchunsaurus parvus. The exact position of this clade within Ornithischia is unresolved in the strict consensus trees of both Butler et al. (2011) and Makovicky et al. (2011), though the maximum agreement subtrees presented by both authors place this clade near the base of Ornithopoda and the strict reduced consensus tree of Han et al. (2012) produces a similar result. However, Makovicky et al. (2011) cautioned that character support for these relationships was weak and that the large number of homoplastic characters displayed by these three taxa hinted at their possibly paraphyly. A different set of relationships is recovered for these taxa in this analysis. A clade consisting of C. parvus and H. griva is situated at the base of Thescelosaurinae within Parksosauridae (Fig. 2), while J. shangyuanensis is positioned outside of Parksosauridae near the base of Neornithischia. Given the incongruence between the results presented here and those of prior studies (e.g., Butler et al., 2011; Makovicky et al., 2011; Han et al., 2012), a brief discussion of the characters supporting the placement of these taxa in the present analysis is warranted.

This study incorporates new character data for J. shangyuanensis based on personal examination of multiple articulated and nearly complete specimens in the collections at Peking University, allowing much of the postcranial skeleton to be analyzed for the first time and for a clearer understanding of the cranial anatomy to be achieved (see Table 2 for a list of specimens examined). As a result, a set of key differences between J. shangyuanensis and the Asian taxa C. parvus and H. griva were noted that are crucial to determining the position of these taxa within Neornithischia. The ventral process of the predentary of J. shangyuanensis is unilobate (72:0), while in C. parvus and H. griva it is bifurcate (72:1). Six premaxillary teeth are present in J. shangyuanensis (112:0) in contrast to C. parvus and H. griva that display five premaxillary teeth (112:1). The morphology of the dentary of J. shangyuanensis is distinctly different than that of C. parvus and H. griva. In J. shangyuanensis, the anterior tip of the dentary is positioned close to the ventral margin (74:2), the ventral and dorsal margins of the dentary converge anteriorly (75:0), and the dorsoventral height of the dentary just anterior to the coronoid process is less than 20% of the total length of the dentary (77:0). Alternatively, in C. parvus and H. griva the anterior tip of the dentary is positioned at midheight (74:1), the ventral and dorsal margins of the dentary are subparallel (75:1), and the dorsoventral height of the dentary just anterior to the rising coronoid process is greater than 20% of the total length of the dentary (77:1). Additionally, the crowns of the dentary teeth in J. shangyuanensis lack a prominent primary ridge (139:0), while a primary ridge is present on the dentary crowns of both C. parvus and H. griva (139:1). In the postcranial skeleton, the lateral surface of the greater trochanter of the femur is convex in J. shangyuanensis (213:0), while the lateral surface of the greater trochanter of the femur is flattened in C. parvus and H. griva (213:1). That character is an unambiguous synapomorphy of Parksosauridae, clearly indicating C. parvus and H. griva are parksosaurids, while J. shangyuanensis is positioned outside of this clade. Overall, the character evidence outlined above strongly argues against a close relationship between J. shangyuanensis and either C. parvus or H. griva.

Prior investigations into the systematic relationships of South American taxa previously referred to either Hypsilophodontidae (e.g., Notohypsilophodon comodorensis) or Iguanodontia (e.g., Anabisetia saldiviai) tended to be relatively restricted in scope, focusing largely on South American taxa (e.g., Coria, 1999; Novas, Cambiaso & Ambrosio, 2004; Calvo, Porfiri & Novas, 2007; Rozadilla et al., 2016). Those investigations often recovered South American taxa in an endemic clade (Coria, 1999; Calvo, Porfiri & Novas, 2007; Rozadilla et al., 2016), or closely situated to one another as part of a South American ‘grade’ of taxa (Novas, Cambiaso & Ambrosio, 2004). Several studies discussed tentative character support for some or all of these taxa forming a clade of strictly South American or Gondwanan taxa (Coria & Calvo, 2002; Novas, Cambiaso & Ambrosio, 2004; Calvo, Porfiri & Novas, 2007; Ibiricu et al., 2010). Thus, the recovery in this study of an iguanodontian clade comprised entirely of Gondwanan taxa is not unexpected. In fact, Coria (1999) previously recovered a clade composed of Gasparinisaura cincosaltensis + Anabisetia saldiviai, the same two South American taxa recovered as a part of this Gondwanan clade. Coria (1999) also suggested that G. cincosaltensis and A. saldiviai had evolved from other Gondwanan taxa and likely dispersed into South America via Antarctica, possibly from Australia, prior to the Cretaceous (Coria, 1999:57). The structure of the strict consensus tree obtained by this study (Fig. 2), the results of the biogeographic reconstructions (Figs. 3 and 4), the inferred distribution of ghost lineages for the taxa recovered within the Gondwanan clade (Fig. 5), and the recent discovery of iguanodontian taxa from Antarctica (Morrosaurus antarcticus and Trinisaura santamartaensis: Coria et al., 2013; Rozadilla et al., 2016) that appear to share a close relationship with several South American iguanodontian taxa (e.g., Rozadilla et al., 2016) all support this interpretation.

No prior analysis recovered a close relationship between any South American and Laurasian taxa, though some authors have suggested certain South American taxa more closely resembled Laurasian taxa than other Gondwanan taxa (e.g., Talenkauen santacrucensis; Novas, Cambiaso & Ambrosio, 2004). The placement of the South American taxa Macrogryphosaurus gondwanicus, Notohypsilophodon comodorensis, and Talenkauen santacrucensis within Thescelosaurinae amongst the North American taxa Thescelosaurus and Parksosaurus warreni provide insight into the evolution of the ornithischian fauna of South America. The South American ornithischian taxa treated in this study are supported as parts of two distinct radiations that dispersed into South America at different times via separate geographic paths. Based on the results presented in this study, basal iguanodontian taxa dispersed into South America from Australia (possibly via Antarctica) during the Late Jurassic or the beginning of the Early Cretaceous (Figs. 4 and 5). Alternatively, thescelosaurine taxa most likely dispersed into South America from Asia (via North America) sometime during the latter portion of the Early Cretaceous, and then diversified, giving rise to Elasmaria (Figs. 4 and 5).

The close relationship between the South American members of Elasmaria and the North American taxa Thescelosaurus and Parksosaurus warreni may also answer some questions regarding the known stratigraphic distribution of parksosaurid taxa in North America during the Cretaceous. During most of the Cretaceous, orodromine taxa were the dominant basal neornithischian taxa present in North American faunas (Sues, 1980; Scheetz, 1999; Weishampel et al., 2004; Varricchio, Martin & Katsura, 2007; Krumenacker, 2010; Brown et al., 2013; Gates et al., 2013). At the end of the Campanian, all orodromine taxa disappear from the North American fossil record. In the Maastrichtian the thescelosaurine taxa Thescelosaurus and Parksosaurus warreni appear in the North American fossil record, which may be an example of faunal replacement (Boyd et al., 2009; Brown, Boyd & Russell, 2011). The results of the LFR biogeographic analysis suggest that the lineages leading to these latter two taxa may have originated in South America, and then dispersed into North America during the Maastrichtian. This observation strengthens the paleontological support for the presence of a land bridge and associated faunal interchange between North and South America during the latest Cretaceous (Brett-Surman & Paul, 1985; Rage, 1986; Hutchinson & Chiappe, 1998; Ezcurra & Agnolin, 2011).

Conclusion

Our understanding of the basal relationships within Ornithischia and amongst the major ornithischian sub-clades is improving, aided by the discovery and description of new taxa, redescripion and re-evaluation of previously named taxa, and the implementation of improved phylogenetic methods and practices. While the phylogenetic hypothesis presented in Fig. 2 is among the best resolved, most stratigraphically consistent hypotheses yet proposed for this clade, additional improvements are needed. Clarification of the possible synonymization of Lesothosaurus diagnosticus and Stormbergia dangershoeki will help resolve questions regarding the split between and early evolution within the clades Thyreophora and Neornithischia. While support is building for the recognition of a clade of ‘hypsilophodontids’ that does not include Hypsilophodon (=Parksosauridae), the exact position of that clade relative to Marginocephalia and its taxonomic contents still fluctuates between analyses. While repositioning Parksosauridae within Ornithopoda substantially increases the length of the strict consensus tree obtained in this study, the future inclusion of additional character information, in terms of new characters, character state observations, and/or terminal taxa, will more robustly test this placement. Additionally, while the current placement of Parksosauridae does decrease the ghost lineages for several parksosaurid taxa (e.g., Thescelosaurus), the ghost lineage for the base of Parksosauridae is still the longest within Ornithischia, indicating there is much regarding the early evolution within Neornithischia yet to learn. Finally, the complete absence of thescelosaurines and the dominance of orodromines in North America during most of the Cretaceous, followed by the abrupt disappearance of orodromines and wide geographic distribution of thescelosaurines in North America during the Maastrichtian remains an interesting puzzle. Increased North American sampling will help to further clarify if this temporal segregation is a real pattern or an artifact of sampling/preservation and additional work on non-cerapodan neornithischians from South America will provide additional information regarding biogeographic dispersals of thescelosaurines during the Cretaceous.

Given the differences between the phylogenetic hypothesis presented in this study and those presented by Butler, Upchurch & Norman (2008), the natural next step is to work towards combining the character state observations presented in those datasets together (i.e., combining congruent characters, reaching a consensus on the description and scoring of character states, and removing or modifying characters and character states that are demonstrated to be influenced by tokogenetic and ontogenetic processes). The resulting dataset, possibly supplemented with additional new characters, would allow for these alternate hypotheses to be compared against a new phylogenetic hypothesis that considers all of the evidence put forth by those two studies. Such a dataset would not only be the most comprehensive test yet of basal ornithischian relationships, but would also form an excellent starting point for constructing a single dataset aimed at testing the evolutionary relationship of all ornithischian dinosaurs at the species level. Achieving such a goal will require the cooperation of numerous researchers, but would allow ornithischian relationships to be tested in such a way that all available character evidence is brought to bear on every evaluation of systematic relationships and ensure that homology statements within Ornithischia remain consistent across analyses. It would also guarantee that analyses of taxa positioned higher up the ornithischian tree will be properly rooted and character states correctly polarized. In the current technological age where data sharing and collaboration are easier than ever and analysis of large datasets is easily accommodated by a variety of search methods, the development of separate, competing phylogenetic datasets for the same taxonomic group should be abandoned and the construction of large scale datasets that include all available character observations is a goal that should be embraced by the entire research community.

Supplemental Information