I need to calculate the occurrences of a motif (including overlaps) in sequences (motif is passed in the first line of standard input and FASTA sequences in subsequent lines). The sequence name starts with >, and after whitespace is just a comment about the sequence that needs to be neglected. The input of program is like:
AT
>seq1 Comment......
AGGTATA
TGGCGCC
>seq2 Comment.....
GGCCGGCGC
GCCGGCGGG
The output should be:
seq1: 2
seq2: 0
I decided to save the first line as a motif, strip the comment from sequence name, join lines of sequence in one line and save sequence names (keys) and sequences (values) in a dictionary. I also wrote a function for motif_count and want to call it on dictionary values and then save it in a separate dictionary for final output. However, it does not work the way I want it as the method motif_count does not produce any result when called on my d values.
Here is my code:
#!/usr/bin/env python3
import sys
sequence = sys.stdin.readlines()
d = {}
temp_genename = None
temp_sequence = None
motif = sequence[0]
def motif_count(motif, temp_sequence):
count = 0
start = -1
while True:
start = temp_sequence.find(motif, start +1)
if start <0:
return count
count += 1
if sequence[1][0] != '>':
print("ERROR")
exit(1)
for line in sequence[1:]:
if line[0] == '>':
temp_genename = line.split(' ')[0].strip()
tmp_genename = temp_genename.replace('>', '')
temp_sequence = ""
else:
temp_sequence += line.strip()
d[tmp_genename] = temp_sequence
for key in d:
motif_count(motif, key)
motif_count
, thewhile
loop does not contain the test that can trigger thereturn
, nor the incrementation ofcount
.exit(1)
should also probably be re-indented. I tried to fix this. Hopefully this corresponds to your actual code now. $\endgroup$motif_count
. To create a new dictionary with values modified by a function, thevalmap
function from thecytoolz
module could be useful. $\endgroup$