dataesr/scientific-paragraphs-categorization huggingface dataset:
Results: We have developed RAPPAS (Rapid Alignment-free Phylogenetic Placement via Ancestral Sequences) which uses an alignment-free approach, removing the hurdle of query sequence alignment as a preliminary step to PP. Our approach relies on the precomputation of a database of k-mers that may be present with non-negligible probability in relatives of the reference sequences. The placement is performed by inspecting the stored phylogenetic origins of the k-mers in the query, and their probabilities. The database can be reused for the analysis of several different metagenomes. Experiments show that the first implementation of RAPPAS is already faster than competing likelihood-based PP algorithms, while keeping similar accuracy for short reads. RAPPAS scales PP for the era of routine metagenomic diagnostics.
Contrary to genomic approaches studying a single genome, metagenomes are generally built from bulk DNA extractions representing complex biological communities that include known and yet unknown species (e. g. protists, bacteria, viruses). Because of the unexplored and diverse origin of the DNA recovered, and because of the sheer volume of sequence data, the bioinformatic analysis of metagenomes is very challenging. Metagenomic sequences are often anonymous, and identifying a posteriori the species from which they originated remains an important bottleneck for analyses based on environmental DNA. Several classes of approaches are possible for taxonomic classification of metagenomic sequences, the choice of which will generally depend on the availability of prior knowledge (e. g. a reference database of known sequences and their relationships) and a trade-off between speed and precision.
For instance, without prior knowledge, the exploration of relatively unknown clades can rely on unsupervised clustering  (Mahé et al. , 2014; Ondov et al. , 2016; Sedlar et al. , 2017). The precision of this approach is limited, and clusters are generally assigned to taxonomic levels via posterior analyses involving local alignment searches against reference databases. An alternative is to contextualise local alignments with taxonomic metadata. For instance, MEGAN  (Huson et al. , 2007 (Huson et al. , , 2016) )  assigns a query metagenomic sequence (QS for simplicity) to the last common ancestor (LCA) of all sequences to which it could be aligned with a statistically-significant similarity score (e. g. using BLAST or other similarity search tools). A faster alternative is to build a database of taxonomically-informed k-mers which are extracted from a set of reference sequences.
, 2016). In the case of viruses (under studied and characterized by rapid evolution), advanced evolutionary modelling is necessary to perform detailed classifications. A possible answer in this direction may be the use of phylogenetic inference methods. However, the standard tools for phylogenetic inference are not applicable to metagenomic datasets, simply because the sheer number of sequences cannot be handled by current implementations  (Izquierdo-Carrasco et al. , 2011). A solution to this bottleneck is phylogenetic placement (PP)  (Matsen et al. , 2010; Berger et al. , 2011). The central idea is to ignore the evolutionary relationships between the QSs and to focus only on elucidating the relationships between a QS and the known reference sequences. A number of reference genomes (or genetic markers) are aligned and a "reference" phylogenetic tree is inferred from the resulting alignment.
But no fundamental algorithmic changes were introduced, and the method still relies on a preliminary alignment. To avoid these bottlenecks, we have developed the phylogenetic placement program RAPPAS (Rapid Alignment-free Phylogenetic Placement via Ancestral Sequences). The novelty of our approach rests on two key ideas: (i) a database of phylogenetically-informed k-mers ("phylo-kmers") is built from the reference tree and reference alignment, (ii) QSs are then placed on the tree by matching their k-mers against the database of phylo-kmers. For a specific group of organisms, the maintenance of the RAPPAS database can be seen as a periodic "build-andupdate" cycle, while the phylogenetic placement of new metagenomes becomes a scalable, alignment-free process compatible with day-to-day routine diagnostics.
Analysing metagenomes with RAPPAS In order to perform phylogenetic placement, one must define a reference dataset (alignment and tree) onto which future QSs will be aligned and placed. Such references are generally based on well-established genetic markers (e. g. 16S rRNA, rbcl, etc. ) or short genomes (e. g. organelles, viruses, etc. ). One of the main differences between our approach and the existing PP methods is in the workflows that they entail. Existing PP algorithms require the alignment of every QS to the reference alignment prior to the placement itself, a complex step that relies on techniques such as profile Hidden Markov Models that in some cases may dominate computation time. The HMMER package  (Eddy, 2011) , and PaPaRa  (Berger and Stamatakis, 2011)  are popular tools for this step.
As we show in sec. 3. 3, the alignment step is particularly onerous in the presence of long reference sequences, in which case this step may actually take longer than placement (Fig. 5A , cf. hmmalign vs. EPA-ng for D155). For each sample analyzed, the alignment phase is then followed by the placement phase, which relies on phylogenetic inference (Figure  1A ). RAPPAS uses a different approach where all computationallyheavy phylogenetic analyses are performed as a preprocessing step, which is only executed once for a given reference dataset, before the analysis of any metagenomic sample. This phase builds a data structure (named pkDB; see next section) which will be the only information needed to place QSs. Many samples can then be routinely placed on the tree without the need for any alignment.
Since the same pkDB can be reused for several query datasets (Figure  1B ), the use of RAPPAS is particularly adept when QSs are placed on standardized, long-lifespan reference trees. On the other hand, when references must be very frequently updated, the overhead for rebuilding the pkDB reduces the advantage of our approach. To follow the standards defined by previous PP methods, the placements are reported in a JSON file following the "jplace" specification  (Matsen et al. , 2012)  (see Suppl. File S1, sec. 3, for more details). The results of RAPPAS can then be exploited in all jplace-compatible software, such as iTOL (visualization), Guppy (diversity analyses) or BoSSA (visualisation and diversity analysis)  (Letunic and Bork, 2016; Matsen and Evans, 2013; Lefeuvre, 2018). The pipelines from query sequence datasets to placement results (jplace files) are depicted.
A. Likelihoodbased software requires, for each new query dataset, the step of aligning the query sequences to the reference alignment via an external tool (red box). The resulting extended alignment is the input for the phylogenetic placement itself (blue box). B. RAPPAS builds a database of k-mers (the pkDB) once for a given reference tree and alignment. Many query datasets can then be placed without alignment, matching their k-mer content to this database. Each operation is run with a separate call to RAPPAS (blue boxes).
Construction of the phylo-kmers DB (pkDB) Before performing placements, a database of phylogenetically informative k-mers must be built; these "phylo-kmers" are those that we would expect to see in a QS with non-negligible probability. The input for this phase is a reference alignment (denoted "refA") and a reference tree (denoted "refT") reconstructed from refA. The refT is interpreted as a rooted tree (as specified by the Newick file) and has lengths associated to its edges. This process is divided in three steps (see Figure  2 ): A) the injection of "ghost nodes" on every edge of refT, B) a process of ancestral state reconstruction based on this extended tree, and C) the generation of the phylo-kmers and their storage in a data structure. The first step enables us to simulate evolution that has occurred and resulted in sequences not present in the reference alignment.
,  +  -1} is given by: 𝑃𝑃 = ∏ 𝑝𝑝 𝜎𝜎 𝑖𝑖 ,𝑛𝑛+𝑖𝑖-1 𝑘𝑘 𝑖𝑖=1 (1) Each k-mer with a probability  larger than a (parameterized) threshold (by default 1/  ) is stored together with (i) the edge corresponding to the ghost node of origin, and (ii) the probability  associated to this k-mer at the ghost node of origin (Figure  2, step C ). In the end, all the k-mers reconstructed above are stored in a data structure that allows fast retrieval of all relevant information associated to a k-mer, namely the list of edges for which it has high-enough probability, and the probabilities associated to each of these edges. This structure will be referred to as the phylo-kmers database (pkDB). We note that when the same k-mer is generated at different positions  of refA, or at different ghost nodes for the same edge, only the highest probability is retained in the pkDB.
We consider this approach realistic in the light of current metagenomic challenges, because it is indeed common that large sections of the tree of life are absent from the reference phylogenies. The procedure was repeated for read lengths of 150, 300, 600, 1200 bp. As in previous studies, the simulated reads are contiguous substrings of the pruned sequences obtained without introducing sequencing errors  (Matsen et al. , 2010; Berger et al. , 2011). A detailed description of the procedure can be found in Suppl. File S1, sec. 4. The impact of sequencing errors was evaluated independently. We repeated the procedure above on a single dataset (D652), and introduced errors in the simulated reads with the Mason read simulator  (Holtgrewe, 2010)  using the Illumina error model. Different rates (0.
Placement consistency on real-world data As a complement to the experiments above, we also evaluated the placements of real-world reads on relatively large reference trees by reproducing the protocol described in  Barbera et al. , 2018. Because no proxy for the correct placement can be defined on these datasets, we evaluate the consistency between placements produced by RAPPAS and by the other methods. For each dataset (see Sec. 2. 7 below) between 10,000 and 15,060 reads were placed on the reference tree using each method (EPA, EPA-ng, PPlacer and RAPPAS). For each pair of methods, the resulting placements were compared using the Phylogenetic Kantorovich-Rubinstein (PKR) metric  (Evans and Matsen, 2012). Naturally, the lower the distance between two methods, the more similar the placements they produce are.
Runtime and memory experiments CPU times were measured for each software with the Unix command time, using a single CPU core and on the same desktop PC (64GB of RAM, Intel i7 3. 6GHz, SSD storage). The fast aligner program hmmalign from HMMER  (Eddy, 2011)  was used to align the QSs to the refA, as in previous placement studies  (Matsen et al. , 2010). Placement runtimes were measured for two of the reference datasets that we describe in the next section, chosen so as to represent different combinations of the two parameters with the largest impact on runtimes (number of taxa and alignment length). Prior to starting the time measurement, for each dataset, GTR+Γ model parameters were estimated using RAxML. Then, dedicated options were used in EPA-ng and PPlacer to load these parameters and avoid optimization overheads (Suppl File S1, sec.
9). Each software was launched independently 3 times and the median CPU time was retained. When available, masking options were kept as defaults. For EPA and EPA-ng, a single thread was explicitly assigned (option -T 1). PPlacer only uses a single thread and required no intervention. Default values were used for all other parameters. RAPPAS is based on Java 8+, which uses a minimum of 2 threads. We ensured fair comparison by using the Unix controller cpulimit, which aligns CPU consumption to the equivalent of a single thread (options --pid process_id -- limit 100). For each RAPPAS process, memory was monitored by programmatically reading its usage at regular intervals.
Software and datasets used in experiments The performance of RAPPAS (v1. 00), both in terms of accuracy and computational efficiency, was compared against that of EPA (from RAxML v8. 2. 11)  (Berger et al. , 2011) , PPlacer (v1. 1alpha19)  (Matsen et al. , 2010)  and  for the experiments on real-world data and v0. 2. 0-beta for all other experiments)  (Barbera et al. , 2018). Reference datasets (refA, refT) used for evaluating accuracy were taken from  (Matsen et al. , 2010; Berger et al. , 2011)  with the exception of D155, which was built for this study and corresponds to a whole genome alignment and tree for the hepatitis C virus. Each dataset is identified by the number of taxa present in the tree (e. g. D150 has 150 taxa). In total, 10 reference datasets were tested. Descriptive metadata reporting the nature of the locus, the length of the alignment and the gap content are displayed in Figure  3  (left headers) and Suppl.
File S2. One of the reference datasets (D140) is based on a reference alignment of protein sequences. Globally, the reference datasets represent a panel of loci commonly used in metabarcoding and metagenomic studies. PKR distances between alternative placement methods were measured using five real world datasets: the three datasets available in  Barbera et al. , 2018 (named bv, neotrop, tara)  and two large datasets of 5088 and 13903 taxa (named green85 and LTP). The latter two datasets were retrieved from the Greengenes database  (DeSantis et al. , 2006 ; release gg_13_8_otus, 85% OTUs tree) and the "All-Species Living Tree Project" (Yilmaz al. , 2014; ssuRNA tree, release 132, june 2018), respectively. On the reference trees for green85 and LTP, we placed the QSs provided with the bv dataset, as they correspond to the same genetic region.
The runtime experiments are based on two of the above reference datasets, D652 (1. 7 kbp) and D155 (9. 5 kbp), and the QSs for placement were obtained in the following way. For D652, real-world 16S rRNA bacterial amplicons of about 150 bp were retrieved from the Earth Microbiome Project  (Thompson et al. , 2017)  via the European Nucleotide Archive  (Silvester et al. , 2018)  using custom scripts from github. com/biocore/emp/. Only sequences longer than 100 bp were retained. A collection of 3. 1×10 8 QSs were thus obtained. For D155, 10 7 Illumina short reads of 150 bp were simulated from the 155 viral genomes in this reference dataset using the Mason read simulator (Holtgrewe, 2010). A substitution rate of 0. 1% and very low indel rate of 0. 0001% were chosen to avoid impacting negatively the runtime of the alignment phase for EPA, EPA-ng and PPlacer.
Accuracy The simulation procedure (see Methods) was repeated for 10 different reference datasets representing genes commonly used for taxonomic identification. Accuracy is measured on the basis of the node distance between observed placement and expected placement (i. e. the number of nodes separating the two edges). The closer to zero the node distances are, the more accurate the method is. Figure  3  summarizes the results of the experiments for each choice of dataset, read length, and software. Each point in the plot shows the mean node distance for the reads generated by one pruned subtree of the reference tree. Each violin plot is based on as many points as there are pruned subtrees (i. e. 100). Horizontal lines link the per software node distance means, calculated as the mean of the 100 means above.
Detailed mean values are reported in Suppl. File S3. As expected, given that they are based on the same general likelihood-based approach, EPA, EPA-ng, and PPlacer show very similar performance for all datasets and read lengths (Figure  3 , blue lines almost perfectly horizontal). As expected, the longer the reads are, the more accurate the placement is (Figure  3 , left to right boxes). The worst performances were observed for the D628 and D855 datasets, which agrees with previous studies  (Matsen et al. , 2010; Berger et al. , 2011). RAPPAS was subjected to the same tests using four k-mer sizes. The corresponding results can be classified into two groups. The first group contains the datasets for which using different values of k has a limited impact (Figure  3 , all except the last three).
For these, RAPPAS shows globally an accuracy comparable to likelihood-based software. In simulations based on the shortest read lengths (Figure  3 , leftmost boxes) similar or even slightly better accuracies (leftmost boxes of D1604, D652) are measured. On the other hand, the longest read lengths generally resulted in lower accuracy compared to other methods (rightmost boxes of D140, D150, D500, D855, D1604). The second group contains datasets where different values of k resulted in variable accuracy (Figure  3 , D218, D628, D155). Two out of these three datasets are also characterized by a high proportion of gaps in their reference alignments (Figure  3 , left headers), which may be an adverse factor for the current version of RAPPAS (discussed below). Finally, D155, which corresponds to the longest reference alignment (9.
5 kbp), is the only dataset where shorter k-mer sizes clearly lowered the accuracy of RAPPAS. Each square is a software tested under 16 conditions corresponding to varying rates of errors introduced into 300 bp reads (substitutions on the horizontal axis, indels on the vertical axis). Values represent the accuracy loss, calculated as the difference (Δ ND) between the mean node distance with indel/substitutions, and the mean node distance in the absence of sequencing errors (the bottom left condition with 0% indels, 0% substitutions). Figure  4  shows the accuracy loss measured on the D652 dataset after introducing different rates of sequencing errors in the simulation procedure (300 bp reads, k=8 for RAPPAS). For all tested software, a very limited accuracy loss (<=0. 2 node distance units) was measured for 0.
1% to 1% substitutions/indel rates. Such error rates are typical for short read technologies  (Glenn et al. , 2011) , suggesting that RAPPAS and all tested software are equally suitable for the placement of sequences produced by current technologies. When higher error rates are considered (5% substitutions/indels), RAPPAS appears to lose more accuracy than competing tools (Figure  4 , top rows and rightmost columns). Interestingly, indels appear to induce a higher accuracy loss than substitutions for every placement method tested. Different read lengths and different values of k for RAPPAS were also tested but showed similar trends (Suppl. File S1, sec. 5). Table  1. Median PKR distances over the placements of real-world reads on five different datasets. The lower the value, the closer the placements on the reference tree are between the two methods.
Placement of real-world reads The PKR distances measured for real-world reads placed on five medium-to-large reference trees are reported in Table  1. The question here is whether RAPPAS produces placements that can be viewed as strong outliers to those of the alignment-based methods. Note that some degree of difference is to be expected, given that RAPPAS has a radically different approach than the other methods. Another prediction is that PPlacer and EPA-ng should produce similar placements, as their default behaviour includes similar heuristics aiming to reduce runtimes. Consistent with the latter prediction, we observe that the PKR distances between PPlacer and EPA-ng are the smallest across all tested datasets. As for the former prediction (RAPPAS vs. the other methods), for the tara and neotrop datasets, PKR distances between RAPPAS and the alignment-based methods are about two to three times higher than the average distance between the alignment-based methods.
dataset #taxa pkDB PPlacer EPA-ng RAPPAS  In each box, the horizontal axis represents the tested software (RAPPAS was tested for several k-mer sizes), while the vertical axis corresponds to the measured node distance (ND) (the closer to 0, the more accurate the placement). For each software, each black dot represents the mean ND for one of the 100 independent pruning experiment (see Methods) and white violins show the approximate distribution of the sample containing these 100 measures. The means of the four samples for RAPPAS are linked by a red line. The means of the three samples for the other methods are linked by a blue line.
Runtime Figure  5  summarises runtimes measured for the placement of 10 3 to 10 7 QSs on two reference datasets corresponding to short and long reference sequences (a marker gene for bacteria, and a complete virus genome, respectively). Figure  5A  shows the times measured independently for each software. To consider the actual time needed in a scenario where a long-lifespan reference tree is used, in Figure  5B , the runtime reported for EPA, EPA-ng, and PPlacer includes both read alignment and placement, while the runtime of RAPPAS does not include pkDB construction (as the pkDB is already available). As expected, EPA, PPlacer, and EPA-ng show a runtime scaling linearly with the number of placed reads for both datasets (Figure  5B ). With the D155 dataset, EPA and PPlacer could not be run with a set of 10 7 reads due to high memory requirements.
Implementation RAPPAS is distributed as a standalone command-line application. It is compatible with any system supporting Java 1. 8+. Ancestral sequence reconstruction is performed via an external program. Currently PhyML 3. 3+  (Guindon et al. , 2010)  and PAML 4. 9  (Yang, 2007)  binaries are supported by RAPPAS. The main limiting factor for RAPPAS is the amount of memory required to build the pkDB (and subsequently, to load the pkDB for the placement phase). Figure  6  reports the memory and drive storage requirements related to building and storing the pkDB for three representative datasets and different values of k. D150 (few taxa, short alignment), D155 (few taxa, long alignment) and D1604 (many taxa, short alignment) show memory usage below 4GB for k = 10, and 12GB for k = 12, which indicates that RAPPAS can be launched on most standard laptop computers, if k is set to a low value.
Advantages of alignment-free placements The novel algorithm introduced in RAPPAS removes a major limitation imposed by previous PP software: these methods rely on the computation of a multiple alignment merging the reference alignment and query sequences, a step which must be repeated for each new set of query sequences (Figure  2-A ). The idea of attaching k-mers to pre-computed metadata to accelerate sequence classification has been explored by algorithms combining taxonomically-informed k-mers and LCA-based classification  (Liu et al. , 2018; Müller et al. , 2017; Wood and Salzberg, 2014; Ounit et al. , 2015). However, none of these methods have the goal of producing evolutionarily precise classifications based on probabilistic models of sequence evolution. RAPPAS adapts this idea to PP by reconstructing sets of probable, phylogenetically-informed k-mers, called phylo-kmers.
Compared to previous PP methods, RAPPAS has two critical algorithmic advantages (see Suppl. File S1, sec. 2 for details). First, its runtime does not depend on the length of the reference alignment, since the sequences in the reference tree, and those that are related to them, are summarized by their phylo-kmers. Note that the alignment phase required by the other PP methods has a runtime that scales linearly with the length of refA (if profile alignment is used  (Eddy, 1998) ). Second, Fig. 5. Placement runtime comparison. Runtimes were measured for a long (D155) and a short (D652) reference alignment, for a number of placement pipelines. A. The execution time for each pipeline is decomposed into "preparation" and "placement" steps. The preparation step for alignment-based methods is the alignment of each read to the reference (denoted "align"), while RAPPAS prepares a database once (noted "DB").