dataesr/scientific-paragraphs-categorization huggingface dataset:
A set of UCSC genome browser tools  (Kuhn et al. , 2013) , including mafsInRegion, maf-Filter, maf2fasta, and several in-house PERL scripts were applied to divide the whole genome alignment into sub-regions. Then, for each sub-region, we included the genome sequences of the non-de novo sequenced rhinoceros samples, via extracting the corresponding sequences based on the region coordinate information on their reference genomes, and conducted multiple sequence alignment using MAFFT  (Katoh and Standley, 2013). We also removed the alignments with effective length < 1,000 bp and effective ratio < 0. 5, of which the effective sites represent nucleotide sites that do not include missing information (Ns) for any taxa in the multiple sequence alignment. After that, the GTR+CAT model of RaxML  (Stamatakis, 2014)  was used to build Maximum Likelihood (ML) phylogenetic trees for each window.
Following genome alignment, 82 genome regions of length > 20,000 bp were selected to infer phylogenetic trees (Figure  S1A ). Those regions were randomly scattered across the genome, thus having little chance of being found in the same recombination blocks. We calculated the frequency of three topologies around each focal internal branch of the species tree (Figure  2B ) using DiscoVista  (Sayyari et al. , 2018)  and NW_utils  (Junier and Zdobnov, 2010) , for the 100 kb windows based gene trees. A multi-species coalescent simulation was also applied to determine the expected gene tree distributions on the basis of the dated species tree using the Phybase package following  Wang al. (2018). Then, we inspected the congruency between the frequency of the three topologies inferred from empirical genomic data and that generated from the simulation for the Rhinocerotinae lineage.
The addition of these data was followed by filtering loci for biased rate estimates as described above, which led to a dataset of similar size (3,163 orthologs with 7,325,631 nucleotides). The species tree topology in these augmented analyses included the tapir as the sister to Rhinocerotidae, while the horse remained as the outgroup. We included a sixth calibration in addition to those described above, this time on the timing of the split between Rhinocerotidae and Tapiridae. The calibration had a hard minimum bound and soft maximum bound (0. 025 prior probability on older ages), ranging between 54 Ma and 64 Ma. The minimum age is based on the earliest appearance of Tapiroidea in the fossil record shortly after the Paleocene-Eocene transition, around 54 Ma, with Vastanolophus holbrooki and Cambaylophus vastanensis in Asia  (Rose et al.
, 2012)  and differences in lineage specific mutation rates  (Zheng and Janke, 2018). We independently mapped all raw reads from each rhinoceros species and the horse to three different reference genomes (greater one-horned rhinoceros, white rhinoceros and Sumatran rhinoceros) and then calculated D-statistics using the doAbbababa function module in ANGSD  (Korneliussen et al. , 2014)  specifying the horse as the outgroup ancestral allele and the following filtering parameters: 1) use only sites where the reads of the outgroup all have the same base; 2) minimum base quality and mapping quality of 20; 3) use only the transversion sites, all individuals have a minimum depth of 3, and maximum depth of 70; 4) a block size of 5 Mb to estimate standard errors using jackknife procedure. We further assessed whether less efficient mapping of the horse outgroup to the rhinoceros reference genomes may be driving some of our conflicting D-statistics results that differed based on the reference genome choice.
Heterozygosity estimation We estimated genome-wide heterozygosity for all the rhinoceros species based on a Site Frequency Spectrum (SFS), of which one diploid individual will generate three allelic states: homozygous ancestral allele state (AA), heterozygotic state (AB) and homozygous derived allele state (BB). Therefore, the whole genome level heterozygosity rate can be calculated as AB/(AA+AB+BB). We applied the doSaf function (-doSaf 1) module in ANGSD (version 0. 924)  (Korneliussen et al. , 2014)  to calculate genome-wide heterozygosity for each rhinoceros species with parameters of: unique mapping, a minimum mapping quality of 20, a minimum base quality of 20, a minimum and maximum depth value of 1/3 and 2 times of the average genome depths, respectively, for each species. For the ancient and historic samples (Javan, Merck's, woolly rhinoceros and Siberian unicorn), transition sites were removed (-noTrans 1) to eliminate the potential influences of DNA damage derived from cytosine deamination.
, 2012)  to create a gene function database. Then, we obtained the heterozygous sites for each rhinoceros species using BCFTOOLS  (Li, 2011)  requiring: minimum and maximum depth of 1/3 and 2x of the average depth; mapping and base quality R 20; variance quality R 30; read supporting number of minor allele account for R 20% of the total depth. After that, functional classes of the variances were estimated using SnpEFF  (Cingolani et al. , 2012)  with default settings. We trimmed 10 bp from the 5 0 and 3 0 ends of each read for the Merck's rhinoceros and Siberian unicorn to alleviate the impact of DNA damages, because both the samples showed abnormal transition/transversion ratios -as high as 5. 367 without end trimming (Figure  S4B ). However, as the transition/transversion ratio of the Merck's rhinoceros stayed abnormally high after trimming, which is consistent with previous reports of its DNA damage pattern (Figure  S4  in  Cappellini et al.
To investigate the putative genetic background of unique rhinoceros biological adaptations, we explored frameshift mutations shared across the rhinoceros lineages. We inspected the frameshift mutations in the alignment results generated by the EXONERATE analysis mentioned above. First, all the frameshift mutations with their location information (Gene IDs and Exon IDs) were extracted out for rhinoceros species with de novo assemblies being available (black, white, Sumatran and greater one-horned rhinoceroses). Then, we filtered out the frameshift mutations that share the same locations and exist in all four rhinoceros species, and recorded their gene information for further examination. Test of heterozygosity, PSMC and ROH estimation for the non-modern rhino samples Since four of the rhinoceros genomes were generated from historic and ancient samples, data could only be generated from them by short read re-sequencing.
Furthermore, given they represent historical or ancient samples, they are expected to contain damaged DNA (principally cytosine deamination (Pa ¨a ¨bo, 1989)) manifested as elevated C > T/G > A transitions at the 5 0 and 3 0 ends of sequence reads  (Briggs et al. , 2007; Brotherton et al. , 2007). Therefore, we tested the viability of restricting a number of the analyses performed to only the transversions, including heterozygosity estimation, historical demographical inference, and ROH estimation (Figures  S3  and S4 ). Also, the influences of the other factors such as reference selection and sequencing depth were evaluated. For the ROH estimation, we applied the ROH inference method to the four rhinoceros species with de novo genome assemblies, using the datasets both with and without transitions.
ADDITIONAL RESOURCES The results of some additional tests could be useful in understanding the sequence dataset for the historic/ancient rhinoceros samples and the comparative genome analyses, and can be found at: https://github. com/liushanlin/rhinoceros-comparative-genome/ blob/main/additional_resources. md   2  GOR: greater one-horned rhinoceros; SUR: Sumatran rhinoceros; WHR: white rhinoceros; SIU: Siberian unicorn; BLR: black rhinoceros; JAR: Javan rhinoceros; MER: Merck's rhinoceros; WOR: woolly rhinoceros. Each dot represents a test using the topology and species in the subtitles. For example, (((GOR, SUR), BLR), Horse) represents one of the topologies in Panel A-III and is shown as one of the dots with the x axis of BLR. For a preset topology (((H1,H2),H3),H4), a negative Dstatistics value means a closer relationship between H1 and H3 compared to that of H2 and H3, while a positive value means that H2 is closer to H3 than H1 is.
15% similar to that of the rhinoceroses. In the bottom panel we only include families in Artiodactyla, Carnivora, (legend continued on next page) Primates and Rodentia that have multiple species of low heterozygosity level (upper panel), and families that have > 2 species in our dataset (Table  S3  for all the details). (C): Factors influencing heterozygosity ratio estimation. Numbers represent the median values for each group. Influence of genome coverage (Depth of 303, 153 and 103 that indicated by different shades of blue) was tested by mapping reads from the black rhinoceros onto the reference genome of white rhinoceros. We noticed that different reference genomes contributed little to the variation of heterozygosity estimation, although when compared to the heterozygosity rate estimated using its own genome as a reference, using other closely related species tends to augment the estimation to a small degree and the inflation worsen when mapped to a higher divergent reference genome.
ACKNOWLEDGMENTS The authors acknowledge support from the  Science for Life Laboratory , the  Garvan Institute of Medical Research , the  Knut and Alice Wallenberg Foundation , and the  National Genomics Infrastructure  funded by the  Swedish Research Council and Uppsala Multidisciplinary Center for Advanced Computational Science  for assistance with massively parallel sequencing and access to the  UPPMAX computational infrastructure. We thank the  Natural History Museum at the University of Oslo  for providing the Javan rhinoceros sample. We thank the  Museum of the Institute of Plant and Animal Ecology (UB RAS, Ekaterinburg)  for providing the sample of Siberian unicorn. M. T. P. G. was supported by  European Research Council (ERC) Consolidator  grant  681396  ( Extinction Genomics ). E. D.
L. was supported by  Independent Research Fund Denmark  grant  8021-00218B. A. C. was supported by an  Australian Research Council   Laureate Fellowship  ( FL140100260 ). T. M. B. is supported by funding from the  ERC  under the European Union's  Horizon 2020 research and innovation program  (grant agreement  864203 ), grant  BFU2017-86471-P (MINECO/FEDER ,  UE), ''Unidad de Excelencia Marı ´a de Maeztu '' funded by the  AEI  ( CEX2018-000792-M ),  Howard Hughes International Early Career , and  Secretaria d'Universitats i Recerca and CERCA Programme del Departament d'Economia i Coneixement de la Generalitat de Catalunya  ( GRC 2017 SGR 880 ). L. D. was supported by the  Swedish Research Council  ( 2017-04647 ) and  Formas  ( 2018-01640 ). We thank  Dmitry Bogdanov  and  Roger Hall  for giving us permission to use their rhinoceros artwork.
REAGENT or RESOURCE SOURCE IDENTIFIER Biological samples Javan rhinoceros This study NHM, UiO 734; NCBI:txid102233  Merck 's rhinoceros This study F-4160; NCBI:txid2003782 Siberian Unicorn This study IPAE 915/2804; NCBI:txid2491732 Woolly rhinoceros This study ND036; NCBI:txid222863 Greater one-horned rhinoceros This study KB14498/SB137; NCBI:txid9809 Black rhinoceros Moodley et al. , 2020 GCA_013634535. 1;NCBI:txid9805 Sumatran rhinoceros Lord al. , 2020 GCA_014189135. 1;NCBI:txid89632 White rhinoceros  SAMN00778988 GCF_000283155. 1;NCBI:txid73337 Deposited data Raw sequence data and assemblies This study  NCBI  Project Number  PRJNA687817  Custom scripts for genome alignment analysis This study https://github. com/liushanlin/rhinoceros-.
Taxonomic distribution of the capybara gut microbiota diverges from typical hindgut fermenters The taxonomic structure of the capybara gut microbiota is mainly constituted by Bacteria, with only a small fraction of reads corresponding to Archaea (16S: 1. 07%, MG: 0. 40% and MT: 1. 55%) and Fungi (MG 0. 12% and MT: 0. 32%). 16S rRNA gene-based taxonomy analysis indicates that the most abundant bacteria found in capybara gut microbiota are members from the phyla Firmicutes (mean + sd: 35. 8 + 12. 4%), Bacteroidetes (31. 5 + 9. 8%), followed by Fusobacteria (15. 3 + 5. 4%) and Proteobacteria (8. 4 + 5. 4%) (Fig. 1a ). A similar taxonomic distribution was observed for 16S rRNA reads recovered from metagenome (16S_MG), MG and MT datasets (correlation coe cient r=0. 96 for MG:MT, r=0. 74 for 16S_MG:MG and r=0.
71 for 16S_MG:MT, P < 0. 05). Bacteroidetes, Fusobacteria and Proteobacteria phyla were signi cantly more abundant in MG than in MT considering cecal samples (Fig. 1b ), whereas Euryarchaeota (MT/MG ratio expressed as mean + SD: 2. 3 ± 0. 6), Fibrobacteres (1. 7 ± 0. 9) and Spirochaetes (1. 6 +/-0. 4) were more prevalent in MT than in MG data of cecal samples (Fig. 1b ). Although the protein-to-RNA ratios are rather dynamic in complex microbial communities according to absolute multi-omics studies  [15] , the higher MT/MG ratios for some phyla, such as Euryarchaeota known to include methanogens species, might correlate to the fact that capybara is a pronounced methane producer among the hindgut fermenters  [12, 16]. Binning of MG assembled contigs based on tetranucleotide frequency and coverage pro le resulted in the reconstruction of 79 unique Metagenome-Assembled Genomes (MAGs) (Suppl.
Fibrobacteres and Bacteroidetes are main degraders of dietary bers in capybara gut In order to understand the ability of capybara to convert plant polysaccharides into free sugars, MG and MT data were investigated to determine the genomic potential of its gut microbiota associated with Carbohydrate-Active enZymes (CAZymes). A total of 6,132 putative CAZymes genes encoding for 105 Glycoside Hydrolases (GH), 11 Carbohydrate Esterases (CE), and 10 Polysaccharide Lyases (PL) families were identi ed, of which 456 genes presented a modular architecture (Suppl. Table  S2 ). The most abundant CAZymes identi ed are members of the families GH3, GH2 and GH1 (by decreasing abundance) (Fig. 2 ), which is in agreement with that reported for other host-associated gut microbiomes such as human, swine and cattle rumen  [25].
These enzymes encompass diversi ed activities including β-glucosidase, β-xylosidase, β-galactosidase and β-mannosidase, and are often associated with the nal steps in the depolymerization cascade of several plant polysaccharides such as cellulose, heteroxylans, mixed-linkage β-glucans and β-mannans. In the CAZyme repertoire of this microbiota neither cellulases from families GH6, GH7 and GH48, nor cellulosomes, assessed by the presence of cohesin and dockerin domains associated with cellulases, could be identi ed in MG or MT datasets. In ruminal anaerobic fungi, these families are found in high abundance, possibly targeting recalcitrant cellulose structures  [26]. However, in the capybara gut microbiota assayed here, fungi were detected only at very low abundance (Fig. 1a ). This suggests that cellulose degradation in the capybara gut might be mainly accomplished by endo-β-1,4-glucanases (EC 3.
bioinformatics. babraham. ac. uk/projects/fastqc/) and ltered using Trimmomatic v. 0. 36  [45]  to remove adapters and low-quality reads, using the following parameters: ILLUMINACLIP:NexteraPE-PE. fa:2:30:10:8:true, LEADING:4, TRAILING:4, MINLEN:60 and SLIDINGWINDOW:4:16. Filtered paired-end reads were merged using fastq_mergepairs function from Usearch v. 10 package [46] (parameters: fastq_maxee 0. 5, fastq_minovlen 50, and fastq_minmergelen 250). Detailed number of reads per sample are summarized in Suppl. Table  S13. Prior to OTUs clustering, primer sequences were removed, and singletons were discarded. Filtered amplicon reads were denoised (error-corrected) using the UPARSE unoise3 function (parameters: -minsize 8 and alpha 0. 2), to likely recover true biological sequences. Prokaryotic taxonomy assignment was performed using sintax function, as implemented in Usearsh v.
9 and markerset 107)  [51]  to recover putative genomes from the metagenomic data. The binned genomes were dereplicated to remove redundancies using dRep v. 2. 0. 5 (parameters: -comp 80 -con 10 -str 100 and -p 10) and analyzed using CheckM v1. 0. 6 with the lineage_wf work ow  [52]  to determine the completeness and contamination ratios of these genomes. Long-reads sequencing were used for MAGs scaffolding using SSPACE-long-reads v1. 1 (parameters: -k 5, -a 0. 7, -x 1, -m 50, -o 20 and -n 1000). Genomes with completeness smaller than 55% and more than 15% contamination rate were discarded. To assign taxonomy to the recovered genomes GTDB-tk tool v. 1. 4 was used with the release 95 of the GTDB database  [53]. Gene prediction and annotation of both the recovered genomes and the co-assembly were performed using Prokka v.
1. 11 with the meta parameter  [54] , and annotation statistics are described in Suppl. Table  S14. KEGG pathways and Kegg Orthologous (KOs) annotation were performed using KOFAM (e-value < 1e-5)  [55]  and Functional Ontology Assignments for Metagenomes (FOAM) database (e-value < 1e-5) [56]. CAZymes annotation was performed according to CAZy database pipelines  [57]. Furthermore, MG and MT reads were mapped to the whole set of genes recovered from the co-assembled metagenome and the set of genes recovered from the MAGs using Kallisto v. 0. 46. 1 with quant function  [58]  to estimate the coverage/abundance of protein coding genes in cecal and rectal samples. Normalized abundance was estimated based on the count/number of reads per kilobase per million mapped reads expressed as TPM.
Phylogenetic analysis and metabolic reconstruction Phylogenetic analysis of the MAG57, reference Bacteroidetes type strains and Prevotellaceae uncultured genomes recovered from UBA project [38] was performed using concatenated 92 single copy core genes according to UBCG method  [59]. CAZymes phylogenetic analysis was carried out using the catalytic domain of each family aligned with MAFFT [60], and using maximum likelihood methods employing the RAxML software [61] with 1,000 rapid bootstrap inferences and LG as the substitution model. To perform the metabolic pathways reconstruction of each recovered MAG, their annotation obtained from KOFAM database were ltered to keep only the top 5 hits of each protein with e-value below the 1e-5 threshold. These ltered annotations were then supplied to the Annotation of Metabolite Origins (AMON) tool  [62] , which based on the KOs annotated in each MAGs predicts the putative metabolites that it can generate.
The nal volume of ltrate ranged from 500 to 650 μL. 1H NMR spectra of samples were acquired using a Varian Inova NMR spectrometer (Agilent Technologies Inc. ) equipped with a 5 mm triple resonance cold probe and operating at a 1H resonance frequency of 599. 84 MHz and constant temperature of 298 K (25 °C). A total of 1024 free induction decays were collected with 32-k data points over a spectral width of 16 ppm. A 1. 5-s relaxation delay was incorporated between scans, during which a continual water presaturation radio frequency (RF) eld was applied. Spectral phase and baseline corrections, as well as the identi cation and quanti cation of metabolites present in samples, were performed using Chenomx NMR Suite 7. 6 software (Chenomx Inc. ).
Kinetic data were calculated from initial velocities and expressed as mean + SD from three independent experiments (n = 3). Binding capacity of the carbohydrate-binding module from theGH10 enzyme and its mutants were evaluated by a nity gel electrophoresis (AGE) according to  [65]. Small Angle X-ray Scattering (SAXS) SAXS data of the CapGH10 β-helix domain were collected at the SAXS1 beamline (Brazilian Synchrotron Light Laboratory, Campinas, Brazil) at protein concentration of 8. 4 mg. mL -1 in 20 mM Hepes buffer pH 7. 5. Buffer scattering were recorded and subtracted from the protein scattering. SAXS patterns were integrated using Fit2D [66] and GNOM [67] was used to evaluate the pair-distance distribution functions p(r). Ab initio molecular envelope was calculated from SAXS data with DAMMIN [68] and an averaged model was generated from several runs using DAMAVER  [69].
The crystallographic coordinates were tted into the SAXS low-resolution model using SUPCOMB  [70]. Crystallization, diffraction data collection and structure determination experiments of the CapGH10 β-helix domain were carried by the vapordiffusion method at 18 °C. Native crystals were grown in 20% (w. v -1 ) PEG 6000, 0. 1 M sodium acetate (pH 5. 0) and 0. 2 M sodium chloride. SeMet crystals were obtained under the same condition with the addition of 20 mM betaine hydrochloride. For cryoprotection, crystals were soaked in the reservoir solution added with glycerol or PEG 400 (20% (w. v -1 )) prior to ash cooling. Diffraction datasets of native and SeMet crystals were collected at the PROXIMA-2A and MX2 beamlines from SOLEIL (Gif-sur-Yvette Cedex, France) and LNLS (Brazilian Synchrotron Light Laboratory, Campinas, Brazil), respectively.
Datasets were indexed, integrated, merged, and scaled using XDS package  [71]. The structure was solved by single anomalous dispersion (SAD) using the programs SHELXC, SHELXD and SHELXE  [72]  for data preparation, anomalous scatters location and phase calculation, respectively. An initial model was built with the AutoBuild Wizard  [73]  from the Phenix package  [74]. The structure was re ned with the programs PHENIX. REFINE  [75]  and REFMAC5 [76], and the models were inspected and manually adjusted according to the computed σ A -weighted (2F o -F c ) and (F o -F c ) electron density maps using COOT  [77]. TLS groups were calculated by TLSMD [78] and applied during the crystallographic re nement. The nal structure was evaluated by MolProbity  [79]  and PDBRedo server  [80]. Structure factors and atomic coordinates were deposited at the Protein Data Bank (PDB) under the accession codes 7JVI.
Acknowledgements acknowledge the  Brazilian Synchrotron Light Laboratory  for the provision of time on the MX2 and SAXS1 beamlines, the  Brazilian Biosciences National Laboratory (LNBio)  for the use of the crystallization (Robolab),  NMR  and spectroscopy facilities, and the  Brazilian Biorenewables National Laboratory (LNBR)  for the use of the characterization of macromolecules and next generation sequencing facilities. LNLS, LNBio and LNBR are operated by the  Brazilian Center for Research in Energy and Materials for the Brazilian Ministry for Science, Technology, Innovations and Communications. We acknowledge SOLEIL for provision of synchrotron radiation facilities at PROXIMA-2A (proposal 20181915) and we would like to thank  William Shepard  and  Martin Savko  for assistance in using the beamline.
Availability of data and material data for this study can be found under the BioProject ID PRJNA563062. The 16S, metagenomic and metatranscriptome reads for cecal and rectal samples are available at SRA under the accession numbers SRR11852069-SRR11852086, SRR11852046-SRR11852057 and SRR11852097-SRR11852108, respectively (Suppl. Table  S15 ). The MAGs can be found at GenBank under the accession numbers JABUSA000000000 -JABUVA000000000. Structural data have been deposited in the Protein Data Bank (https://www. rcsb. org/) under accession code 7JVI (CapCBMXX). All other data generated or analyzed during this study are included in this published article (and its Supplementary information les). Data and code used for microbiome analyses are publicly available at https://github. com/gpersinoti/capybara_microbiome.