DGEN2.jar reinstantiates GT probability class for case found of bug when used for many taxa

Found internally that gene tree probabilities were not being calculated correctly in a case with many taxa. This was due to not reinstantiating the class that calculates the gene tree probabilities before passing in a network with new branch lengths. The class requires reinstantiating because of its typical use case where it saves things in an internal state.

Author: chilleo

CommandLineFiles/DoWholePaper.py

@@ -1,12 +1,10 @@
 from RunDGEN import *
 
-i = 0.01
-j = 'cat'
-
-h = j/i
-
 k = 0
 
+#8 tax - sites ignored? (was just using this line for tsting, prolly shouldnt rerun)
+#run_saved_dgen('/Users/leo/rice/res/data/dgen/tmp/figBig8/stat_8tax.txt', ['/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim1/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim2/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim3/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim4/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim5/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim6/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim7/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim8/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim9/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim10/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim11/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim12/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim13/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim14/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim15/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim16/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim17/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim18/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim19/seqfileNamed','/Users/leo/rice/res/data/dgen/simulations/12taxa/bigNetwork/firstTry/sim20/seqfileNamed'], verbose=True, plot='/Users/leo/rice/res/data/dgen/tmp/figBig8/plot_big_8', meta='Dgen')
+
 # cichlid alignment - i put the commands for these in the java code
 #run_saved_dgen('/Users/leo/rice/res/data/dgen/tmp/figC/stat_cichlid.txt', ['/Users/leo/rice/res/data/cichlid/alignment/cichlid6tax.phylip-sequential.txt'], verbose=False, plot='/Users/leo/rice/res/data/dgen/tmp/figC/plot_figC', meta='Dgen')
 run_saved_dgen('/Users/leo/rice/res/data/dgen/tmp/figC/stat_cichlid.txt', ['/Users/leo/rice/res/data/cichlid/alignment/cichlid6tax.phylip-sequential.txt'], verbose=True, plot='/Users/leo/rice/res/data/dgen/tmp/figC/plot_figCVerbose', meta='Dgen')

CommandLineFiles/RunDGEN.py

@@ -283,6 +283,13 @@ def calculate_windows_to_DGEN(alignments, taxa_order, outgroup, list_of_tree_and
                 elif bases.issubset(possibleBases) == False:
                     num_ignored += 1
 
+                #final catch all else for ignored sites (forgot to add back in check for non biallelic (they were ignored but not counted in verbose output. this line now properly counts them for verbose mode))
+                #includes also sites that only have all A's or whatever, basically anything non strictly biallelic
+                else:
+                    num_ignored += 1
+
+                #should i add count here for sites that specifically violate biallelic as in 3 or 4 dif letters? (len(bases)>2)?
+
                 # Increment the site index
                 site_idx += 1
 
@@ -297,6 +304,9 @@ def calculate_windows_to_DGEN(alignments, taxa_order, outgroup, list_of_tree_and
                     lineAdd.append(groupCount)
                 list_of_tree_and_net_invariants_counts.append(lineAdd)
 
+            #for debugging, how many sites are actually observed that go into the observed calculations?
+
+
             # calculate new version of chi squared test
             # chi squared = sum of (obs - exp)^2 / exp
             # obs is Na and exp is Nt * |Na| / |Nt|