test/unit/bio/sequence/test_common.rb in bio-1.4.2 vs test/unit/bio/sequence/test_common.rb in bio-1.4.3

@@ -237,10 +237,14 @@
       assert(rseqs[0] != rseqs[1])
       assert(rseqs[0] != rseqs[2])
       assert(rseqs[1] != rseqs[2])
     end
 
+  end #class TestSequenceCommonRandomize
+
+  class TestSequenceCommonRandomizeChi2 < Test::Unit::TestCase
+
     def chi2(hist, f, k)
       chi2 = 0
       (0...k).each do |i|
         chi2 += ((hist[i] - f) ** 2).quo(f)
       end
@@ -248,11 +252,11 @@
     end
     private :chi2
 
     # chi-square test for distribution of chi2 values from
     # distribution of index('a')
-    def randomize_equiprobability
+    def randomize_equiprobability_chi2
       # Reference: http://www.geocities.jp/m_hiroi/light/pystat04.html
       seq = TSequence.new('ccccgggtta') # length must be 10
       k = 10
       hist = Array.new(k, 0)
       iter = 200
@@ -278,14 +282,42 @@
 
       chi2_f = chi2_iter.quo(k).to_f
       chi2_chi2 = chi2(chi2_hist, chi2_f, k)
       #$stderr.puts chi2_chi2
 
+      chi2_chi2
+    end
+    private :randomize_equiprobability_chi2
+
+    def randomize_equiprobability(&block)
       ## chi-square test, freedom 9, significance level 5%
-      #assert_operator(16.919, :>, chi2_chi2, "test of chi2 < 16.919 failed (#{chi2_chi2})")
+      #critical_value = 16.919
+      #significance_level_message = "5%"
+      #
       # chi-square test, freedom 9, significance level 1%
-      assert_operator(21.666, :>, chi2_chi2, "test of chi2 < 21.666 failed (#{chi2_chi2})")
+      critical_value = 21.666
+      significance_level_message = "1%"
+
+      # max trial times till the test sucess
+      max_trial = 10
+
+      values =[]
+      max_trial.times do |i|
+        chi2_chi2 = randomize_equiprobability_chi2(&block)
+        values.push chi2_chi2
+        # immediately breaks if the test succeeds
+        break if chi2_chi2 < critical_value
+        $stderr.print "Bio::Sequence::Common#randomize test of chi2 (=#{chi2_chi2}) < #{critical_value} failed (expected #{significance_level_message} by chance)"
+        if values.size < max_trial then
+          $stderr.puts ", retrying."
+        else
+          $stderr.puts " #{values.size} consecutive times!"
+        end
+      end
+
+      assert_operator(values[-1], :<, critical_value,
+                      "test of chi2 < #{critical_value} failed #{values.size} times consecutively (#{values.inspect})")
     end
     private :randomize_equiprobability
 
     def test_randomize_equiprobability
       randomize_equiprobability { |seq| seq.randomize }
@@ -317,10 +349,10 @@
     #    end
     #    newseq
     #  end
     #end
 
-  end
+  end #class TestSequenceCommonRandomizeChi2
 
 
   class TestSequenceCommonSubseq < Test::Unit::TestCase
     #def subseq(s = 1, e = self.length)