-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathPcirc_unmapped_reads_blast.py
39 lines (31 loc) · 1.24 KB
/
Pcirc_unmapped_reads_blast.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
import argparse
import os
parser = argparse.ArgumentParser(description='Unmapped reads alignment')
parser.add_argument("-g",
"--genome",
help='reference genome files are Fasta')
parser.add_argument("-q",
"--query",
help='reference query files are Fasta')
parser.add_argument("-o",
"--output_file",
help='Output_file saved reasult')
parser.add_argument("-e",
"--evalue",
default=1e-5,
help='Expectation value (E) threshold for saving hits')
parser.add_argument("-p",
"--num_threads",
default=4)
args = parser.parse_args()
evalue = args.evalue
num_threads = args.num_threads
output_file = args.output_file
query_filename = args.query
genome_filename = args.genome
genome_path = os.path.realpath(genome_filename)
genome_database = genome_path[:genome_path.rfind('/')] + '/genome_database'
os.system("makeblastdb -in %s -dbtype nucl -out %s" %
(genome_filename, genome_database))
os.system("blastn -db %s -query %s -out %s -outfmt 6 -evalue %f -num_thread %d" %
(genome_database, query_filename, output_file, evalue, num_thread))