-
Notifications
You must be signed in to change notification settings - Fork 0
/
GC.py
71 lines (50 loc) · 1.89 KB
/
GC.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
from sys import argv
from sys import exit
from itertools import islice
from re import findall
# 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
read(data)
# Keeps track of highest GC count
hfreq = 0
hseq = ""
# Counts amount of GC bases per sequence and records highest
for sequence in dataset:
freq = len(findall(r'[G C]', dataset[sequence]))
freq = freq * 100 / len(dataset[sequence])
if freq > hfreq:
hfreq = freq
hseq = sequence
# Formats highest sequence for writing
result = hseq + "\n" + str(hfreq) + "%"
# Writes highest sequence to file
with open(argv[2], 'w') as file:
file.writelines(result)
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
record(beginning, i + 1, label, file)
return
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
return
if __name__ == "__main__":
main()