-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathmksequence.py
More file actions
95 lines (89 loc) · 3.2 KB
/
mksequence.py
File metadata and controls
95 lines (89 loc) · 3.2 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
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
import numpy as np
import random as rm
import os
from progressbar import ProgressBar
pbar = ProgressBar()
os.system('cls')
states = ['T', 'C', 'A', 'G']
transitionName = [['TT', 'TC', 'TA', 'TG'], ['CT', 'CC', 'CA', 'CG'],
['AT', 'AC', 'AA', 'AG'], ['GT', 'GC', 'GA', 'GG']]
years = int(input('Enter the number of years to simulate the sequence over:'))
T = []
C = []
A = []
G = []
transitionMatrix = [T, C, A, G]
lam = 7.3333 * 10 ** -10
for base in states:
if base == 'T':
T.append(1 - (3 * lam))
else:
T.append(lam)
for base in states:
if base == 'C':
C.append(1 - (3 * lam))
else:
C.append(lam)
for base in states:
if base == 'A':
A.append(1 - (3 * lam))
else:
A.append(lam)
for base in states:
if base == 'G':
G.append(1 - (3 * lam))
else:
G.append(lam)
testSequence = 'CTTCCAGCCCGCGGACCGATGCTGCCGGCGGCCGCTCGCCCCCTGTGGGGGCCTTGCCTTGGGCTTCGGGCCGCTGCGTTCCGCCTTGCCAGGCGACAGGTGCCATGTGTCTGTGCCGTGCGACATATGAGGAGCAGCGGCCATCAGAGGTGTGAGGCCCTCGCTGGTGCAACGCCCCCAAGGAGTACCCCCCCAAGATACAGCAGCTGGTCCAGGACATACTCTCTTGGAAATCTCAGACCTCAACGAGCTCCTGAAGAAAACGTTGAAGTCGGGCTTGTGCCGATGGGTGGTGTGATGTCTGGGGCTGTCCCTGCTGCGAGGCGGTGGAAGAAGATATCCCCATAGCGAAAGAACGGACACATTTCACACCGAGGCGAAGCCCGTGGACAAAGTGAAGCTGATCAAGGAAATCAAGAAGGCATCAACCTCGTCCAGGCAAAGAAGCTGGTGGAGTCCCTGCCCCAGGAAATGTCGCCAAAGCTGAGGCGGAGAAGATCAAGGCGGCCCTGGAGGCGGTGTGGTTCTGGAGTAGCCTCCAGCTCGGAGGACTTGTGTTCAGGGGTCCTGGGCCCCGGG'
# genbank ascesion #E01991
newSequence = []
for base in pbar(testSequence):
currentBase = base
i = 0
while i != years:
if currentBase == 'T':
transition = np.random.choice(transitionName[0], replace=True, p=transitionMatrix[0])
if transition == 'TT':
pass
elif transition == 'TC':
currentBase = 'C'
elif transition == 'TA':
currentBase = 'A'
else:
currentBase = 'G'
elif currentBase == 'C':
transition = np.random.choice(transitionName[1], replace=True, p=transitionMatrix[1])
if transition == 'CT':
currentBase = 'T'
elif transition == 'CC':
pass
elif transition == 'CA':
currentBase = 'A'
else:
currentBase = 'G'
elif currentBase == 'A':
transition = np.random.choice(transitionName[2], replace=True, p=transitionMatrix[2])
if transition == 'AT':
currentBase = 'T'
elif transition == 'AC':
currentBase = 'C'
elif transition == 'AA':
pass
else:
currentBase = 'G'
else:
transition = np.random.choice(transitionName[3], replace=True, p=transitionMatrix[3])
if transition == 'GT':
currentBase = 'T'
elif transition == 'GC':
currentBase = 'C'
elif transition == 'GA':
currentBase = 'A'
else:
pass
i += 1
newSequence.append(currentBase)
x = sum(i != j for i, j in zip(testSequence, newSequence))
print(x)
print(testSequence)
print(''.join(newSequence))