1 """
   2 nussinov.py 
   3 NussinovAlgorithm code python 
   4 
   5 input.txt:
   6 GGGAAAUCG
   7 
   8 c:\workplace> python nussinov.py <input.txt >output.txt
   9 
  10 output.txt (which will be created in the same directory):
  11 (2,9)
  12 (3,8)
  13 (4,7)
  14 """
  15 """nussinov"""
  16 
  17 class NussinovMatrix:
  18     def __init__(self, rnaseq):
  19         self.rnaseq = rnaseq
  20         self.L = len(rnaseq)
  21         self.matrix = []
  22         for i in range(self.L):
  23             self.matrix.append(range(self.L))
  24         for i in range(self.L):
  25 	    for j in range(self.L):
  26 	        self.matrix[i][j] = 'n'
  27 	self.fillInitialMatrix()
  28   	
  29     def fillInitialMatrix(self):   	       	
  30 	for i in range(1, self.L):
  31 	    self.matrix[i][i-1] = 0
  32 	for i in range(self.L):
  33 	    self.matrix[i][i] = 0
  34     	
  35     def fillMatrix(self):
  36         endOfSeq = self.L + 1
  37         for oneDiagonal in range(1, endOfSeq):
  38             i = 1
  39 	    for oneCell in range(1, endOfSeq):
  40 	        j = i + oneDiagonal
  41 		if j > self.L: break
  42 		else:
  43 		    self.writeValue(i, j)
  44 		    i += 1
  45 	    endOfSeq -= 1
  46 
  47     def writeValue(self, i, j):
  48     	unpairi = self.getValue(i+1, j)
  49 	unpairj = self.getValue(i, j-1)
  50 	pair = self.getValue(i+1, j-1) + self.score(i, j)
  51 	
  52 	ks = range(i+1, j)
  53 	if ks == []: bifurcation = 0
  54 	else:
  55 	    bifurcations = []
  56 	    for k in ks:
  57 	        case = self.getValue(i, k) + self.getValue(k+1, j)
  58 	        bifurcations.append(case)
  59     	    bifurcation = max(bifurcations) 
  60 	result = max(pair, unpairi, unpairj, bifurcation)
  61 	self.setValue(i, j, result)
  62 
  63     def score(self, ai, aj):
  64         i = self.rnaseq[ai-1]
  65 	j = self.rnaseq[aj-1]
  66         RNAComplementary = {'A':'U', 'U':'A', 'G':'C', 'C':'G'}
  67         return RNAComplementary[i] == j
  68 
  69     def getInitialMatrix(self):
  70     	return self.matrix
  71     def getFinalMatrix(self):
  72     	self.fillMatrix()
  73 	return self.matrix
  74     def getValue(self, i, j):
  75         return self.matrix[i-1][j-1]
  76     def setValue(self, i, j, aValue):
  77     	self.matrix[i-1][j-1] = aValue
  78 
  79 class Nussinov(NussinovMatrix):
  80     def __init__(self, aRnaseq):
  81         NussinovMatrix.__init__(self, aRnaseq)
  82         matrix = self.getFinalMatrix()
  83 	self.stack = []
  84 	self.stack.append((1, self.L))
  85 	self.basepair = []
  86 	self.iterateTraceback()
  87 
  88     def iterateTraceback(self):
  89         while 1:
  90             try: popedItem = self.stack.pop()
  91 	    except IndexError: break
  92 	    self.doRecursion(popedItem)
  93 
  94     def doRecursion(self, popedItem):
  95 	i, j = popedItem
  96 	if i >= j: pass
  97 	elif self.getValue(i+1, j) == self.getValue(i,j):
  98 	    self.stack.append((i+1, j))
  99 	elif self.getValue(i, j-1) == self.getValue(i,j):
 100 	    self.stack.append((i, j-1))
 101 	elif self.getValue(i+1, j-1) + self.score(i,j) == self.getValue(i,j):
 102 	    self.basepair.append((i,j))
 103 	    self.stack.append((i+1, j-1))
 104 	else:
 105 	    ks = range(i+1, j)
 106 	    for k in ks:
 107 	        if self.getValue(i,k) + self.getValue(k+1,j) == self.getValue(i,j):
 108 		    self.stack.append((k+1,j))
 109 		    self.stack.append((i,k))
 110 		    break
 111     	
 112     def getBasePair(self):
 113         return self.basepair
 114     	
 115 def writeResult():
 116     inputRnaTotal = ''
 117     while 1:
 118         try: inputRna=raw_input()
 119         except EOFError: break
 120         if not inputRna: break
 121         inputRnaTotal += inputRna.strip() 
 122 
 123     a = Nussinov(inputRnaTotal)
 124     for each in a.getBasePair():
 125         print each
 126     	
 127 if __name__=='__main__':
 128     writeResult()
web biohackers.net