Conditional Probability- Python

user_11981324 0
Viewed 889

I'm working on this python problem:

Given a sequence of the DNA bases {A, C, G, T}, stored as a string, returns a conditional probability table in a data structure such that one base (b1) can be looked up, and then a second (b2), to get the probability p(b2 | b1) of the second base occurring immediately after the first. (Assumes the length of seq is >= 3, and that the probability of any b1 and b2 which have never been seen together is 0. Ignores the probability that b1 will be followed by the end of the string.)

You may use the collections module, but no other libraries.

However I'm running into a roadblock:

word = 'ATCGATTGAGCTCTAGCG'

def dna_prob2(seq):
    tbl = dict()
    levels = set(word)
    freq = dict.fromkeys(levels, 0)
    for i in seq:
        freq[i] += 1
    for i in levels:
        tbl[i] = {x:0 for x in levels}
    lastlevel = ''
    for i in tbl:
        if lastlevel != '':
             tbl[lastlevel][i] += 1
        lastlevel = i
    for i in tbl:
        print(i,tbl[i][i] / freq[i])
    return tbl

tbl['T']['T'] / freq[i]

Basically, the end result is supposed to be the final line tbl you see above. However, when I try to do that in print(i,tbl[i][i] /freq[i), and run dna_prob2(word), I get 0.0s for everything.

Wondering if anyone here can help out.

Thanks!

probability python
2 Answers

I am not sure what it is your code is doing, but this works:

def makeprobs(word):
  singles = {}
  probs = {}
  thedict={}
  ll = len(word)
  for i in range(ll-1):
    x1 = word[i]
    x2 = word[i+1]
    singles[x1] = singles.get(x1, 0)+1.0
    thedict[(x1, x2)] = thedict.get((x1, x2), 0)+1.0
  for i in thedict:
    probs[i] = thedict[i]/singles[i[0]]
  return probs
user_29772560

I finally got back to my professor. This is what it was trying to accomplish:

word = 'ATCGATTGAGCTCTAGCG'

def dna_prob2(seq):
    tbl = dict()
    levels = set(seq)     
    freq = dict.fromkeys(levels, 0)
    for i in seq:
        freq[i] += 1
    for i in levels:
        tbl[i] = {x:0 for x in levels}
    lastlevel = ''
    for i in seq:            
        if lastlevel != '':
             tbl[lastlevel][i] += 1
        lastlevel = i
    return tbl, freq          

condfreq, freq = dna_prob2(word)  
print(condfreq['T']['T']/freq['T']) 
print(condfreq['G']['A']/freq['A']) 
print(condfreq['C']['G']/freq['G'])

Hope this helps.

user_119813240
Related