Running blast command line inside python script

Question

[flogin]

Hello guys, I'm writing a python script to automatize several analysis that I need to do, one of those are BLAST analyses, so, in this way I create these structure, where the user put the file with queries and the number of threads, but when I try to run, an error message appears: [blastall] ERROR: Number of processors to use [threads] is bad or out of range [? to ?]

import os, argparse

parser = argparse.ArgumentParser('Automatizating BLAST analaysis...')

parser.add_argument("-in", "--input", help="genome in fasta file format",  required=True)
parser.add_argument("-th", "--threads", help="number of threads",  required=True)
args = parser.parse_args()

read_file = args.input
threads = int(args.threads)
blast_output = args.input+'.tblastn'

first_blast = 'blastall -p tblastn -d ./database.fasta -i read_file -e 0.01 -o blast_output -a threads -m 8 -M BLOSUM45 -W 2 -t 30'
os.system(first_blast)

Can anyone explain to me this error? I'm wondering if the variable threads are parsing in an incorrect manner, but I doesn't have much experience with chmod commands inside python scripts...

Answer 1

[Mensur Dlakic]

The way you wrote this, threads is not interpreted as a variable. It is just a string inside quotation marks that is interpreted literally.

Try this:

first_blast = 'blastall -p tblastn -d ./database.fasta -i read_file -e 0.01 -o blast_output -a ' + str(threads) + ' -m 8 -M BLOSUM45 -W 2 -t 30'

[finswimmer]

I think it's the same for read_file and blast_output. I would also recommend to use f-strings (available since python 3.6) to format the string:

first_blast = f"blastall -p tblastn -d ./database.fasta -i {read_file} -e 0.01 -o {blast_output} -a {threads} -m 8 -M BLOSUM45 -W 2 -t 30"

[flogin]

Yeah, this is a good way, but I was using an institutional server with python 3.5 :/

[flogin]

Thanks for the advice, I rewrite my code, but now this message error appears: [blastall] ERROR: Arguments must start with '-' (the offending argument #6 was: '0.01')

[Mensur Dlakic]

That means the error is somewhere around -e 0.01. Maybe, as finswimmer suggested, you also need to take the read_file and blast_output variables out of quotation marks?

first_blast = 'blastall -p tblastn -d ./database.fasta -i ' + str(read_file) + ' -e 0.01 -o ' + str(blast_output) + ' -a ' + str(threads) + ' -m 8 -M BLOSUM45 -W 2 -t 30'