dataesr/scientific-paragraphs-categorization huggingface dataset:
To support confidence on this identified group of sequences, we scanned them to see if we could find motifs that are known to characterise the ABC transporter domain (https://www. ebi. ac. uk/interpro/entry/IPR000515). For this, we considered all 1109 environmental sequences annotated as ABC transporter type 1 by MetaCLADE (Additional file 1: Table  S6 ) and selected the ones with an E value smaller than 1e-4. There are 945 sequences with an average length of 36. 9 residues. We ran MEME  [68]  on them to find the 10 most significant motifs. Among these motifs, we identified the EAA motif, a 20 amino acid-conserved sequence known to occur in ABC transporter type 1 (https:// www. ebi. ac. uk/interpro/entry/IPR000515) (See motif 9 in Additional file 1: Table  S7 ; strictly speaking, we found a portion of the known EAA motif, where the submotif EAA----G occurs  [69].
). The consensus sequences of the 10 motifs (for example, FNLLGDGLRDALDPR for motif 1 and GAILTEAALSFLGLG for motif 9 in Additional file 1: Table  S7 ) were used as query to search the NR database. For all consensus sequences, most of the hits found (> 95%) matched ABC transporters. In rare cases, BLAST  [70]  retrieved, in addition to the great majority of ABC transporters, other transport systems permeases or hypothetical proteins. The presence of known motifs favourably supports the finding, and MetaCLADE proves to be able to extract useful functional information even from very short reads. In this respect, Additional file 1: Table  S6  shows that MetaCLADE annotations can substantially change the estimations of domain abundance in MG samples compared to estimations realised with HMMer.
Sensitivity of MetaCLADE on the distribution of species generating the models To analyse MetaCLADE's sensitivity on the distribution of species from which models have been generated compared to species where reads come from, we verified the distribution of species generating models used for the annotation of nine simulated datasets of reads. These datasets contain short fragments coming from species belonging to bacteria, viruses, archaea and eukaryotes. Specifically, they are constructed by gradually incrementing (by 10%) the number of eukaryotic sequences in them (see the "Methods" section). The nine resulting datasets have been annotated with MetaCLADE, and the origin of the CCMs used is reported in Fig. 5. As expected, in the annotation process, MetaCLADE tends to use models close to the communities represented in the dataset; namely, the number of models generated from eukaryotic sequences used for annotation is proportional to the quantity of eukaryotic reads in the datasets.
Comparison with UProC UProC is a very fast protein classification tool designed for large-scale sequence analysis. It is much faster than profile-based methods, like MetaCLADE, and in MG datasets was demonstrated to achieve high sensitivity  [38]. We tested whether MetaCLADE, on real datasets, is identifying more domains than UProC or not. The answer depends on the average length of the reads in the dataset. We took the five environmental samples considered in the section "Improved annotation of MG/MT datasets compared to InterProScan" (see Additional file 1: Table  S4 ), ran UProC on them and compared to Meta-CLADE results. For the Rainforest MG dataset, the analysis is reported in Fig. 6a. As illustrated in Fig. 6b , the 5439 reads that are uniquely annotated by UProC are either of very small size, < 50 aa, or much larger, > 150 aa.
In contrast, the 22,059 MetaCLADE exclusive annotations do Fig. 5  Distributions of species whose sequences generated models for MetaCLADE annotation. Analysis of nine simulated datasets named "Eukx" containing x% of reads coming from eukaryotic sequences. For increasing values of x, one observes a proportionally higher number of CCMs coming from eukaryotes (red) that have been used for the annotation of the dataset. The proportion of bacteria (violet), archaea (yellow) and viruses (green) is reported for each dataset not concern very small reads but rather reads with larger size > 50 aa, and particularly > 150 aa. The E value density distribution curve of UProC annotations (see Fig. 6c  and "Methods" section) highlights reasonably low E values showing a high confidence in most UProC domain annotations.
The second best curve is MetaCLADE's curve, placed on its left, followed by the InterProScan curves. A quite large number of reads exclusively predicted by UProC is also predicted by MetaCLADE but not selected by MetaCLADE because of its probability threshold set at 0. 9. More precisely, they cover less than 50% of the UProC exclusively predicted reads. By looking at the confidence of these domain predictions, a large number of these predicted reads have very low probability (see Fig. 6d ). This observation suggests that, for these domains, the CLADE's library does not cover properly the spread of evolutionary variability of the domains. In Additional file 1: Figures  S12-S15 , we report the analysis of the other four MG and MT datasets. All these datasets confirm that most uniquely annotated UProC reads are very short (< 50 aa) compared to reads uniquely annotated by MetaCLADE (see Additional file 1: Figures  S12B-S15B ).
The density curves of the UProC's associated E values occupy the right-hand side of the plots (see Additional file 1: Figures  S12C-S15C ) with a peak that indicates an optimal average E value for UProC reads. Moreover, from the distribution of domain lengths for uniquely annotated reads in Additional file 1: Figure  S16 , we observe that, depending on the dataset, for a given domain, the read multiplicity can be much larger for MetaCLADE than for UProC. The reason is anchored on the way UProC and MetaCLADE handle sequence information. Indeed, if the space of sequences in the FULL Pfam dataset associated to a domain spans well over diversified homologous sequences, then UProC, which exploits all FULL Pfam sequences, can produce highly confident predictions. If the number of Meta-CLADE models is small compared to the number of diversified FULL Pfam sequences, it might represent well only a part of the diversification and fail in predicting homology on the unrepresented pool of sequences.
A first comparison with HMM-GRASPx using an assembly-based approach Assembly algorithms can be used prior to functional annotation of a given MG dataset. We evaluated MetaCLADE's performance against HMM-GRASPx [37], a state-of-the-art assembly-based annotation method. HMM-GRASPx is characterised by a profile-guided assembly phase in which assembled protein contigs are verified through a HMMER realignment of Pfam profiles. Hence, a domain family is assigned to the reads by mapping them back to verified contigs. In order to provide a fair comparison of the two methods, we decided to annotate the assembly generated by HMM-GRASPx with MetaCLADE and transfer this annotation on the reads by mapping them back to the assembly as done by HMM-GRASPx. Unmapped reads were also separately annotated with MetaCLADE.
In parallel, we considered another assembly approach, based on the construction and annotation of a gene catalog  [71, 72] , where input reads are assembled using a canonical assembly pipeline, ORFs are predicted and putative gene sequences are clustered in order to create a non-redundant set (i. e. the gene catalog). The latter was annotated with MetaCLADE. Reads were finally mapped back to the catalog and the annotation transferred accordingly. Unmapped reads were also separately annotated with MetaCLADE. We compared MetaCLADE based on both assembly approaches against HMM-GRASPx. To assess annotation performances, the experiment was run on the dataset used in  [37]  which was augmented in order to better reflect the size of modern metagenomic datasets. More precisely, a set of 20 million paired-end reads was generated from a simulated marine dataset with uneven coverage and read length of 100 bp.
A total of 303 Pfam domain families that are involved in some important metabolic pathways were selected as input for HMM-GRASPx and MetaCLADE as in  [37]  (details are reported in the "Methods" section). The results of our comparison are reported in Table  2  showing that MetaCLADE always outperforms HMM-GRASPx regardless of the assembly method considered. MetaCLADE, however, performed better considering the gene catalog as a reference. The same observation still holds when clans are evaluated, that is when domain hits of the same clan (with respect to the gold standard) are counted as true positives. The same comparison had also been carried out on a simulated dataset of about 10 million 200-bp paired-end reads (keeping the same average coverage and composition of the 100-bp read dataset).
Here again, Meta-CLADE outperformed HMM-GRASPx (Additional file 1: Table  S13 ), especially when considering HMM-GRASPx's assembly as a reference for MetaCLADE. As a matter of fact, the lower number of reads led to a more fragmented gene catalog and a lower performance of MetaCLADE on these assembled sequences. However, it was not possible to increase the read number (e. g. to 20 million) due to the excessive amount of computational resources (more than 128 GB of RAM) demanded by HMM-GRASPx. In fact, the application of HMM-GRASPx seems limited to datasets of modest size (in terms of both read and profile number). Overall, even though MetaCLADE achieved very good performances in this assembly-based scenario, we should emphasise that it had been specifically tailored to work with relatively short fragments.
Comparison with HMM-GRASPx and UProC on datasets of 100-and 200-bp reads To investigate further MetaCLADE performance with very short reads, characterising the growing number of Illumina sequencing MG datasets available, we generated two sets of sequences, 100 bp (1,226,882 reads) and 200 bp (682,380 reads) long, from the Guerrero Negro Hypersaline Microbial Mats dataset (GNHM; see Additional file 1). GNHM was used in  [38]  to demonstrate UProC performance on short reads versus profile-based methods and, indeed, in  [42] , it was shown that profile hidden Markov models substantially lose sensitivity on shorter reads. The two datasets have been evaluated by considering the annotations obtained either with Pfam (version 27) or with CLADE domain library (also based on Pfam27). Table  3  shows a slightly better performance of UProC compared to HMM-GRASPx and MetaCLADE on GNHM for 100-bp reads with a gold standard set by Pfam annotation.
The two plots in Additional file 1: Figure  S17  display the relation between precision and recall, in order to evaluate the performance of the two tools on the two datasets, when the gold standard is CLADE and clans are considered. MetaCLADE displays a slightly better behaviour than UProC for fixed recall values. Notice that for very small recall values and, hence, very high hit scores, UProC detects a higher number of false positives (yet quite small) compared to MetaCLADE. This is seen with the behaviour of the curves in the zoomed plots in Additional file 1: Figure  S17.
A measure of the improvement obtained by combining MetaCLADE and UProC By combining MetaCLADE and UProC, we obtain an improved quality of the functional assignments of sequences of GNHM on both the 100-bp and the 200-bp datasets. Indeed, for the two datasets, we considered the annotation realised by MetaCLADE augmented with UProC annotation on those reads left unannotated by MetaCLADE, and called this approach Meta-CLADE+UProC (see the "Methods" section). Vice versa, we considered the annotation realised by UProC augmented with MetaCLADE annotation on those reads left unannotated by UProC, and called this approach UProC+MetaCLADE (see the "Methods" section). On both GNHM datasets, for the 100-and 200-bp reads, UProC+MetaCLADE outperformed on exact domain annotation and MetaCLADE+UProC on clan annotation.
This is expected because UProC, based on word matching, is intuitively similar to very conserved CCMs that are very close to known sequences. Most importantly, the performance of their combination is bringing a clear improvement in terms of number of correctly predicted domains independently on the gold standard and on the dataset (Tables  3  and 4 ). In Additional file 1: Figure  S17 , the precision-recall curve for Meta-CLADE+UProC shows that the addition of unique UProC annotations to MetaCLADE increases the number of false positives and therefore decreases precisions, as also seen in Table  3. The advantage in using a combined approach like MetaCLADE+UProC relies on the increase of the recall (since false negatives decrease) counterbalanced by a small decrease in precision (since false positives slightly increase).
The combination MetaCLADE+UProC has been tested also on the Rainforest MG dataset, where it produced a much larger number of read predictions as illustrated in Fig. 6e , accompanied by a distribution of E values showing high confidence (see Fig. 6f  and compare it to Fig. 4a, c ). MetaCLADE+UProC behaviour on the other four MG/MT samples in Table  1  is reported in Additional file 1: Figures  S12E-S15E  and shows a high improvement in performance associated.
Note that MetaCLADE package provides the program to pre-compute gathering thresholds for all domain models. This allows the user to compute appropriate thresholds based on new CCMs. Constructing a library of conserved small motifs. The search for sequence motifs in an environmental sample might be realised with a computationally costly "all against all" read comparison. Alternatively, starting from the most conserved patterns comprised in CCMs, we can generate a repertoire of significant motifs specific of each domain in order to improve hit selection criteria. A systematic classification of these motifs might lead to datasets of motifs that could be used as environmental signatures of metabolic activities. These "environmental patterns" could be also used to find new domains in environmental samples with Meta-CLADE.
New criteria to filter overlapping hits in MetaCLADE. Different domain hits could be selected by exploiting further the characteristics of the two-dimensional space of sequences pre-computed for the domains. For instance, one could privilege the domain hits with larger bitscore/mean-bit-score distance from the closest negative in the space. These filtering conditions could improve the annotation and need to be tested at large scale. MetaCLADE differences with CLADE. MetaCLADE has been designed with the purpose of annotating MG and MT reads. It exploits the multi-source annotation strategy introduced in CLADE and the CLADE model library, but it handles the models and their output in a different manner. Indeed, the CLADE pipeline combines the output of its rich database of probabilistic models with a machine-learning strategy in order to determine a set of best predictions for each domain sequence.
Conclusion MG and MT datasets have been explored mostly to learn about which and in what abundance species are present in the community. Learning about the functional activity of the community and its subcommunities is a crucial step to understand species interactions and large-scale environmental impact. Ecological questions, such as how limited availability of abiotic factors in an ocean shape most abundant genes in a community, or how temperature affects eukaryotic phytoplankton growth strategies, for example, can be approached with an accurate domain annotation and a precise functional mapping. In this respect, one might need to zoom into functional activities and metabolic pathways employed by the environmental communities that might involve non-highly expressed genes. This means searching for lowly abundant domains that, through cooperation, might imply important functional effects.
Methods This section explains MetaCLADE's methodology and the datasets used in the analyses. The differences between MetaCLADE and CLADE  [43]  are presented in the "Discussion" section. The time complexity is explained in Additional file 1. The testing datasets were designed to fit current technological characteristics. We considered that Illumina represents nowadays the dominant technology in most sequencing projects. Currently, the HiSeq and MiSeq platforms are able to produce pair-end reads of 150 bp and 300 bp, respectively (with a much higher put for the first one). Such fragments in practice might be even shorter after the required low-quality-base trimming. For these reasons, we considered testing datasets of reads of increasing lengths. More precisely, we chose 100 bp and 200 bp as read lengths.
The multi-source annotation strategy and the CLADE library Widely used search methods based on sequence-profile and profile-profile comparison, such as PSI-BLAST  [45] , HMMer  [47]  and HHblits  [49] , are based on a monosource annotation strategy, where a single probabilistic model, generated from the consensus of a set of homologous sequences, is used to represent a protein domain. The mono-source strategy typically performs well when sequences are highly conserved. In this case, the consensus model captures the most conserved features in domain sequences and it can be successfully used to find new domains in databases of sequences, sharing the same features as the original sequence. However, when sequences have highly diverged, consensus signals become too weak to generate a useful probabilistic representation and models constructed by global consensus do not characterise domain features properly.
To overcome this fundamental bottleneck, CLADE  [43] , a domain annotation tool tailored to full genomes, introduced a multi-source annotation strategy, where several probabilistic models are used to represent a protein domain. For each Pfam domain D i , CLADE considers the FULL set of homologous sequences S i in Pfam  [61]  associated to D i , and for some representative sequences s j in S i (see below), it constructs a model by retrieving with PSI-BLAST  [45]  a set of sequences similar to s j from the NCBI NR database. The probabilistic model generated in this way displays features that are characteristic of the sequence s j and that might be very different for other sequences s k in S i. The more divergent the homologous domain sequences s j and s k are, the more models constructed from these sequences are expected to display different features.
Identification of domain hits MetaCLADE takes as input a set of MG/MT sequences to be annotated and the CLADE model library. More specifically, the input sequences coming from a dataset of reads are expected to be (subsequences of) open reading frames (ORFs). Alternatively, one can use MetaCLADE on the six reading frame translations of the reads. Each sequence is scanned with the model library in order to identify all domain hits. Each hit is defined by a bit-score, that is the PSI-BLAST/HMMer score associated to the match, and by a mean-bit-score, that is the bitscore of the hit divided by its length. These two scores are used to evaluate the likelihood of the hit to represent a true annotation (see the "Selection of domain hits" section; for the computation of the likelihood see the "Pre-computed two-dimensional gathering thresholds for domain identification" section).
In particular, only one probability space is estimated for all CCMs of a domain. All positive and negative sequences, generated for all CCMs, are considered in the same sequence space and the associated probability space is computed. Figure  7a  illustrates an example of separation of the spaces of positive and negative hits for the CCMs and SCM of a Pfam domain (PF01036), analysed with the naive Bayes classifier. Note that short fragments have small bit-scores with a possibly large mean-bit-score and that negative sequences are characterised by small bitscores and small mean-bit-scores. Also, the identification of fragmented coding regions (especially important for the annotation of MG/MT datasets, where ORFs that are present in MG/MT reads are fragmented) will only be realised through very low bit-scores because only parts of domains are present.
Also, the usage of a two-dimensional sequence space, determined by bit-scores and mean-bit-scores, improves the separation of positive and negative sequences in MetaCLADE compared to HMMer (hmmscan). In the plot describing a b Fig. 7  Naive Bayes classifier analysis of Pfam domains. a Naive Bayes classifier analysis of the training set evaluating thresholds for positive (cyan dots) and negative (orange dots) hits discrimination of the Pfam domain "bacteriorhodopsin-like protein" (PF01036). The space of positive hits is coloured blue, and the space of negative hits is the complementary one, coloured from green to dark red (see colour scaling). From dark red to blue, each coloured area is characterised by a different probability for a sequence to be a positive sequence. The space of solutions for all CCMs confounded (left) and the SCM (pHMM; right) are described.
The red dots correspond to the MG sequences identified by the models in the Equatorial Pacific (EPAC) MT dataset discussed below  [57]. The green vertical line on the SCM plot (right) corresponds to hmmscan GA threshold (= 24, for this domain). Note that the grid is discrete and that rectangular regions are coloured with respect to probability intervals. This means that two adjacent regions with the same colour might have associated two different probabilities. b The probability space of CCMs and SCMs generated for MetaCLADE; all domains included. Each point of the plot is the average of the corresponding points in all spaces of solutions computed for CCMs and SCMs, respectively. Examples of such spaces of solutions are reported in a for the Pfam domain PF01036. The colour scale is defined with respect to probability values associated to the regions.
The green vertical line in the SCM plot (right) corresponds to the average hmmscan GA threshold; all domain confounded. Note that if we take the bit-score only as a threshold in our database, we obtain that an average bit-score at 25 in a is statistically meaningful, independently of the mean-bit-score. This value is used as a default threshold in http://www. ebi. ac. uk/Tools/hmmer/search/phmmer the SCM space (Fig. 7a , right), the Pfam GA threshold excludes most of the Bacteriorhodopsin-like protein sequences detected by CCMs. More generally, for all Pfam domains, we computed the difference between the GA threshold associated to the SCM and the mean of the bitscores for the five best negative sequences identified by the SCM. The distribution of the differences, displayed in Additional file 1: Figure  S5 , shows a small standard deviation suggesting that MetaCLADE and Pfam/HMMer estimation of the (one dimensional) cutoff is similar.
Comparison to InterProScan To compare MetaCLADE to InterProScan, we used several tools and domain model libraries: Pfam  [61] , TIGRFAM  [64] , Gene3D  [63]  and PRINTS & ProSite  [65, 66]. Annotations produced by these tools were downloaded from available annotation files provided by the EBI metagenomics pipeline http://www. ebi. ac. uk/ metagenomics/pipelines/ (versions 1. 0 or 2. 0 depending on the dataset). The EBI metagenomics pipeline uses InterProScan (v5. 0 for pipeline 1. 0 and 5. 9-50 for pipeline 2. 0) as tool to annotate predicted ORF (using FragGe-neScan 1. 15) using the domain libraries above. HMMer v3. 0 was downloaded from http://hmmer. org and run with default parameters and curated inclusion thresholds (option -cut_ga).