Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 RESEARCH ARTICLE A genomic exploration of the early evolution of extant cats and their sabre-toothed relatives [version 2; peer review: 2 approved] Michael V Westbury 1, Ross Barnett 1, Marcela Sandoval-Velasco1, Graham Gower 1, Filipe Garrett Vieira 1, Marc de Manuel2, Anders J Hansen 1, Nobuyuki Yamaguchi 3, Lars Werdelin4, Tomas Marques-Bonet2,5-7, M Thomas P Gilbert 1,8*, Eline D Lorenzen1* 1The GLOBE Institute, University of Copenhagen, Øster Voldgade 5-7, Copenhagen, Denmark 2Institute of Evolutionary Biology (UPF-CSIC), PRBB, Dr. Aiguader 88, 08003 Barcelona, Spain 3Institute of Tropical Biodiversity and Sustainable Development, Universiti Malaysia Terengganu, 21030 Kuala Nerus, Terengganu, Malaysia 4Department of Palaeobiology, Swedish Museum of Natural History, Box 50007, SE-104 05 Stockholm, Sweden 5Catalan Institution of Research and Advanced Studies (ICREA), Passeig de Lluís Companys, 23, 08010, Barcelona, Spain 6CNAG-CRG, Centre for Genomic Regulation (CRG), Barcelona Institute of Science and Technology (BIST), Baldiri i Reixac 4, 08028 Barcelona, Spain 7Institut Català de Paleontologia Miquel Crusafont, Universitat Autònoma de Barcelona, Edifici ICTA-ICP, c/ Columnes s/n, 08193 Cerdanyola del Vallès, Barcelona, Spain 8Department of Natural History, NTNU University Museum, Norwegian University of Science and Technology (NTNU), NO-7491 Trondheim, Norway * Equal contributors v2 First published: 24 Mar 2021, 1:25 Open Peer Review https://doi.org/10.12688/openreseurope.13104.1 Latest published: 21 Jun 2021, 1:25 https://doi.org/10.12688/openreseurope.13104.2 Reviewer Status Invited Reviewers Abstract Background: The evolutionary relationships of Felidae during their 1 2 Early–Middle Miocene radiation is contentious. Although the early common ancestors have been subsumed under the grade-group version 2 Pseudaelurus, this group is thought to be paraphyletic, including the (revision) early ancestors of both modern cats and extinct sabretooths. 21 Jun 2021 Methods: Here, we sequenced a draft nuclear genome of Smilodon populator, dated to 13,182 ± 90 cal BP, making this the oldest version 1 palaeogenome from South America to date, a region known to be 24 Mar 2021 report report problematic for ancient DNA preservation. We analysed this genome, together with genomes from other extinct and extant cats to 1. Daniel Förster , Leibniz Institute for Zoo investigate their phylogenetic relationships. Results: We confirm a deep divergence (~20.65 Ma) within sabre- and Wildlife Research (IZW), Berlin, Germany toothed cats. Through the analysis of both simulated and empirical data, we show a lack of gene flow between Smilodon and 2. Manuel Salesa , Museo Nacional de contemporary Felidae. Page 1 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Conclusions: Given that some species traditionally assigned to Pseudaelurus originated in the Early Miocene ~20 Ma, this indicates Ciencias Naturales-CSIC, Madrid, Spain that some species of Pseudaelurus may be younger than the lineages Any reports and responses or comments on the they purportedly gave rise to, further supporting the hypothesis that Pseudaelurus was paraphyletic. article can be found at the end of the article. Keywords Smilodon, ancient DNA, genomics, gene flow, Felidae, palaeogenome, phylogeny This article is included in the Excellent Science gateway. Corresponding author: Michael V Westbury (m.westbury@sund.ku.dk) Author roles: Westbury MV: Conceptualization, Formal Analysis, Visualization, Writing – Original Draft Preparation, Writing – Review & Editing; Barnett R: Conceptualization, Investigation; Sandoval-Velasco M: Investigation; Gower G: Formal Analysis, Software; Vieira FG: Data Curation, Formal Analysis; de Manuel M: Formal Analysis; Hansen AJ: Funding Acquisition; Yamaguchi N: Writing – Review & Editing; Werdelin L: Writing – Review & Editing; Marques-Bonet T: Supervision; Gilbert MTP: Conceptualization, Funding Acquisition, Project Administration, Supervision, Writing – Review & Editing; Lorenzen ED: Funding Acquisition, Supervision, Writing – Review & Editing Competing interests: No competing interests were disclosed. Grant information: This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No [681396]), (project Extinction Genomics). This project has also received funding from the European Union’s Seventh Framework Programme for research, technological development and demonstration under grant agreement no. FP7-PEOPLE-2011-IEF-298820. The work was also supported by the Independent Research Fund Denmark | Natural Sciences, Forskningsprojekt 1, grant no. 8021-00218B to EDL. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. Copyright: © 2021 Westbury MV et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. How to cite this article: Westbury MV, Barnett R, Sandoval-Velasco M et al. A genomic exploration of the early evolution of extant cats and their sabre-toothed relatives [version 2; peer review: 2 approved] Open Research Europe 2021, 1:25 https://doi.org/10.12688/openreseurope.13104.2 First published: 24 Mar 2021, 1:25 https://doi.org/10.12688/openreseurope.13104.1 Page 2 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Apart from a few anomalous mass occurrences in tar pits (e.g. REVISED Amendments from Version 1       Rancho La Brea, USA and Talara, Peru), Smilodon is rela- We have made some minor updates to the manuscript based on tively rare in the fossil record, which is not uncommon for an the reviewers’ suggestions. This includes fixing some typos and apex carnivore. Moreover, when it is present, there are only adding two new tables to the main text from the extended data. limited remains5. Any further responses from the reviewers can be found at the end of the article The evolutionary relationships of Smilodon to other extinct and extant large cats are not resolved. An early ancient DNA study using both mitochondrial and nuclear DNA placed Smilodon within Felinae6. However, this was later shown to likely reflect Plain language summary contaminant DNA from a domestic cat7. A later study using Here we sequenced the genome of the extinct sabre-toothed complete mitochondrial genomes found Smilodon to be highly cat Smilodon populator. By comparing this genome to those of divergent from both Felinae (with an estimated divergence living cat species as well as another extinct sabre-toothed cat of ~20 Ma), and another extinct Machairodontinae lineage Homotherium latidens, we were able to assess their evolution- (Homotherium (tribe Homotheriini)), diverging ~18 Ma8. ary relationships to one another. We show that not only are However, conclusions based exclusively on a single maternally sabre-toothed and extant cats very divergent from one another, inherited locus can be biased by interspecific hybridisation and but there were also highly divergent species within sabre-toothed incomplete lineage sorting, which has been well-documented cats. This high level of divergence within sabre-toothed cats in living cats9. makes it difficult to place the common ancestor of all cat spe- cies based on just the fossil record. Moreover, we were able To date, only three studies present authentic DNA from to show that it was very unlikely that sabre-toothed cats Smilodon8,10,11, all of which is mitochondrial DNA. Although hybridised with any living cat species. there is a large collection from Rancho La Brea, retrieving endogenous DNA from the material is not considered feasible7, Introduction greatly restricting the number of specimens available for Based on the fossil record, the origin of Felidae, colloqui- DNA study. The relative absence of studies of Smilodon likely ally referred to as cats, is well resolved to the Middle-Late in part reflects their rarity in the fossil record outside of tar Oligocene (~30–27 million years ago (Ma))1. However, the pits. In addition, as Smilodon are assumed to have been mixed subsequent radiation of Felidae in the Early–Middle Miocene wood-edge ambush predators12, the majority of fossils are (~23–15 Ma) is less well understood, encompassing the evo- found in open-air sites with poor DNA preservation; cave sites lution and radiation from the common ancestor of a number (which afford better preservation) are rare. Finally, the detri- of species subsumed under the paraphyletic grade-group mental effects of temperature at equatorial and near-equatorial Pseudaelurus, which was widely distributed across Europe, latitudes further exacerbate DNA preservation. Asia, and North America1. The putative paraphyly of Pseu- daelurus is exemplified by the fact that it likely includes not To elucidate the evolutionary relationship of Smilodon to only the early ancestors of modern cats (Felinae) but also those extant and extinct Felidae, and to gain insights into the early of the extinct sabretooths (Machairodontinae). radiation of Felidae, we sequenced a draft nuclear genome of a single Smilodon populator individual from the Ultima In contrast to the early radiation of Felidae, the phylogenetic Esperanza region of Chile. The specimen was radiocarbon relationships among extant cats are relatively well under- dated to 13,182 ± 90 calibrated years before present (cal BP), stood; extant cats are believed to share a most recent common and represents the oldest palaeogenome from South America, ancestor in the Late Miocene (~11 Ma)2. However, despite a region in which ancient DNA preservation is expected to be this, there is uncertainty surrounding the constituents of the limited13. long-stem lineage of Felinae after its split from Machairodon- tinae (which extends through the Early and Middle Miocene) Results and discussion until its radiation in the Late Miocene. The evolutionary Using an African lion assembly as the reference genome, we relationships within Machairodontinae are even less well under- successfully mapped a draft nuclear genome of a single Smilodon stood, and increased knowledge of the evolutionary relation- populator individual to an average genome-wide coverage of ships both between Felinae and Machairodontinae, as well ~0.7x, with an average read depth of ~2x. We achieved this as within Machairodontinae itself, would provide important using a combination of both target enrichment and shotgun insights to help resolve the complex early radiation of Felidae. sequencing. We performed target enrichment using baits created from either DNA extracted from a modern lion (whole-genome The last surviving members of the Machairodontinae belonged capture), or based on the exome annotations of the lion reference to the genus Smilodon (tribe Smilodontini). While once genome (exome capture). The discrepancy between genome- widespread across the continents of North America (S. fatalis) wide coverage and average read-depth likely reflects the use of and South America (S. populator), the genus went extinct ~10 captured data, and lack of a closely related reference genome. thousand years ago (kya)3,4. Smilodon are not known further If we were to have a conspecific reference genome, we would north than southernmost Canada (~42 °N), and their range expect a more even genome-wide coverage, more comparable extended until the tip of the South American continent (~53 °S). to the read depth. Page 3 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 We sequenced approximately 800 million reads from four clearly see a much higher D-statistic (D= 0.54-0.67, Z= 286.6- independently constructed libraries (Table 1). Investiga- 442.5) when placing the sabre-toothed cats parapyletically, tions into mapped reads showed high levels of fragmenta- compared to when they are placed as sister taxa (D= 0.12-0.20, tion (average read length <90bp for all libraries) and high Z= 39.2-67.7). A high D-statistic could be interpreted as rates of C-T transitions across the reads (Figure S1, Extended gene flow between either H1 and H3, or between H2 and H3. data14), both typical of authentic ancient DNA. Comparisons However, another possible explanation could be the incorrect among libraries showed that the capture experiments did yield input topology. Although paraphyly of the sabre-toothed cats significant levels of endogenous DNA, but were not highly dif- followed by gene flow could explain the observed pattern, a ferent from the shotgun data but had lower levels of complexity. more likely explanation for the higher D-statistic when placing Moreover, although the exome capture was able to increase the sabre-tooths paraphyletically, would be the monophyly of the relative coverage of the exome, it also included a lot of the sabre-tooths as seen in our NJ tree. This result further sup- non-coding whole genomic data (Table 1). ports Smilodon and Homotherium as more closely related to one another, than either are to any extant cat species. To investigate the topological placement of Smilodon, we computed a neighbour-joining (NJ) tree using genome-wide To further contextualise this relationship, we built a dated phy- transversion pairwise distances, rooted using Crocuta crocuta logenetic tree (Figure 1). However, due to the low coverage (spotted hyena). Similar to the mitochondrial genome, we of our Smilodon genome, traditional methods for phylogeny found Smilodon and Homotherium to be sister taxa, and dating are likely unsuitable. Therefore, we devised a method to Machairodontinae to be sister taxon to all extant cats (Figure S2, overcome this by computing pairwise genetic drift distances Extended data14). Support for this topology was high; we using F-statistics (in our case F2) (Table S3, Extended data14)16. repeated this analysis independently for the 100 longest scaf- The pairwise F2 statistics were built into an unrooted NJ tree folds (27.2Mb - 5.3Mb), and found identical placement of (Figure S3, Extended data14). The relationships recov- Smilodon for each scaffold. ered in this tree were the same as those obtained using the transversional genetic distances. We further tested the closer affinity of Smilodon and Homoth- erium to each other, relative to any living Felinae, using To calibrate our phylogenetic tree, we estimated the average D-statistics15 to test for topology. For this, we took advantage rate of genetic drift between all pairs of genera within Felinae of the topological input requirement for the D-statistics test based on the F2 results and previously calculated divergence ([[[H1, H2], H3], Outgroup]). We computed the D-statistic dates17. We made a number of assumptions about our data, with Smilodon and Homotherium as sister taxa (i.e. in the H1 including (i) a strict molecular clock, (ii) constant population and H2 positions) (Table S1, Extended data14), each species of sizes through time, and (iii) that Felinae drift rates are similar Felinae as H3, and the spotted hyena as Outgroup, and com- to those in Machairodontinae. pared this to the D-statistic and significance from 0 (Z-score) when placing Smilodon and Homotherium paraphyletically (i.e. Using the F2 statistics and average drift rate, we estimated in the H1 and H3 positions) (Table S2, Extended data14). We a divergence date between Machairodontinae and Felinae of Table 1. Smilodon mapping results. WGC - whole genome capture, SE - single end, PE - paired end. Library name BU19_5 BU23_1 Smi745_1 Smi745_2 Library type Exome capture WGC WGC/Shotgun mix Shotgun # SE reads 145,712,209 7,856,010 125,060,258 64,850,185 # PE reads 164,381,221 309,399,989 - - # PE reads merged 146,380,424 278,659,340 - - # Retained reads 313,124,333 342,949,720 121,988,548 60,427,680 # Mapped reads 6,070,092 6,474,904 5,191,598 7,136,656 Proportion duplicates 0.893 0.830 0.731 0.100 Endogenous proportion 0.019 0.019 0.043 0.118 Exome coverage 1.237 0.217 0.605 0.272 Genome coverage 0.171 0.234 0.139 0.179 Average read length (bp) 68.7 88.3 65.4 61.3 Page 4 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Figure 1. Dated Felidae phylogenetic tree based on genome-wide pairwise F2 statistics, calibrated using an average genetic drift rate estimated from within Felinae (yellow shading). Blue bar shows 95% confidence interval for the divergence between the Smilodon (red) Homotherium lineages, calculated from the F2 standard error. Smilodon illustration by Binia De Cahsan and included with permission. ~22.1 Ma (Figure 1), similar to the estimate of ~22.5 Ma cal- to. This is strong support for the hypothesis that Pseudaelurus as culated using a high-coverage exome dataset17. Although the a grade-group is paraphyletic. Furthermore, current phylog- correlation between these two results was not unexpected, as enies using well-known characters (such as length of upper our tree was calibrated using the average within Felinae drift canines and presence of serrations or crenulations) resolve the calculated using the divergence times of Barnett et al.17, their Smilodontini-Homotheriini divergence to only ~11-10 Ma, similar divergence estimates suggest extrapolating Felinae while older machairodonts constitute a stem lineage of drift rates to other subfamilies of Felidae is valid. uncertain relationship18. Therefore, given the problematic nature of machairodont phylogenetics, the new molecular We further tested for the robustness of this method by both information analysed here provides an important reference recalculating the divergence times within Felinae, and by down- point for identifying morphologically Smilodontini or Homoth- sampling two Felinae species (Caracal caracal - caracal and eriini characters in specimens from >11 Ma, to help resolve Prionailurus bengalensis - leopard cat) to ~0.7x. We recovered machairodont phylogeny back to the early Miocene. similar divergence estimates to those previously reported, and to those produced without downsampling, providing confidence We next investigated the evolutionary relationships between in the use of this methodology for low-coverage genomes Smilodon and Felinae by searching for signatures of gene flow, (Tables S3 and S4, Extended data14). Using this methodology, using two independent methods (f316 and D319). The results we estimated the divergence date of Smilodon/Homotherium were assessed using simulated data (Tables S6, S7, and S8, to be ~20.65 Ma (95% CI 26.07-15.25 Ma) (Figure 1; Table S5, Extended data14). As lineages leading to Smilodon and Extended data14). Furthermore, although based on low- Homotherium diverged relatively soon after the Machairodontinae/ coverage data and a number of data assumptions, the congru- Felinae split, we were able to test for ancient gene flow (up to ence of our results with previous mitogenome-based estimates 20 Ma) between Smilodon or Homotheirum and stem Felinae. from the same individuals8, despite the use of different data Ancient gene flow may have prevented the diverging and calibration methods, further adds confidence to our method lineages in early Felidae from accumulating obvious morpho- of phylogenetic assessment. logical differences, thus preventing the confident phylogenetic placement of lineages during this time period. Species traditionally assigned to Pseudaelurus originated in the Early Miocene, ~20 Ma. Given the deep divergence Moreover, we tested for more recent gene flow events between between Machairodontinae and Felinae, as well as within either Smilodon or Homotherium and the entire Panthera big Machairodontinae, this indicates that some species of Pseudae- cat lineage (all Panthera species grouped together), due to lurus are younger than the lineages they purportedly gave rise their potential size overlap, as well as with single species that Page 5 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 may have more recently met Smilodon due to spatial proximity using the int13.14c calibration curve. It has been identified (i.e. from South America: Panthera onca (jaguar), Puma concolor as a left tibia of Smilodon populator and is part of the Kruimel (puma), and Leopardus pardalis (ocelot)). collection, an assortment of megafaunal remains recovered from the Última Esperanza region of Chile, most likely from the site We did not find any indication of gene flow between any of of Cueva del Milodón. This and other specimens of Smilodon the lineages tested, regardless of method (Tables S6 and S7, from the Kruimel Collection were described and figured by Extended data14). The lack of more recent gene flow events McDonald and Werdelin (2018); specimen ZMA20.042 is is not surprising, due to the relatively ancient divergence of presented in Figure 4.6D5. We additionally included genomic Machairodontinae and Felinae. However, the lack of ancient gene information from another extinct sabre-toothed cat, 12 extant flow signatures do not necessarily mean that gene flow was cat species, and the spotted hyena (Table 2). not present during the early divergence of these lineages. Rather, it may result from the inadequacy of the methods in Ancient DNA extraction, library preparation, and uncovering such events. sequencing Samples of cortical bone were taken from the long bone We assessed the power of the D3 statistic to detect gene flow, element (approx. 1 cm3) using a Dremel drill, and reduced to which occurred at different times (20 Ma-50 kya), using sim- powder in a Mikrodismembrator. Two independent DNA extrac- ulated data. When using a simple demographic model of tions were performed as described in Orlando et al.21 in a constant population size, mutation rate, and recombination rate dedicated ancient DNA laboratory, with negative controls. We on error-free simulated data, we detect significant levels of gene built each DNA extract and negative control into genomic librar- flow back to ~16 Ma (Table S8, Extended data14). However, ies using the NEB E6070 kit following a modified version of empirical data do not always behave in such a simple manner, the protocol used by Vilstrup et al.22. Briefly, the extract (30 µl) and a variety of factors may influence the results of the D3 was end-repaired and cleaned using a MinElute column, the statistic. These include, but are not restricted to: violations of collected flow-through was adapter-ligated and cleaned using a the infinite sites model, different mutation rates across lineages, QiaQuick column, and the adapter fill-in reaction was per- ancestral populations structure, and introgression with unsam- formed on the flowthrough. To complete the library build, we pled lineages19. However, our results suggest that D3 may be performed a final incubation at 37°C (30 min) followed by suitable for highly divergent lineages with ancient gene flow inactivation overnight at -20°C. events and therefore violations of the infinite sites model may not be as problematic. Thus, although we are unable to exclude the possibility of Table 2. List of the additional species ancient gene flow during the early radiation of Felidae, we are included in this study. Original sources and somewhat more confident that there was no more recent gene accession codes for all raw reads used in this flow between either Smilodon or Homotherium and Felinae study can be found in Table S9 (see Extended data14). within the last 16 Ma years. If very ancient gene flow had occurred prior to 16 Ma, it is likely that recombination and genetic drift would have either highly fragmented the intro- Species name: Common name: gressed regions, or completely removed them from contempo- Felis catus Domestic cat rary Felinae individuals. This would render the regions too small to be detected with current methods, which use genome-wide Prionailurus bengalensis Leopard cat summary statistics or regions of phylogenetic incongruence Lynx pardinus Iberian lynx with known evolutionary relationships, to infer gene flow. Acinonyx jubatus Cheetah Our study exemplifies how even a draft palaeogenome from Caracal caracal Caracal an extinct species can provide important information into their evolutionary history. Through the sequencing of a single Smilo- Neofilis nebulosa Clouded leopard don populator genome, we provide insights into Felidae’s early Panthera uncia Snow leopard radiation in the Early–Middle Miocene (~23 - 15 Ma), which could not be uncovered using genetic data from extant species Panthera onca Jaguar alone. Panthera pardus Leopard Panthera leo Lion Methods Sample information Crocuta crocuta Spotted hyena Specimen ZMA20.042 is from the Naturalis museum, Leiden, Homotherium latidens Scimitar-toothed cat the Netherlands, and was radiocarbon dated to 11,335±30 uncalibrated years before present (2-sigma range of 13,269- Puma concolor Puma 13,095 calibrated years before present) (Stafford: UCIAMS- 142836). We calibrated the radiocarbon date using Calib v7.0420 Leopardus pardalis Ocelot Page 6 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 For each library, we performed two independent indexing realignment step to adjust for increased error rates at the PCR amplifications (Veriti thermal cycler, Applied Biosys- end of short reads in the presence of indels. In the absence tems) in a 50 µl volume reaction, using 25 µl of library, 25 µl of a curated dataset of indels, this step relied on a set of indels PCR master mix, and 12 cycles of PCR reactions. The final identified in the specific sample being processed. Read concentrations in the PCR master mix were 1.25 U Accu- damage patterns were assessed and base quality scores recali- Prime™ Pfx DNA Polymerase (Invitrogen, Cat # 12344-024), 1x brated around read terminal damage patterns using mapDamage AccuPrime™ Pfx reaction mix (Invitrogen, Cat # 12344-024), v2.0.530,31. 0.4 mg/ml BSA, 120 nM primer InPE, 120 nM of a multiplex- ing indexing primer containing a unique six-nucleotide index Data processing of the extant Felidae species, excluding Puma code (Illumina – sequences TGCAGG, CGATGA, GCGAGA, concolor and Leopardus pardalis, followed the same pipeline or CAGCAC). PCR cycling conditions consisted of an initial with the following minor adjustments; no minimum read length denaturation step at 95°C for 2 min, followed by 12 cycles of cut-off or missing data cut-off was applied during the 95°C denaturation for 15 s, 60°C annealing for 30 s, and 68°C adapter trimming step, and bases were not recalibrated using extension for 30 s, and a final extension step at 68°C for 7 min. mapDamage. These steps were removed as the data from the Indexed libraries were checked for presence of DNA on a 2% extant species would not display the highly fragmented and dam- agarose gel before purification using the QIAquick column age patterns found in ancient DNA. The Puma concolor and system (Qiagen, Cat # 28104) and quantification on an Agilent Leopardus pardalis samples had Illumina adapter and short 2100 BioAnalyzer. Quantified libraries were communally sequences trimmed using skewer 0.2.232 but followed the same pooled in equimolar ratios and sequenced as 100 bp single-end protocol as the other extant species for the rest of the processing reads on an Illumina HiSeq2000 platform at the Danish National steps. High-Throughput Sequencing Centre and 100 bp paired-end reads on an Illumina Hiseq2000 at BGI Copenhagen. Neighbour-joining tree To build a NJ phylogenetic tree, we computed an identity by Genome capture state distance matrix considering only transversion differ- We assessed the shotgun-sequenced libraries for endogenous ences using ANGSD v0.92133, and specifying the following content and selected the libraries with the highest levels of parameters; call the consensus base (-doIBS 2), minimum endogenous DNA for two sets of capture experiments. The mapping and quality of 30 (-minmapq 30, -minq 30), only first set used biotinylated RNA probes transcribed from fresh include a site if all individuals are covered (-minInd), remove DNA extracted from a modern lion, for the purpose of enrich- secondary alignments (-remove_bads 1), only include reads that ing the entire nuclear genome (whole-genome capture). The sec- map to a single location (-uniqueonly 1), compute major and ond method used biotinylated RNA baits, assembled based on the minor alleles based on genotype likelihoods (-domajorminor 1), exonic annotations of lion genomic data23 (exome capture). remove transitions (-rmtrans 1), print all sites (-minminor 0), Both types of baits were generated by Arbor Biosciences (Ann and use the GATK algorithm to compute genotype likelihoods Arbor, MI, USA) and carried out using the myBaits target (-GL 2) and only include scaffolds over 1Mb in length (-rf). enrichment kit and the instructions described in manual V3. After filtering, 221,350,529 sites remained. After capture and cleanup, enriched libraries were re-amplified We further checked for support of the genome-wide for further sequencing using either a Phusion polymerase topology by building a distance matrix for the 100 (New England Biolabs, Cat # M0530S) or a KAPA HiFi longest scaffolds independently, resulting in 100 independ- HotStart polymerase (Roche, Cat # KK2801 07959052001) with ent distance matrices based on 171,281 - 3,344,878 sites. The primers IS5_reamp.P5 and IS6_reamp.P7 over 14 cycles24. resultant distance matrices were converted into NJ trees We quantified the resultant enriched libraries on a TapeSta- using fastME v2.1.534. tion 2200 instrument and sequenced them on an Illumina Hiseq2000 at the Danish National High-throughput Sequencing D-statistics topology test Centre. To investigate the closer relationship of Smilodon and Homotherium to each other relative to other extant cat species, Data processing pipeline we implemented a D-statistics test for topology35. Although Post-sequencing read processing of the Smilodon populator D-statistics is most commonly used to find evidence of gene was performed using the PALEOMIX v1.2.5 pipeline25. flow, it can also be used to test for phylogenetic relationships. Adaptor sequence removal and trimming of low-quality bases This test takes advantage of the fact that D-statistics (BaseQ < 5 or Ns) was done with AdapterRemoval v2.0.026. This relies on a predefined four-taxon topology as input step also removed all reads shorter than 30 bp in length or [[[H1,H2],H3],outgroup]. A high D-score is most commonly with more than 50 bp of missing data. Trimmed reads were used to infer post-divergence gene flow, but it can also be mapped against an African lion reference genome (NCBI caused by more recent common ancestry brought about by an Accession JAAVKH00000000027) with BWA-MEM v0.7.5a28, incorrect predefined topology. utilising default parameters. PCR duplicates were identi- fied and filtered based on the 5’-end mapping coordinate using Taking the latter into account, we performed a number of Picard v2.18.0. GATK v3.8.029 was used to perform an indel D-statistics tests including the topologies: [[[Smilodon, Page 7 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Homotherium], extant cat], Crocuta crocuta], [[[Smilodon, extant this analysis we assumed a known species tree and divergences cat], Homotherium], Crocuta crocuta], and [[[Homotherium, within Felinae previously found in Barnett et al.17, a genome- extant cat], Smilodon], Crocuta crocuta]. In this test, ‘extant wide constant mutation/drift rate, and no variation in drift cat’ was replaced with each living Felidae species included in rates between lineages. Uncertainty around the divergence this study. We used ANGSD v0.92133 to perform the D-statistics between Smilodon and Homotherium was calculated using with the following parameters; -minmapq 30 -minq 30, -minind 1.96 x the standard deviation on either side of the F2 statistic. 15, -remove_bads 1, -uniqueonly 1, -domajorminor 1, -rmtrans The standard deviation was calculated using the standard 1, -GL 2, calculate D in a block size of 1Mb (-blocksize), error x √ N (where N = number of blocks using for jackknifing). and the spotted hyena as the ancestral/outgroup sequence (-anc). Assessing gene flow We implemented two independent analyses to test for the Dated phylogeny presence of signs of gene flow between Machairodontinae and We computed the consensus haploid base calls (-dohaplocall Felinae. We computed F3 statistics16 to assess whether there 2) for all scaffolds greater than 1MB in length using ANGSD, were any signals of gene flow between a predefined triplet specifying the following parameters; minimum mapping and [[A,B],C] using the same PLINK file computed above for quality of 30 (-minmapq 30, -minq 30), only include a site if all the F2 statistics. We used popstats.py36 for five independent individuals are covered (minInd), remove secondary align- triplet combinations, which we present as ((A,B),C) in the ments (-remove_bads 1), only include reads that map to a sin- subsequent text. We placed Smilodon and Homotherium in gle location (-uniqueonly 1), compute major and minor alleles the A and B positions, while alternating C. based on genotype likelihoods (-domajorminor 1), remove transitions (-rmtrans 1), print all sites (-minminor 0), and use First, we investigated signs of very ancient gene flow between the GATK algorithm to compute genotype likelihoods (-GL 2). Machairodontinae and Felinae by specifying all extant Felinae After filtering, 268,152,250 sites remained. as one population. Next, we looked for gene flow with the Panthera spp. big cats by specifying all Panthera as a single We converted the haploid call file into a PLINK file format population. Finally, we computed F3 three times independ- using the haplo2plink command from the ANGSD toolsuite. ently to test for signs of very recent gene flow between Using the resultant PLINK file, we calculated F2 statistics15 Machairodontinae and either of three extant cat species occu- for each pairwise combination of our Smilodon populator pying South America, which may have come into contact with individual, Homotherium latidens, and 12 Felinae species using Smilodon based on geography (Panthera onca (jaguar), Puma the popstats.py script36. We used genetic drift calculated via concolor (puma), and Leopardus pardalis (ocelot)). F-statistics as opposed to absolute pairwise distance as it should be more suitable for ancient DNA data37. From these To complement this analysis, we also computed D3 statistics18 pairwise F2 comparisons, we built a distance matrix that on the same five triplet comparisons. For this, we also was converted into a newick tree file using PHYLIP v3.696 produced a haplocall file in ANGSD using the same filtering neighbor38 for visualisation. parameters as above, but using a random base call (-dohaplo- call 1), as opposed to the consensus base call, while specifying From the distance matrix, we computed the Felinae average the same additional parameters specified above. The resultant rate of genetic drift (F2) to be 0.000305 per million years. For haploid output was then converted into a geno file and run this, we calculated the average F2 between pairwise through the popgenWindows.py python script. We ran the comparisons of all genera within Felinae that had divergence popgenWindows.py script using default parameters, speci- time estimates available in Barnett et al.17. These included fying a window size of 1MB, and the minimum number of Acinonyx, Caracal, Felis, Lynx, Neofelis, Panthera, and sites per window as 1kb. From this output, we calculated D3 Prionailurus, and resulted in 21 comparisons (Table S3, for each window independently by applying the equation Extended data14). We used the formula of (F2 × 0.5)/divergence [[BC-AC]/[BC+AC]]. Mean values, standard deviations, and time to estimate the rate of genetic drift per million years for all significance from 0 were measured in R v3.6.0 using the pnorm 21 comparisons. We calculated the mean of these 21 com- function40. parisons, giving us the Felinae mean rate of genetic drift per million years. This mean rate was used in conjunction with the Evaluating D3 for detecting ancient gene flow previously calculated average F2 to calibrate the tree includ- To evaluate the adequacy of the D3 method for detecting very ing the Smilodon populator, Homotherium latidens, puma, ancient gene flow events (up to 20 Ma), we simulated three and ocelot. diploid sequences representing Homotherium (A), Smilodon (B), and a Felinae species (C), using msprime41. As input, we To test for the robustness of our method to low-coverage estimated the average transversion mutation rate within Felidae data, we downsampled both Caracal and Prionailurus to using the average pairwise distance between Homotherium comparable coverage to Smilodon populator (~0.7x) using and Felinae and divided this by two before multiplying it by SAMtools v1.639, and recomputed the F2 statistics. We used the Machairodontinae/Felinae divergence time calculated the same average rate of genetic drift to estimate the diver- above (22.1 Ma). We computed pairwise comparisons between gence of these genera from their closest relatives (i.e. Homotherium and each Felinae species included in the study Prionailurus - Felis, and Caracal - Acinonyx/Felis/Lynx). For using ANGSD and took the average pairwise distance from Page 8 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 these comparisons. We excluded the Smilodon from this Data availability calculation due to the lower quality of the genome. Underlying data NCBI BioProject: Raw sequencing reads for the Smilodon To calculate the pairwise distances we used the consensus individual. Accession number PRJNA691254; https://www.ncbi. identity by state parameter (-doIBS 2) in ANGSD and applied the nlm.nih.gov/bioproject/691254. following filters; -minmapq 30, -minq 30, -minind 15, -unique- only 1, -remove_bads 1, -domajorminor 1, -rmtrans 1, -minminor Extended data 0, -makematrix 1, and only included scaffolds >1MB in length. Zenodo: A genomic exploration of the early evolution of This gave us an average pairwise distance of 0.01207, which extant cats and their sabre-toothed relatives - extended data. we converted into an average transversion mutation rate of https://doi.org/10.5281/zenodo.492245014. 2.730769e-10 per year. We converted this to a generational mutation rate using a generation time calculated for lions of This project contains the following extended data: 6.5 years42; if the generation time of the investigated species - Westbury et al extended data.pdf (Supplementary differed from that of lion, the mutation rate would adjust Figures S1-S3 and Tables S1-S9) accordingly, and we therefore expect minimal impact of this parameter on the final results, especially as we ran the simulations specifying years as opposed to number of generations. For the Data are available under the terms of the Creative Commons recombination rate, we used a previously published sex-averaged Attribution 4.0 International license (CC-BY 4.0). recombination rate of 1.9 cM/Mb (1.9e-8)43. Using the above information, we ran five independent simula- Acknowledgements tions, each consisting of 2,000 1 Mb windows with a constant We thank the laboratory technicians of the Centre for effective population size of 40,000 individuals, a genera- GeoGenetics and the staff of the Danish National High- tion time of 6.5 years, the above-mentioned mutation and Throughput DNA Sequencing Centre for technical assistance. recombination rates, a 5% pulse of migration (m=0.05), and We thank Tom Stafford Jr and Stafford Research LLC for radio- a different timing of the migration pulse for each of the five carbon dating and discussion. We thank Jean-Marie Rouillard runs (20 Ma, 18 Ma, 17 Ma, 16 Ma, 15 Ma, 10 Ma, 5 Ma, (Arbor Biosciences) for help with bait design. We would like 50 kya). The scripts for these simulation runs can be found to thank Binia De Cahsan for the illustration in Figure 1. on GitHub. We calculated the D3 statistic from this output Lastly, we would like to thank Reinier van Zelst and Caroline using the tskit toolkit44 using the d3.py script. Significance was Pepermans Naturalis, Leiden, the Netherlands, for access calculated as it was done for the empirical data above. to Smilodon populator samples. References 1. Werdelin L, Yamaguchi N, Johnson WE, et al.: Phylogeny and evolution of cats 3330–3336.e5. (Felidae). Biology and conservation of wild felids. Oxford, 2010; 59–82. PubMed Abstract | Publisher Full Text Reference Source 9. Li G, Figueiró HV, Eizirik E, et al.: Recombination-Aware Phylogenomics 2. Johnson WE, Eizirik E, Pecon-Slattery J, et al.: The late Miocene radiation of Reveals the Structured Genomic Landscape of Hybridizing Cat Species. Mol modern Felidae: a genetic assessment. Science. 2006; 311(5757): 73–77. Biol Evol. 2019; 36(10): 2111–2126. PubMed Abstract | Publisher Full Text PubMed Abstract | Publisher Full Text | Free Full Text 3. Werdelin L, McDonald HG, Shaw CA: Smilodon: The Iconic Sabertooth (JHU 10. Metcalf JL, Turney C, Barnett R, et al.: Synergistic roles of climate warming Press). 2018. and human occupation in Patagonian megafaunal extinctions during the Reference Source Last Deglaciation. Sci Adv. 2016; 2(6): e1501682. 4. Kurtén B, Werdelin L: Relationships between North and South American PubMed Abstract | Publisher Full Text | Free Full Text Smilodon. J Vert Paleontol. 1990; 10(2): 158–169. 11. Barnett R, Barnes I, Phillips MJ, et al.: Evolution of the extinct Sabretooths Publisher Full Text and the American cheetah-like cat. Curr Biol. 2005; 15(15): R589–90. 5. McDonald HG, Werdelin L: The sabertooth cat Smilodon populator (Carnivora: PubMed Abstract | Publisher Full Text Felidae), from Cueva del Milodón, Chile. In Smilodon: The Iconic Sabertooth, L. 12. Bocherens H, Cotte M, Bonini R, et al.: Paleobiology of sabretooth cat Werdelin, H. G. McDonald, and C. A. and Shaw, eds. (Baltimore: John Hopkins Smilodon populator in the Pampean Region (Buenos Aires Province, University Press), 2018; 53–75. Argentina) around the Last Glacial Maximum: Insights from carbon and Reference Source nitrogen stable isotopes in bone collagen. Palaeogeogr Palaeoclimatol 6. Janczewski DN, Yuhki N, Gilbert DA, et al.: Molecular phylogenetic inference Palaeoecol. 2016; 449: 463–474. from saber-toothed cat fossils of Rancho La Brea. Proc Natl Acad Sci U S A. Publisher Full Text 1992; 89(20): 9769–9773. 13. Hofreiter M, Paijmans JLA, Goodchild H, et al.: The future of ancient DNA: PubMed Abstract | Publisher Full Text | Free Full Text Technical advances and conceptual shifts. Bioessays. 2015; 37(3): 284–293. 7. Gold DA, Robinson J, Farrell AB, et al.: Attempted DNA extraction from a PubMed Abstract | Publisher Full Text Rancho La Brea Columbian mammoth (Mammuthus columbi): prospects for 14. Westbury M, Barnett R, Sandoval-Velasco M, et al.: A genomic exploration ancient DNA from asphalt deposits. Ecol Evol. 2014; 4(4): 329–336. of the early evolution of extant cats and their sabre-toothed relatives PubMed Abstract | Publisher Full Text | Free Full Text - extended data. [Data set]. Zenodo. 2021. 8. Paijmans JLA, Barnett R, Gilbert MTP, et al.: Evolutionary History of http://www.doi.org/10.5281/zenodo.4434076 Saber-Toothed Cats Based on Ancient Mitogenomics. Curr Biol. 2017; 27(21): 15. Green RE, Krause J, Briggs AW, et al.: A draft sequence of the Neandertal Page 9 of 16 Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 genome. Science. 2010; 328(5979): 710–722. damage patterns in ancient DNA sequences. Bioinformatics. 2011; 27(15): PubMed Abstract | Publisher Full Text | Free Full Text 2153–2155. 16. Reich D, Thangaraj K, Patterson N, et al.: Reconstructing Indian population PubMed Abstract | Publisher Full Text history. Nature. 2009; 461(7263): 489–494. 31. Jónsson H, Ginolhac A, Schubert M, et al.: mapDamage2.0: fast approximate PubMed Abstract | Publisher Full Text | Free Full Text Bayesian estimates of ancient DNA damage parameters. Bioinformatics. 17. Barnett R, Westbury MV, Sandoval-Velasco M, et al.: Genomic adaptations 2013; 29(13): 1682–1684. and evolutionary history of the extinct scimitar-toothed cat, Homotherium PubMed Abstract | Publisher Full Text | Free Full Text latidens. Curr Biol. 2020; 30(24): 5018–5025.e5. 32. Jiang H, Lei R, Ding SW, et al.: Skewer: a fast and accurate adapter trimmer PubMed Abstract | Publisher Full Text | Free Full Text for next-generation sequencing paired-end reads. BMC Bioinformatics. 2014; 15: 182. 18. Werdelin L, Flink T: The phylogenetic context of Smilodon. In Smilodon: PubMed Abstract | Publisher Full Text | Free Full Text The Iconic Sabertooth, L. Werdelin, H. G. McDonald, and C. A. and Shaw, eds. (Baltimore: John Hopkins University Press), 2018; 14–29. 33. Korneliussen TS, Albrechtsen A, Nielsen R: ANGSD: Analysis of Next Reference Source Generation Sequencing Data. BMC Bioinformatics. 2014; 15(1): 356. PubMed Abstract | Publisher Full Text | Free Full Text 19. Hahn MW, Hibbins MS: A Three-Sample Test for Introgression. Mol Biol Evol. 2019; 36(12): 2878–2882. 34. Lefort V, Desper R, Gascuel O: FastME 2.0: A Comprehensive, Accurate, and PubMed Abstract | Publisher Full Text Fast Distance-Based Phylogeny Inference Program. Mol Biol Evol. 2015; 32(10): 2798–2800. 20. Reimer PJ, Bard E, Bayliss A, et al.: IntCal13 and Marine13 Radiocarbon Age PubMed Abstract | Publisher Full Text | Free Full Text Calibration Curves 0–50,000 Years cal BP. Radiocarbon. 2013; 55(4): 1869–1887. Publisher Full Text 35. Westbury MV, Hartmann S, Barlow A, et al.: Extended and Continuous Decline in Effective Population Size Results in Low Genomic Diversity in the World’s 21. Orlando L, Ginolhac A, Zhang G, et al.: Recalibrating Equus evolution using Rarest Hyena Species, the Brown Hyena. Mol Biol Evol. 2018; 35(5): 1225–1237. the genome sequence of an early Middle Pleistocene horse. Nature. 2013; PubMed Abstract | Publisher Full Text | Free Full Text 499(7456): 74–78. PubMed Abstract | Publisher Full Text 36. Skoglund P, Mallick S, Bortolini MC, et al.: Genetic evidence for two founding populations of the Americas. Nature. 2015; 525(7567): 104–108. 22. Vilstrup JT, Seguin-Orlando A, Stiller M, et al.: Mitochondrial phylogenomics of PubMed Abstract | Publisher Full Text | Free Full Text modern and ancient equids. PLoS One. 2013; 8(2): e55950. PubMed Abstract | Publisher Full Text | Free Full Text 37. Harris AM, DeGiorgio M: Admixture and Ancestry Inference from Ancient and Modern Samples through Measures of Population Genetic Drift. Hum 23. Cho YS, Hu L, Hou H, et al.: The tiger genome and comparative analysis with Biol. 2017; 89(1): 21–46. lion and snow leopard genomes. Nat Commun. 2013; 4: 2433. PubMed Abstract | Publisher Full Text PubMed Abstract | Publisher Full Text | Free Full Text 38. Felsenstein J: PHYLIP (Phylogeny Inference Package) version 3.6. 2005. 24. Meyer M, Kircher M: Illumina sequencing library preparation for highly Reference Source multiplexed target capture and sequencing. Cold Spring Harb Protoc. 2010; 39. Li H, Handsaker B, Wysoker A, et al.: The Sequence Alignment/Map format 2010(6): pdb.prot5448. and SAMtools. Bioinformatics. 2009; 25(16): 2078–2079. PubMed Abstract | Publisher Full Text PubMed Abstract | Publisher Full Text | Free Full Text 25. Schubert M, Ermini L, Der Sarkissian C, et al.: Characterization of ancient and 40. R Core Team: R: A language and environment for statistical computing. modern genomes by SNP detection and phylogenomic and metagenomic 2013. analysis using PALEOMIX. Nat Protoc. 2014; 9(5): 1056–1082. Reference Source PubMed Abstract | Publisher Full Text 41. Kelleher J, Etheridge AM, McVean G: Efficient Coalescent Simulation and 26. Schubert M, Lindgreen S, Orlando L: AdapterRemoval v2: rapid adapter Genealogical Analysis for Large Sample Sizes. PLoS Comput Biol. 2016; 12(5): trimming, identification, and read merging. BMC Res Notes. 2016; 9: 88. e1004842. PubMed Abstract | Publisher Full Text | Free Full Text PubMed Abstract | Publisher Full Text | Free Full Text 27. de Manuel M, Barnett R, Sandoval-Velasco M, et al.: The evolutionary history 42. Black SA, Fellous A, Yamaguchi N, et al.: Examining the extinction of the of extinct and living lions. Proc Natl Acad Sci U S A. 2020; 117(20): Barbary lion and its implications for felid conservation. PLoS One. 2013; 8(4): 10927–10934. e60174. PubMed Abstract | Publisher Full Text | Free Full Text PubMed Abstract | Publisher Full Text | Free Full Text 28. Li H, Durbin R: Fast and accurate short read alignment with 43. Li G, Hillier LW, Grahn RA, et al.: A High-Resolution SNP Array-Based Linkage Burrows-Wheeler transform. Bioinformatics. 2009; 25(14): 1754–1760. Map Anchors a New Domestic Cat Draft Genome Assembly and Provides PubMed Abstract | Publisher Full Text | Free Full Text Detailed Patterns of Recombination. G3 (Bethesda). 2016; 6(6): 1607–1616. 29. McKenna A, Hanna M, Banks E, et al.: The Genome Analysis Toolkit: a PubMed Abstract | Publisher Full Text | Free Full Text MapReduce framework for analyzing next-generation DNA sequencing 44. Ralph P, Thornton K, Kelleher J: Efficiently Summarizing Relationships in data. Genome Res. 2010; 20(9): 1297–1303. Large Samples: A General Duality Between Statistics of Genealogies and PubMed Abstract | Publisher Full Text | Free Full Text Genomes. Genetics. 2020; 215(3): 779–797. 30. Ginolhac A, Rasmussen M, Gilbert MTP, et al.: mapDamage: testing for PubMed Abstract | Publisher Full Text | Free Full Text Page 10 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Open Peer Review Current Peer Review Status: Version 1 Reviewer Report 19 April 2021 https://doi.org/10.21956/openreseurope.14174.r26653 © 2021 Salesa M. This is an open access peer review report distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Manuel Salesa Departamento de Paleobiología, Museo Nacional de Ciencias Naturales-CSIC, Madrid, Spain In the present work, the authors show the results of the sequencing of a draft nuclear genome of the largest sabre-toothed felid, Smilodon populator, with an age of around 13,000 years ago, making this the oldest palaeogenome from South America to date, a region known to be problematic for ancient DNA preservation. The analysis of this genome, together with genomes from other extinct and extant cats, allows the authors to suggest a series of very interesting phylogenetic proposals. Among them, they confirm a deep divergence (around 20.65 Ma) within Machairodontinae, something that was evident from the fossil record, but never stated from palaeogenomics. I haven’t found major concerns about this work. The methodology seems impeccable, the figure is informative, the taxonomic citations are correct, and the discussion is well supported by the data. The style of the text is clear, and the manuscript is correctly written. I just detected a grammatical mistake when citing S. populator on page 3, which is cited as “popular”. Finally, I find the subject treated in the manuscript suitable for its indexing on Open Research Europe. MINOR COMMENTS Although fossils of Early and Middle Miocene felids are known since the XIX century, our understanding of their phylogenetic relationships has not changed greatly over this time. This was in part due to the fact that those fossils of primitive felids, mostly fragmentary dental and postcranial material, show a very similar morphology to that of extant species, and only sabre- toothed taxa exhibit evident differences in their elongated and flattened upper canines or their huge carnassial teeth. In the present paper, Westbury et al. have sequenced the genome of one of the most iconic sabre-toothed felids, the giant Smilodon populator, and then they compared it to those of extant species of Felidae, as well as to that of another extinct sabre-toothed species, the relatively slender and cursorial Homotherium latidens. The results were expectable but very interesting: they not only show the strong divergence between Machairodontinae (sabre-toothed felids) and Felinae (extant felids), but also between the two analysed species of sabre-toothed felids, S. populator and H. latidens. These two derived forms were the last sabre-toothed felids that Page 11 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 were present in the faunas, and both were extinct as recently as around 28,000–11,000 years. Smilodon and Homotherium illustrate two very different ways to be a sabre-toothed felid: S. populator, included in a lineage known as Smilodontini, was a strong-built animal, with short and powerful limbs and very elongated upper canines, whereas H. latidens was a more slender predator, with relatively long limbs (the hind limbs much shorter than the forelimbs) and shorter canines. Thus, the divergence found by Westbury et al. fits well with the idea that the sabre- toothed felids developed at least two different morphotypes very early in their evolution, a robust, forest-dweller type illustrated by the Middle Miocene Pseudaelurus quadridentatus and the Late Miocene Promegantereon ogygia, and a slenderer model, probably occupying more open environments, whose first well-known representative is the Late Miocene Machairodus aphanistus. These two lineages evolved in parallel, refining each of the morphotypes, occupying different habitats thanks to different locomotor adaptations, thus allowing the presence in the faunas of two large predators through habitat segregation. Concerning the early divergence between sabre-toothed felids and felines, the present work by Westbury et al. points towards a date of around 22.1–22.5 Ma. This is interesting, as we have little data on the felids of that moment, but all the available fossils suggest that the split between felines and macairodontines occurred at that time. There are several European localities with fossils slightly younger than that age (around 18 Ma) that can be identified as early felines, as their dentitions exhibit short canines of a rounded section, although they also show some very primitive features such as the presence of m2; and besides these forms, there are other fossils that show an incipient development of the macairodontine model illustrated by the laterally flattened upper canines of Pseudaelurus quadridentatus. Now seems clear that this latter genus should be considered as the first known sabre-toothed felid, as it has laterally flattened upper canines, a relatively robust first metacarpal, and an overall postcranial anatomy very similar to that of the younger Promegantereon ogygia (Beaumont, 19751; Salesa et al., 20102). Nevertheless, as the present work by Westbury et al. remarks, the genus Pseudaelurus has traditionally included felids with evident differences in their dental morphology, both from North America and Eurasia; the paraphyly of this genus was already suggested by Werdelin et al. (2010)3 but this new work reinforces this idea. And finally, it is very interesting that, although Westbury et al. remark that they were unable to exclude the possibility of ancient gene flow during the early radiation of Felidae, they are somewhat more confident that there was no more recent gene flow between either Smilodon or Homotherium and Felinae within the last 16 Ma. Felinae and Machairodontinae developed a body plan devoted to preying upon different kinds of prey, and their dental and skeletal morphology suggests an early independent evolution of each model. This absence of gene flow also suggests the existence of several ethological and physiological differences between felines and sabre- toothed felids whose trace unfortunately vanished when the latter disappeared. References 1. de Beaumont G: Recherches sur les Félidés (Mammifères, Carnivores) du Pliocène inférieur des sables à Dinotherium des environs d’Eppelsheim (Rheinhessen). Archives des Sciences Physiques et Naturelles, Genève. 28: 369-405 2. SALESA M, ANTÓN M, TURNER A, ALCALÁ L, et al.: Systematic revision of the Late Miocene sabre- toothed felid Paramachaerodus in Spain. Palaeontology. 2010; 53 (6): 1369-1391 Publisher Full Text 3. Werdelin L, Yamaguchi N, Johnson WE, O’Brien SJ: Phylogeny and evolution of cats (Felidae). Biology and Conservation of Wild Felids - Oxford University Press. 2010. 59-82 Page 12 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Is the work clearly and accurately presented and does it cite the current literature? Yes Is the study design appropriate and does the work have academic merit? Yes Are sufficient details of methods and analysis provided to allow replication by others? Yes If applicable, is the statistical analysis and its interpretation appropriate? Yes Are all the source data underlying the results available to ensure full reproducibility? Yes Are the conclusions drawn adequately supported by the results? Yes Competing Interests: No competing interests were disclosed. Reviewer Expertise: Palaeobiology, Systematics, Functional anatomy I confirm that I have read this submission and believe that I have an appropriate level of expertise to confirm that it is of an acceptable scientific standard. Author Response 10 Jun 2021 Michael Westbury, The GLOBE Institute, University of Copenhagen, Øster Voldgade 5-7, Copenhagen, Denmark We would like to thank Manuel Salesa for praising our work and giving a very interesting summary of the results within the context of prior literature. We believe this additional information published alongside our manuscript to be helpful to the reader to better understand the background of the study. We have also corrected the mistake when citing S. populator on page 3, which was cited as “popular”. Competing Interests: No competing interests were disclosed. Reviewer Report 08 April 2021 https://doi.org/10.21956/openreseurope.14174.r26645 © 2021 Förster D. This is an open access peer review report distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Page 13 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Daniel Förster Evolutionary Genetics Department, Leibniz Institute for Zoo and Wildlife Research (IZW), Berlin, Germany Westbury et al. use nuclear data of a Smilodon populator individual to resolve the phylogenetic relationship of this species with another extinct sabre-toothed cat (Homotherium latidens), and of the Machairodontinae (here, Smilodon and Homotherium) with the extant cats (Felinae). They illuminate the deeper divergences within this Carnivoran family, i.e. the Felidae. The study confirms that Smilodon and Homotherium are sister lineages, previously suggested by mtDNA data, and, as the authors point out, these results have implications for resolving the phylogeny of the machairodonts based on morphological characters. Gene flow between lineages during the Early Miocene radiation of the Felidae was not detected, within machairodonts or between machairodonts and stem Felinae, potentially owing to available data and analytical approaches that were possible. An absence of gene flow between machairodonts and extant cats in the more recent past is not too surprising. An interesting study illuminating the history of a charismatic group. Unfortunately, perhaps owing to genome coverage, the authors did not examine genomic adaptations as Barnett et al. 2020 did for Homotherium latidens. Tables & Figures It seems unnecessary to avoid any Tables in the main text. For example, having the table of species analysed (Table S10) in the Methods section would benefit clarity, as these data suddenly find mention in the data processing section, as well as being featured in multiple subsequent sections. The mapping results (Table S1) are summarized in the Results & Discussion section, but the table could feature in the main text as well. Plain language summary “assess” rather than “access”? Introduction The first sentence reads rather clunky. Results & Discussion The start of this section would benefit from more details about (1) baits used for capture (based on data available for the lion; exome capture missing here?), and (2) mapping approach used (i.e. to which species: lion). Can their approach to detect gene flow be applied to some Felinae (e.g. Panthera) in their dataset to confirm that it does indeed function for more recent events? Methods Exome capture used DNA baits? (typo?) Why was the lion genome used as reference? Why not a genome with higher contiguity? Did this Page 14 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 not potentially result in unnecessary loss of data during filtering (e.g. due to scaffold lengths)? Is the work clearly and accurately presented and does it cite the current literature? Yes Is the study design appropriate and does the work have academic merit? Yes Are sufficient details of methods and analysis provided to allow replication by others? Yes If applicable, is the statistical analysis and its interpretation appropriate? Yes Are all the source data underlying the results available to ensure full reproducibility? Yes Are the conclusions drawn adequately supported by the results? Yes Competing Interests: No competing interests were disclosed. Reviewer Expertise: Conservation genetics (incl felids), evolutionary genetics (incl felids), archival and ancient DNA (incl felids). I confirm that I have read this submission and believe that I have an appropriate level of expertise to confirm that it is of an acceptable scientific standard. Author Response 10 Jun 2021 Michael Westbury, The GLOBE Institute, University of Copenhagen, Øster Voldgade 5-7, Copenhagen, Denmark We would like to thank Daniel Förster for taking the time to review our manuscript and for providing useful feedback that we have now incorporated into our manuscript. We have included our responses to each comment below. Tables & Figures It seems unnecessary to avoid any Tables in the main text. For example, having the table of species analysed (Table S10) in the Methods section would benefit clarity, as these data suddenly find mention in the data processing section, as well as being featured in multiple subsequent sections. Answer: We have included the table listing the names of the additional species included in this study to the methods section. The mapping results (Table S1) are summarized in the Results & Discussion section, but the table could feature in the main text as well. Answer: We have now shifted the mapping results table to the main text. Plain language summary “assess” rather than “access”? Answer: Fixed Page 15 of 16 Open Research Europe Open Research Europe 2021, 1:25 Last updated: 30 JUN 2021 Introduction The first sentence reads rather clunky. Answer: We have now removed “a family within Carnivora” to make it smoother to read. Results & Discussion The start of this section would benefit from more details about (1) baits used for capture (based on data available for the lion; exome capture missing here?), and (2) mapping approach used (i.e. to which species: lion). Answer: We have now added more information about the baits and mapping approach to this section. Can their approach to detect gene flow be applied to some Felinae (e.g. Panthera) in their dataset to confirm that it does indeed function for more recent events? Answer: We have now tested the D3 method on some Felinae species. We find that our method does not detect gene flow between Panthera species. However, we do find significant signs of different levels of gene flow between Panthera (more specifically the snow leopard) and its sister genus, Neofelis (clouded leopard). Although the lack of detection between Panthera is unexpected based on the literature (e.g. Figueiro et al 2017 DOI: 10.1126/sciadv.1700299), this may be due to our strict filtering, artificially making the closely related species look more similar to each other. We note that the results between more distantly related Panthera species are only just not significant (p=0.066). We have included a table with the results here. Methods Exome capture used DNA baits? (typo?) Answer: Thank you for finding this. Yes it is a typo and has now been corrected to RNA baits. Why was the lion genome used as reference? Why not a genome with higher contiguity? Did this not potentially result in unnecessary loss of data during filtering (e.g. due to scaffold lengths)? Answer: We used the lion as reference genome as the baits were designed with lion sequences and to make results more comparable with our previous study (Barnett et al 2020) where we investigated the Homotherium genome within the context of extant Felinae species which were also all mapped to the same reference genome. Although there are currently other Felinae chromosome level assemblies available, they were not available at the time of the bait design. Despite this, even after filtering of shorter scaffolds, we were still left with 613 scaffolds >1Mb in length equating to 2,107,406,870bp which was more than sufficient for our purposes especially as we either relied on the complete dataset (trees) or 1Mb sliding windows (gene flow analyses). Competing Interests: No competing interests were disclosed. Page 16 of 16