Filtering BLAST results.

Question

[bpz]

Hello. I need some advice about how I can handle some BLAST results.

I have a set of Illumina reads obtained from a metagenomic sample. I need to obtain from these reads the ones belonging to a certain organism, while eliminating others. I prepared a BLAST database from the genomes (mitogenomes) of the organism I am interested in, plus 3 others I need to "screen-out". The query for the BLAST were the Illumina reads.

I have the results (in tabular output) and I need to keep the reads that are a match to my problem organism, and eliminate the ones that are also a (better) match to the other 3.

Any ideas on how can I achieve this? Thanks for your help.

[Biomonika (Noolean)]

Show a header of your table.

[bpz]

HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       NC_008332.1                              100.00  89      0       0       1       89      146573  146485  1e-41    165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       gi|40794996|gb|AY506529.1|      100.00  89      0       0       1       89      137429  137341  1e-41    165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       gi|93116077|gb|DQ490952.1|      100.00  89      0       0       1       89      642229  642141  1e-41    165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       gi|93116122|gb|DQ490953.1|      100.00  89      0       0       1       89      53028   52940   1e-41      165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       NC_007982.1                                100.00  89      0       0       1       89      137429  137341  1e-41    165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       gi|102567891|gb|DQ645538.1|     100.00  89      0       0       1       89      71502   71590   1e-41      165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       gi|102567891|gb|DQ645538.1|     100.00  89      0       0       1       89      246653  246565  1e-41    165
HWI-ST1379:164:c5k28acxx:3:1101:2770:2209       gi|110287586|gb|DQ645537.1|     100.00  89      0       0       1       89      164556  164644  1e-41    165

[bpz]

Thanks a lot for your answers. I think BBsplit is the best tool for my problem. I will give it a try.

Answer 1

[felipehcoutinho]

Assuming you have BLAST output files in the standard order, the identifier of the queries (names of illumina reads) that matched the sequences in your database should be in the first column. What I would do is write a script (e.g. using perl or python) to run through your BLAST table(s) and keep track of the query identifiers that match your mitogenomes, check if these queries have matches with a higher bitscore (last column of the BLAST output) to the the genomes you wanna screen out, and finally go through the file with all the illumina reads and filter only the ones you are interested in. You also should decide what you wanna do with the reads that match none of the genomes in the database. Hope this helps.

Answer 2

[h.mon]

It seems BBSplit does exactly what you want, binning reads according to which reference they map.

If you intend to bin (possibly millions of) raw Illumina reads, Blast is not the best tool, it will be too slow.

Answer 3

[Josh Herr]

This should be pretty straight forward as you should have your output in tabular format. I do this quite frequently.

1. Sort on hits to your target genomes and output a list of the FASTA sequence headers.
2. Then I use a tool (the one I use is seqtk, but there are others) to input a list of these headers to subset the sequences out of your FASTA file into a new smaller FASTA file.

From there, you should be good -- let us know if you're still having problems selecting your BLAST hits.