Changeset 203

Timestamp:

Nov 24, 2006, 1:27:16 PM (16 years ago)

Author:

Johan Enell

Message:

Fixes #33

File:

• trunk/se/lu/onk/ZTest/src/ztest/ZTest.java (modified) (3 diffs)

trunk/se/lu/onk/ZTest/src/ztest/ZTest.java

@@ -24 +24 @@
 package ztest;
 
+import org.apache.commons.math.MathException;
+import org.apache.commons.math.special.Erf;
+
 import java.util.HashMap;
 import java.util.List;
@@ -31 +34 @@
 public class ZTest
 {
-  private int precision = (int) Math.pow(10, 6);
-  
-  private int depth = 6;
-  
-  private double[] normsdistV;
-  
   public void calculate(BASEFileSpotSection<Reporter, Spot> bfss)
   {
-    normsdistV = new double[(precision * depth) + 1];
     List<Integer> assays = bfss.findFieldIntList("assays");
     HashMap<Integer, Double> means = getAssayMeans(bfss, assays);
     HashMap<Integer, Double> sds = getAssaySDS(bfss, assays, means);
-    
-    fillNormsdist();
     
     for (int i = 0; i < bfss.getReporterSize(); i++)
@@ -87 +81 @@
   private double normsdist(double z)
   {
-    int zz = (int) Math.abs(Math.round(z * precision));
-    double value;
-    if (z < 0) throw new IllegalArgumentException("normsdist cant handle a negative z score");
-    else if (z <= (depth - 1)) value = normsdistV[zz];
-    else 
+    double value = 0;
+    try
     {
-      value = 1;
+      if (z < -26) value = 0;
+      else if (z > 26) value = 1;
+      else value = 0.5 + Erf.erf(z/Math.sqrt(2))/2;
+    }
+    catch (MathException e)
+    {
+      System.err.println("This should never happen but if it do... The plug-in is crap.");
+      System.err.println(z);
+      e.printStackTrace();
+      System.exit(0);
     }
     return value;
   }
   
-  private double distribution(double t)
-  {
-    return Math.pow(Math.E, -0.5 * Math.pow(t, 2)) / Math.sqrt(2 * Math.PI);
-  }
-
   private double z(double x, double popMean, double popSD)
   {
     return (x - popMean) / popSD;
-  }
-  
-  private void fillNormsdist()
-  {
-    normsdistV[0] = 0.5;
-    for (int i = 0; i < normsdistV.length-1; i++)
-    {
-      double w = (double)i/precision;
-      normsdistV[i+1] = normsdistV[i] + distribution(w)/precision;
-    }
   }