-
Notifications
You must be signed in to change notification settings - Fork 0
/
CONS.py
83 lines (60 loc) · 2.35 KB
/
CONS.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
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
from sys import argv
from sys import exit
from itertools import islice
import numpy as np
import time
start_time = time.time()
# Stores labels and sequences
dataset = {}
def main():
# Ensures proper command line arguments
if len(argv) != 3:
exit("Please use proper input")
# Reads file into string
with open(argv[1], 'r') as file:
data = file.readlines()
# Organizes data into dict dataset
length = read(data)
# Creates profile matrix
profile = np.zeros([4, length])
# Dict for matching bases and their matrix row
letters = {'A': 0, 'C': 1, 'G': 2, 'T': 3}
# Adds base frequency to profile matrix
for sequence in dataset:
for i in range(length):
base = dataset[sequence][i]
profile[letters[base]][i] += 1
# Finds matrix row containing the max in each column
intcstring = np.argmax(profile, 0)
# For converting row from number to letter
strcstring = ""
# Reverses dictionary for converting numbers to letters
letters = dict(map(reversed, letters.items()))
# Converts numbers in intcstring to letters
for number in intcstring:
strcstring += letters[number]
# Writes result to file
with open(argv[2], 'w') as file:
file.write(strcstring)
print("--- %s seconds ---" % (time.time() - start_time))
def read(file):
beginning = 0
label = None # For first iteration
for i, line in enumerate(file):
if line[0] == '>': # If new sequence is found
record(beginning, i, label, file) # Add the previous sequence to dataset
label = line[1:].rstrip() # Record label of current sequence
beginning = i + 1 # Record beginning of current sequence
elif i == len(file) - 1: # If end of file is reached
length = record(beginning, i + 1, label, file)
return length
def record(beginning, end, label, file):
if label is None: # First iteration contains no data
return
dataset[label] = "" # So future values can be concatenated
for line in islice(file, beginning, end, 1): # Read only lines in sequence
dataset[label] += line.rstrip() # Add these lines to the dataset
length = len(dataset[label])
return length
if __name__ == "__main__":
main()