-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathgenmut
More file actions
executable file
·49 lines (42 loc) · 1.53 KB
/
genmut
File metadata and controls
executable file
·49 lines (42 loc) · 1.53 KB
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
#!/usr/bin/env python3
import argparse
import os
from readfasta import read_record
parser = argparse.ArgumentParser(description='Simulate splice-site-generating point mutaitons')
parser.add_argument('fasta', type=str, metavar='<fasta>', help='input path to fasta file')
parser.add_argument('-o', type=str, metavar='<out_dir>', required=False, help='output dir')
arg = parser.parse_args()
for idn, seq in read_record(arg.fasta):
file_name = os.path.basename(arg.fasta)
base_name = os.path.splitext(file_name)[0]
OUT_DIR='.'
if arg.o: OUT_DIR = arg.o
os.makedirs(f'{OUT_DIR}/{base_name}', exist_ok=True)
count = 1
for i in range(0, len(seq)-1, 2):
# donor site
cur_id = None
cur_seq = None
if seq[i] == 'G' and seq[i+1] != 'T':
cur_id = f'{idn} {i+1} d {seq[i+1]}~T'
cur_seq = f'{seq[:i+1]}T{seq[i+2:]}'
if seq[i] != 'G' and seq[i+1] == 'T':
cur_id = f'{idn} {i} d {seq[i]}~G'
cur_seq = f'{seq[:i]}G{seq[i+1:]}'
# acceptor site
if seq[i] == 'A' and seq[i+1] != 'G':
cur_id = f'{idn} {i+1} a {seq[i+1]}~G'
cur_seq = f'{seq[:i+1]}G{seq[i+2:]}'
if seq[i] != 'A' and seq[i+1] == 'G':
cur_id = f'{idn} {i} a {seq[i]}~A'
cur_seq = f'{seq[:i]}A{seq[i+1:]}'
# out
if cur_id is not None:
with open(f'{OUT_DIR}/{base_name}/mut_{count}_{base_name}.fa', 'w') as fh:
fh.write(f'>{cur_id}\n')
for i in range(0, len(cur_seq), 80):
fh.write(cur_seq[i:i+80])
count+=1
#print(f'>{cur_id}')
#for i in range(0, len(cur_seq), 80):
#print(cur_seq[i:i+80])