Answers.Body: <p>One thing that hasn't been suggested so far would be phylogenetic placement. First option there would be placing against one larger tree roughly touching the domains you mentioned to get an overview. Second option would be having three trees (bacteria, plants, viruses) made up of relevant taxa for the type of sample you have.</p>
  
  <p>This would give you a detailed and accurate view of your sample, but it does involve some work. Sadly the computational throughput for placement is still a bit behind BLAST, and not many complete pipelines exist, but in my opinion it is worth the extra mile (I might be biased. I'm working to resolve as many of these issues as possible).</p>
  
  <p>If you want to look into it, I can reccomend <a href="https://academic.oup.com/sysbio/article/64/1/e26/2847641" rel="nofollow noreferrer">this excellent paper by Erick Matsen</a>. Since your samples sound vaguely soil-ish, I will also link <a href="https://www.nature.com/articles/s41559-017-0091" rel="nofollow noreferrer">this recent Nature paper</a> using placement to analyze rainforest soil samples.</p>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-09-15T14:10:56.140Z
Answers.LastActivityDate: 2017-09-15T14:10:56.140Z
Answers.UserDisplayName: Pierre Barbera
Answers.Body: <p>I think you're just looking for <code>nr</code>. It is absolutely not limited to prokaryotes, far from it. It is, according to the blurb on the NCBI's blast page:</p>
  
  <blockquote>
    <p>The nucleotide collection consists of GenBank+EMBL+DDBJ+PDB+RefSeq sequences, but excludes EST, STS, GSS, WGS, TSA, patent sequences as well as phase 0, 1, and 2 HTGS sequences. The database is non-redundant. Identical sequences have been merged into one entry, while preserving the accession, GI, title and taxonomy information for each entry.</p>
  </blockquote>
  
  <p>The link you gave in your question doesn't describe the <code>nr</code> database, but the collection of RefSeq non-redundant sequences. The two are not the same thing at all. The nr databases contains sequences from all available species so it should be perfect for what you're after. </p>
  
  <p>Now, I admit I haven't tried this with millions of redas, but I did used to run many blast searches with several thousand sequences using the command line blast tools and the remote nr database. I just checked the <a href="https://www.ncbi.nlm.nih.gov/books/NBK279690/" rel="nofollow noreferrer">blast+ docs</a> now and see no mention of a rate limit. Which isn't to say there isn't one, for sure, but since each query is run sequentially when you give a multifasta file, there seems to reason to limit since the servers should be able to queue the jobs as needed. </p>
  
  <p>So, don't use the API or roll your own, just <a href="https://blast.ncbi.nlm.nih.gov/Blast.cgi?CMD=Web&amp;PAGE_TYPE=BlastDocs&amp;DOC_TYPE=Download" rel="nofollow noreferrer">download NCBI's blast client</a> and try that:</p>
  
  <pre><code>blastn -query query.fa -db nr -task blastn 
  </code></pre>
  
  <p>If that complains, you'll have to download the <code>nr</code> DB (well, since you want nucleotide blast, you probably want its nucleotide cousin, <code>nt</code>) which you can find here: <a href="ftp://ftp.ncbi.nlm.nih.gov/blast/db/" rel="nofollow noreferrer">ftp://ftp.ncbi.nlm.nih.gov/blast/db/</a>. Just get all the <code>foo.nt*</code> files.</p>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-01T17:58:48.210Z
Answers.Text: According to the [BLAST API Usage guidelines](https://blast.ncbi.nlm.nih.gov/Blast.cgi?CMD=Web&PAGE_TYPE=BlastDocs&DOC_TYPE=DeveloperInfo), this isn’t really supported, even if there is no hard limit. But you seem to be right about NR, this seems to be what I’m after, after all.
Answers.UserDisplayName: Konrad Rudolph
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-01T18:27:24.883Z
Answers.Text: @KonradRudolph I'm hoping that refers to using the API in your own tools and not to using their own standalone blast+ program. I may very well be wrong, but I'm not sure if the instructions for the API apply to them as well. I would expect them to have built their tool in such a way as to not overload the servers. I figure it's at least worth a try.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-01T17:36:43.390Z
Answers.LastActivityDate: 2017-08-01T17:36:43.390Z
Answers.UserDisplayName: terdon
Answers.Body: <p>You might want to look into <a href="https://en.wikipedia.org/wiki/MinHash" rel="nofollow noreferrer">minihashes</a> as a more efficient way of clustering and comparing your sequences to known species and samples. <a href="https://genomebiology.biomedcentral.com/articles/10.1186/s13059-016-0997-x" rel="nofollow noreferrer">In the mash paper</a> the authors show how it would be feasible to create minihash sketches that you could search against from large metagenomics datasets (check the "Clustering massive metagenomic datasets" section in the results):</p>
  
  <blockquote>
    <p>For a large-scale test, samples from the Human Microbiome Project [36]
    (HMP) and Metagenomics of the Human Intestinal Tract [37] (MetaHIT)
    were combined to create a ~10 TB 888-sample dataset. Importantly, the
    size of a Mash sketch is independent of the input size, requiring only
    70 MB to store the combined sketches (s = 10,000, k = 21) for these
    datasets.</p>
  </blockquote>
  
  <p>Moreover, the creators of <a href="https://sourmash.readthedocs.io/en/latest/" rel="nofollow noreferrer">sourmash</a> (a tool using the same approach) provide <a href="https://sourmash.readthedocs.io/en/latest/databases.html" rel="nofollow noreferrer">precomputed sketches for the whole of NCBI's RefSeq database</a> (includes bacterial, viral and fungal genomes). Creating a similar sketch for the whole of ENA/SRA (or a subset of it) could require some time but would allow for very fast and accurate searches.</p>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-03T21:38:09.977Z
Answers.LastActivityDate: 2017-08-03T21:38:09.977Z
Answers.UserDisplayName: mgalardini
Answers.Body: <p>As you point out, BLAST is not going to work well for reads from a typical RNASeq experiment.</p>
  
  <p>From an RNASeq dataset, you might be able to get away with just a ribosomal screen using <a href="https://www.arb-silva.de/projects/ssu-ref-nr/" rel="noreferrer">Silva</a>, because the ribosome is the most abundantly transcribed gene. This will usually be the case even if polyA selection or ribosomal depletion was carried out on a sample.</p>
  
  <p>However, bear in mind that there will be a lot of false results due to database incompleteness. I used <a href="http://ccb.jhu.edu/software/kraken/" rel="noreferrer">Kraken</a> when I did something like this a few years ago, and ended up with one sequence annotated as a fruit bat (where my subject was a water-living planarian flatworm). <a href="https://github.com/infphilo/centrifuge/" rel="noreferrer">Centrifuge</a> would be a better option now; the authors even claim that <a href="http://genome.cshlp.org/content/26/12/1721.full" rel="noreferrer">they can store the entirety of nr in a 70 GB index</a> with Centrifuge (which would be better / more comprehensive than just a ribosomal screen).</p>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-01T21:54:51.507Z
Answers.LastActivityDate: 2017-08-01T21:54:51.507Z
Answers.UserDisplayName: gringer
Answers.Body: <p>I have had good success using the MEGAN software to detect and visualize contamination.  My input data was WGS not RNA-seq But I would check it out.  It can create a very nice weight phylogenetic tree for all the detected organisms.</p>
  
  <p><a href="https://ab.inf.uni-tuebingen.de/software/megan6" rel="nofollow noreferrer">https://ab.inf.uni-tuebingen.de/software/megan6</a></p>
  
  <p><a href="https://i.sstatic.net/ttT6l.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/ttT6l.png" alt="enter image description here"></a></p>
  
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-03T22:32:43.370Z
Answers.Text: I've downloaded MEGAN but my (admittedly cursory) reading of the documentation made it seem as if the software only visualised results computed by other tools.
Answers.UserDisplayName: Konrad Rudolph
Answers.ContentLicense: CC BY-SA 3.0
Answers.CreationDate: 2017-08-03T18:55:25.983Z
Answers.LastActivityDate: 2017-08-03T18:55:25.983Z
Answers.UserDisplayName: Bioathlete
Body: <p>I’ve got an RNA-seq dataset with a large proportion of environmental RNA “contamination”. BLASTing random reads reveals that much of the data comes from bacterial, plant and viral RNA. My target organism only accounts for ~5% of the RNA-seq read data.</p>
  
  <p>I would like to get a more or less (more, if possible) comprehensive overview of what species are found in my sample — bacteria, plants, animals (?), but also viruses. How can I perform this?</p>
  
  <p>I have been unable to find a “standard” way of performing a metagenome screen. BLAST online services all seem to be severely rate limited (certainly unable to handle uploading ~80 M reads). Installing BLAST (or a similar tool) locally is of course not an issue, but I can’t find a comprehensive database that spans all phyla — the best I’ve been able to find are databases, such as <a href="https://www.ncbi.nlm.nih.gov/refseq/about/nonredundantproteins/" rel="noreferrer">NCBI-NR</a>, that are restricted to single phyla or classes, or to bacteria.</p>
  
  <p>Ideally I’d like a workflow that I feed with my RNA-seq data, and get out a list of species with % of reads mapped.</p>
  
Comments.ContentLicense: CC BY-SA 3.0
Comments.CreationDate: 2017-08-01T16:46:41.617Z
Comments.Text: NCBI's nr database isn't restricted to any species as far as I know.  Why do you say it's limited to  single phyla or classes or to bacteria? Wouldn't blasting against nr be good enough? And are they really that rate limited? I admit I haven't done this in a few years, but back in the day, I could use NCBI's command-line blast client which would use a remote nr DB and could happily handle several thousand sequences. Have you tried that? Perhaps by splitting your input reads?
Comments.UserDisplayName: terdon
Comments.ContentLicense: CC BY-SA 3.0
Comments.CreationDate: 2017-08-01T16:50:45.040Z
Comments.Text: @terdon From what I’ve read, it’s only for procaryotes. But I may be completely wrong. Regarding the rate limit, the NCBI BLAST FAQ states that they don’t support batch usage, and that large jobs shouldn’t be run through the URL API.
Comments.UserDisplayName: Konrad Rudolph
Comments.ContentLicense: CC BY-SA 3.0
Comments.CreationDate: 2017-08-01T19:16:24.470Z
Comments.Text: Nr covers a wide range of taxa as @terdon says. You can download it and build your own blast db if you have ~100gb disk space and a medium memory server. However, BLAST will be slow, I would recommend DIAMOND instead
Comments.UserDisplayName: Chris_Rands
Tags: <rna-seq><metagenome>
Title: How do I efficiently perform a metagenome screen of “all” species?
ThreadText: # How do I efficiently perform a metagenome screen of “all” species?
  
  (Asked by: username_0 on 2017-08-01)
  
  I’ve got an RNA-seq dataset with a large proportion of environmental RNA “contamination”. BLASTing random reads reveals that much of the data comes from bacterial, plant and viral RNA. My target organism only accounts for ~5% of the RNA-seq read data.
  
  I would like to get a more or less (more, if possible) comprehensive overview of what species are found in my sample — bacteria, plants, animals (?), but also viruses. How can I perform this?
  
  I have been unable to find a “standard” way of performing a metagenome screen. BLAST online services all seem to be severely rate limited (certainly unable to handle uploading ~80 M reads). Installing BLAST (or a similar tool) locally is of course not an issue, but I can’t find a comprehensive database that spans all phyla — the best I’ve been able to find are databases, such as [NCBI-NR](https://www.ncbi.nlm.nih.gov/refseq/about/nonredundantproteins/), that are restricted to single phyla or classes, or to bacteria.
  
  Ideally I’d like a workflow that I feed with my RNA-seq data, and get out a list of species with % of reads mapped.
  
  **Tags:** rna-seq, metagenome
  
  > **Comments:**
  > **username_1**: NCBI's nr database isn't restricted to any species as far as I know.  Why do you say it's limited to  single phyla or classes or to bacteria? Wouldn't blasting against nr be good enough? And are they really that rate limited? I admit I haven't done this in a few years, but back in the day, I could use NCBI's command-line blast client which would use a remote nr DB and could happily handle several thousand sequences. Have you tried that? Perhaps by splitting your input reads?  
  > **username_0**: @username_1 From what I’ve read, it’s only for procaryotes. But I may be completely wrong. Regarding the rate limit, the NCBI BLAST FAQ states that they don’t support batch usage, and that large jobs shouldn’t be run through the URL API.  
  > **username_2**: Nr covers a wide range of taxa as @username_1 says. You can download it and build your own blast db if you have ~100gb disk space and a medium memory server. However, BLAST will be slow, I would recommend DIAMOND instead  
  
  ---
  
  ## 5 answers
  
  ### Answer by username_5 on 2017-08-01. Score: 7
  
  As you point out, BLAST is not going to work well for reads from a typical RNASeq experiment.
  
  From an RNASeq dataset, you might be able to get away with just a ribosomal screen using [Silva](https://www.arb-silva.de/projects/ssu-ref-nr/), because the ribosome is the most abundantly transcribed gene. This will usually be the case even if polyA selection or ribosomal depletion was carried out on a sample.
  
  However, bear in mind that there will be a lot of false results due to database incompleteness. I used [Kraken](http://ccb.jhu.edu/software/kraken/) when I did something like this a few years ago, and ended up with one sequence annotated as a fruit bat (where my subject was a water-living planarian flatworm). [Centrifuge](https://github.com/infphilo/centrifuge/) would be a better option now; the authors even claim that [they can store the entirety of nr in a 70 GB index](http://genome.cshlp.org/content/26/12/1721.full) with Centrifuge (which would be better / more comprehensive than just a ribosomal screen).
  
  
  ---
  
  ### Answer by username_3 on 2017-09-15. Score: 4
  
  One thing that hasn't been suggested so far would be phylogenetic placement. First option there would be placing against one larger tree roughly touching the domains you mentioned to get an overview. Second option would be having three trees (bacteria, plants, viruses) made up of relevant taxa for the type of sample you have.
  
  This would give you a detailed and accurate view of your sample, but it does involve some work. Sadly the computational throughput for placement is still a bit behind BLAST, and not many complete pipelines exist, but in my opinion it is worth the extra mile (I might be biased. I'm working to resolve as many of these issues as possible).
  
  If you want to look into it, I can reccomend [this excellent paper by Erick Matsen](https://academic.oup.com/sysbio/article/64/1/e26/2847641). Since your samples sound vaguely soil-ish, I will also link [this recent Nature paper](https://www.nature.com/articles/s41559-017-0091) using placement to analyze rainforest soil samples.
  
  
  ---
  
  ### Answer by username_1 on 2017-08-01. Score: 4
  
  I think you're just looking for `nr`. It is absolutely not limited to prokaryotes, far from it. It is, according to the blurb on the NCBI's blast page:
  
  > The nucleotide collection consists of GenBank+EMBL+DDBJ+PDB+RefSeq sequences, but excludes EST, STS, GSS, WGS, TSA, patent sequences as well as phase 0, 1, and 2 HTGS sequences. The database is non-redundant. Identical sequences have been merged into one entry, while preserving the accession, GI, title and taxonomy information for each entry.
  
  The link you gave in your question doesn't describe the `nr` database, but the collection of RefSeq non-redundant sequences. The two are not the same thing at all. The nr databases contains sequences from all available species so it should be perfect for what you're after.
  
  Now, I admit I haven't tried this with millions of redas, but I did used to run many blast searches with several thousand sequences using the command line blast tools and the remote nr database. I just checked the [blast+ docs](https://www.ncbi.nlm.nih.gov/books/NBK279690/) now and see no mention of a rate limit. Which isn't to say there isn't one, for sure, but since each query is run sequentially when you give a multifasta file, there seems to reason to limit since the servers should be able to queue the jobs as needed.
  
  So, don't use the API or roll your own, just [download NCBI's blast client](https://blast.ncbi.nlm.nih.gov/Blast.cgi?CMD=Web&PAGE_TYPE=BlastDocs&DOC_TYPE=Download) and try that:
  
  ```
  blastn -query query.fa -db nr -task blastn 
  
  ```
  
  If that complains, you'll have to download the `nr` DB (well, since you want nucleotide blast, you probably want its nucleotide cousin, `nt`) which you can find here: <ftp://ftp.ncbi.nlm.nih.gov/blast/db/>. Just get all the `foo.nt*` files.
  
  > **Comments:**
  > **username_0**: According to the [BLAST API Usage guidelines](https://blast.ncbi.nlm.nih.gov/Blast.cgi?CMD=Web&PAGE_TYPE=BlastDocs&DOC_TYPE=DeveloperInfo), this isn’t really supported, even if there is no hard limit. But you seem to be right about NR, this seems to be what I’m after, after all.  
  > **username_1**: @username_0 I'm hoping that refers to using the API in your own tools and not to using their own standalone blast+ program. I may very well be wrong, but I'm not sure if the instructions for the API apply to them as well. I would expect them to have built their tool in such a way as to not overload the servers. I figure it's at least worth a try.  
  
  ---
  
  ### Answer by username_4 on 2017-08-03. Score: 4
  
  You might want to look into [minihashes](https://en.wikipedia.org/wiki/MinHash) as a more efficient way of clustering and comparing your sequences to known species and samples. [In the mash paper](https://genomebiology.biomedcentral.com/articles/10.1186/s13059-016-0997-x) the authors show how it would be feasible to create minihash sketches that you could search against from large metagenomics datasets (check the "Clustering massive metagenomic datasets" section in the results):
  
  > For a large-scale test, samples from the Human Microbiome Project [36]
  > (HMP) and Metagenomics of the Human Intestinal Tract [37] (MetaHIT)
  > were combined to create a ~10 TB 888-sample dataset. Importantly, the
  > size of a Mash sketch is independent of the input size, requiring only
  > 70 MB to store the combined sketches (s = 10,000, k = 21) for these
  > datasets.
  
  Moreover, the creators of [sourmash](https://sourmash.readthedocs.io/en/latest/) (a tool using the same approach) provide [precomputed sketches for the whole of NCBI's RefSeq database](https://sourmash.readthedocs.io/en/latest/databases.html) (includes bacterial, viral and fungal genomes). Creating a similar sketch for the whole of ENA/SRA (or a subset of it) could require some time but would allow for very fast and accurate searches.
  
  
  ---
  
  ### Answer by username_6 on 2017-08-03. Score: 1
  
  I have had good success using the MEGAN software to detect and visualize contamination. My input data was WGS not RNA-seq But I would check it out. It can create a very nice weight phylogenetic tree for all the detected organisms.
  
  <https://ab.inf.uni-tuebingen.de/software/megan6>
  
  [![enter image description here](https://i.sstatic.net/ttT6l.png)](https://i.sstatic.net/ttT6l.png)
  
  > **Comments:**
  > **username_0**: I've downloaded MEGAN but my (admittedly cursory) reading of the documentation made it seem as if the software only visualised results computed by other tools.  
  
  ---
  
  