-
Notifications
You must be signed in to change notification settings - Fork 4
/
extract_contigs.py
executable file
·66 lines (53 loc) · 1.99 KB
/
extract_contigs.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
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
#! /usr/bin/env python
# -*- coding: utf-8 -*-
import sys
import os
import re
from Bio import SeqIO
from Bio.SeqUtils import GC
from optparse import OptionParser
parser = OptionParser(usage='''%prog [-i] [-o] [-l] [-r] [-h]
*************************************************************************
Extracts contigs from a fasta file using a contig ID list (header part).
Possible to extract all from list or not in list (reverse -r).
-i input (fasta file)
-l list (ID list)
-o output
BioPython needed!
http://biopython.org/wiki/Main_Page
contact: [email protected]
*************************************************************************
''', version="%prog 1.0")
parser.add_option("-i", "--input", dest="InputFile", type='string', help='specify file to be searched', default=None)
parser.add_option("-l", "--list", dest="ListFile", type='string', help='search list',default=None)
parser.add_option("-o", "--outfile", dest="OutputFile", type='string', help='specify output file',default=None)
parser.add_option("-r", "--reverse", dest="reverse", action="store_true", help='reverse hit', default=False)
(options, args) = parser.parse_args()
if options.OutputFile is None:
parser.error("no output file given")
if options.InputFile is None:
parser.error("no input file given")
if os.path.exists(options.OutputFile):
sys.stderr.write('\nFile exists!\n\n')
overwrite=input('overwrite? (y/n): ')
if overwrite != 'y':
sys.exit()
#here, (with open('file', 'w' as handle) you do not need to close the outfile
sys.stderr.write('\nGoing through file...\n\n')
wanted = set()
with open(options.ListFile) as f:
for line in f:
line = line.strip()
if line != "":
wanted.add(line)
fasta_sequences = SeqIO.parse(open(options.InputFile),'fasta')
end = False
with open(options.OutputFile, "w") as f:
for seq in fasta_sequences:
if options.reverse == True:
if seq.id not in wanted:
SeqIO.write([seq], f, "fasta")
else:
if seq.id in wanted:
SeqIO.write([seq], f, "fasta")
sys.stderr.write('\nDone\n\n')