Problems Using Formatdb And Fastacmd

Question

[Sakti]

Hello!!

I am having an extremely difficult time with formatdb and fastcmd from NCBI blast.

I am trying to make a db from single mouse chromosome from NCBI using formatdb. Afterwards I try to obtain sequence from defined positions of the db. I use the commands:

formatdb -p F -i 12.fa -o T -t C57Bl6_chr12 -n C57Bl6_chr12
fastacmd -d C57Bl6_chr12 -p T -L 1000000,10001000 -o test.txt

[fastacmd] ERROR: ERROR: Cannot initialize readdb for C57Bl6_chr12 database

It is very annoying that the formatdb is not making the readdb for my file, do you know what I can do about it??

Thanks!!

[Sakti]

Oh god, you are right, I am so sorry to have spammed the site with such a silly thing... very much appreciated Neilws!!!! And about retrieving sequences with fastacmd, well, it is at least fastest than my previous experiences with SeqIO when managing huge data volumes, but once again, it might be just the way I progam it... thanks !!!!

[Istvan Albert]

we all make silly mistakes all the time - in fact the simpler the problem the more stupefying it is

[Neilfws]

Yes, we've all experienced hours of frustration only to realise that the error is a simple typo or wrong switch.

Answer 1

[Neilfws]

The error has a simple cause: you are running fastacmd with the option "-p T", which means "assume sequence is protein." I note that you used the correct option, "-p F" when you ran formatdb.

When I run formatdb on a chromosome 12 sequence downloaded from NCBI and I look at formatdb.log, I see that 7 sequences were formatted, corresponding to 7 contigs. If this is the case for you, do not be surprised when fastacmd generates an empty output file.

Finally - there are easier ways to extract sub-sequences from larger sequences (using the original FASTA file) than dumping from a BLAST database.