Relative performance of Bayesian morphological clock and parsimony methods for phylogenetic reconstructions: Insights from the case of Myomiminae and Dryomyinae glirid rodents

Extinct organisms provide crucial information about the origin and time of origination of extant groups. The importance of morphological phylogenetics for rigorously dating the tree of life is now widely recognized and has been revitalized by methodological developments such as the application of tip-dating Bayesian approaches. Traditionally, molecular clocks have been node calibrated. However, node-calibrations are often unsatisfactory because they do not allow the fossil age to inform about phylogenetic hypothesis. The introduction of tip-calibrations allow fossil species to be included alongside their living relatives, and the absence of molecular sequence data for these taxa remedied by supplementing the sequence alignments for living taxa with phenotype character matrices for both living and fossil taxa. So, only phylogenetic analyses that take into account morphological characters can incorporate both fossil and extant species. Herein we present an unprecedented morphological dataset for a vast group of glirid rodents, to which different phylogenetic methodologies have been applied. We have compared the tree topologies resulting from


Introduction
Establishing an evolutionary timescale is fundamental yet elusive goal of the earth and life sciences.Morphological data are a crucial complement of phylogenetic diversity and the only available information to draw phylogenetic scenario and to reconcile the fossil record with molecular trees (Hunt and Slater 2016).Maximum parsimony is the most widely applied method for analysing morphological data.However, since the development of Bayesian methods using fossil taxa as tips (Pyron, 2011;Ronquist et al., 2012), fossil species can be included alongside their living relatives into the tree building process instead of being merely used to indirectly calibrate the divergences between living lineages (node-dating).This has boosted the importance of morphological phylogenetics for dating the tree of life, particularly through the development of birth-death tree models, such as the fossilized birth-death (FBD) models (Heath et al., 2014;Stadler et al., 2018).Recent numerical methods such as STRAP (Stratigraphic Tree Analysis for Palaeontology) (Bell and Lloyd, 2014) have bolstered the incorporation of stratigraphic data into phylogenetic analyses.They allow to date phylogenies of fossil taxa resulting from parsimony and to compare average stratigraphic congruence metrics of parsimony and Bayesian trees (Bell and Lloyd, 2014;O'Connor and Wills, 2016;King and Beck, 2019;King, 2021).So, both in, developments on tree models and comparative methods make now possible to compare phylogenies of extinct taxa obtained by means of evolutionary models with those resulting from maximum parsimony.
Gliridae is one of the oldest families of living rodents.They are known since Eocene times, which suggests a late Paleocene-early Eocene origin, but it was not until the late Early Miocene when they experienced their greatest diversification and widest distribution (Daams and de Bruijn, 1995).Extant dormice are represented by 9 genera and 29 species (Holden-Musser et al., 2016) but they were much more diverse in the past.Inside this family of rodents, three subfamilies are recognized: Graphiurinae, Leithiinae and Glirinae (Montgelard et al. 2003).
Several palaeontologists attempted to arrange fossil genera and species of glirids in a bunch of subfamilies (Daams and de Bruijn, 1995;Berger, 2008;Martín Suárez, 2013).However, their classifications were not constructed on the basis of phylogenetic analyses and they differ from one another as well as show discrepancies with the neontological classification based on molecular results.For instance, Daams and de Bruijn (1995) recognized five subfamilies of glirids, based on the morphology of the cheek-teeth of their extinct and extant representatives.
However, two of the subfamilies (Dryomyinae and Myomiminae) proposed by them belong in fact to the Leithiinae (for discussion about the validity of the subfamilies and their synonymies see Wilson andReeder's, 2005 andMac Kenna andBell, 1997) and the fossil taxa that should be assigned to these groups are still matter of discussion (Freudenthal and Martín Suárez, 2013).Here, we use as empirical data a morphological dataset of Myomiminae and Dryomyinae species, two controversial groups of glirid rodents and apply parsimony and model-based approaches to test the performance of these approaches in a morphological dataset of understudied fossils.Moreover, this work allow providing the first divergence date estimates for different fossil lineages belonging to the largest groups of glirids.The results allow to discuss the different propositions for the taxonomic arrangement of the group based on dental morphology as well as to those resulting from molecular phylogenies (Mongerald et al., 2003).

Material
We have chosen rodents as study material because their teeth have a high morphological diversity and outstanding preservation in the fossil record.Moreover, rodent teeth can be readily characterized by discrete characters, providing an ideal morphological dataset to be analysed.Amongst rodents, we have turned our attention to Gliridae because it is one of the oldest families of living rodents.The systematic study presented below is based on the examination of original specimens and data from the literature listed in (Appendix I).

Methods
Eogliravus moltzeri, a basal Gliridae, is selected as outgroup.A total of 57 phylogenetically informative dental characters from 65 taxa have been coded.The list of characters has been modified from Sinitsa and Nesin (2018) (Appendix 2).The modified data matrix has been entirely coded by us.

Parsimony analysis:
The data matrix has been built using Mesquite 3.04 (Maddison and Maddison, 2009) and the analysis ran in TNT (Goloboff et al., 2008) v.1.135with the 'traditional search' option (using TBR).The analysis has been performed applying implied weighting as there is evidence that it improves the results of the phylogenetical analyses when working with morphological dataset (Goloboff et al., 2008).K values larger than the default (3.0) are more accurate to perform this kind of analyses (Goloboff et al., 2017), so K = 10 (weak implied weighting) and collapsing all branches with support = 0 has been used.Branch support has been estimated through the complementary index Relative Bremer Support (Goloboff and Farris, 2001).Owing to the lack of a priori information, all characters were unordered.
Bayesian analyses: Bayesian undated and dated analyses have been carried out with Mr. Bayes v.3.2.6 (Ronquist et al., 2012) using the CIPRES Science Gateway v.3.3 (Miller et al., 2010).The morphological dataset has been analysed with the MkV model (Lewis, 2001) under the γ model, 30 million generations and four independent runs.The prior on the gamma shape parameter is an exponential distribution.Following Simões et al. (2020), the fossilized birthdeath tree model with sampled ancestors (FBD-SA) (Stadler, 2010), under relaxed-clock models has been applied.In order to put the clock rate prior the median value for tree height Median (TL) from previous non-clock analysis (TL=9.77)divided by the age of the tree has been used.The age of the root is calculated by dividing the age of Eogliravus (45.25 Ma) by the maximum divergence age of Eogliravus and the species of the ingroup (50.5 Ma, bottom of MP10).So, the prior is: 9.77/47.87=0.204.The fossil ages used in this work for tip dating correspond to the uniform prior distributions on the age range of the stratigraphic occurrence of the fossils.Fossil recovery rate has been set up in 0.06.
To time-scale, plot the phylogenies against stratigraphy and to assess the stratigraphic congruence of the parsimony and Bayesian-derived phylogenies resulting from our analyses we have used the package STRAP for R with the default number of permutations (1000) (Bell and Lloyd, 2014).Following the methodology detailed by Bell and Lloyd (2014), first and last appearances of each taxa have been compiled as the lowermost and uppermost stratigraphic occurrences.Time-scaling has been performed with the DatePhylo function and plotting the trees using the geoscalePhylo function.

Parsimony analysis
A total of two most parsimonious trees with 353 steps long, a consistency index (CI) of 0.252 and a retention index (RI) of 0.723 have been obtained and have been used to compute the strict consensus tree (357 steps) (Fig. 1).Relative Bremer Support values for each node are shown in Appendix III.The cladogram is, on the whole, pectinate (unbalance).The most basal taxon is Peridyromys murinus.One node further (node 86) from the base we find P. soondari and one node up Pseudodryomys granatensis (node 85).The next node (node 84) gives rise to Peridyromys turbatus and Peridyromys darocensis, as sister species to one another, on one side and the rest of the ingroup on the other side.At the following node (node 83), on one side Margaritamys llulli groups with Carbomys sacaresi (node 124) and together they form the sister group of Miodyromys praecox plus Vasseuromys rambliensis (node 108).Together the two twosomes constitute the sister group to the rest of the ingroup.The latter clade divides into two large clades.The smaller of these two clades originating from node 116, includes on one hand the species belonging to Tempestia (node 122) and on the other hand the species Altomiramys, Simplomys, Armantomys, Quercomys and Prearmantomys.The topology of the tree shows that Altomiramys daamsi+more derived species form a monophyletic group (Armantomyini).Our results evidence that Quercomys and Praearmantomys should be considered as synonyms of Armantomys.This is in agreement with Daams (1990), who already considered Quercomys as a synonym of Armantomys.Most of the species of Simplomys form a polytomy (except for S. aljaphi, which is basal to all the other ones).According to this analysis the origin of the Armantomyini is set up circa 24 Ma.Pseudodryomys ibericus is the most basal species of the larger clade (originated from node 81).One node crownward (node 80), a small clade (Myodyromyini, whose origin is set up circa 19 Ma) comprising species of Miodyromys and Myomimus qafzensis (species originated from node 97) is sister group to the remaining, more derived species of the ingroup.At the following node (node 79), there is on one side all species originating from node 95 (Myomimini).According to these results the Myomimini would have originated circa 14 Ma.This clade is represented by all species of Myomimus (except for My.qafzensis) and it is sister group of a large clade that includes the remaining species of the ingroup (originating from node 78).The most basal taxon of this latter is Pseudodryomys rex followed by Miodyromys prosper.The next node (node 123) gives rise to the sister species Prodryomys gregarius and Prodryomys brailloni, and the rest of the ingroup, of which Prodryomys satus is the most basal taxon.The rest of the ingroup is symmetrically distributed into two lineages.One (all taxa originating from node 102) consists of most of the species belonging to Vasseuromys (except for V. elegans and V. duplex) and the two species of Ramys.Our results evidence that this monophyletic group (tribe Vasseuromyini) would have originated circa 24 Ma.The other one splits into two sister clades: V. elegans plus more derived taxa on one side and the species belonging to Microdyromys and Dryomys on the other.Vasseuromys elegans form a monophyletic group together with V. duplex, Gliruloides and Glirulus (node 89) that we name here as Glirulini (originating from node 89).The results obtained by Sinitsa and Nesin (2018) did not support the monophyly of most taxa belonging to Vasseuromys but fully agree with our results in the phylogenetic position of the taxa that we include into the Glirulini.Glirulini is sister clade to Microdyromys spp.plus Dryomys spp.This latter form a monophyletic group originating from node 69 that we call here the tribe Dryomyini.The origin of both Glirulini and Dryomyini is set up circa 29.8 Ma.According to the topology of the tree, the subfamily Dryomyinae coined by Daams and de Bruijn (1995) would correspond, according to this analysis, with the species that originate from node 73, which forms a monophyletic group.In contrast, our results suggest that, as currently understood, the subfamily Myomiminae, which was first created by Daams (1981) to include Myomimus and Peridyromys and, later on, enlarged by Daams and de Bruijn (1995) to include the genera Miodyromys, Peridyromys, Dryomimus, Vasseuromys, Pseudodryomys, Praearmantomys, Armantomys, Nievella, Tempestia, Altomiramys, Prodryomys, Carbomys, Margaritamys and Ramys, is polyphyletic.The subfamily Peridyromyinae (Freudenthal and Martín Suárez, 2013) is likewise invalidated.In fact, our results, (Fig. 1) provide evidence of the successive developments of various clades, which could be considered as tribes: Armantomyini, Miodyromyini, Vasseuromyini, Myomimini, Glirulini and Dryomyini.The synapomorphies characterizing each of these groups are listed in Appendix IV.Interestingly, Dryomys is phylogenetically closer to Glirulus than to Myomimus (Fig. 1), which is in disagreement with the results obtained by molecular phylogenies (Montgelard et al., 2003).

Dated Bayesian analysis
The results of the topology produced by tip-dating Bayesian methods (Fig. 2), evidence that some clades merit the status of tribe on their own.
Tribe Armantomyini: The Armantomyini (Armantomys bijmai + Quercomys daamsi + Prearmantomys crusafonti + more derived species of Armantomys) + (Altomiramys daamsi + Simplomys spp.)) is the most basal clade of the tree.All species belonging to these genera are characterized by having a very simplified dental pattern.The results of this analysis allow us to re-evaluate the status of the genera Quercomys and Praearmantomys.Quercomys was coined by Lacomba and Martinez-Salanova (1988) with Quercomys bijmai as type species.
Later on, Daams (1990) synonymized Quercomys with Armantomys (De Bruijn, 1966), on the basis of many similarities.The results of this work agree with the suggestion of Daams (1990) in considering Quercomys bijmai and Quercomys daamsi as belonging to the genus Armantomys.Moreover, according to our results, Praearmantomys crusafonti De Bruijn, 1966, the only species of the genus, should also be considered as a synonym of Armantomys.
Armantomys group is sister clade of Altomiramys daamsi plus Simplomys spp., which are sister clades of each other.The three genera are here included in the new tribe Armantomyini, the species of which are characterized by having a simplified teeth pattern.These results are in agreement with those of Daams and De Bruijn (1995), who suggested that Altomiramys and Pseudodryomys (from which most of their species have been reallocated to the genus Simplomys) were phylogenetically closely-related.The tribe Armantomyini diverges early (approximately 33.16 Ma) from the remaining species of the ingroup.Its origin is set up about 26.10 Ma in Spain, then spread to France and also entered Central Europe (Germany, Swiss) (Gómez Cano et al., 2011 ;Gómez Cano et al., 2014 ;Prieto et al., 2019).These results support the conclusions of Lu et al. (2021) concerning the monophyly of Simplomys and Armantomys and the early divergence of this clade of glirids.Our analysis includes all described species of Simplomys, which are characterized by having a very simplified dental pattern with four main ridges, one or two accessory ridges for the upper molars and one for the lower molars.The topology of the tree shows Simplomys aljaphi as the most basal species of Simplomys.
Interestingly, S. aljaphi has a less simplified dental pattern than the remaining species of the genus.This character can be considered as less derived and justify the basal position of this species inside Simplomys.The next node gives rise to the sister taxa Simplomys simplicidens and Simplomys robustus on one hand, which are sister group of Simplomys meulenorum + (Simplomys hugi + Simplomys julii) on the other one.This is not surprising due to the similarity between the dental pattern of Simplomys simplicidens and Simplomys robustus, except for the size, the width of the lingual part of the anteroloph and the size variation of the centrolophid.
The same is true for Simplomys meulenorum, Simplomys hugi and Simplomys julii, which are characterized by having a reduced centrolophid or no centrolophid.In addition, Simplomys meulenorum and Simplomys hugi show a reduced posteroloph on the M2.
Tribe Peridyromyini: One node further from the base there is a small clade with Vasseuromys rambliensis and Peridyromys murinus as sister species to one another, on one side and the rest of the ingroup on the other side.Vasseuromys rambliensis should be reallocated to the genus Peridyromys, from which P. murinus is the type species.Both taxa are here included in the new tribe Peridyromyini, whose divergence from the rest of the ingroup is set up at about 31.42Ma.This tribu has originated approximately 25.37 Ma.

Tribes Glirulini, Microdyromyini and Vasseuromyini:
The results of the phylogenetic analysis presented in this work provide evidence that the genus Vasseuromys is not monophyletic.This is in agreement with the results of the phylogenetic analyses carried out by Wu et al. (2016) and Sinitsa and Nesin (2018).Wu et al. (2016) considered Vasseuromys duplex and Glirulus ekremi as being similar to Gliruloides zhoui.Our results show that Gliruloides zhoui is the sister species of Vasseuromys duplex plus Glirulus spp.Vasseuromys duplex could be reallocated to the genus Glirulus, as Sinitsa and Nesin (2018) proposed previously.In fact, "Vasseuromys" duplex plus more derived taxa are here considered as members of the new tribe Glirulini, which diverged from Gliruloides approximately 24.73 Ma.In addition, according to our analysis, Vasseuromys bergasensis and Vasseuromys priscus should also be excluded from Vasseuromys.In fact, V. bergasensis and one node up, V. priscus, are the most basal taxa of the Microdyromys clade and should be reallocated to this latter genus that we include in the new tribe Microdyromyni.Our results agree with those of Sinitsa and Nesin (2018) who pointed out that these two species do not belong to Vasseuromys.However, they found that V. priscus and V. bergasensis were more distantly related than what we find in our analysis.
However, according to our analysis, the clade Vasseuromys contains all species of Vasseuromys excluding V. duplex, V. bergasensis, V. priscus and V. rambliensis.This clade also includes as terminal taxa the species belonging to the genus Ramys.We can therefore suggest that Ramys is a synonym of Vasseuromys.Moreover, we consider the clade Vasseuromys bacchius plus more derived taxa as the new tribe Vasseuromyini.This tribe diverged from the Mycrodyromyini circa 26 Ma.
Tribe Prodryomyini: According to our analysis, Prodryomys forms a clade with Miodyromys prosper in a basalmost position.Miodyromys prosper was first described by Thaler (1966) as Peridyromys prosper.The attribution to the genus Miodyromys was made later by Daams and De Bruijn (1995).Interestingly, Aguilar et al. (1999) andSigé et al. (1997) mentioned this species as Prodryomys prosper, which is in agreement with our results.Wu (1990) discussed the difficulty of telling apart Peridyromys from Prodryomys.However, we agree with García-Paredes et al. (2009) in considering Prodryomys as a genus clearly distinct from Pseudodryomys.We here consider the group "Miodyromys" prosper plus more derived taxa as the members of the new tribe Prodryomyini, whose split from the rest of the ingroup is set up at about 25.9 Ma.According to our results this tribe originated circa 22.77 Ma.
Tribe Miodyromyini: Our results show a clade composed of all the species belonging to the genus Peridyromys, except for Peridyromys murinus (type species of the genus), Miodyromys hamadryas, Miodyromys aegercii, Miodyromys vagus and Pseudodryomys granatensis.All these species share characters such as the presence of 2 to 4 extra ridges on the lower molars or the "moderately chaotic" pattern of the upper molars.So, they may be considered to be all allocated to the genus Miodyromys, the type species of which is M. hamadryas ( Major, 1899) and could form the tribe Miodyromyini, whose divergence from more derived glirids took place approximately 24.1 Ma.The origin of this tribu is set up at 24 Ma.
Miodyromys praecox plus more derived taxa: According to the topology of the tree, Miodyromys praecox is the most basal taxon of a clade, in which Tempestia is sister group of Margaritamys plus Carbomys.All the species belonging to this clade share some characters such as the presence of centrolophs not connected to each other, sometimes only one centroloph is present, and the lack of accessory ridges.The above-mentioned characters that are present in Miodyromys praecox but not in other species of Miodyromys contributes to its phylogenetical situation at the base of this clade.
Tribe Myomimini and Dryomyini: Our work suggests that the genera Myomimus and Dryomys are both monophyletic and sister taxa, which is in full agreement with the results of molecular analyses (Montgelard et al., 2003).The latter indeed provided evidence that Myomimus is the sister group of Dryomys plus Eliomys.All these taxa would belong to the subfamily Leithiinae, to which the genus Muscardinus should be allocated.According to Montgelard et al. (2003), the subfamily Leithiinae is the sister group of the Glirinae, which are represented by Glirulus and Glis in that work.Our results are in line with a closer relationship of Myomimus to Dryomys than to Glirulus.The divergence between the tribes Myomimini and Dryomyini has been estimated circa 16Ma.However Myomyimini originated much earlier (15.2Ma) than the Dryomyini (10.53 Ma).
Interestingly, most of the tribes above mentioned originated and diversified around the boundary Oligocene-Miocene, which can be correlated to Mi-1 glacial event detected at global scale.It would be, therefore, interesting to further investigate if this diversification event observed on glirids could have a direct relationship with this global climatic event.
Moreover, the results of this analysis suggests that neither the subfamilies Myomiminae nor the Peridyromyinae are monophyletic and the classic classifications of the Gliridae (Daams, 1981;Daams and de Bruijn, 1995;Freudenthal and Martín Suárez, 2013) should be questioned.

Bayesian versus Parsimony
The main tribes (such as Armantomyini, Miodyromyini, Vasseuromyini, Glirulini, Myomimini and Dryomyini.)are recovered in both the parsimony and Bayesian analyses (Figs 1 and 2).However, the relationships of these tribes with each other differ significantly depending on the methodology employed.For instance, the topology of the tree resulting from parsimony (Fig. 1) shows the tribes Armantomyini, Vasseuromyini and Glirulini placed as less basal groups than what is found in the topology of the Bayesian tree (Fig. 2).Another striking difference between the two topologies is that in parsimony, Dryomys and Microdyromys form a monophyletic group whereas in the Bayesian tip-dating analysis, Dryomys is monophyletic and sister group of the monophyletic Myomimus.Stratigraphic congruence indices show high levels of stratigraphic fit for all metrics for Bayesian (undated and tip-dating) and parsimony methods (Supplementary Tab.S1 available on Dryad).The strict consensus tree recovered by tip-dating has significantly higher stratigraphic congruence than that resulting from undated Bayesian and parsimony analyses.Taking into account that the topology obtained by Parsimony shows some anomalous results (e.g.basal position of younger taxa such as Tempestia) and that time-scaling these topologies (Fig. 1) results in improbable long ghost lineages for a large number taxa (e.g.Dryomys spp.), the topology obtained via the tip-dating Bayesian analyses (Fig. 2) is most likely to be more accurate with this kind of morphological dataset (fossil rodent teeth).This is in agreement with the hypothesis of López-Antoñanzas and Peláez-Campomanes (2022).Moreover, the topology obtained by applying Bayesian methodology is congruent with the results of molecular analyses (Montgelard et al., 2003).
This suggests that the Bayesian method has yielded more accurate results than Parsimony in this study.

Conclusion
Our analyses evidence that tip-dating Bayesian analyses (morphological clock) deal better than traditional Parsimony analyses when working with high homoplastic characters (such as rodent teeth).The calibrated Bayesian analysis invalidates the classic classification, particularly in that concerning the subfamilies Dryomyinae and Myomiminae.This results are in agreement with those obtained by means of molecular studies.The topology of the strict consensus tree resulting from tip-dating Bayesian methods evidence that the following clades: Armantomyini, Peridyromyini, Glirulini, Microdyromyini, Vasseuromyini, Prodryomyini, Miodyromyini, Myomimyini and Dryomyini merit the status of tribe on their own.