dataesr/scientific-paragraphs-categorization huggingface dataset:
Benchmarking BIGSI We benchmarked the empirical scaling properties of BIGSI against the Sequence Bloom Tree (SBT)13 and the Split Sequence Bloom Tree15 (SSBT) using a dataset of 10,000 random microbial sequence datasets from the ENA. We compared build and query times, and peak storage requirements on 21 subsets ranging from 100 to 10,000 samples. For each subset, we built a BIGSI, an SBT and an SSBT using two different parameter settings: one optimized for speed (recommended by the authors in Lemma 1, and subsequent text, of Solomon et al13 and the SBT/SSBT user manual) (SBT-fast/SSBT-fast), and another optimized for compression (SBT-small/SSBT-small) (Online methods). For database sizes of more than 2,000 samples, peak storage requirements of SBT-fast and SSBT-fast exceeded our available disk space (1 Tb; by comparison the total storage required for the input data was 403 Gb; Online methods).
We queried the resulting indexes for 2,157 antimicrobial resistance genes with a mean length of 937 bp and total query length 2,021,655 bp. SBT-fast, SBT-small, SSBT-fast, SSBT-small and BIGSI returned near-identical hits from the exact match search, with only 1 difference across all queries (measured on N = 1,000 database size). With the inexact search, concordance was >99% with the differences likely due to differences between methods in the construction of the underlying bloom filters. Inexact (T = 40%) query times versus peak storage requirements for the various methods are shown in Figure  3 ; BIGSI maintains good query performance in small space for all input data sizes, whereas (since they need to build an uncompressed tree of bloom filters before compressing it) both SBT and SSBT trade storage requirements for performance (we quantify performance below).
SSBT-fast and SBT-fast have query time comparable to or better than BIGSI but require orders of magnitude more storage to build the uncompressed tree. SBT-small and SSBT-small have almost equivalent storage as BIGSI, but slower query performance, due to saturation of internal bloom filters within the sequence bloom tree, a result of the total number of unique k-mers being significantly larger than the number of unique k-mers in any individual dataset. With 2,000 samples in the index, the SBT-fast query was only slightly slower than BIGSI's (292s vs 274s) but required 132X (878GB/6. 6GB) more storage during construction, and 9X (57. 5GB/6. 6GB) more storage after compression. Further details in Online Methods. Similar results were found for exact queries (T = 100%, Supplementary Figure  4 ); exact match searches with BIGSI were up to 2-3X faster than inexact queries.
To ensure our SBT comparison was fair, we also benchmarked BIGSI and SBT on the same human RNA-seq dataset used in 13,15,30, using the prebuilt SBT index provided by the authors. We measured the query time of 1000 RNA transcripts randomly selected from the 214,294 known transcripts (reported in 13), and compressed index size. BIGSI has faster query time than SBT in smaller space for these datasets (360 s and 144 Gb for BIGSI, and 2221 s and 200 Gb for SBT (T = 70%); details in Online methods and full results, including for exact search, in Supplementary Table  2 ). Finally, we simulated the scaling of storage requirements required to build SBT-fast, and BIGSI for dataset of up to 1 million genomes of two types: genomes with high proportions of k-mer sharing (e. g. human), and genomes with lower proportions of k-mer sharing (e.
An index of all bacterial and viral sequence in the ENA We downloaded the entire set of 469,654 bacterial, viral and parasitic WGS datasets in the ENA as of December 2016. After a step to exclude eukaryotic genomes on the basis of size (Online Methods), we were left with 447,833 datasets, from which we created a BIGSI index. This index, named the 'all-microbial index' required 1. 5 TB of storage, <1% of the original data size (170 TB) and contained more than 60 billion unique k-mers. Data download took 6 weeks, and we built bloom filters 'on the fly'. Combining the complete set of bloom filters took approximately 2 days (Online Methods). We estimate that the intermediate storage required to build an SBT-fast/SSBT-fast of the same data with recommended bloom filter size equal to the number of unique k-mers in the collection, would be >6.
7PB. We used our 'all-microbial index' for several analyses. First, we estimated the species and abundances present in each dataset that within the all-microbial-index using the Bayesian abundance estimator Bracken31, which parses output from the read-classifier Kraken32 (Online Methods). Using this method, we found that more than 90% of the input datasets were from only 20 genera, and 65% of the input datasets were from the 5 most common bacterial genera (Salmonella, Streptococcus, Staphylococcus, Escherichia and Mycobacterium); counts for the most prevalent bacterial genera in the all-microbial index are shown in Supplementary Figure  6.
Ultrafast gene search using the all-microbial index We searched for exact matches to the mobilized colistin resistance genes MCR-1/2/333-37 in the all-microbial-index. Searching for all 3 genes in 10x more genomes than has been previously reported took just 1. 73 s seconds in total. We did not detect MCR-2. MCR-1 was present in 169 datasets from 3 species (Escherichia coli, Salmonella enterica, Enterobacter aerogenes) and MCR-3 was present in 34 datasets from 3 species (E. coli, S. enterica, Klebsiella pneumoniae) (Supplementary Data 1). Our all microbial index therefore enables almost instant detection of drug resistance genes, or indeed any other gene studied, across the global corpus.
Host-range of plasmids in the all-microbial index We downloaded a set of 2,827 plasmids from the ENA (Online Methods, Supplementary Data 2) and ran an inexact (T = 40%) search for them in the all-microbial-index. We restricted analysis to hits with T > 90%. The total length of query sequences was 227 Mbp, and the query took 2,120 CPU hours (11 days) on a single server using 8 cores and 1. 5 Gb RAM per process. The search returned 665,619 hits with 121,758 unique accessions across 258 genera. Since contamination could confound observations of a plasmid in a genus, we excluded all datasets containing evidence of more than one genus at abundance above 0. 1%. Only 41% (184,652) of datasets and 62% of search hits passed this filter. In the filtered output we identified plasmids shared by closely related genera, for example Escherichia, Shigella and Salmonella or Enterococcus, Streptococcus and Staphylococcus.
We found 37 plasmids present in at least five datasets from at least two genera (Figure  4 , Supplementary Data 2, 3); and 5 plasmids that were present in multiple orders and families. The plasmid pETHIS-1 (GenBank identifier: AF012911) was found in 5 phyla, 10 taxonomic classes and 17 genera. This plasmid is used as an expression vector and its identification in so many species serves as a positive control, confirming that BIGSI can detect common plasmids. Of more biological interest, the Tn916 conjugative transposon encoding tetracycline resistance that was first found in Enterococcus faecium and is known Sampling biases in the ENA preclude inference about plasmid distribution, but they do allow us to ask if plasmids bearing antibiotic resistance genes are more widely distributed than those plasmids that do not have antibiotic resistance genes.
We undertook a similar analysis of the whole ENA, by using an exact match search of the allmicrobial-index for nucleotide alleles corresponding to previously identified MOB and T4SS "types" (alleles defined in amino acid space). We applied the same contamination filter mentioned above (restricting to 184,652 datasets), and required exact matches to at least one MOB and one T4SS be present to be considered a putative conjugative system. 36,030 datasets met these criteria, giving an estimate of 19. 5% of bacteria in the ENA containing a conjugative system, which is consistent with the previous estimate of 18% in Guglielmini et al39. This proportion varied by phylum (minimum 0. 5% in Spirochaetes and maximum 31. 7% in Firmicutes, full data in Supplementary Table  3 , and MOB type distribution in Supplementary Figure  8 ).
These data could be mined further to explore the potential spread of antibiotic resistance genes. For example, by analysing datasets shown to contain MOB T we observed genetic flux between Staphylococcus and Streptococcus, but not between either of them and Salmonella. We also found MOB Q in Salmonella and Streptococcus but not Staphylococcus, indicating a different probability of cross-genus (and cross-phylum) transfer by conjugation (Supplementary Data 4).
Antibiotic resistance genes in the ENA To understand whether antibiotic resistance gene prevalence in the ENA has changed over time, we downloaded the set of 2,157 sequences associated with antibiotic resistance from the the Comprehensive Antibiotic Resistance Database (CARD) database (v1. 1. 7)40 searched for this set of genes in the all-microbial-index with thresholds of 100% and 70% (full results for both searches are in Supplementary Data 5). An exact search for a single gene took an average of 1. 1s and returned an average of 438 hits, whereas the corresponding figures for an inexact search (T = 70%) were 34. 4s and 5,320 hits. We show the prevalence over time in 3 genera, and in total, in Figure  6. Focusing only on Staphylococcus, we found that the proportion of datasets containing mecA gene (methicillin resistance), fell from 70% in 2013 to 40% in 2016, and that prevalence of tet and aac genes also decreased (Figure  6b ).
For Klebsiella however, the prevalence of all resistance genes that we screened for in the allmicrobial dataset increased over time (Figure  6c ). In M. tuberculosis, antibiotic resistance mainly arises through amino acid mutations in specific genes 1,41. We genotyped the all-microbial-index, which contains 30,226 M. tuberculosis genome datasets, for the 206 resistance mutations enumerated by Walker et al. 27. This exercise took only 103 minutes on a single-core, which is approximately 10,000x faster than typing each dataset individually with the fast resistance prediction software Mykrobe predictor1. Our results indicate an increase in prevalence of MDR-TB in the ENA since 2011 (Figure  6d ). It should be noted that one cannot infer anything about global prevalence from this, as the contents of the ENA are subject to sampling biases (that is, are dependent on specific research studies rather than containing an unbiased sampling of isolates from tuberculosis patients).
Importantly, new datasets can be readily and quickly added to our index, allowing it to grow as new datasets are sequenced. Our indexing method works for both raw data and assembled genomes, so it is poised for a future in which finished reference genomes are routine. BIGSI was designed with SNP or indel genotyping and allele search queries in mind. It allows a user to identify datasets worthy of in-depth study, as well as enabling global monitoring of specific alleles. Nevertheless, BIGSI has limitations, which it shares with previous k-mer tools such as the SBT; we give the main ones here, but also see Online Methods. It is a k-mer index and is therefore as lossy as all de Bruijn methods. Reconstruction of stored genomes is impossible and repeat regions cannot be resolved. Also BIGSI does not store coverage information which rules out queries in which copy number is important.
One example of such a query is detecting azithromycin resistance in Neisseria gonorrhoeae, in which the level of resistance is mediated by the number of rRNA genes containing a particular SNP45. BIGSI only supports nucleotide k-mers although we would note that an extension to amino acid search would be straightforward, because bloom filters are agnostic to what they are storing. Finally, BIGSI is optimized for very diverse datasets, in which the combined unique k-mer count is much higher than that of any individual sample e. g. the entire microbial ENA/SRA. Where samples are less diverse (eg. human), we expect SSBT or Mantis to be a better choice, as they take advantage of compression from sample similarity. In the future we envisage that sequencing data volumes will continue to grow and the userbase that wants to access and exploit those data will increase to include clinical and public health practitioners.
BIGSI construction and querying BIGSI indexes a set of N (number of datasets) bloom filters by position in the bloom filter. Each bloom filter must be constructed with the same parameters (m, η), where m is the bloom filter's length in bits and η is the number of hash functions applied to each k-mer. The same hash functions must also be used to construct each bloom filter. To construct a BIGSI, the N bloom filters are column-wise concatenate into a matrix. The row index and row bit-vectors are then inserted into a hash table or key-value store as key-value pairs so that row lookups can be done in O(1) time. This set of key-value pairs can be stored on disk, in memory or distributed across several machines and is indexed via a hash index (a b-tree would be an alternative option). To insert a new bloom filter we simply append it as a column to the existing bitmatrix.
Parameter choices The choice of BIGSI parameters (m, η), depends on: the maximum number of k-mers expected in any dataset (K max ), the number of datasets (N) expected, the smallest number of unique k-mers in a query sequence (L min ) to be supported, the k-mer size (k) and the maximum number of acceptable false discoveries per query (q max ). The expected number of false discoveries (V) for any query can be calculated as q = E[V] = Np L where p is the false positive rate of the bloom filter and L is the number of unique k-mers in the query, assuming independence of the k-mers and bloom filters. Parameters m and η determine false positive rate for a bloom filter with K max elements -if there are fewer elements, then the false positive rate will be lower. We assume below that all bloom filters have the maximum number of unique k-mers inserted to give an upper bound on error rate.
However, the independence assumption mentioned above is strictly speaking false for real data archives (such as the ENA/SRA) which have biased distribution of datasets across the phylogeny. If there is an enrichment of datasets from some genus, then those columns (bloom filters) will be more similar, and conditional on a false positive in one column, the probability of a false positive in the similar columns will be elevated. To keep q below a chosen threshold q < q max for a given N and k, p must be chosen to satisfy q for the smallest allowable number of unique k-mers in a query L min : p L min = q max N. Therefore, the desired bloom filter false positive rate is p = q max N 1 L min. Since, for a given number of inserted k-mers (n), and desired false positive rate (p), optimal bloom filter parameters can be determined by the following formula48 which becomes: m = - K max ln(q max /N) L min ln(2) 2 , η = - ln(q max /N).
BIGSI parameters for all-microbial-index We assume initially that a bacterial dataset contains at most 10 million k-mers since bacterial genomes are generally under 6 Mb in length, leaving 4 million k-mers available for some sequencing errors which escape de-noising, and plasmid variation. Unless otherwise specified we use BIGSI parameters m=25,000,000, η = 3 for all analyses. For these parameters, the upper bound on number of false discoveries, q max for an SNP allele from a probe (flanks plus allele) of length (L min 2k-1=61) is 10 -9 (if K max = 10 7 , N = 10 6 ).
K-mer identity and scoring If, for example, the k-mer size is 31, each SNP difference between the query and the nearest sequence in the index causes a window of 31 absent k-mers. Therefore k-mer identity drops more rapidly than sequence identity; e. g for k=31, matches with sequence identity above 80% can be detected (Supplementary Figure  2 ). Therefore, BIGSI most naturally enables exact or close match searches, or situations where combinatorial searching is feasible, such as querying all sequences 2 SNP different from a given allele, or all single amino acid changes in a gene. Datasets containing matches for a more divergent allele can be sought by searching for subsequences (seeds). Although BIGSI does not carry out an alignment, an approximation to a megaBLAST alignment score can be inferred from the presence/absence pattern of k-mers in the query.
Benchmarking query time and storage requirements of BIGSI, SBT, and SSBT We randomly chose 10,000 microbial cleaned de Bruijn graphs from the all-microbial-index accessions and we then further randomly sub-sampled these into collections of 10, 100, 200, 300, 400, 500, 600, 700, 800, 900, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000 datasets. A BIGSI of each set of datasets was built with parameters (m = 2. 5 × 10 7 ;η = 3). A SBT and SSBT were built for each set of datasets with η = 1 and bloom filter size (m) equal to the count of the total number of k-mers in the collections of graphs, as recommend in the text following Lemma 1 of Solomon et al13, called "SBT-fast" and "SSBT-fast" respectively. Redis hyperloglog was used to count the unique k-mers in the set of cleaned graphs for each increment.
A SBT and a SSBT were built for each dataset with m = 2. 5 × 10 7 ;η = 3 (the same bloom filter parameters as BIGSI), called "SBT-small" and "SSBT-small". Construction and query time analyses were run on an Amazon Web Service i3. 8xlarge instance with 32vCPUs, 224 GiB of memory and 4 x 1. 9 TB non-volatile SSD-backed instance storage. SBT-fast construction exceeded 1TB for 3000 datasets, the practical limit we set, and as a result SBT-fast was not built for increments above 2,000 datasets. The search time comparison was run with 'bt query' and 'bigsi search -seqfile $f', using kmer thresholds 40% and 100%. A full table of results can be found in Supplementary Data 6. As mentioned above, SBT-fast and SSBT-fast construction failed for databases > 2,000 due to excessive disk requirements, but we were able to calculate a lower bound for the peak disk-use for database sizes > 2,000; query times were extrapolated linearly (these points are shown as triangles on Figure  3 ; details of lower bound below).
There were approximately 4. 4 × 10 9 unique k-mers in the union of all 10,000 datasets, almost 1000 times more than in a typical individual dataset. It would require ~11TB of storage to build SBT-fast, 350X more space than is required by BIGSI. By contrast, querying SBT-small took 65X longer and SSBT-small 448X longer than a BIGSI of the same database size (N=10,000). SSBT-small was nearly 7X slower than SBT-small, perhaps due to the lower proportion of shared k-mers between datasets. Even for this small benchmark dataset, constructing unsaturated (i. e. fast) SBTs and SSBTs for larger numbers of datasets is prohibitive in terms of storage requirements; we address scaling further below. Query times for SBT/SSBT-slow were unacceptable for this benchmark, and so we excluded them as candidates for storing the ENA/SRA, which is 50X bigger than this benchmark dataset.
Simulation of storage requirements for a BIGSI for N datasets is given by: BIGSI storage [bytes] = mN 8 Although it is possible to append to an SBT incrementally, as new microbial datasets will keep adding new k-mers, this will lead to saturation of the root-level bloom filter in the SBT, and a collapse in query performance. This can be avoided by reconstructing the SBT, ensuring the bloom filters are large enough to support the full set of k-mers. This was borne out by our benchmarking. In simulating scaling to a million genomes, we therefore focussed on SBT-fast rather than SBT-small. As the best case for a binary tree with N leaves is 2N-1 nodes, we estimate: SBT storage [bytes] = (2N -1)N k 8 where N k is the total number of k-mers in the combined set of datasets and also equal to the size of the bloom filter required.
This is a lower bound for the peak storage use. See Supplementary Data 6, for the close correspondence between this theoretical estimate and the empirical measurement in our benchmark datasets. As a result, we can calculate lower bounds for the peak disk usage for SBT/SSBT; when the explosion of disk usage made construction of indexes for the benchmark datasets with size>2000 datasets, we used this lower bound to plot (suitably labelled) extrapolated datapoints.
Benchmarking against Mantis We attempted to benchmark against Mantis but despite assistance from the authors we were unable to resolve a number of issues: large numbers of false positive hits on some data, and a bug causing segmentation faults when trying to build more than 3000 datasets. We reluctantly excluded Mantis from our benchmarking due to limitations of time. As currently implemented, the Mantis data structure does not support incremental insertion, as it needs up-front the full set of k-mers and for each k-mer, the list of datasets containing it ("colour class", stored as a bit-vector). There is currently an intermediate stage where the uncompressed colour-class matrix is held in RAM, scaling quadratically in number of datasets.
Benchmarking query time and storage requirements of BIGSI on RNA-seq data De Bruijn graphs (k=31) were constructed and cleaned from the downloaded RNA-seq fastq files listed in13. A BIGSI was built with bloom filter parameters m = 4000000000; η = 1 by chunking into 1,600 batches, building in parallel and then combining into a final index. The SBT and the 214,294 transcripts were provided by the authors (personal correspondence). 1,000 transcripts were randomly selected from the full set and queried with 'bt query' and 'bigsi search -seqfile $f' respectively. Since the SBT was pre-built by the authors, and therefore compressed, we used a RocksDB back-end store for BIGSI, which used Snappy (https://google. github. io/snappy/) for bitslice compression. This gave us a compressed BIGSI index, allowing us to compare compressed index size.
Using an array of BIGSI to support variable dataset size A limitation of BIGSI is that K max , the maximum number of k-mers per dataset, must be set in advance. One way to extend BIGSI to datasets with varying k-mer cardinality (e. g. metagenomic data) is to build multiple BIGSIs with different K max , e. g. K max = 10 5 , 10 6 , 10 7 , etc …, and insert each sequence into the appropriate level by k-mer counting before insertion. Queries can then be sent to all the indexes in parallel and then results combined.
Genotyping accuracy measurement on TB Conservative SNP calls were made using Cortex17 (independent workflow, k=31) on 3480 Mycobacterium tuberculosis datasets from Walker et al27. Singleton variants were discarded, and a de-duplicated list of 68,695 SNPs was constructed. We generated "probe sets" consisting of a reference and alternate alleles of these variants from the NC_000962. 3 reference. An index of the 3,480 datasets was built and 100 random datasets were genotyped at the 68,695 sites as follows: Each allele of the probe-set is searched for in the BIGSI resulting in Boolean presence/absence of each allele. This requires querying for multiple probes for each variant. If only a reference allele is present the genotype is returned as 0/0, if only an alternate allele 1/1, if both 0/1 and if neither -/-.
We compared the concordance of genotypes of the 100 random datasets with those generated with the samtools pipeline from Walker et al27, excluding filtered positions. As described in the main text, the concordance between methods was 99. 997% with a total of 286/682,690 discrepancies. The majority of these discrepancies (203/286) were mixed (heterozygous) calls from BIGSI; samtools had been run with a haploid model it did not make any mixed calls, so we expect some of these were correct.