BioInformaticsAlgorithms/RandomizedMotifSearch.py at master · Andrewwu73/BioInformaticsAlgorithms · GitHub

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
def scoreMotifs(scores):
    consensus = ""
    count = 0
    for l in range(len(scores[0])):
        tempString= ""
        for a in scores:
            tempString = tempString + a[l]
        a = tempString.count("A")
        c = tempString.count("C")
        g = tempString.count("G")
        t = tempString.count("T")
        y = max(a, c, g, t)
        if(y==a):
            consensus = consensus + "A"
        elif(y==c):
            consensus = consensus + "C"
        elif(y==g):
            consensus = consensus + "G"
        elif(y==t):
            consensus = consensus + "T"
    for m in range(len(scores[0])):
        for xd in scores:
            if(xd[m]!=consensus[m]):
                count = count + 1

    return count
bestMotifs = []
length, num = map(int, input("").split(" "))
dna = []
for k in range(num):
    dna.append(input(""))
import random
for holy in range(2000):
    motifs = []
    for m in range(num):
        xd = random.randint(0, len(dna[m])-length)
        motifs.append(dna[m][xd:xd+length])
    if(len(bestMotifs)==0):
        bestMotifs = motifs
    while True:
        profile = [[], [], [], []]
        for xd in range(length):
            profile[0].append(1)
            profile[1].append(1)
            profile[2].append(1)
            profile[3].append(1)
        for motif in motifs:
            for l in range(length):
                if(motif[l]=="A"):
                    profile[0][l] = profile[0][l]+1
                elif(motif[l]=="C"):
                    profile[1][l]=profile[1][l] + 1
                elif(motif[l]=="G"):
                    profile[2][l] = profile[2][l] + 1
                elif(motif[l]=="T"):
                    profile[3][l] = profile[3][l] + 1
        motifs = []
        for j in range(num):
            mostProb = ""
            prob = -1
            for haha in range(len(dna[j])-length+1):
                temp = 1
                for y in range(length):
                    if(dna[j][haha+y]=="A"):
                        temp = temp * profile[0][y]
                    elif(dna[j][haha+y]=="C"):
                        temp = temp * profile[1][y]
                    elif(dna[j][haha+y]=="G"):
                        temp = temp * profile[2][y]
                    elif(dna[j][haha+y]=="T"):
                        temp = temp * profile[3][y]
                if(temp>prob):
                    prob = temp
                    mostProb = dna[j][haha:haha+length]
            motifs.append(mostProb)
        if(scoreMotifs(motifs)<scoreMotifs(bestMotifs)):
            bestMotifs = motifs
        else:
            break

for m in bestMotifs:
    print(m)