tree_script that takes into account all multiple sites of the same types

sun13005 · Apr 19, 2017 · 3173c3e · 3173c3e
1 parent 12f0bdf
commit 3173c3e
Showing 1 changed file with 29 additions and 15 deletions.
diff --git a/tree_test.py b/tree_test.py
@@ -25,19 +25,23 @@ def __init__(self, ndnum):              # initialization function
 
     def allocatePartial(self, patterns):
         # skip internals for now
-        print 'self.number, self.edgelen =',self.number, self.edgelen
+#         print 'self.number, self.edgelen =',self.number, self.edgelen
 
         if self.number > 0:
         # patterns['ACGTT'] = 51
             npatterns = len(patterns)
             self.partial = [0.0]*(4*npatterns)
 
-            print 'self.partial=', self.partial
+#             print 'self.partial=', self.partial
 
             for i,pattern in enumerate(patterns.keys()):
+#                 print 'patterns=', patterns
+#                 print 'pattern=', pattern
+#                 print patterns[pattern]
+
 
                 base = pattern[self.number-1]
-                print 'i = %d, pattern = %s, base = %s' % (i,pattern,base)
+#                 print 'i = %d, pattern = %s, base = %s' % (i,pattern,base)
                 if base == 'A':
                     self.partial[i*4 + 0] = 1.0
                 elif base == 'C':
@@ -50,15 +54,18 @@ def allocatePartial(self, patterns):
                     assert(False), 'oops, something went horribly wrong!'
 
 
-            print self.partial
-            print '****************************'    
+#             print self.partial
+#             print '****************************************************************************************************************'    
 
         else:
             npatterns = len(patterns)
-            print 'npatterns=', npatterns
+#             print 'npatterns=', npatterns
             self.partial = [0.0]*(4*npatterns)
-
+            like_list = []
             for i,pattern in enumerate(patterns.keys()):
+#                 print 'patterns=', patterns
+#                 print 'pattern=', pattern
+                num_pattern = patterns[pattern]
 #                 print 'edge_length==', self.lchild.edgelen  
 #                 print 'test1=', self.lchild.partial[i*4 + 0], self.lchild.rsib.partial[i*4 + 0]
 #                 print 'test2=', self.lchild.partial[i*4 + 1], self.lchild.rsib.partial[i*4 + 1]
@@ -68,7 +75,6 @@ def allocatePartial(self, patterns):
 #                 print 'test5=', self.lchild.edgelen
 #                 print 'test6=', self.lchild.rsib.edgelen
 #                 print'____________________'
-
                 pAA = (0.25+0.75*exp(-4.0*(self.lchild.edgelen)/3.0))*(self.lchild.partial[i*4 + 0])
                 pAC = (0.25-0.25*exp(-4.0*(self.lchild.edgelen)/3.0))*(self.lchild.partial[i*4 + 1])
                 pAG = (0.25-0.25*exp(-4.0*(self.lchild.edgelen)/3.0))*(self.lchild.partial[i*4 + 2])
@@ -128,16 +134,24 @@ def allocatePartial(self, patterns):
                 pfromT_rchild = pTA2+pTC2+pTG2+pTT2
                 self.partial[i*4 + 3] = pfromT_lchild*pfromT_rchild
 ##############################################################################################################################
-            site1= (sum(self.partial[0:4]))*0.25
-            site2= (sum(self.partial[4:8]))*0.25
-            Like = log(site1)+log(site2)
+
+#                 print '~~~~~~~~~~~~~~~~~~~~'
+#                 print 'num_pattern=', num_pattern
+#                 print self.partial
+#                 print self.partial[i*4:i*4+4]
+#                 site_like2 = sum(self.partial[i:i+4])*num_pattern
+                site_like = (log((sum(self.partial[i*4:i*4+4]))*0.25))*num_pattern
+#                 print 'site_like=', site_like
+                like_list.append(site_like)
+#                 print '~~~~~~~~~~~~~~~~~~~~'
 
-            print self.partial
-            print 'site1=', site1
-            print 'site2=', site2
+
+
+            Like = sum(like_list)
             print 'Like=', Like
 
-            print '****************************'    
+            print '**************************************************'    
+
     def __str__(self):
         # __str__ is a built-in function that is used by print to show an object
         descendants_as_string = ','.join(['%d' % d for d in self.descendants])