dataesr/scientific-paragraphs-categorization huggingface dataset:
Landmark analysis under parsimony Several methods can be applied to morphometric data for estimating phylogenetic relationships (e. g. neighbour joining (NJ), linear parsimony, squared-change parsimony (SCP), maximum likelihood (ML)). Thanks to the nature of shape data, these methods can be applied in the shape space to reduce the length of the tree. While NJ works directly on distances without explicitly reconstructing ancestral character states, the three other methods estimate ancestral shapes from landmark data in order to minimise the morphological evolution of the tree according to an optimisation logic. LAUP differs from these methods because it considers individual landmarks as a unit for the optimisation procedure and introduces several new concepts to parsimony analysis. We will provide a brief account of the method to facilitate the comprehension of this paper, but interested readers are referred to the original publications  (Catalano, Goloboff & Giannini, 2010; Goloboff & Catalano, 2011 , 2016; Catalano & Goloboff, 2012).
Reference phylogeny We constructed a phylogeny of geoemydids from molecular data to use as a reference against which the analyses of morphometric data could be compared. We collected DNA sequences of 64 geoemydid species and two testudinid outgroups from seven loci: 12s, C-MOS, COI-5p, CYTB, R35, RAG1, and RAG2. COI-5p sequences were obtained from DNA barcode sequences available on the BOLD workbench  (Ratnasingham & Hebert, 2007). Sequences of the other loci were collected from GenBank. We only selected sequences that were used in published studies. Our least inclusive analytical unit is the species. In some cases, we used sequences from different subspecies in order to maximise locus sampling per species. We aligned all sequences with MUSCLE v. 3. 8. 31  (Edgar, 2004) , and performed manual corrections as necessary (Dataset S1).
Selection of partition scheme and transition models was performed using the greedy merging strategy  (Lanfear et al. , 2014)  implemented in IQ-TREE version 1. 5. 3 (multi-thread)  (Nguyen et al. , 2015). For the ML phylogenetic estimation, we ran 15 searches on IQ-TREE with an edge-proportional partition model  (Chernomor, Von Haeseler & Minh, 2016)  and the default search heuristic settings. Branch support was assessed with ultra-fast bootstrapping  (Minh, Nguyen & Von Haeseler, 2013). Finally, the resulting ML tree was adapted to match the species sampling of each morphometric dataset (see below) by dropping superfluous species and replacing the  testudinid outgroups by the emydid Malaclemys terrapin for the reasons explained in the following section.
Specimen and landmark sampling We examined specimens from collections in the following institutions in the USA: the Field Museum of Natural History (FMNH) in Chicago, Illinois, the Chelonian Research Institute (PCHP) in Oviedo, Florida, the Museum of Comparative Zoology (MCZ) in Cambridge, Massachusetts, and the Yale Peabody Museum (YPM) in New Haven, Connecticut. We restricted our sampling to specimens with visible bony sutures and without supernumerary bone plates or scutes. The full details of all the specimens examined is given in Dataset S2. In total, 103 carapaces and 98 plastra were used in this study (Table  S1 ), all of adult specimens, representing 47 out of the 71 currently recognised geoemydid species. On average, each species is represented by 2. 3 carapace specimens and 2. 1 plastron specimens (Table  S2 ).
We consider, in contrast, that the shell of emydids is more conservative and morphologically closer to the shape of ancestral testugurians. Therefore, we selected the emydid M. terrapin to be the outgroup for our analyses. We defined 152 landmarks for the carapace (Table  1 ; Dataset S3) and 48 (Table  2 ; Dataset S4) landmarks for the plastron (Fig. 1 ). Most of these landmarks (150 and 46, respectively) are laterally symmetric, i. e. they have a homologous counterpart on the opposite side of the shell. Because of this redundancy, we reduced the number of effective landmarks for the phylogenetic analyses to 77 landmarks on the carapace and 25 on the plastron (see Data filtering and Procrustes superimposition). Our landmark definition scheme corresponds mostly to the landmarks of Claude et al.
(  2003 ), but we added landmarks for neural bones 1-5 (the number of neural bones in the Note: The landmarks marked with asterisks (49-52) were not included in the analyses for this study. We show them in this list because they are present in the raw data we provide. These landmarks actually belong to the carapace, and were acquired together with the plastron landmarks in order to make it possible to realign the plastra with their respective carapaces. That alignment was not done in this study, because it's not applicable to species with mobile plastra and to collection specimens in which the carapace and the plastron were preserved separately. posterior half of the carapace is often variable within a species), the marginal series, and the anterior peripheral series. We also adapted the definitions of the plastral landmarks to preserve all information related to right/left variation in our raw data, although side variation was not addressed in the present study.
Data filtering and Procrustes superposition We excluded right/left variation and developmental abnormalities from the data. For this, we first filtered our data to exclude specific landmarks from regions of 10 carapace specimens with obvious developmental abnormalities, as inferred from extreme right/left asymmetry and by comparison to conspecifics. No plastra had abnormalities requiring this treatment. Then, we performed generalised procrustes analyses (GPA) and generated a consensus of the two digitisation replicates of each specimen. The specimen consensus of carapaces and plastra were symmetrised estimating the missing landmarks by reflection using the R package StereoMorph v. 1. 5. 1  (Olsen & Westneat, 2015). As side variation was removed from the data, we dropped the left side landmarks for further analyses.
Therefore, the effective number of landmarks for the analyses became 77 for the carapace and 25 for the plastron. From the symmetrised specimens, we computed a generalised Procrustes consensus configuration for each species. The code used to process the landmark data is given in Dataset S6. In box turtles (Cuora), leaf turtles (Cyclemys), and Notochelys, the plastron is not always a single structure as it is usually divided along the transverse hinge that develops in these groups forming an anterior and a posterior plastral lobe. If there is a 'natural resting angle' between the two lobes, it would seldom be observed in osteological specimens. Instead of imposing an arbitrary angle between the lobes for the generalised Procrustes superimposition with the plastra of unhinged turtles, we artificially divided the unhinged plastra into anterior and posterior lobes.
Datasets In order to assess the performance of individual superimpositions created in the previous step, and some of their combinations, we defined five datasets for the analyses (Table  3 ). For the sake of clarity, we use italics to refer to datasets. When we do not use italics, we refer to the anatomical structures after which the datasets are named. The species of the shell dataset correspond to the set intersection of the species of the carapace and plastron datasets. There is no missing data in any of the datasets. The superposition steps were performed on the anterior lobe, posterior lobe, and carapace datasets independently. The shell and plastron datasets combine those independently superposed datasets, and therefore are not representations of the intact shell and plastron shapes.
The number of grid divisions and grid-nesting iterations have been shown to have great impact on score estimation  (Goloboff & Catalano, 2011) , with a significant trade-off in computation time. In order to estimate reasonable grid subdivision and nesting parameters, we performed a series of trial analyses on the shell dataset with incremental values. We tried values of 6, 8, and 10 for the number of grid divisions, and values of 1 and 2 iterations of grid nesting with a cell window size of 1. These ranges of values were determined by our computational resources and are not atypical for the LAUP literature  (Goloboff & Catalano, 2011; Catalano, Ercoli & Prevosti, 2015; Perrard, Lopez-Osorio & Carpenter, 2016). Combined, these are six parameter settings, for each of which we ran three LAUP searches of two random addition sequences followed by tree-bisection-and-reconnection (TBR) branch swapping.
LAUP search strategies We performed LAUP searches on all the datasets using GPA alignments and the dynamic alignment based on parsimony proposed by  Catalano & Goloboff (2012). The heuristic search strategy for the GPA alignments consisted of 20 replicates of random addition sequences followed by TBR branch swapping using the default settings of TNT's MULT command. For the dynamic alignments, we used a search strategy conceived by Santiago Catalano (E. Ascarrunz, 2017, personal communication). It first performs a random pairwise alignment of the configurations followed by a tree search with a random addition sequence and a round of TBR. The tree obtained from that initial search is then fed into a loop of dynamic alignments followed by a tree search with TBR, until the same topology is obtained in two successive iterations.
Squared-change parsimony To our knowledge, there is no published software capable of performing heuristic tree searches using SCP as an optimality criterion (but see  Klingenberg & Gidaszewski, 2010  for a branch-and-bound implementation). Therefore, we wrote an R script for this purpose. The script computes the unweighted sum of squared changes on unrooted trees using the method described by  Rohlf (2001;  see also  McArdle & Rodrigo, 1994; Claude, 2008). This implementation yields practically identical results to the anc. recon function for estimating ML ancestral states  (Goolsby, 2017)  included in the Rphylopars package  (Goolsby, Bruggeman & Ané, 2017) , when all branch lengths are set to unity (a punctuational model corresponding to unweighted SCP). This implementation is also slightly faster for datasets of the size used in this study.
The code is written in pure R (Dataset S7), except for a single instruction that calls the C++ Eigen library  (Bates & Eddelbuettel, 2013; Jacob, Guennebaud & Eigen Project Contributors, 2017)  to perform matrix multiplication. The search strategy consisted of performing 100 random addition sequences, each followed by hill-climbing optimisation using first random SPRs, and secondly with nearest-neighbour interchanges (NNI). The 25 trees with the best sum of squared changes scores were then used as starting trees for further optimisation with an implementation of the stochastic NNI perturbation search strategy of  Nguyen et al. (2015) , which was run 10 times. All the branch swapping operations made use of functions from the Phangorn package v. 2. 1. 1  (Schliep, 2011). Bootstrapping was performed on 200 resampled matrices, each analysed with a single random addition sequence and the hill-climbing optimisation described above, but with less thorough settings.
In addition to comparing the trees obtained from landmark data against the reference tree to assess error in topology estimation, we compared the landmark trees against each other to explore possible associations between topologies estimated from different datasets and different methods. For this, we computed pairwise distance matrices of all the trees using the quartet distance  (Estabrook, McMorris & Meacham, 1985)  and the Robinson-Foulds distance, which are closely related to the dQ and CD dissimilarity measures, respectively, but satisfy the conditions of true metrics (dQ and CD fail to satisfy the identity of indiscernibles and the triangle inequality for trees with polytomies). We used metric multidimensional scaling (MDS) (see  Hillis, Heath & John, 2005  for a similar application) on the distance matrices to produce visualisations of the tree space projected on two dimensions and applied the partitioning around medoids algorithm (PAM)  (Kaufman & Rousseeuw, 1990 ) in order to identify clusters of topologically similar trees.
We determined the optimal number of clusters to be found by PAM, the k parameter, by means of the elbow method  (Thorndike, 1953). In addition, we also applied PAM with k set to 5 and 6, to see if trees would cluster following the five datasets and/or the six methods of phylogenetic analysis. Comparisons between trees estimated from landmark data required pruning out all the species not common to all datasets. After the pruning step, all the trees used for comparisons had 42 species.
Utility software Computations were performed on a desktop computer using GNU Parallel build 20161222  (Tange, 2011)  to make efficient use of all the CPU cores. A few utility functions in our scripts make use of the the phytools package v. 0. 6  (Revell, 2012). Tree plots were produced with FigTree v. 1. 4. 2  (Rambaut, 2014)  and the APE package v. 5. 1  (Paradis, Claude & Strimmer, 2004).
Parameters for LAUP spatial optimisation We identify a single topology as the most parsimonious in replications of four parameter combinations, including the most exhaustive one (10 grid divisions, two nesting iterations). We consider that topology to be the best estimate of the MPT for the GPA alignment of the shell dataset. Scaled quartet distances between the best tree and the trees obtained in other trials were never greater than 0. 072. As expected, parsimony scores improve with more thorough approximation settings (Table  4 ). Using two levels of grid nesting produced greater improvements in the tree score than increasing the number of grid subdivisions. The two-level nesting setting was also more effective in increasing the chance of finding the optimal topology. More exhaustive spatial optimisation settings slow down tree searches dramatically, raising search times from between 53 min and 1 h 30 min to over 26 h on mid-2014 MacBook Pro with a 2.
Phylogenetic analyses of morphometric data For the optimality criterion methods (GPA and dynamic LAUP, linear parsimony, SCP, and ML) the tree scores obtained from the morphometric data analyses were always better than the scores of the reference tree (Table  5 , trees available in Dataset S8). That suggests that any limited ability to recover relationships from morphometric data was not the result of suboptimal tree searches. This conclusion is also supported by the fact that we used the reference tree as the starting tree in LAUP, ensuring that the topological vicinity of the reference tree was explored. Searches initiated with the reference tree found better parsimony scores with dynamic LAUP three out of five times, whereas all the GPA LAUP searches initiated from random addition sequences outperformed the searches initiated from the reference tree.
Note: The reference tree scores for dynamic LAUP were computed by realigning the configuration in relation to the reference tree. Note that log-likelihood scores are positive. Quantitatively, the overall performance of all methods with all datasets was poor (Table  5 ). The optimal trees obtained with each method had low accuracy, with TA dQ values ranging in 0. 179-0. 445. Accuracy estimated with the CD (TA CD ) was even lower, in the range of 0. 108-0. 319, in some cases about half the value of TA dQ. This large difference is likely due to the greater sensitivity of TA CD to displacements of a few outlier species, an issue that also affects other split-based measures like the Robinson-Foulds distance  (Böcker, Canzar & Klau, 2013). Nonetheless, TA dQ and TA CD broadly coincide on the accuracy of the trees relative to each other; the Pearson correlation coefficient between TA dQ and TA CD is 0.
841. No analytical method or type of alignment showed a significant advantage or disadvantage over the rest. However, there were more consistent differences in performance between the different datasets (Table  5 , see also Figs. 3 4 5 ). Analyses of the anterior lobe landmarks tended to perform the worst (median TA dQ = 0. 198, median TA CD = 0. 135), closely followed by the posterior lobe and plastron datasets with similar performance (median TA dQ = 0. 229 and median TA CD = 0. 149 for posterior lobe, and median TA dQ = 0. 217 and median TA CD = 0. 149 for plastron). The analyses of the shell and carapace datasets tended to be the most accurate (carapace median TA dQ = 0. 342, median TA CD = 0. 254; shell median TA dQ = 0. 296, median TA CD = 0. 2). The 70%-rule consensus of bootstrap trees improved accuracy (Table  6 ).
Predictably, this enhanced accuracy came at the expense of major losses of resolution; the most resolved consensus tree was obtained from the carapace dataset with neighbour-joining, with 52% of resolved nodes. Even that best-resolved tree is disappointing, as most of the resolved nodes only join pairs of species, failing to find major clades. The clades recovered more often in the 70%-rule bootstrap consensus included Batagur, Pangshura, and Orlitia (Fig. 3 ). We found no notable differences in total resolution or artefactual resolution (incorrectly resolved nodes) across the different methods. Qualitatively, some trees present interesting features (Fig. 4 ). The reference tree has a large and well-supported clade of Asian turtles including Batagur, Hardella, Pangshura, Morenia, Orlitia, and Malayemys (Fig.
2 ). The trees from morphometric data often recover close relationships between the species in that Batagur-Malayemys clade, although some of those 'genera' are only occasionally recovered as monophyletic. The trees obtained from the carapace dataset tended to show all the Mauremys species clustering together either as a clade or as a paraphyletic assemblage, except for the linear parsimony tree. This is surprising, as Mauremys seems difficult to recover as a clade from its carapace alone (never recovered as monophyletic or paraphyletic in Garbin, Ascarrunz & Joyce, 2018), which does not seem very distinctive among geoemydids from our visual inspection of collection specimens. The anterior lobe, posterior lobe, and plastron datasets yielded trees that consistently found Cuora as monophyletic, but with the two Cyclemys species often arranged as paraphyletic to Cuora (Figs.
In our results, although the clustering of Cuora and Cyclemys clearly reflects an evolutionary convergence, it is encouraging that Cuora is always recovered as a clade, possibly because of features such as the rounded margins of its plastral lobes and the particular geometry of the entoplastral region, among others. The strength of that plastral signal is such, that all shell trees also recover Cuora as monophyletic. It is noteworthy that many of the trees inferred from datasets that include plastral lobes also grouped together as a clade or paraphyletic assemblage all or almost all Rhinoclemmys species (Fig. 4 ), but this result is not as consistent across methods. Quantitative comparison between trees is concordant with some of the previous observations (Fig. 5 ). The elbow method indicates that the tree spaces produced by quartet and Robinson-Foulds distances can be optimally partitioned into three clusters with the PAM algorithm (k = 3).
Note: The resolution is given as the number of nodes observed in the tree divided by the number of nodes in a fully resolved tree with the same number of species. Trees were treated as unrooted. distance tree space). The trees inferred from the plastron dataset are found to cluster either with the posterior lobe trees (quartet distance tree space) or included in the cluster with the shell and posterior lobe trees (Robinson-Foulds distance tree space). These results are indicative of the strong influence of the plastral configurations, the anterior lobe in particular, in the overall signal detected by all methods. This influence must be due to specific geometric features, as the carapace dataset contains more than twice as many landmarks as the plastron dataset. The GPA superposition of the anterior lobe and posterior lobe have a between-species variance of 0.
012 and 0. 015, respectively, and the whole plastron has a variance of 0. 013. In contrast, the carapace superposition has a variance of 0. 002. Trees from the shell dataset still retain some of the better properties seen in the carapace trees, such as the better grouping of Mauremys and finding Notochelys close to Cyclemys. We also applied PAM with k values of 5 and 6, to see if the trees would cluster reflecting the five datasets or the six methods of phylogenetic analysis used. We found five clusters that for the most part reflect the five landmark datasets. This corroborates that the datasets are the most important factor shaping the tree space even within the three main clusters indicated by the elbow method. PAM with k = 6 does not yield clusters that reflect the six phylogenetic methods.
Fossil placement performance The fossil placement analyses consist in removing a single species from the reference tree and inserting it back into the tree based on the landmark data. This is done for all species in the tree, one by one. The results show some differences in performance across datasets, with the plastron dataset yielding the greatest scaled nodal errors, and carapace the least (Fig. 6 ). Additionally, ML performed consistently better than GPA LAUP and SCP, as measured by the median and the integrated cumulative distribution of the median scaled error. The differences are small and appear most pronounced for the carapace dataset, where the median scaled placement error with ML (0. 077) was less than half the median error with GPA LAUP and SCP (both 0. 167). For the carapace dataset, scaled placement errors in the range 0.
067-0. 09 correspond to placements at only one node away from the original position in the reference tree. With ML, the 70th percentile of the scaled placement errors reflects up to three nodes of separation between the reinsertion point and the original ('correct') position on the reference tree. The fossil placement analyses were performed with the fully resolved reference tree that includes many branches that are poorly supported by molecular data (Fig. 2 ). Therefore, at least a small amount of the placement error can be reasonably ascribed to the incorrect branching patterns in the reference tree. Considering this, the performance of fossil placement analyses is much better than the performance of full phylogenetic analyses, at least for the carapace and shell datasets.
Performance of methods of phylogenetic analysis No method consistently yielded more accurate and precise trees. In this respect, our results are in agreement with the study of  Catalano & Torres (2017) , who compiled 41 landmark datasets and analysed them with five phylogenetic methods including neighbourjoining of Procrustes distances, linear parsimony, and LAUP (GPA). They found little difference in the performance of all methods. Their results and ours suggest that the current methods of phylogenetic inference and landmark datasets available do not directly address the major obstacles in making reliable phylogenetic inferences from geometric information. In our case the results suggest that the poor performance is particularly driven by low phylogenetic information and/or high noise (convergence) in the data, and the differences between trees from different methods are attributable to comparatively minor effects of the search strategies and optimality criteria.