1
$\begingroup$

I have a set of 100 amino acid sequences and I want to perform a BLASTP sesrch against the refseq_protein database. Accordingly I had set up the standalone version of BLAST (Version 2.11.0+) and downloaded the refseq_protein database from NCBI using the following code

wget ftp://ftp.ncbi.nlm.nih.gov/refseq/release/complete/*.faa.gz

The database gets downloaded as 3027 zipped files containing FASTA sequences. I unzipped all these files and concatenated them into a single file refseq_protein.faa (which is around 95 GB in size). Now when I run the following Python code

from Bio.Blast.Applications import NcbimakeblastdbCommandline
from Bio.Blast.Applications import NcbiblastpCommandline

cline = NcbimakeblastdbCommandline(dbtype = "prot", input_file = "D:\\refseq_protein.faa", out 
= "refseq_protein")

blastp_cline = NcbiblastpCommandline(query = "D:\\DEP_sequences.fasta', db = 
"refseq_protein", evalue = 0.01, outfmt = "7 sseqid evalue qcovs pident")

cline()

response = blastp_cline()

the NcbimakeblastdbCommandline function keeps creating multiple .phr, .pin, .psq etc files which take up a lot of space (In a demo run it had created ~30GB of these files and was still running). I'm afraid this will exhaust the entire space available on my internal hard drive. I'm wondering if there is a way to estimate the total size of the files which NcbimakeblastdbCommandline would create. This will help me in deciding whether or not to switch to an external storage to perform the BLASTP search.

I am aware of the fact that pre-formatted refseq_protein database exists but I'm not sure what value is to be passed in the db parameter of the NcbiblastpCommandline function, because it asks for the name of the database against which the BLASTP search is to be performed. In the approach that I chose, I had the liberty to set the name of the database.

Any suggestions on how to solve this issue would be appreciated.

$\endgroup$

1 Answer 1

2
$\begingroup$

Welcome to the Bioinformatics site! Hopefully you find the below helpful.

Regarding the creation of a RefSeq protein database locally...as you mentioned, there are pre-made RefSeq BLAST databases available, and looking at the website, the current DB (as of 12/27/2021) is roughly 62.9GB compressed - and is probably ~2x larger decompressed, so I’d guess that the database you are constructing manually will be at least 125GB in size - and likely more.

The above is regarding the technical issue you are currently facing - but I'd back up and look closer at what you are trying to do - align the 100 amino acid sequences to RefSeq. Here are some alternatives that may be -much- easier.

  1. BLAST your sequences on the NCBI website. This may not be as satisfying as building out your local tools and databases, but may get the job done more efficiently. This of course depends on whether this is a one-off thing, or something you plan to do frequently, in which case a local database may make sense.
  2. Download the premade RefSeq BLAST database. You mentioned you were having issues with how to refer to this database (which may be a useful question in and of itself), but generally an alias file exists within the premade database (entitled refseq_protein.pal), and calling this file will allow you to BLAST all the separate chunks. So, the command line might be something like this:
blastp -db refseq_protein -query amino_acid_100.fa

There are a ton more options you may want to use (try blastp -help) , but the above should work.

  1. Use the -remote option in blastp to run from the command line against the NCBI server. The below line should perform the BLAST alignment against the NCBI refseq protein database.
blastp -remote -db refseq_protein -query amino_acid_100.fa
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.