Entering edit mode
8.8 years ago
Nano
▴
20
I have made a Local Blast of one query file (containing thousands of sequences) against many individual database formatted files (each containing thousands of sequences). The result is good. Script below
#!/bin/bash
for file in *.fasta
do
pb tblastn -query queryFile.fasta -dbnuc $file -out $file\.table -evalue 0.00001 -outfmt 6
done
The results in each file look like this:
BGIBMGA000034 c9147_g1_i1 33.95 162 97 2 42 193 120 605 1e-22 93.2
BGIBMGA000034 c13162_g1_i2 51.35 74 36 0 42 115 300 79 4e-20 84.0
BGIBMGA000034 c15905_g1_i1 36.00 125 79 1 44 167 493 119 4e-19 82.0
BGIBMGA000034 c15528_g1_i1 29.70 165 104 1 42 194 241 735 1e-15 75.1
BGIBMGA000034 c13162_g1_i3 42.65 68 39 0 129 196 2104 1901 2e-08 52.8
BGIBMGA000060 c15177_g1_i1 44.85 689 269 10 51 715 1423 3228 2e-177 561
BGIBMGA000060 c15177_g1_i1 37.74 106 61 1 286 386 1369 1686 2e-12 70.9
BGIBMGA000060 c15177_g1_i3 44.70 689 266 11 55 715 1435 3240 8e-176 556
BGIBMGA000060 c15177_g1_i3 37.74 106 61 1 286 386 1369 1686 2e-12 71.2
BGIBMGA000060 c15177_g1_i2 55.99 384 160 2 341 715 502 1653 3e-136 440
BGIBMGA000060 c15821_g1_i2 30.14 282 174 8 454 714 1955 2794 8e-24 108
The first columns are query IDs while the second column are matched IDs from database files, the values are normal Blast values.
Question: How can I get the code to also output the full sequence (including the header >xxx|yyy
) of the matching database sequence as part of the result.
For example (output):
BGIBMGA000034 c9147_g1_i1 33.95 162 97 2 42 193 120 605 1e-22 93.2
>c9147_g1_i1
AGAGGACCATTGCACCGGGAAATTCCTCCCAATTGGGGTTCCAAAGGGGCCCCTTAAGGGAAGGAAGGTTTCCRCCGGCCTTTTTTCCCTTGGGGAAAGGGGTTTCCCCCAAGGTTCCAACCAATTGGGGGCCTTCTTTAAGGG........
Thanks
You can't.
You need to run blastdbcmd to retrieve your sequences. Write the id's and the corresponding start and stop positions to a file and send that file as a query to blastdbcmd.
If I'm not mistaken I think the format of the query should be:
Accession_id_1 50 410
Accession_id_2 113 460