Changeset 197


Ignore:
Timestamp:
Nov 21, 2006, 1:53:52 PM (15 years ago)
Author:
Johan Enell
Message:

callculates with good acc

Location:
trunk/se/lu/onk/ZTest/src/ztest
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/se/lu/onk/ZTest/src/ztest/Reporter.java

    r187 r197  
    3030  private double m;
    3131  private double p;
     32  private int size;
    3233 
    3334 
     
    6162    return m;
    6263  }
     64
     65  public void setSize(int size)
     66  {
     67    this.size = size;
     68  }
     69
     70  public int getSize()
     71  {
     72    return size;
     73  }
    6374}
  • trunk/se/lu/onk/ZTest/src/ztest/Spot.java

    r187 r197  
    2626public class Spot
    2727{
    28   float m;
     28  double m;
    2929 
    30   public Spot(float m)
     30  public Spot(Double m)
    3131  {
    3232    this.m = m;
    3333  }
    3434
    35   public final float getM()
     35  public final double getM()
    3636  {
    3737    return m;
    3838  }
    3939
     40  @Override
     41  public String toString()
     42  {
     43    return "["+m+"]";
     44  }
     45
    4046 
    4147}
  • trunk/se/lu/onk/ZTest/src/ztest/Start.java

    r187 r197  
    7373        {
    7474          int index = mCol + i * assayFields.size();
    75           float m = Float.NaN;
     75          Double m = Double.NaN;
    7676          if (!data[index].equals(""))
    7777          {
    78             m = new Float(data[index]);
     78            m = new Double(data[index]);
    7979          }
    8080          spots[i] = new Spot(m);
     
    130130      out.println("<meta http-equiv=\"Content-Type\" content=\"text/html; charset=ISO-8859-1\" />");
    131131      out.println("<title>BASE 1.2.16dev - Experiment Enell:Enell Plugin</title>");
     132      out.println("<link rel=\"stylesheet\" type=\"text/css\" href=\"base.css\" />");
    132133      out.println("</head>");
    133134      out.println("<table>");
    134       out.println("<tr><th>Reporter</th><th>M</th><th>p</th><th>Rank</th><th>Number of expected</th><th>Ratio between expected and observed</th></tr>");
     135      out.println("<tr><th>Reporter</th><th>Average M</th><th>P</th><th>Number of values</th><th>Rank</th<th>Expected</th><th>FDR</th></tr>");
     136     
    135137      for (int i = 0; i < bfss.getReporterSize() && i < 100; i++)
    136138      {
     139        int rank = i + 1;
    137140        Reporter r = bfss.getReporter(i);
    138         int numExpected = (int) (r.getP() * numOfReporters);
    139         double ratio = numExpected / (i + 1.0);
    140 
    141         out.print("<tr>");
     141        double numExpected = r.getP() * numOfReporters;
     142        double fdr = (numExpected / rank);
     143       
     144        out.print("<tr class=\"row"+(i%2+1)+"\">");
    142145        out.printf(
    143           "<td>%1$d</td><td>%2$g</td><td>%3$g</td><td>%4$d</td><td>%5$d</td><td>%6$g</td>",
     146          "<td>%d</td>" +
     147          "<td>%f</td>" +
     148          "<td>%f</td>" +
     149          "<td>%d</td>" +
     150          "<td>%d</td>" +
     151          "<td>%f</td>" +
     152          "<td>%.2f</td>",
    144153          r.getReporter(),
    145154          r.getM(),
    146155          r.getP(),
    147           i + 1,
     156          r.getSize(),
     157          rank,
    148158          numExpected,
    149           ratio);
     159          fdr);
    150160        out.println("</tr>");
    151161      }
  • trunk/se/lu/onk/ZTest/src/ztest/ZTest.java

    r187 r197  
    2424package ztest;
    2525
    26 import java.util.Arrays;
    2726import java.util.HashMap;
    2827import java.util.List;
     
    3231public class ZTest
    3332{
    34 
     33  private int precision = (int) Math.pow(10, 6);
     34 
     35  private int depth = 6;
     36 
     37  private double[] normsdistV;
     38 
    3539  public void calculate(BASEFileSpotSection<Reporter, Spot> bfss)
    3640  {
     41    normsdistV = new double[(precision * depth) + 1];
    3742    List<Integer> assays = bfss.findFieldIntList("assays");
    38     HashMap<Integer, Float> means = getAssayMeans(bfss, assays);
    39     HashMap<Integer, Float> sds = getAssaySDS(bfss, assays, means);
    40 
     43    HashMap<Integer, Double> means = getAssayMeans(bfss, assays);
     44    HashMap<Integer, Double> sds = getAssaySDS(bfss, assays, means);
     45   
     46    fillNormsdist();
     47   
    4148    for (int i = 0; i < bfss.getReporterSize(); i++)
    4249    {
     
    4855      for (int j = 0; j < assays.size(); j++)
    4956      {
    50         double m = bfss.getSpot(i + j).getM();
     57        double m = bfss.getSpot((i*assays.size()) + j).getM();
    5158        if (!Double.isNaN(m))
    5259        {
     
    5764        }
    5865      }
    59       popSD = popSD / popSize;
     66      popSD = Math.sqrt(popSD) / popSize;
    6067      popMean = popMean / popSize;
    6168      meanM = meanM / popSize;
    6269
    6370      double z = z(meanM, popMean, popSD);
     71     
    6472      double normsdist = normsdist(Math.abs(z));
    6573      double p = p(normsdist);
    66      
    67       bfss.getReporter(i).setP(p);
    68       bfss.getReporter(i).setM(meanM);
     74
     75      Reporter r = bfss.getReporter(i);
     76      r.setP(p);
     77      r.setM(meanM);
     78      r.setSize(popSize);
    6979    }
    7080  }
     
    7282  private double p(double normsdist)
    7383  {
    74     return 2 * (1 / normsdist);
     84    return 2 * (1 - normsdist);
    7585  }
    7686
    7787  private double normsdist(double z)
    7888  {
    79     return 1 / Math.sqrt(2 * Math.PI) * Math.pow(Math.E, -0.5 * Math.pow(z, 2));
     89    int zz = (int) Math.abs(Math.round(z * precision));
     90    double value;
     91    if (z < 0) throw new IllegalArgumentException("normsdist cant handle a negative z score");
     92    else if (z <= (depth - 1)) value = normsdistV[zz];
     93    else
     94    {
     95      value = 1;
     96    }
     97    return value;
     98  }
     99 
     100  private double distribution(double t)
     101  {
     102    return Math.pow(Math.E, -0.5 * Math.pow(t, 2)) / Math.sqrt(2 * Math.PI);
    80103  }
    81104
     
    84107    return (x - popMean) / popSD;
    85108  }
     109 
     110  private void fillNormsdist()
     111  {
     112    normsdistV[0] = 0.5;
     113    for (int i = 0; i < normsdistV.length-1; i++)
     114    {
     115      double w = (double)i/precision;
     116      normsdistV[i+1] = normsdistV[i] + distribution(w)/precision;
     117    }
     118  }
    86119
    87   private HashMap<Integer, Float> getAssaySDS(BASEFileSpotSection<Reporter, Spot> bfss, List<Integer> assays, HashMap<Integer, Float> means)
     120  private HashMap<Integer, Double> getAssaySDS(BASEFileSpotSection<Reporter, Spot> bfss, List<Integer> assays, HashMap<Integer, Double> means)
    88121  {
    89     int[] size = new int[assays.size()];
    90     float[] sqrSum = new float[assays.size()];
    91     Arrays.fill(sqrSum, 0);
    92 
    93     int count = bfss.findIntOpt("count");
    94     for (int i = 0; i < count; i += assays.size())
     122    HashMap<Integer, Double> sds = new HashMap<Integer, Double>();
     123    for (int i = 0; i < assays.size(); i++)
    95124    {
    96       for (int j = 0; j < assays.size(); j++)
     125      int size = 0;
     126      double sqrSum = 0;
     127      for (int j = 0; j < bfss.getSpotSize(); j += assays.size())
    97128      {
    98         float m = bfss.getSpot(i + j).getM();
    99         if (!Float.isNaN(m))
     129        double m = bfss.getSpot(i + j).getM();
     130        if (!Double.isNaN(m))
    100131        {
    101           sqrSum[j] += Math.pow(m - means.get(assays.get(j)), 2);
    102           size[j]++;
     132          sqrSum += Math.pow(m - means.get(assays.get(i)), 2);
     133          size++;
    103134        }
    104135      }
    105     }
    106     HashMap<Integer, Float> sds = new HashMap<Integer, Float>();
    107     for (int i = 0; i < assays.size(); i++)
    108     {
    109       sds.put(assays.get(i), sqrSum[i] / size[i]);
     136      sds.put(assays.get(i), Math.sqrt(sqrSum / size));
    110137    }
    111138
     
    113140  }
    114141
    115   private HashMap<Integer, Float> getAssayMeans(BASEFileSpotSection<Reporter, Spot> bfss, List<Integer> assays)
     142  private HashMap<Integer, Double> getAssayMeans(BASEFileSpotSection<Reporter, Spot> bfss, List<Integer> assays)
    116143  {
    117     int[] size = new int[assays.size()];
    118     float[] sum = new float[assays.size()];
    119     Arrays.fill(sum, 0);
    120 
    121     int count = bfss.findIntOpt("count");
    122     for (int i = 0; i < count; i += assays.size())
     144    HashMap<Integer, Double> means = new HashMap<Integer, Double>();
     145    for (int i = 0; i < assays.size(); i++)
    123146    {
    124       for (int j = 0; j < assays.size(); j++)
     147      int size = 0;
     148      double sum = 0;
     149      for (int j = 0; j < bfss.getSpotSize(); j += assays.size())
    125150      {
    126         float m = bfss.getSpot(i + j).getM();
    127         if (!Float.isNaN(m))
     151        double m = bfss.getSpot(i + j).getM();
     152        if (!Double.isNaN(m))
    128153        {
    129           sum[j] += m;
    130           size[j]++;
     154          sum += m;
     155          size++;
    131156        }
    132157      }
     158      means.put(assays.get(i), sum / size);
    133159    }
    134     HashMap<Integer, Float> means = new HashMap<Integer, Float>();
    135     for (int i = 0; i < assays.size(); i++)
    136     {
    137       means.put(assays.get(i), sum[i] / size[i]);
    138     }
    139 
    140160    return means;
    141161  }
Note: See TracChangeset for help on using the changeset viewer.