#! /usr/bin/env python # numSaturated.py # P.Clote # Computes the number of saturated sec str, assuming # every base can base pair with every other base. Used for # combinatorics study of saturated secondary structures. import sys,math THRESHOLD = 1 def maxNumBasePairs(n): return (n-THRESHOLD)/2 def saturated(n): #D[(m,k)] is number of loc opt str on m having k base pairs D = {} for m in range(THRESHOLD+2): D[(m,0)] = 1L for k in range(1,maxNumBasePairs(m)+1): D[(m,k)] = 0L for m in range(THRESHOLD+2,n+1): D[(m,0)] = 0L #once m>threshold+1, there are no loc opt with 0 bp for k in range(1,maxNumBasePairs(m)+1): sum = 0L if k-1 <= maxNumBasePairs(m-2): sum += D[(m-2,k-1)] #Case 1: (1,m) is base pair for i in range(m-THRESHOLD-1,m): #Case 2: (1,m-1) OR (1,m-2) OR ... OR (1,m-THRESHOLD-1) is a base pair #Note that for imaxNumBasePairs(x): continue sum += D[(x,k-1)] for i in range(2,m-THRESHOLD): #Case 3: (i,m) is base pair for k0 in range(0,maxNumBasePairs(i-1)+1): #k0 many base pairs in sequence 1,...,i-1 k1 = k-k0-1 #want k0+k1+1 = k many base pairs if k1<0 or k1>maxNumBasePairs(m-i-1): continue sum += D[(i-1,k0)]*D[(m-i-1,k1)] D[(m,k)] = sum return D def main(n): D = saturated(n) for k in range(maxNumBasePairs(n)+1): if D[(n,k)]>0: print "%d\t%s" % (k,D[(n,k)]) print for k in range(maxNumBasePairs(n),-1,-1): k0 = maxNumBasePairs(n)-k if D[(n,k)]>0: print "%d\t%s" % (k0,D[(n,k)]) if __name__ == '__main__': if len(sys.argv)!= 2: print "Usage: %s integer" % sys.argv[0] sys.exit(1) main(int(sys.argv[1]))