mev-4.0.01/source/org/tigr/microarray/mev/cluster/algorithm/impl/OneWayANOVA.java

Code
Comments
Other
Rev Date Author Line
2 26 Feb 07 jari 1 /*
2 26 Feb 07 jari 2 Copyright @ 1999-2003, The Institute for Genomic Research (TIGR).
2 26 Feb 07 jari 3 All rights reserved.
2 26 Feb 07 jari 4 */
2 26 Feb 07 jari 5 /*
2 26 Feb 07 jari 6  * $RCSfile: OneWayANOVA.java,v $
2 26 Feb 07 jari 7  * $Revision: 1.6 $
2 26 Feb 07 jari 8  * $Date: 2005/03/10 15:45:21 $
2 26 Feb 07 jari 9  * $Author: braistedj $
2 26 Feb 07 jari 10  * $State: Exp $
2 26 Feb 07 jari 11  */
2 26 Feb 07 jari 12
2 26 Feb 07 jari 13 package org.tigr.microarray.mev.cluster.algorithm.impl;
2 26 Feb 07 jari 14
2 26 Feb 07 jari 15 import java.util.Random;
2 26 Feb 07 jari 16 import java.util.Vector;
2 26 Feb 07 jari 17
2 26 Feb 07 jari 18 import org.tigr.microarray.mev.cluster.Cluster;
2 26 Feb 07 jari 19 import org.tigr.microarray.mev.cluster.Node;
2 26 Feb 07 jari 20 import org.tigr.microarray.mev.cluster.NodeList;
2 26 Feb 07 jari 21 import org.tigr.microarray.mev.cluster.NodeValue;
2 26 Feb 07 jari 22 import org.tigr.microarray.mev.cluster.NodeValueList;
2 26 Feb 07 jari 23 import org.tigr.microarray.mev.cluster.algorithm.AbortException;
2 26 Feb 07 jari 24 import org.tigr.microarray.mev.cluster.algorithm.AbstractAlgorithm;
2 26 Feb 07 jari 25 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmData;
2 26 Feb 07 jari 26 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmEvent;
2 26 Feb 07 jari 27 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmException;
2 26 Feb 07 jari 28 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmParameters;
2 26 Feb 07 jari 29 import org.tigr.microarray.mev.cluster.gui.impl.owa.OneWayANOVAInitBox;
2 26 Feb 07 jari 30 import org.tigr.util.FloatMatrix;
2 26 Feb 07 jari 31 import org.tigr.util.QSort;
2 26 Feb 07 jari 32
2 26 Feb 07 jari 33 import JSci.maths.statistics.FDistribution;
2 26 Feb 07 jari 34
2 26 Feb 07 jari 35 /**
2 26 Feb 07 jari 36  *
2 26 Feb 07 jari 37  * @author  nbhagaba
2 26 Feb 07 jari 38  * @version 
2 26 Feb 07 jari 39  */
2 26 Feb 07 jari 40 public class OneWayANOVA extends AbstractAlgorithm {
2 26 Feb 07 jari 41     
2 26 Feb 07 jari 42     public static final int FALSE_NUM = 12;
2 26 Feb 07 jari 43     public static final int FALSE_PROP = 13;    
2 26 Feb 07 jari 44     
2 26 Feb 07 jari 45     private boolean stop = false;
2 26 Feb 07 jari 46     
2 26 Feb 07 jari 47     private int function;
2 26 Feb 07 jari 48     private float factor;
2 26 Feb 07 jari 49     private boolean absolute, calculateAdjFDRPVals;
2 26 Feb 07 jari 50     private FloatMatrix expMatrix;
2 26 Feb 07 jari 51     
2 26 Feb 07 jari 52     private Vector[] clusters;
2 26 Feb 07 jari 53     private int k; // # of clusters
2 26 Feb 07 jari 54     
2 26 Feb 07 jari 55     private int numGenes, numExps, numGroups;
2 26 Feb 07 jari 56     private float alpha, falseProp;
2 26 Feb 07 jari 57     private boolean usePerms, drawSigTreesOnly;
2 26 Feb 07 jari 58     private int numPerms, falseNum;
2 26 Feb 07 jari 59     private int correctionMethod;    
2 26 Feb 07 jari 60     int[] groupAssignments; 
2 26 Feb 07 jari 61     private double[] origPVals;
2 26 Feb 07 jari 62     
2 26 Feb 07 jari 63     float currentP = 0.0f;
2 26 Feb 07 jari 64     //float currentF = 0.0f;
2 26 Feb 07 jari 65     int currentIndex = 0; 
2 26 Feb 07 jari 66     double constant;
2 26 Feb 07 jari 67     
2 26 Feb 07 jari 68     AlgorithmEvent event;
2 26 Feb 07 jari 69     
2 26 Feb 07 jari 70     Vector fValuesVector = new Vector();
2 26 Feb 07 jari 71     Vector rawPValuesVector = new Vector();
2 26 Feb 07 jari 72     Vector adjPValuesVector = new Vector();
2 26 Feb 07 jari 73     //Vector pValuesVector = new Vector(); 
2 26 Feb 07 jari 74     Vector dfNumVector = new Vector();
2 26 Feb 07 jari 75     Vector dfDenomVector = new Vector();
2 26 Feb 07 jari 76     Vector ssGroupsVector = new Vector();
2 26 Feb 07 jari 77     Vector ssErrorVector = new Vector();
2 26 Feb 07 jari 78     private boolean[] isSig;
2 26 Feb 07 jari 79     
2 26 Feb 07 jari 80     private int hcl_function;
2 26 Feb 07 jari 81     private boolean hcl_absolute;
2 26 Feb 07 jari 82     
2 26 Feb 07 jari 83     private boolean useFastFDRApprox = true;
2 26 Feb 07 jari 84     
2 26 Feb 07 jari 85     /**
2 26 Feb 07 jari 86      * This method should interrupt the calculation.
2 26 Feb 07 jari 87      */
2 26 Feb 07 jari 88     public void abort() {
2 26 Feb 07 jari 89         stop = true;        
2 26 Feb 07 jari 90     }
2 26 Feb 07 jari 91     
2 26 Feb 07 jari 92     /**
2 26 Feb 07 jari 93      * This method execute calculation and return result,
2 26 Feb 07 jari 94      * stored in <code>AlgorithmData</code> class.
2 26 Feb 07 jari 95      *
2 26 Feb 07 jari 96      * @param data the data to be calculated.
2 26 Feb 07 jari 97      */
2 26 Feb 07 jari 98     public AlgorithmData execute(AlgorithmData data) throws AlgorithmException {
2 26 Feb 07 jari 99   groupAssignments = data.getIntArray("group-assignments");
2 26 Feb 07 jari 100   
2 26 Feb 07 jari 101   AlgorithmParameters map = data.getParams();
2 26 Feb 07 jari 102   function = map.getInt("distance-function", EUCLIDEAN);
2 26 Feb 07 jari 103   factor   = map.getFloat("distance-factor", 1.0f);
2 26 Feb 07 jari 104   absolute = map.getBoolean("distance-absolute", false);
2 26 Feb 07 jari 105         usePerms = map.getBoolean("usePerms", false);
2 26 Feb 07 jari 106         numPerms = map.getInt("numPerms", 0);
2 26 Feb 07 jari 107
2 26 Feb 07 jari 108         hcl_function = map.getInt("hcl-distance-function", EUCLIDEAN);
2 26 Feb 07 jari 109         hcl_absolute = map.getBoolean("hcl-distance-absolute", false);        
2 26 Feb 07 jari 110         
2 26 Feb 07 jari 111   boolean hierarchical_tree = map.getBoolean("hierarchical-tree", false);
2 26 Feb 07 jari 112         if (hierarchical_tree) {
2 26 Feb 07 jari 113             drawSigTreesOnly = map.getBoolean("draw-sig-trees-only");
2 26 Feb 07 jari 114         }        
2 26 Feb 07 jari 115   int method_linkage = map.getInt("method-linkage", 0);
2 26 Feb 07 jari 116   boolean calculate_genes = map.getBoolean("calculate-genes", false);
2 26 Feb 07 jari 117   boolean calculate_experiments = map.getBoolean("calculate-experiments", false);
2 26 Feb 07 jari 118   
2 26 Feb 07 jari 119   this.expMatrix = data.getMatrix("experiment");
2 26 Feb 07 jari 120   
2 26 Feb 07 jari 121   numGenes = this.expMatrix.getRowDimension();
2 26 Feb 07 jari 122   numExps = this.expMatrix.getColumnDimension();
2 26 Feb 07 jari 123   alpha = map.getFloat("alpha", 0.01f);
2 26 Feb 07 jari 124   correctionMethod = map.getInt("correction-method", OneWayANOVAInitBox.JUST_ALPHA);        
2 26 Feb 07 jari 125         numGroups = map.getInt("numGroups", 3);
2 26 Feb 07 jari 126         
2 26 Feb 07 jari 127         calculateAdjFDRPVals = false;
2 26 Feb 07 jari 128         
2 26 Feb 07 jari 129         if (correctionMethod == FALSE_NUM) {
2 26 Feb 07 jari 130             falseNum = map.getInt("falseNum", 10);
2 26 Feb 07 jari 131         }
2 26 Feb 07 jari 132         
2 26 Feb 07 jari 133         if (correctionMethod == FALSE_PROP) {
2 26 Feb 07 jari 134             falseProp = map.getFloat("falseProp", 0.05f);
2 26 Feb 07 jari 135         }        
2 26 Feb 07 jari 136         
2 26 Feb 07 jari 137         getFDfSSValues();
2 26 Feb 07 jari 138         if ((correctionMethod == FALSE_NUM) || (correctionMethod == FALSE_PROP)) {
2 26 Feb 07 jari 139             rawPValuesVector = getRawPValuesFromFDist();
2 26 Feb 07 jari 140             origPVals = new double[rawPValuesVector.size()];
2 26 Feb 07 jari 141             for (int i = 0; i < origPVals.length; i++) {
2 26 Feb 07 jari 142                 origPVals[i] = ((Float)(rawPValuesVector.get(i))).doubleValue();
2 26 Feb 07 jari 143             }
2 26 Feb 07 jari 144             adjPValuesVector = new Vector();
2 26 Feb 07 jari 145             for (int i = 0; i < origPVals.length; i++) {
2 26 Feb 07 jari 146                 adjPValuesVector.add(new Float(0f));
2 26 Feb 07 jari 147             }
2 26 Feb 07 jari 148         }  else {        
2 26 Feb 07 jari 149             if (!usePerms) {
2 26 Feb 07 jari 150                 rawPValuesVector = getRawPValuesFromFDist();
2 26 Feb 07 jari 151             } else {
2 26 Feb 07 jari 152                 rawPValuesVector = getRawPValsFromPerms();
2 26 Feb 07 jari 153             }
2 26 Feb 07 jari 154             
2 26 Feb 07 jari 155         }
2 26 Feb 07 jari 156         
2 26 Feb 07 jari 157         
2 26 Feb 07 jari 158         //Vector clusterVector = sortGenesBySignificance();
2 26 Feb 07 jari 159         Vector clusterVector = new Vector();
2 26 Feb 07 jari 160         Vector sigGenes = new Vector();
2 26 Feb 07 jari 161         Vector nonSigGenes = new Vector();
2 26 Feb 07 jari 162         
2 26 Feb 07 jari 163         if ( (correctionMethod == OneWayANOVA.FALSE_NUM) || (correctionMethod == OneWayANOVA.FALSE_PROP) ) {
2 26 Feb 07 jari 164             boolean[] isGeneSig = new boolean[1];
2 26 Feb 07 jari 165             if (correctionMethod == OneWayANOVA.FALSE_NUM) {
2 26 Feb 07 jari 166                 isGeneSig = isGeneSigByFDRNum();
2 26 Feb 07 jari 167             } else {
2 26 Feb 07 jari 168                 isGeneSig = isGeneSigByFDRPropNew2();
2 26 Feb 07 jari 169             }
2 26 Feb 07 jari 170             /*            
2 26 Feb 07 jari 171             if (!calculateAdjFDRPVals) {
2 26 Feb 07 jari 172                 adjustedPVals = new double[numGenes];
2 26 Feb 07 jari 173             }
2 26 Feb 07 jari 174              */
2 26 Feb 07 jari 175             //System.out.println("isGeneSig.length = " + isGeneSig.length);
2 26 Feb 07 jari 176             for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 177                 if (isGeneSig[i]) {
2 26 Feb 07 jari 178                     sigGenes.add(new Integer(i));
2 26 Feb 07 jari 179                 } else {
2 26 Feb 07 jari 180                     nonSigGenes.add(new Integer(i));
2 26 Feb 07 jari 181                 }
2 26 Feb 07 jari 182             }
2 26 Feb 07 jari 183             
2 26 Feb 07 jari 184         }  else {      
2 26 Feb 07 jari 185             
2 26 Feb 07 jari 186             adjPValuesVector = getAdjPVals(rawPValuesVector, correctionMethod);
2 26 Feb 07 jari 187             event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numGenes);
2 26 Feb 07 jari 188             fireValueChanged(event);
2 26 Feb 07 jari 189             event.setId(AlgorithmEvent.PROGRESS_VALUE);
2 26 Feb 07 jari 190             
2 26 Feb 07 jari 191             for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 192                 if (stop) {
2 26 Feb 07 jari 193                     throw new AbortException();
2 26 Feb 07 jari 194                 }
2 26 Feb 07 jari 195                 event.setIntValue(i);
2 26 Feb 07 jari 196                 event.setDescription("Finding significant genes: Current gene = " + (i + 1));
2 26 Feb 07 jari 197                 fireValueChanged(event);
2 26 Feb 07 jari 198                 
2 26 Feb 07 jari 199                 float currAdjP = ((Float)(adjPValuesVector.get(i))).floatValue();
2 26 Feb 07 jari 200                 //System.out.println("currAdjP: gene" + i + " = " + currAdjP);
2 26 Feb 07 jari 201                 if (correctionMethod == OneWayANOVAInitBox.ADJ_BONFERRONI) {// because we have to break out of the loop if a non-sig value is encountered
2 26 Feb 07 jari 202                     if (isSig[i]) {
2 26 Feb 07 jari 203                         sigGenes.add(new Integer(i));
2 26 Feb 07 jari 204                     } else {
2 26 Feb 07 jari 205                         nonSigGenes.add(new Integer(i));
2 26 Feb 07 jari 206                     }
2 26 Feb 07 jari 207                 } else {
2 26 Feb 07 jari 208                     if (currAdjP <= alpha) {
2 26 Feb 07 jari 209                         sigGenes.add(new Integer(i));
2 26 Feb 07 jari 210                     } else {
2 26 Feb 07 jari 211                         nonSigGenes.add(new Integer(i));
2 26 Feb 07 jari 212                     }
2 26 Feb 07 jari 213                 }
2 26 Feb 07 jari 214             }
2 26 Feb 07 jari 215         }
2 26 Feb 07 jari 216         
2 26 Feb 07 jari 217         clusterVector.add(sigGenes);
2 26 Feb 07 jari 218         clusterVector.add(nonSigGenes);
2 26 Feb 07 jari 219         k = clusterVector.size();
2 26 Feb 07 jari 220
2 26 Feb 07 jari 221         FloatMatrix fValuesMatrix = new FloatMatrix(fValuesVector.size(), 1);  
2 26 Feb 07 jari 222         //FloatMatrix pValuesMatrix = new FloatMatrix(pValuesVector.size(), 1);   
2 26 Feb 07 jari 223         FloatMatrix rawPValuesMatrix = new FloatMatrix(rawPValuesVector.size(), 1);
2 26 Feb 07 jari 224         FloatMatrix adjPValuesMatrix = new FloatMatrix(adjPValuesVector.size(), 1);
2 26 Feb 07 jari 225         FloatMatrix dfNumMatrix = new FloatMatrix(numGenes, 1);
2 26 Feb 07 jari 226         FloatMatrix dfDenomMatrix = new FloatMatrix(numGenes, 1);
2 26 Feb 07 jari 227
2 26 Feb 07 jari 228         for (int i = 0; i < fValuesVector.size(); i++) {
2 26 Feb 07 jari 229             fValuesMatrix.A[i][0] = ((Float)(fValuesVector.get(i))).floatValue();
2 26 Feb 07 jari 230         }  
2 26 Feb 07 jari 231         
2 26 Feb 07 jari 232         for (int i = 0; i < rawPValuesVector.size(); i++) {
2 26 Feb 07 jari 233             rawPValuesMatrix.A[i][0] = ((Float)(rawPValuesVector.get(i))).floatValue();
2 26 Feb 07 jari 234             adjPValuesMatrix.A[i][0] = ((Float)(adjPValuesVector.get(i))).floatValue();
2 26 Feb 07 jari 235         } 
2 26 Feb 07 jari 236         
2 26 Feb 07 jari 237         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 238             //dfNumMatrix.A[i][0] = (float)(getDfNum(i));
2 26 Feb 07 jari 239             //dfDenomMatrix.A[i][0] = (float)(getDfDenom(i));            
2 26 Feb 07 jari 240             dfNumMatrix.A[i][0] = ((Integer)(dfNumVector.get(i))).floatValue();
2 26 Feb 07 jari 241             dfDenomMatrix.A[i][0] = ((Integer)(dfDenomVector.get(i))).floatValue();
2 26 Feb 07 jari 242             if (dfNumMatrix.A[i][0] <= 0) dfNumMatrix.A[i][0] = Float.NaN;
2 26 Feb 07 jari 243             if (dfDenomMatrix.A[i][0] <= 0) dfDenomMatrix.A[i][0] = Float.NaN;
2 26 Feb 07 jari 244         }
2 26 Feb 07 jari 245         
2 26 Feb 07 jari 246         FloatMatrix ssGroupsMatrix = new FloatMatrix(numGenes, 1);
2 26 Feb 07 jari 247         FloatMatrix ssErrorMatrix = new FloatMatrix(numGenes, 1);
2 26 Feb 07 jari 248         
2 26 Feb 07 jari 249         for (int i = 0; i < ssGroupsMatrix.getRowDimension(); i++) {
2 26 Feb 07 jari 250             //float[] currentGene = getGene(i);
2 26 Feb 07 jari 251             //constant = getConstant(currentGene);
2 26 Feb 07 jari 252             /*
2 26 Feb 07 jari 253             if (i == 11) {
2 26 Feb 07 jari 254                 System.out.println("Gene A3: groupsSS = " + getGroupsSS(currentGene));
2 26 Feb 07 jari 255             }
2 26 Feb 07 jari 256              */
2 26 Feb 07 jari 257             ssGroupsMatrix.A[i][0] = ((Double)(ssGroupsVector.get(i))).floatValue();//(float)(getGroupsSS(currentGene));
2 26 Feb 07 jari 258             /*
2 26 Feb 07 jari 259             if (i == 11) {
2 26 Feb 07 jari 260                 System.out.println("Gene A3: ssGroupsMatrix.A[11][0] = " + ssGroupsMatrix.A[i][0]);
2 26 Feb 07 jari 261             }
2 26 Feb 07 jari 262              */
2 26 Feb 07 jari 263             ssErrorMatrix.A[i][0] = ((Double)(ssErrorVector.get(i))).floatValue();//(float)(getTotalSS(currentGene) - getGroupsSS(currentGene));
2 26 Feb 07 jari 264         }
2 26 Feb 07 jari 265         
2 26 Feb 07 jari 266   clusters = new Vector[k];
2 26 Feb 07 jari 267   
2 26 Feb 07 jari 268   for (int i = 0; i < k; i++) {
2 26 Feb 07 jari 269       clusters[i] = (Vector)(clusterVector.get(i));
2 26 Feb 07 jari 270   }
2 26 Feb 07 jari 271   
2 26 Feb 07 jari 272   FloatMatrix means = getMeans(clusters);
2 26 Feb 07 jari 273   FloatMatrix variances = getVariances(clusters, means);        
2 26 Feb 07 jari 274
2 26 Feb 07 jari 275   AlgorithmEvent event = null;
2 26 Feb 07 jari 276   if (hierarchical_tree) {
2 26 Feb 07 jari 277       event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, clusters.length, "Calculate Hierarchical Trees");
2 26 Feb 07 jari 278       fireValueChanged(event);
2 26 Feb 07 jari 279       event.setIntValue(0);
2 26 Feb 07 jari 280       event.setId(AlgorithmEvent.PROGRESS_VALUE);
2 26 Feb 07 jari 281       fireValueChanged(event);
2 26 Feb 07 jari 282   }        
2 26 Feb 07 jari 283
2 26 Feb 07 jari 284   Cluster result_cluster = new Cluster();
2 26 Feb 07 jari 285   NodeList nodeList = result_cluster.getNodeList();
2 26 Feb 07 jari 286   int[] features;        
2 26 Feb 07 jari 287   for (int i=0; i<clusters.length; i++) {
2 26 Feb 07 jari 288       if (stop) {
2 26 Feb 07 jari 289     throw new AbortException();
2 26 Feb 07 jari 290       }
2 26 Feb 07 jari 291       features = convert2int(clusters[i]);
2 26 Feb 07 jari 292       Node node = new Node(features);
2 26 Feb 07 jari 293       nodeList.addNode(node);
2 26 Feb 07 jari 294       if (hierarchical_tree) {
2 26 Feb 07 jari 295                 if (drawSigTreesOnly) {
2 26 Feb 07 jari 296                     if (i == 0) {
2 26 Feb 07 jari 297                         node.setValues(calculateHierarchicalTree(features, method_linkage, calculate_genes, calculate_experiments));
2 26 Feb 07 jari 298                         event.setIntValue(i+1);
2 26 Feb 07 jari 299                         fireValueChanged(event);                       
2 26 Feb 07 jari 300                     }                    
2 26 Feb 07 jari 301                 } else {                
2 26 Feb 07 jari 302                     node.setValues(calculateHierarchicalTree(features, method_linkage, calculate_genes, calculate_experiments));
2 26 Feb 07 jari 303                     event.setIntValue(i+1);
2 26 Feb 07 jari 304                     fireValueChanged(event);
2 26 Feb 07 jari 305                 }                
2 26 Feb 07 jari 306       }
2 26 Feb 07 jari 307   }
2 26 Feb 07 jari 308         
2 26 Feb 07 jari 309   // prepare the result
2 26 Feb 07 jari 310   AlgorithmData result = new AlgorithmData();
2 26 Feb 07 jari 311   result.addCluster("cluster", result_cluster);
2 26 Feb 07 jari 312   result.addParam("number-of-clusters", String.valueOf(clusters.length));    
2 26 Feb 07 jari 313   result.addMatrix("clusters_means", means);
2 26 Feb 07 jari 314   result.addMatrix("clusters_variances", variances); 
2 26 Feb 07 jari 315         //result.addMatrix("pValues", pValuesMatrix);
2 26 Feb 07 jari 316         result.addMatrix("rawPValues", rawPValuesMatrix);
2 26 Feb 07 jari 317         result.addMatrix("adjPValues", adjPValuesMatrix);
2 26 Feb 07 jari 318         result.addMatrix("fValues", fValuesMatrix);
2 26 Feb 07 jari 319         result.addMatrix("dfNumMatrix", dfNumMatrix);
2 26 Feb 07 jari 320         result.addMatrix("dfDenomMatrix", dfDenomMatrix);
2 26 Feb 07 jari 321         result.addMatrix("ssGroupsMatrix", ssGroupsMatrix);
2 26 Feb 07 jari 322         result.addMatrix("ssErrorMatrix", ssErrorMatrix);
2 26 Feb 07 jari 323         result.addMatrix("geneGroupMeansMatrix", getAllGeneGroupMeans());
2 26 Feb 07 jari 324         result.addMatrix("geneGroupSDsMatrix", getAllGeneGroupSDs());
2 26 Feb 07 jari 325   return result;        
2 26 Feb 07 jari 326         //return null; //for now
2 26 Feb 07 jari 327     }
2 26 Feb 07 jari 328
2 26 Feb 07 jari 329     private NodeValueList calculateHierarchicalTree(int[] features, int method, boolean genes, boolean experiments) throws AlgorithmException {
2 26 Feb 07 jari 330   NodeValueList nodeList = new NodeValueList();
2 26 Feb 07 jari 331   AlgorithmData data = new AlgorithmData();
2 26 Feb 07 jari 332   FloatMatrix experiment = getSubExperiment(this.expMatrix, features);
2 26 Feb 07 jari 333   data.addMatrix("experiment", experiment);
2 26 Feb 07 jari 334         data.addParam("hcl-distance-function", String.valueOf(this.hcl_function));
2 26 Feb 07 jari 335         data.addParam("hcl-distance-absolute", String.valueOf(this.hcl_absolute));
2 26 Feb 07 jari 336   data.addParam("method-linkage", String.valueOf(method));
2 26 Feb 07 jari 337   HCL hcl = new HCL();
2 26 Feb 07 jari 338   AlgorithmData result;
2 26 Feb 07 jari 339   if (genes) {
2 26 Feb 07 jari 340       data.addParam("calculate-genes", String.valueOf(true));
2 26 Feb 07 jari 341       result = hcl.execute(data);
2 26 Feb 07 jari 342       validate(result);
2 26 Feb 07 jari 343       addNodeValues(nodeList, result);
2 26 Feb 07 jari 344   }
2 26 Feb 07 jari 345   if (experiments) {
2 26 Feb 07 jari 346       data.addParam("calculate-genes", String.valueOf(false));
2 26 Feb 07 jari 347       result = hcl.execute(data);
2 26 Feb 07 jari 348       validate(result);
2 26 Feb 07 jari 349       addNodeValues(nodeList, result);
2 26 Feb 07 jari 350   }
2 26 Feb 07 jari 351   return nodeList;
2 26 Feb 07 jari 352     }
2 26 Feb 07 jari 353     
2 26 Feb 07 jari 354     private void addNodeValues(NodeValueList target_list, AlgorithmData source_result) {
2 26 Feb 07 jari 355   target_list.addNodeValue(new NodeValue("child-1-array", source_result.getIntArray("child-1-array")));
2 26 Feb 07 jari 356   target_list.addNodeValue(new NodeValue("child-2-array", source_result.getIntArray("child-2-array")));
2 26 Feb 07 jari 357   target_list.addNodeValue(new NodeValue("node-order", source_result.getIntArray("node-order")));
2 26 Feb 07 jari 358   target_list.addNodeValue(new NodeValue("height", source_result.getMatrix("height").getRowPackedCopy()));
2 26 Feb 07 jari 359     }
2 26 Feb 07 jari 360     
2 26 Feb 07 jari 361     private FloatMatrix getSubExperiment(FloatMatrix experiment, int[] features) {
2 26 Feb 07 jari 362   FloatMatrix subExperiment = new FloatMatrix(features.length, experiment.getColumnDimension());
2 26 Feb 07 jari 363   for (int i=0; i<features.length; i++) {
2 26 Feb 07 jari 364       subExperiment.A[i] = experiment.A[features[i]];
2 26 Feb 07 jari 365   }
2 26 Feb 07 jari 366   return subExperiment;
2 26 Feb 07 jari 367     }
2 26 Feb 07 jari 368     
2 26 Feb 07 jari 369     /**
2 26 Feb 07 jari 370      * Checking the result of hcl algorithm calculation.
2 26 Feb 07 jari 371      * @throws AlgorithmException, if the result is incorrect.
2 26 Feb 07 jari 372      */
2 26 Feb 07 jari 373     private void validate(AlgorithmData result) throws AlgorithmException {
2 26 Feb 07 jari 374   if (result.getIntArray("child-1-array") == null) {
2 26 Feb 07 jari 375       throw new AlgorithmException("parameter 'child-1-array' is null");
2 26 Feb 07 jari 376   }
2 26 Feb 07 jari 377   if (result.getIntArray("child-2-array") == null) {
2 26 Feb 07 jari 378       throw new AlgorithmException("parameter 'child-2-array' is null");
2 26 Feb 07 jari 379   }
2 26 Feb 07 jari 380   if (result.getIntArray("node-order") == null) {
2 26 Feb 07 jari 381       throw new AlgorithmException("parameter 'node-order' is null");
2 26 Feb 07 jari 382   }
2 26 Feb 07 jari 383   if (result.getMatrix("height") == null) {
2 26 Feb 07 jari 384       throw new AlgorithmException("parameter 'height' is null");
2 26 Feb 07 jari 385   }
2 26 Feb 07 jari 386     }
2 26 Feb 07 jari 387     
2 26 Feb 07 jari 388     private int[] convert2int(Vector source) {
2 26 Feb 07 jari 389   int[] int_matrix = new int[source.size()];
2 26 Feb 07 jari 390   for (int i=0; i<int_matrix.length; i++) {
2 26 Feb 07 jari 391       int_matrix[i] = (int)((Integer)source.get(i)).intValue();
2 26 Feb 07 jari 392   }
2 26 Feb 07 jari 393   return int_matrix;
2 26 Feb 07 jari 394     }
2 26 Feb 07 jari 395     
2 26 Feb 07 jari 396     private FloatMatrix getMeans(Vector[] clusters) {
2 26 Feb 07 jari 397   FloatMatrix means = new FloatMatrix(clusters.length, numExps);
2 26 Feb 07 jari 398   FloatMatrix mean;
2 26 Feb 07 jari 399   for (int i=0; i<clusters.length; i++) {
2 26 Feb 07 jari 400       mean = getMean(clusters[i]);
2 26 Feb 07 jari 401       means.A[i] = mean.A[0];
2 26 Feb 07 jari 402   }
2 26 Feb 07 jari 403   return means;
2 26 Feb 07 jari 404     }
2 26 Feb 07 jari 405     
2 26 Feb 07 jari 406     private FloatMatrix getMean(Vector cluster) {
2 26 Feb 07 jari 407   FloatMatrix mean = new FloatMatrix(1, numExps);
2 26 Feb 07 jari 408   float currentMean;
2 26 Feb 07 jari 409   int n = cluster.size();
2 26 Feb 07 jari 410   int denom = 0;
2 26 Feb 07 jari 411   float value;
2 26 Feb 07 jari 412   for (int i=0; i<numExps; i++) {
2 26 Feb 07 jari 413       currentMean = 0f;
2 26 Feb 07 jari 414       denom = 0;
2 26 Feb 07 jari 415       for (int j=0; j<n; j++) {
2 26 Feb 07 jari 416     value = expMatrix.get(((Integer) cluster.get(j)).intValue(), i);
2 26 Feb 07 jari 417     if (!Float.isNaN(value)) {
2 26 Feb 07 jari 418         currentMean += value;
2 26 Feb 07 jari 419         denom++;
2 26 Feb 07 jari 420     }
2 26 Feb 07 jari 421       }
2 26 Feb 07 jari 422       mean.set(0, i, currentMean/(float)denom);
2 26 Feb 07 jari 423   }
2 26 Feb 07 jari 424   
2 26 Feb 07 jari 425   return mean;
2 26 Feb 07 jari 426     }
2 26 Feb 07 jari 427     
2 26 Feb 07 jari 428     private FloatMatrix getVariances(Vector[] clusters, FloatMatrix means) {
2 26 Feb 07 jari 429   final int rows = means.getRowDimension();
2 26 Feb 07 jari 430   final int columns = means.getColumnDimension();
2 26 Feb 07 jari 431   FloatMatrix variances = new FloatMatrix(rows, columns);
2 26 Feb 07 jari 432   for (int row=0; row<rows; row++) {
2 26 Feb 07 jari 433       for (int column=0; column<columns; column++) {
2 26 Feb 07 jari 434     variances.set(row, column, getSampleVariance(clusters[row], column, means.get(row, column)));
2 26 Feb 07 jari 435       }
2 26 Feb 07 jari 436   }
2 26 Feb 07 jari 437   return variances;
2 26 Feb 07 jari 438     }
2 26 Feb 07 jari 439     
2 26 Feb 07 jari 440     int validN;
2 26 Feb 07 jari 441     
2 26 Feb 07 jari 442     private float getSampleNormalizedSum(Vector cluster, int column, float mean) {
2 26 Feb 07 jari 443   final int size = cluster.size();
2 26 Feb 07 jari 444   float sum = 0f;
2 26 Feb 07 jari 445   float value;
2 26 Feb 07 jari 446   validN = 0;
2 26 Feb 07 jari 447   for (int i=0; i<size; i++) {
2 26 Feb 07 jari 448       value = expMatrix.get(((Integer) cluster.get(i)).intValue(), column);
2 26 Feb 07 jari 449       if (!Float.isNaN(value)) {
2 26 Feb 07 jari 450     sum += Math.pow(value-mean, 2);
2 26 Feb 07 jari 451     validN++;
2 26 Feb 07 jari 452       }
2 26 Feb 07 jari 453   }
2 26 Feb 07 jari 454   return sum;
2 26 Feb 07 jari 455     }
2 26 Feb 07 jari 456     
2 26 Feb 07 jari 457     private float getSampleVariance(Vector cluster, int column, float mean) {
2 26 Feb 07 jari 458   return(float)Math.sqrt(getSampleNormalizedSum(cluster, column, mean)/(float)(validN-1));
2 26 Feb 07 jari 459   
2 26 Feb 07 jari 460     }  
2 26 Feb 07 jari 461     
2 26 Feb 07 jari 462     private boolean[] isGeneSigByFDRPropNew2() throws AlgorithmException {
2 26 Feb 07 jari 463         double[] nonNanPVals = new double[origPVals.length];
2 26 Feb 07 jari 464         for (int i = 0; i < origPVals.length; i++) { //gets rid of NaN's for sorting
2 26 Feb 07 jari 465             if (Double.isNaN(origPVals[i])) {
2 26 Feb 07 jari 466                 nonNanPVals[i] = Double.POSITIVE_INFINITY;
2 26 Feb 07 jari 467             } else {
2 26 Feb 07 jari 468                 nonNanPVals[i] = origPVals[i];
2 26 Feb 07 jari 469             }
2 26 Feb 07 jari 470         } 
2 26 Feb 07 jari 471         
2 26 Feb 07 jari 472         QSort sortOrigPVals = new QSort(nonNanPVals, QSort.ASCENDING);
2 26 Feb 07 jari 473         double[] sortedOrigPVals = sortOrigPVals.getSortedDouble();
2 26 Feb 07 jari 474         int[] sortedIndices = sortOrigPVals.getOrigIndx();  
2 26 Feb 07 jari 475         boolean[] isGeneSig = new boolean[numGenes];
2 26 Feb 07 jari 476         for (int i = 0; i < isGeneSig.length; i++) {
2 26 Feb 07 jari 477             isGeneSig[i] = false;
2 26 Feb 07 jari 478         } 
2 26 Feb 07 jari 479         
2 26 Feb 07 jari 480         double[] yKArray = getYKArray();
2 26 Feb 07 jari 481         
2 26 Feb 07 jari 482         if (sortedOrigPVals[0] >= yKArray[0]) {
2 26 Feb 07 jari 483             return isGeneSig;
2 26 Feb 07 jari 484             
2 26 Feb 07 jari 485         } else {
2 26 Feb 07 jari 486             isGeneSig[sortedIndices[0]] = true;
2 26 Feb 07 jari 487             if (useFastFDRApprox) {
2 26 Feb 07 jari 488                 
2 26 Feb 07 jari 489                 for (int i = 1; i < sortedOrigPVals.length; i++) {
2 26 Feb 07 jari 490                     int rGamma = (int)(Math.floor((i+1)*falseProp));
2 26 Feb 07 jari 491                     int rMinusOneGamma = (int)(Math.floor(i*falseProp));
2 26 Feb 07 jari 492                     double yKRGamma = yKArray[rGamma];
2 26 Feb 07 jari 493                     //System.out.println("rGamma = " + rGamma + ", (r - 1)Gamma = " + rMinusOneGamma + ", yKRGamma = " + yKRGamma);
2 26 Feb 07 jari 494                     if ((rGamma > rMinusOneGamma) || (sortedOrigPVals[i] < yKRGamma)) {
2 26 Feb 07 jari 495                        isGeneSig[sortedIndices[i]] = true; 
2 26 Feb 07 jari 496                     } else {
2 26 Feb 07 jari 497                         break;
2 26 Feb 07 jari 498                     }
2 26 Feb 07 jari 499                 }
2 26 Feb 07 jari 500                 
2 26 Feb 07 jari 501             }
2 26 Feb 07 jari 502         }
2 26 Feb 07 jari 503         
2 26 Feb 07 jari 504         return isGeneSig;
2 26 Feb 07 jari 505     }    
2 26 Feb 07 jari 506     
2 26 Feb 07 jari 507     private double[] getYKArray() throws AlgorithmException {
2 26 Feb 07 jari 508         AlgorithmEvent event2 = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numPerms);
2 26 Feb 07 jari 509         fireValueChanged(event2);
2 26 Feb 07 jari 510         event2.setId(AlgorithmEvent.PROGRESS_VALUE);  
2 26 Feb 07 jari 511         int maxRGamma = (int)(Math.floor(numGenes*falseProp));
2 26 Feb 07 jari 512         double[][] pValArray =new double[maxRGamma + 1][numPerms];     
2 26 Feb 07 jari 513         
2 26 Feb 07 jari 514         Vector validExpts = new Vector();
2 26 Feb 07 jari 515         for (int i = 0; i < groupAssignments.length; i++) {
2 26 Feb 07 jari 516             if (groupAssignments[i] != 0) validExpts.add(new Integer(i));
2 26 Feb 07 jari 517         }
2 26 Feb 07 jari 518         
2 26 Feb 07 jari 519         int[] validArray = new int[validExpts.size()];
2 26 Feb 07 jari 520         for (int j = 0; j < validArray.length; j++) {
2 26 Feb 07 jari 521             validArray[j] = ((Integer)(validExpts.get(j))).intValue();
2 26 Feb 07 jari 522         }          
2 26 Feb 07 jari 523         
2 26 Feb 07 jari 524         for (int i = 0; i < numPerms; i++) {
2 26 Feb 07 jari 525             if (stop) {
2 26 Feb 07 jari 526                 throw new AbortException();
2 26 Feb 07 jari 527             }
2 26 Feb 07 jari 528             event2.setIntValue(i);
2 26 Feb 07 jari 529             event2.setDescription("Permuting matrix: Current permutation = " + (i + 1));
2 26 Feb 07 jari 530             fireValueChanged(event2);            
2 26 Feb 07 jari 531             int[] permutedExpts = getPermutedValues(numExps, validArray);
2 26 Feb 07 jari 532             FloatMatrix permutedMatrix = getPermutedMatrix(expMatrix, permutedExpts);
2 26 Feb 07 jari 533             float[] currPermFVals = getPermutedFVals(permutedMatrix);
2 26 Feb 07 jari 534             int[][] dfs = getDfs(permutedMatrix);
2 26 Feb 07 jari 535             double[] currPermPVals = getParametricPVals(currPermFVals, dfs[0], dfs[1]);
2 26 Feb 07 jari 536             
2 26 Feb 07 jari 537             for (int j = 0; j < currPermPVals.length; j++) {
2 26 Feb 07 jari 538                 if (Double.isNaN(currPermPVals[j])) currPermPVals[j] = Double.POSITIVE_INFINITY;
2 26 Feb 07 jari 539                 // this is to push NaNs to the end of the sorted array, since otherwise they would mess up quantile calculations
2 26 Feb 07 jari 540             }
2 26 Feb 07 jari 541             
2 26 Feb 07 jari 542             QSort sortCurrPVals = new QSort(currPermPVals, QSort.ASCENDING);
2 26 Feb 07 jari 543             double[] sortedCurrPVals = sortCurrPVals.getSortedDouble();
2 26 Feb 07 jari 544             //uPlusOneSmallestPVector.add(new Double(sortedCurrPVals[u]));    
2 26 Feb 07 jari 545             for (int j = 0; j < pValArray.length; j++) {
2 26 Feb 07 jari 546                 pValArray[j][i] = sortedCurrPVals[j];
2 26 Feb 07 jari 547             }         
2 26 Feb 07 jari 548         }  
2 26 Feb 07 jari 549         
2 26 Feb 07 jari 550         double[] yKArray = new double[pValArray.length];
2 26 Feb 07 jari 551         
2 26 Feb 07 jari 552         for (int i = 0; i < pValArray.length; i++) {
2 26 Feb 07 jari 553             double[] currRow = new double[pValArray[i].length];
2 26 Feb 07 jari 554             
2 26 Feb 07 jari 555             for (int j = 0; j < currRow.length; j++) {
2 26 Feb 07 jari 556                 currRow[j] = pValArray[i][j];
2 26 Feb 07 jari 557             }
2 26 Feb 07 jari 558             
2 26 Feb 07 jari 559             for (int j = 0; j < currRow.length; j++) {
2 26 Feb 07 jari 560                 if (Double.isNaN(currRow[j])) currRow[j] = Double.POSITIVE_INFINITY;
2 26 Feb 07 jari 561                 // this is to push NaNs to the end of the sorted array, since otherwise they would mess up quantile calculations 
2 26 Feb 07 jari 562             }
2 26 Feb 07 jari 563             
2 26 Feb 07 jari 564             QSort sortCurrRow = new QSort(currRow, QSort.ASCENDING);
2 26 Feb 07 jari 565             double[] sortedCurrRow = sortCurrRow.getSortedDouble();
2 26 Feb 07 jari 566             int selectedIndex = (int)Math.floor(sortedCurrRow.length*alpha) - 1;
2 26 Feb 07 jari 567             if (selectedIndex < 0) selectedIndex = 0;
2 26 Feb 07 jari 568             yKArray[i] = sortedCurrRow[selectedIndex]; 
2 26 Feb 07 jari 569             //System.out.println("i= " + i + ", selectedIndex = " + selectedIndex + ", yKArray[i] = " + yKArray[i]);
2 26 Feb 07 jari 570         }
2 26 Feb 07 jari 571         
2 26 Feb 07 jari 572         return yKArray;        
2 26 Feb 07 jari 573     }
2 26 Feb 07 jari 574     
2 26 Feb 07 jari 575     private boolean[] isGeneSigByFDRNum() throws AlgorithmException {
2 26 Feb 07 jari 576         double[] nonNanPVals = new double[origPVals.length];
2 26 Feb 07 jari 577         for (int i = 0; i < origPVals.length; i++) { //gets rid of NaN's for sorting
2 26 Feb 07 jari 578             if (Double.isNaN(origPVals[i])) {
2 26 Feb 07 jari 579                 nonNanPVals[i] = Double.POSITIVE_INFINITY;
2 26 Feb 07 jari 580             } else {
2 26 Feb 07 jari 581                 nonNanPVals[i] = origPVals[i];
2 26 Feb 07 jari 582             }
2 26 Feb 07 jari 583         }
2 26 Feb 07 jari 584         QSort sortOrigPVals = new QSort(nonNanPVals, QSort.ASCENDING);
2 26 Feb 07 jari 585         double[] sortedOrigPVals = sortOrigPVals.getSortedDouble();
2 26 Feb 07 jari 586         int[] sortedIndices = sortOrigPVals.getOrigIndx();        
2 26 Feb 07 jari 587         boolean[] isGeneSig = new boolean[numGenes];
2 26 Feb 07 jari 588         for (int i = 0; i < isGeneSig.length; i++) {
2 26 Feb 07 jari 589             isGeneSig[i] = false;
2 26 Feb 07 jari 590         }
2 26 Feb 07 jari 591         for (int i = 0; i < falseNum; i++) {
2 26 Feb 07 jari 592             isGeneSig[sortedIndices[i]] = true;          
2 26 Feb 07 jari 593         }
2 26 Feb 07 jari 594         if (useFastFDRApprox) {
2 26 Feb 07 jari 595             double yK = getYConservative(alpha, falseNum); 
2 26 Feb 07 jari 596             //System.out.println("yK = " + yK);
2 26 Feb 07 jari 597             for (int i = falseNum; i < sortedOrigPVals.length; i++) {
2 26 Feb 07 jari 598                 if (sortedOrigPVals[i] < yK) {
2 26 Feb 07 jari 599                     isGeneSig[sortedIndices[i]] = true;
2 26 Feb 07 jari 600                 } else {
2 26 Feb 07 jari 601                     break;
2 26 Feb 07 jari 602                 }
2 26 Feb 07 jari 603             }
2 26 Feb 07 jari 604             
2 26 Feb 07 jari 605         } else {// if (!useFastFDRAprpox)
2 26 Feb 07 jari 606             for (int i = falseNum; i < origPVals.length; i++) {
2 26 Feb 07 jari 607                 
2 26 Feb 07 jari 608             }
2 26 Feb 07 jari 609         }
2 26 Feb 07 jari 610         //return null; //for now
2 26 Feb 07 jari 611         return isGeneSig;
2 26 Feb 07 jari 612     }    
2 26 Feb 07 jari 613     
2 26 Feb 07 jari 614
2 26 Feb 07 jari 615     private double getYConservative(double alphaQuantile, int u) throws AlgorithmException {
2 26 Feb 07 jari 616         AlgorithmEvent event2 = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numPerms);
2 26 Feb 07 jari 617         fireValueChanged(event2);
2 26 Feb 07 jari 618         event2.setId(AlgorithmEvent.PROGRESS_VALUE);
2 26 Feb 07 jari 619         Vector uPlusOneSmallestPVector = new Vector(); 
2 26 Feb 07 jari 620         
2 26 Feb 07 jari 621         Vector validExpts = new Vector();
2 26 Feb 07 jari 622         for (int i = 0; i < groupAssignments.length; i++) {
2 26 Feb 07 jari 623             if (groupAssignments[i] != 0) validExpts.add(new Integer(i));
2 26 Feb 07 jari 624         }
2 26 Feb 07 jari 625         
2 26 Feb 07 jari 626         int[] validArray = new int[validExpts.size()];
2 26 Feb 07 jari 627         for (int j = 0; j < validArray.length; j++) {
2 26 Feb 07 jari 628             validArray[j] = ((Integer)(validExpts.get(j))).intValue();
2 26 Feb 07 jari 629         }        
2 26 Feb 07 jari 630         for (int i = 0; i < numPerms; i++) {
2 26 Feb 07 jari 631             if (stop) {
2 26 Feb 07 jari 632                 throw new AbortException();
2 26 Feb 07 jari 633             }
2 26 Feb 07 jari 634             event.setIntValue(i);
2 26 Feb 07 jari 635             event.setDescription("Permuting matrix: Current permutation = " + (i + 1));
2 26 Feb 07 jari 636             fireValueChanged(event);            
2 26 Feb 07 jari 637             int[] permutedExpts = getPermutedValues(numExps, validArray);
2 26 Feb 07 jari 638             FloatMatrix permutedMatrix = getPermutedMatrix(expMatrix, permutedExpts);
2 26 Feb 07 jari 639             float[] currPermFVals = getPermutedFVals(permutedMatrix);
2 26 Feb 07 jari 640             int[][] dfs = getDfs(permutedMatrix);
2 26 Feb 07 jari 641             double[] currPermPVals = getParametricPVals(currPermFVals, dfs[0], dfs[1]);
2 26 Feb 07 jari 642             
2 26 Feb 07 jari 643             for (int j = 0; j < currPermPVals.length; j++) {
2 26 Feb 07 jari 644                 if (Double.isNaN(currPermPVals[j])) currPermPVals[j] = Double.POSITIVE_INFINITY;
2 26 Feb 07 jari 645                 // this is to push NaNs to the end of the sorted array, since otherwise they would mess up quantile calculations
2 26 Feb 07 jari 646             }
2 26 Feb 07 jari 647             
2 26 Feb 07 jari 648             QSort sortCurrPVals = new QSort(currPermPVals, QSort.ASCENDING);
2 26 Feb 07 jari 649             double[] sortedCurrPVals = sortCurrPVals.getSortedDouble();
2 26 Feb 07 jari 650             uPlusOneSmallestPVector.add(new Double(sortedCurrPVals[u]));           
2 26 Feb 07 jari 651         }
2 26 Feb 07 jari 652         
2 26 Feb 07 jari 653         double[] uPlusOneSmallestArray = new double[uPlusOneSmallestPVector.size()];
2 26 Feb 07 jari 654         for(int i = 0; i < uPlusOneSmallestPVector.size(); i++) {
2 26 Feb 07 jari 655             uPlusOneSmallestArray[i] = ((Double)(uPlusOneSmallestPVector.get(i))).doubleValue();
2 26 Feb 07 jari 656         }
2 26 Feb 07 jari 657         
2 26 Feb 07 jari 658         QSort sortUPlusOneArray = new QSort(uPlusOneSmallestArray, QSort.ASCENDING);
2 26 Feb 07 jari 659         double[] sortedUPlusOneArray = sortUPlusOneArray.getSortedDouble();
2 26 Feb 07 jari 660         
2 26 Feb 07 jari 661         int selectedIndex = (int)Math.floor(sortedUPlusOneArray.length*alphaQuantile) - 1;
2 26 Feb 07 jari 662         //System.out.println("Selected index (before setting to zero) = " + selectedIndex);
2 26 Feb 07 jari 663         if (selectedIndex < 0) selectedIndex = 0;        
2 26 Feb 07 jari 664         
2 26 Feb 07 jari 665         return sortedUPlusOneArray[selectedIndex];        
2 26 Feb 07 jari 666         
2 26 Feb 07 jari 667         //return 0d; //for now;
2 26 Feb 07 jari 668     }
2 26 Feb 07 jari 669     
2 26 Feb 07 jari 670     private int[][] getDfs(FloatMatrix permutedMatrix) {
2 26 Feb 07 jari 671         int[][] dfs = new int[2][numGenes];
2 26 Feb 07 jari 672         for (int gene = 0; gene < numGenes; gene++) {
2 26 Feb 07 jari 673             dfs[0][gene] = getDfNum(gene, permutedMatrix); 
2 26 Feb 07 jari 674             dfs[1][gene] = getDfDenom(gene, permutedMatrix);
2 26 Feb 07 jari 675         }
2 26 Feb 07 jari 676         return dfs;
2 26 Feb 07 jari 677     }
2 26 Feb 07 jari 678     
2 26 Feb 07 jari 679     private double[] getParametricPVals(float[] fVals, int[] dfNums, int[] dfDenoms) {
2 26 Feb 07 jari 680         double[] pVals = new double[numGenes];
2 26 Feb 07 jari 681         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 682             double currF = (double)fVals[i];
2 26 Feb 07 jari 683             int currDfNum = dfNums[i];
2 26 Feb 07 jari 684             int currDfDenom = dfDenoms[i];
2 26 Feb 07 jari 685             if (Double.isNaN(currF) || (currDfNum <= 0) || (currDfDenom <= 0) ) {
2 26 Feb 07 jari 686                 pVals[i] = Double.NaN;
2 26 Feb 07 jari 687             } else {
2 26 Feb 07 jari 688                 FDistribution fDist = new FDistribution(currDfNum, currDfDenom);     
2 26 Feb 07 jari 689                 //System.out.println("i (gene) = " + i + ", currF = " + currF + ", currDfNum = " + currDfNum + ", currDfDenom = " + currDfDenom);
2 26 Feb 07 jari 690                 double cumulProb = fDist.cumulative(currF);
2 26 Feb 07 jari 691                 double pValue = 1 - cumulProb; //                
2 26 Feb 07 jari 692                 if (pValue > 1) {
2 26 Feb 07 jari 693                     pValue = 1.0d;
2 26 Feb 07 jari 694                 }   
2 26 Feb 07 jari 695                 pVals[i] = pValue;
2 26 Feb 07 jari 696             }
2 26 Feb 07 jari 697         }
2 26 Feb 07 jari 698         return pVals;
2 26 Feb 07 jari 699     }
2 26 Feb 07 jari 700     
2 26 Feb 07 jari 701     private void getFDfSSValues() throws AlgorithmException {
2 26 Feb 07 jari 702         event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numGenes);
2 26 Feb 07 jari 703         fireValueChanged(event);
2 26 Feb 07 jari 704         event.setId(AlgorithmEvent.PROGRESS_VALUE);   
2 26 Feb 07 jari 705         
2 26 Feb 07 jari 706         for (int gene = 0; gene < numGenes; gene++) {
2 26 Feb 07 jari 707             float currentF = 0.0f;
2 26 Feb 07 jari 708             if (stop) {
2 26 Feb 07 jari 709                 throw new AbortException();
2 26 Feb 07 jari 710             }
2 26 Feb 07 jari 711             event.setIntValue(gene);
2 26 Feb 07 jari 712             event.setDescription("Calculating F and df: Current gene = " + (gene + 1));
2 26 Feb 07 jari 713             fireValueChanged(event);            
2 26 Feb 07 jari 714             //boolean sig = false;
2 26 Feb 07 jari 715             float[] geneValues = new float[numExps];
2 26 Feb 07 jari 716             int n = 0;
2 26 Feb 07 jari 717             for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 718                 geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 719                 if (!Float.isNaN(geneValues[i])) {
2 26 Feb 07 jari 720                     n++;
2 26 Feb 07 jari 721                 }
2 26 Feb 07 jari 722             }
2 26 Feb 07 jari 723             
2 26 Feb 07 jari 724             if (n == 0) {
2 26 Feb 07 jari 725                 currentF = Float.NaN;
2 26 Feb 07 jari 726                 //currentP = Float.NaN;
2 26 Feb 07 jari 727                 //return false;
2 26 Feb 07 jari 728             }
2 26 Feb 07 jari 729             
2 26 Feb 07 jari 730             constant = getConstant(geneValues);
2 26 Feb 07 jari 731             
2 26 Feb 07 jari 732             double totalSS = getTotalSS(geneValues);
2 26 Feb 07 jari 733             double groupsSS = getGroupsSS(geneValues);
2 26 Feb 07 jari 734             double errorSS = totalSS - groupsSS;
2 26 Feb 07 jari 735             
2 26 Feb 07 jari 736             if ((Double.isNaN(totalSS))||(Double.isNaN(groupsSS))||(Double.isNaN(errorSS))) {
2 26 Feb 07 jari 737                 currentF = Float.NaN;
2 26 Feb 07 jari 738                 //currentP = Float.NaN;
2 26 Feb 07 jari 739                 //return false;
2 26 Feb 07 jari 740             }
2 26 Feb 07 jari 741             
2 26 Feb 07 jari 742             //int totalDF = validN - 1;
2 26 Feb 07 jari 743             int groupsDF = getDfNum(gene);
2 26 Feb 07 jari 744             int errorDF = getDfDenom(gene);
2 26 Feb 07 jari 745             
2 26 Feb 07 jari 746             double groupsMS = groupsSS / groupsDF;
2 26 Feb 07 jari 747             double errorMS = errorSS / errorDF;
2 26 Feb 07 jari 748                         
2 26 Feb 07 jari 749             if (!Float.isNaN(currentF)) {
2 26 Feb 07 jari 750                 double fValue = groupsMS/errorMS;
2 26 Feb 07 jari 751                 currentF = (float)(fValue);
2 26 Feb 07 jari 752             }
2 26 Feb 07 jari 753             
2 26 Feb 07 jari 754             if (Float.isInfinite(currentF)) {
2 26 Feb 07 jari 755                 currentF = Float.NaN;
2 26 Feb 07 jari 756             }
2 26 Feb 07 jari 757             //System.out.println("currentF for gene " + gene + " = " + currentF);
2 26 Feb 07 jari 758             
2 26 Feb 07 jari 759             fValuesVector.add(new Float(currentF));
2 26 Feb 07 jari 760             dfNumVector.add(new Integer(groupsDF));
2 26 Feb 07 jari 761             dfDenomVector.add(new Integer(errorDF));
2 26 Feb 07 jari 762             ssGroupsVector.add(new Double(groupsSS));
2 26 Feb 07 jari 763             ssErrorVector.add(new Double(errorSS));
2 26 Feb 07 jari 764         }
2 26 Feb 07 jari 765     }
2 26 Feb 07 jari 766     
2 26 Feb 07 jari 767     private Vector getRawPValuesFromFDist() {
2 26 Feb 07 jari 768         Vector rawPVals = new Vector();
2 26 Feb 07 jari 769         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 770             double currF = ((Float)(fValuesVector.get(i))).doubleValue();
2 26 Feb 07 jari 771             int currDfNum = ((Integer)(dfNumVector.get(i))).intValue(); 
2 26 Feb 07 jari 772             int currDfDenom = ((Integer)(dfDenomVector.get(i))).intValue();
2 26 Feb 07 jari 773             
2 26 Feb 07 jari 774             if (Double.isNaN(currF) || (currDfNum <= 0) || (currDfDenom <= 0) ) {
2 26 Feb 07 jari 775                 rawPVals.add(new Float(Float.NaN));
2 26 Feb 07 jari 776             } else {
2 26 Feb 07 jari 777                 FDistribution fDist = new FDistribution(currDfNum, currDfDenom);     
2 26 Feb 07 jari 778                 //System.out.println("i (gene) = " + i + ", currF = " + currF + ", currDfNum = " + currDfNum + ", currDfDenom = " + currDfDenom);
2 26 Feb 07 jari 779                 double cumulProb = fDist.cumulative(currF);
2 26 Feb 07 jari 780                 double pValue = 1 - cumulProb; //                
2 26 Feb 07 jari 781                 if (pValue > 1) {
2 26 Feb 07 jari 782                     pValue = 1.0d;
2 26 Feb 07 jari 783                 }                
2 26 Feb 07 jari 784                 rawPVals.add(new Float(pValue));                
2 26 Feb 07 jari 785             }
2 26 Feb 07 jari 786         }
2 26 Feb 07 jari 787         return rawPVals;
2 26 Feb 07 jari 788     }
2 26 Feb 07 jari 789     
2 26 Feb 07 jari 790     private Vector getRawPValsFromPerms() throws AlgorithmException {
2 26 Feb 07 jari 791         event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numPerms);
2 26 Feb 07 jari 792         fireValueChanged(event);
2 26 Feb 07 jari 793         event.setId(AlgorithmEvent.PROGRESS_VALUE);
2 26 Feb 07 jari 794         
2 26 Feb 07 jari 795         float[] permPValues = new float[numGenes];
2 26 Feb 07 jari 796         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 797             permPValues[i] = 0f;
2 26 Feb 07 jari 798         }
2 26 Feb 07 jari 799         float[] origFValues = new float[fValuesVector.size()];
2 26 Feb 07 jari 800         for (int i = 0; i < fValuesVector.size(); i++) {
2 26 Feb 07 jari 801             origFValues[i] = ((Float)(fValuesVector.get(i))).floatValue();
2 26 Feb 07 jari 802         }
2 26 Feb 07 jari 803         
2 26 Feb 07 jari 804         Vector validExpts = new Vector();
2 26 Feb 07 jari 805         for (int i = 0; i < groupAssignments.length; i++) {
2 26 Feb 07 jari 806             if (groupAssignments[i] != 0) validExpts.add(new Integer(i));
2 26 Feb 07 jari 807         }
2 26 Feb 07 jari 808         
2 26 Feb 07 jari 809         int[] validArray = new int[validExpts.size()];
2 26 Feb 07 jari 810         for (int j = 0; j < validArray.length; j++) {
2 26 Feb 07 jari 811             validArray[j] = ((Integer)(validExpts.get(j))).intValue();
2 26 Feb 07 jari 812         }  
2 26 Feb 07 jari 813         for (int i = 0; i < numPerms; i++) {
2 26 Feb 07 jari 814             if (stop) {
2 26 Feb 07 jari 815                 throw new AbortException();
2 26 Feb 07 jari 816             }
2 26 Feb 07 jari 817             event.setIntValue(i);
2 26 Feb 07 jari 818             event.setDescription("Permuting matrix: Current permutation = " + (i + 1));
2 26 Feb 07 jari 819             fireValueChanged(event);            
2 26 Feb 07 jari 820             int[] permutedExpts = getPermutedValues(numExps, validArray);
2 26 Feb 07 jari 821             FloatMatrix permutedMatrix = getPermutedMatrix(expMatrix, permutedExpts);
2 26 Feb 07 jari 822             float[] currPermFVals = getPermutedFVals(permutedMatrix);
2 26 Feb 07 jari 823             for (int j = 0; j < numGenes; j++) {
2 26 Feb 07 jari 824                 if (currPermFVals[j] > origFValues[j]) permPValues[j]++;
2 26 Feb 07 jari 825             }
2 26 Feb 07 jari 826         }
2 26 Feb 07 jari 827         
2 26 Feb 07 jari 828         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 829             if (Float.isNaN(origFValues[i])) {
2 26 Feb 07 jari 830                 permPValues[i] = Float.NaN;
2 26 Feb 07 jari 831             } else {
2 26 Feb 07 jari 832                 permPValues[i] = (float)permPValues[i]/(float)numPerms;
2 26 Feb 07 jari 833             }
2 26 Feb 07 jari 834         }
2 26 Feb 07 jari 835         
2 26 Feb 07 jari 836         Vector permPValsVector = new Vector();
2 26 Feb 07 jari 837         
2 26 Feb 07 jari 838         for (int i = 0; i < permPValues.length; i++) {
2 26 Feb 07 jari 839             permPValsVector.add(new Float(permPValues[i]));
2 26 Feb 07 jari 840         }
2 26 Feb 07 jari 841         
2 26 Feb 07 jari 842         return permPValsVector;
2 26 Feb 07 jari 843     }
2 26 Feb 07 jari 844     
2 26 Feb 07 jari 845     private FloatMatrix getPermutedMatrix(FloatMatrix inputMatrix, int[] permExpts) {
2 26 Feb 07 jari 846         FloatMatrix permutedMatrix = new FloatMatrix(inputMatrix.getRowDimension(), inputMatrix.getColumnDimension());
2 26 Feb 07 jari 847         for (int i = 0; i < inputMatrix.getRowDimension(); i++) {
2 26 Feb 07 jari 848             for (int j = 0; j < inputMatrix.getColumnDimension(); j++) {
2 26 Feb 07 jari 849                 permutedMatrix.A[i][j] = inputMatrix.A[i][permExpts[j]];
2 26 Feb 07 jari 850             }
2 26 Feb 07 jari 851         }
2 26 Feb 07 jari 852         return permutedMatrix;
2 26 Feb 07 jari 853     }    
2 26 Feb 07 jari 854     
2 26 Feb 07 jari 855     private float[] getPermutedFVals(FloatMatrix permMatrix) {        
2 26 Feb 07 jari 856         float[] permFVals  = new float[numGenes];
2 26 Feb 07 jari 857         
2 26 Feb 07 jari 858         for (int gene = 0; gene < numGenes; gene++) {
2 26 Feb 07 jari 859             float currentF = 0.0f;
2 26 Feb 07 jari 860             //boolean sig = false;
2 26 Feb 07 jari 861             float[] geneValues = new float[numExps];
2 26 Feb 07 jari 862             int n = 0;
2 26 Feb 07 jari 863             for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 864                 geneValues[i] = permMatrix.A[gene][i];
2 26 Feb 07 jari 865                 if (!Float.isNaN(geneValues[i])) {
2 26 Feb 07 jari 866                     n++;
2 26 Feb 07 jari 867                 }
2 26 Feb 07 jari 868             }
2 26 Feb 07 jari 869             
2 26 Feb 07 jari 870             if (n == 0) {
2 26 Feb 07 jari 871                 currentF = Float.NaN;
2 26 Feb 07 jari 872                 //currentP = Float.NaN;
2 26 Feb 07 jari 873                 //return false;
2 26 Feb 07 jari 874             }
2 26 Feb 07 jari 875             
2 26 Feb 07 jari 876             constant = getConstant(geneValues);
2 26 Feb 07 jari 877             
2 26 Feb 07 jari 878             double totalSS = getTotalSS(geneValues);
2 26 Feb 07 jari 879             double groupsSS = getGroupsSS(geneValues);
2 26 Feb 07 jari 880             double errorSS = totalSS - groupsSS;
2 26 Feb 07 jari 881             
2 26 Feb 07 jari 882             if ((Double.isNaN(totalSS))||(Double.isNaN(groupsSS))||(Double.isNaN(errorSS))) {
2 26 Feb 07 jari 883                 currentF = Float.NaN;
2 26 Feb 07 jari 884                 //currentP = Float.NaN;
2 26 Feb 07 jari 885                 //return false;
2 26 Feb 07 jari 886             } 
2 26 Feb 07 jari 887             
2 26 Feb 07 jari 888             int groupsDF = getDfNum(gene, permMatrix);
2 26 Feb 07 jari 889             int errorDF = getDfDenom(gene, permMatrix);            
2 26 Feb 07 jari 890             
2 26 Feb 07 jari 891             double groupsMS = groupsSS / groupsDF;
2 26 Feb 07 jari 892             double errorMS = errorSS / errorDF;
2 26 Feb 07 jari 893             
2 26 Feb 07 jari 894             if (!Float.isNaN(currentF)) {
2 26 Feb 07 jari 895                 double fValue = groupsMS/errorMS;
2 26 Feb 07 jari 896                 currentF = (float)(fValue);
2 26 Feb 07 jari 897             }
2 26 Feb 07 jari 898             
2 26 Feb 07 jari 899             //System.out.println("currentF for gene " + gene + " = " + currentF);
2 26 Feb 07 jari 900             permFVals[gene] = currentF;
2 26 Feb 07 jari 901         }  
2 26 Feb 07 jari 902         return permFVals;
2 26 Feb 07 jari 903     }
2 26 Feb 07 jari 904     
2 26 Feb 07 jari 905     private int[] getPermutedValues(int arrayLength, int[] validArray) {//returns an integer array of length "arrayLength", with the valid values (the currently included experiments) permuted
2 26 Feb 07 jari 906         int[] permutedValues = new int[arrayLength];
2 26 Feb 07 jari 907         for (int i = 0; i < permutedValues.length; i++) {
2 26 Feb 07 jari 908             permutedValues[i] = i;
2 26 Feb 07 jari 909         }
2 26 Feb 07 jari 910         
2 26 Feb 07 jari 911         int[] permutedValidArray = new int[validArray.length];
2 26 Feb 07 jari 912         for (int i = 0; i < validArray.length; i++) {
2 26 Feb 07 jari 913             permutedValidArray[i] = validArray[i];
2 26 Feb 07 jari 914         }
2 26 Feb 07 jari 915         
2 26 Feb 07 jari 916         for (int i = permutedValidArray.length; i > 1; i--) {
2 26 Feb 07 jari 917             Random generator2 =new Random();
2 26 Feb 07 jari 918             //Random generator2 = new Random(randomSeeds[i - 2]);
2 26 Feb 07 jari 919             int randVal = generator2.nextInt(i - 1);
2 26 Feb 07 jari 920             int temp = permutedValidArray[randVal];
2 26 Feb 07 jari 921             permutedValidArray[randVal] = permutedValidArray[i - 1];
2 26 Feb 07 jari 922             permutedValidArray[i - 1] = temp;
2 26 Feb 07 jari 923         }  
2 26 Feb 07 jari 924         
2 26 Feb 07 jari 925         for (int i = 0; i < validArray.length; i++) {
2 26 Feb 07 jari 926             //permutedValues[validArray[i]] = permutedValues[permutedValidArray[i]];
2 26 Feb 07 jari 927             permutedValues[validArray[i]] = permutedValidArray[i];
2 26 Feb 07 jari 928         }
2 26 Feb 07 jari 929         
2 26 Feb 07 jari 930         try {
2 26 Feb 07 jari 931             Thread.sleep(10);
2 26 Feb 07 jari 932         } catch (Exception exc) {
2 26 Feb 07 jari 933             exc.printStackTrace();
2 26 Feb 07 jari 934         }
2 26 Feb 07 jari 935         
2 26 Feb 07 jari 936         
2 26 Feb 07 jari 937         return permutedValues;        
2 26 Feb 07 jari 938     }    
2 26 Feb 07 jari 939     
2 26 Feb 07 jari 940     private Vector getAdjPVals(Vector rawPVals, int adjMethod) throws AlgorithmException {
2 26 Feb 07 jari 941         event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numGenes);
2 26 Feb 07 jari 942   fireValueChanged(event);
2 26 Feb 07 jari 943   event.setId(AlgorithmEvent.PROGRESS_VALUE); 
2 26 Feb 07 jari 944         
2 26 Feb 07 jari 945         Vector adjPVals = new Vector();
2 26 Feb 07 jari 946         if (adjMethod == OneWayANOVAInitBox.JUST_ALPHA) {
2 26 Feb 07 jari 947             adjPVals = (Vector)(rawPVals.clone());
2 26 Feb 07 jari 948         } 
2 26 Feb 07 jari 949         if (adjMethod == OneWayANOVAInitBox.STD_BONFERRONI) {
2 26 Feb 07 jari 950             for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 951     if (stop) {
2 26 Feb 07 jari 952         throw new AbortException();
2 26 Feb 07 jari 953     }
2 26 Feb 07 jari 954     event.setIntValue(i);
2 26 Feb 07 jari 955     event.setDescription("Computing adjusted p-values: Current gene = " + (i + 1));  
2 26 Feb 07 jari 956     fireValueChanged(event);  
2 26 Feb 07 jari 957                 float currP = ((Float)(rawPVals.get(i))).floatValue();
2 26 Feb 07 jari 958                 float currAdjP = (float)(currP*numGenes);
2 26 Feb 07 jari 959                 if (currAdjP > 1.0f) currAdjP = 1.0f;
2 26 Feb 07 jari 960                 adjPVals.add(new Float(currAdjP));
2 26 Feb 07 jari 961             }
2 26 Feb 07 jari 962         }
2 26 Feb 07 jari 963         if (adjMethod == OneWayANOVAInitBox.ADJ_BONFERRONI) {
2 26 Feb 07 jari 964             adjPVals = getAdjBonfPVals(rawPVals);
2 26 Feb 07 jari 965         }
2 26 Feb 07 jari 966         if (adjMethod == OneWayANOVAInitBox.MAX_T) {
2 26 Feb 07 jari 967             adjPVals = getMaxTPVals();
2 26 Feb 07 jari 968         }
2 26 Feb 07 jari 969         
2 26 Feb 07 jari 970         return adjPVals;
2 26 Feb 07 jari 971     }
2 26 Feb 07 jari 972     
2 26 Feb 07 jari 973     private Vector getMaxTPVals() throws AlgorithmException {
2 26 Feb 07 jari 974         double[] origFValues = new double[numGenes];
2 26 Feb 07 jari 975         double[] descFValues = new double[numGenes];
2 26 Feb 07 jari 976         int[] descGeneIndices = new int[numGenes];
2 26 Feb 07 jari 977         double[] adjPValues = new double[numGenes];
2 26 Feb 07 jari 978         double[][] permutedRankedFValues = new double[numPerms][numGenes];
2 26 Feb 07 jari 979         double[][] uMatrix = new double[numGenes][numPerms];     
2 26 Feb 07 jari 980         
2 26 Feb 07 jari 981         event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numPerms);
2 26 Feb 07 jari 982         fireValueChanged(event);
2 26 Feb 07 jari 983         event.setId(AlgorithmEvent.PROGRESS_VALUE);    
2 26 Feb 07 jari 984         
2 26 Feb 07 jari 985         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 986             origFValues[i] = ((Float)(fValuesVector.get(i))).doubleValue();
2 26 Feb 07 jari 987         }
2 26 Feb 07 jari 988         
2 26 Feb 07 jari 989         QSort sortDescFValues = new QSort(origFValues, QSort.DESCENDING);
2 26 Feb 07 jari 990         descFValues = sortDescFValues.getSortedDouble();
2 26 Feb 07 jari 991         descGeneIndices = sortDescFValues.getOrigIndx(); 
2 26 Feb 07 jari 992         
2 26 Feb 07 jari 993         Vector validExpts = new Vector();
2 26 Feb 07 jari 994         for (int i = 0; i < groupAssignments.length; i++) {
2 26 Feb 07 jari 995             if (groupAssignments[i] != 0) validExpts.add(new Integer(i));
2 26 Feb 07 jari 996         }
2 26 Feb 07 jari 997         
2 26 Feb 07 jari 998         int[] validArray = new int[validExpts.size()];
2 26 Feb 07 jari 999         for (int j = 0; j < validArray.length; j++) {
2 26 Feb 07 jari 1000             validArray[j] = ((Integer)(validExpts.get(j))).intValue();
2 26 Feb 07 jari 1001         }        
2 26 Feb 07 jari 1002         
2 26 Feb 07 jari 1003         for (int i = 0; i < numPerms; i++) {
2 26 Feb 07 jari 1004             if (stop) {
2 26 Feb 07 jari 1005                 throw new AbortException();
2 26 Feb 07 jari 1006             }
2 26 Feb 07 jari 1007             event.setIntValue(i);
2 26 Feb 07 jari 1008             event.setDescription("Permuting matrix: Current permutation = " + (i+1));
2 26 Feb 07 jari 1009             fireValueChanged(event);   
2 26 Feb 07 jari 1010             
2 26 Feb 07 jari 1011             int[] permutedExpts = getPermutedValues(numExps, validArray);
2 26 Feb 07 jari 1012             FloatMatrix permutedMatrix = getPermutedMatrix(expMatrix, permutedExpts);
2 26 Feb 07 jari 1013             float[] currPermFVals = getPermutedFVals(permutedMatrix);  
2 26 Feb 07 jari 1014             
2 26 Feb 07 jari 1015             if (Double.isNaN(currPermFVals[descGeneIndices[numGenes - 1]])) {
2 26 Feb 07 jari 1016                 uMatrix[numGenes - 1][i] = Double.NEGATIVE_INFINITY;
2 26 Feb 07 jari 1017             } else {
2 26 Feb 07 jari 1018                 uMatrix[numGenes - 1][i] = currPermFVals[descGeneIndices[numGenes - 1]];
2 26 Feb 07 jari 1019             }   
2 26 Feb 07 jari 1020             
2 26 Feb 07 jari 1021             for (int j = numGenes - 2; j >= 0; j--) {
2 26 Feb 07 jari 1022                 if (Double.isNaN(currPermFVals[descGeneIndices[j]])) {
2 26 Feb 07 jari 1023                     uMatrix[j][i] = uMatrix[j+1][i];
2 26 Feb 07 jari 1024                 } else {
2 26 Feb 07 jari 1025                     uMatrix[j][i] = Math.max(uMatrix[j+1][i], currPermFVals[descGeneIndices[j]]);
2 26 Feb 07 jari 1026                 }
2 26 Feb 07 jari 1027                 //System.out.println("uMatrix[" + j + "][" + i + "] = " + uMatrix[j][i]);
2 26 Feb 07 jari 1028             }            
2 26 Feb 07 jari 1029             
2 26 Feb 07 jari 1030         }
2 26 Feb 07 jari 1031         
2 26 Feb 07 jari 1032         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 1033             int pCounter = 0;
2 26 Feb 07 jari 1034             for (int j = 0; j < numPerms; j++) {
2 26 Feb 07 jari 1035
2 26 Feb 07 jari 1036                 if (uMatrix[i][j] >= descFValues[i]) {
2 26 Feb 07 jari 1037                     pCounter++;
2 26 Feb 07 jari 1038                 }
2 26 Feb 07 jari 1039             }
2 26 Feb 07 jari 1040             adjPValues[descGeneIndices[i]] = (double)pCounter/(double)numPerms;
2 26 Feb 07 jari 1041         }  
2 26 Feb 07 jari 1042         
2 26 Feb 07 jari 1043         int NaNPCounter = 0;
2 26 Feb 07 jari 1044         for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 1045             if (Double.isNaN(origFValues[i])) {
2 26 Feb 07 jari 1046                 adjPValues[i] = Double.NaN;
2 26 Feb 07 jari 1047                 NaNPCounter++;
2 26 Feb 07 jari 1048                 //System.out.println("NaN index = " + i);
2 26 Feb 07 jari 1049             }             
2 26 Feb 07 jari 1050         } 
2 26 Feb 07 jari 1051         //double[] pStarValues = new double[adjPValues];
2 26 Feb 07 jari 1052         //pStartValues[descGeneIndices[0]]
2 26 Feb 07 jari 1053         for (int i = 1; i < numGenes - NaNPCounter; i++) { // enforcing monotonicity
2 26 Feb 07 jari 1054             adjPValues[descGeneIndices[i]] = Math.max(adjPValues[descGeneIndices[i]], adjPValues[descGeneIndices[i - 1]]); 
2 26 Feb 07 jari 1055         }        
2 26 Feb 07 jari 1056         
2 26 Feb 07 jari 1057         Vector adPVector = new Vector();
2 26 Feb 07 jari 1058         
2 26 Feb 07 jari 1059         for (int i = 0; i < adjPValues.length; i++) {
2 26 Feb 07 jari 1060             adPVector.add (new Float(adjPValues[i]));
2 26 Feb 07 jari 1061         }
2 26 Feb 07 jari 1062         return adPVector;
2 26 Feb 07 jari 1063     }
2 26 Feb 07 jari 1064     
2 26 Feb 07 jari 1065     
2 26 Feb 07 jari 1066     private Vector getAdjBonfPVals(Vector rawPVals) {
2 26 Feb 07 jari 1067         float[] rawPValArray = new float[rawPVals.size()];
2 26 Feb 07 jari 1068         isSig = new boolean[rawPValArray.length];
2 26 Feb 07 jari 1069         for (int i = 0; i < isSig.length; i++) {
2 26 Feb 07 jari 1070             isSig[i] = false;
2 26 Feb 07 jari 1071         }        
2 26 Feb 07 jari 1072         for (int i = 0; i < rawPValArray.length; i++) {
2 26 Feb 07 jari 1073             rawPValArray[i] = ((Float)(rawPVals.get(i))).floatValue();
2 26 Feb 07 jari 1074         }
2 26 Feb 07 jari 1075         float[] adjPValArray = new float[rawPValArray.length];
2 26 Feb 07 jari 1076         
2 26 Feb 07 jari 1077         QSort sortRawPs = new QSort(rawPValArray, QSort.ASCENDING);
2 26 Feb 07 jari 1078         float[] sortedRawPVals = sortRawPs.getSorted();
2 26 Feb 07 jari 1079         int[] origIndices = sortRawPs.getOrigIndx();
2 26 Feb 07 jari 1080         int n = numGenes;
2 26 Feb 07 jari 1081         adjPValArray[origIndices[0]] = (float)(sortedRawPVals[0]*n);        
2 26 Feb 07 jari 1082         for (int i = 1; i < numGenes; i++) {
2 26 Feb 07 jari 1083             if (sortedRawPVals[i - 1] < sortedRawPVals[i]) n--;   
2 26 Feb 07 jari 1084             if (n <= 0) n = 1;
2 26 Feb 07 jari 1085             adjPValArray[origIndices[i]] = (float)(sortedRawPVals[i]*n);
2 26 Feb 07 jari 1086         }
2 26 Feb 07 jari 1087         
2 26 Feb 07 jari 1088         for (int i = 0; i < adjPValArray.length; i++) {
2 26 Feb 07 jari 1089             if (adjPValArray[i] > 1.0f) adjPValArray[i] = 1.0f;        
2 26 Feb 07 jari 1090         }
2 26 Feb 07 jari 1091         
2 26 Feb 07 jari 1092         for (int i = 0; i < origIndices.length; i++) {// break out of loop as soon as non-significant value is encountered 
2 26 Feb 07 jari 1093             if (adjPValArray[origIndices[i]] > (double)alpha) {
2 26 Feb 07 jari 1094                 break;
2 26 Feb 07 jari 1095             } else {
2 26 Feb 07 jari 1096                 if (adjPValArray[origIndices[i]] <= (double)alpha) {
2 26 Feb 07 jari 1097                     isSig[origIndices[i]] = true;
2 26 Feb 07 jari 1098                 }
2 26 Feb 07 jari 1099             }
2 26 Feb 07 jari 1100         }        
2 26 Feb 07 jari 1101         
2 26 Feb 07 jari 1102         Vector adjBonPVals = new Vector();
2 26 Feb 07 jari 1103         for (int i = 0; i < adjPValArray.length; i++) {
2 26 Feb 07 jari 1104             adjBonPVals.add(new Float(adjPValArray[i]));
2 26 Feb 07 jari 1105         }
2 26 Feb 07 jari 1106         return adjBonPVals;
2 26 Feb 07 jari 1107     }
2 26 Feb 07 jari 1108     
2 26 Feb 07 jari 1109     /*
2 26 Feb 07 jari 1110     private Vector sortGenesBySignificance() throws AlgorithmException {
2 26 Feb 07 jari 1111   Vector sigGenes = new Vector();
2 26 Feb 07 jari 1112   Vector nonSigGenes = new Vector();
2 26 Feb 07 jari 1113         
2 26 Feb 07 jari 1114   AlgorithmEvent event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, numGenes);
2 26 Feb 07 jari 1115   fireValueChanged(event);
2 26 Feb 07 jari 1116   event.setId(AlgorithmEvent.PROGRESS_VALUE);    
2 26 Feb 07 jari 1117
2 26 Feb 07 jari 1118   if ((correctionMethod == OneWayANOVAInitBox.JUST_ALPHA)||(correctionMethod == OneWayANOVAInitBox.STD_BONFERRONI)) {
2 26 Feb 07 jari 1119       sigGenes = new Vector();
2 26 Feb 07 jari 1120       nonSigGenes = new Vector();   
2 26 Feb 07 jari 1121       for (int i = 0; i < numGenes; i++) {
2 26 Feb 07 jari 1122     if (stop) {
2 26 Feb 07 jari 1123         throw new AbortException();
2 26 Feb 07 jari 1124     }
2 26 Feb 07 jari 1125     event.setIntValue(i);
2 26 Feb 07 jari 1126     event.setDescription("Current gene = " + (i + 1));
2 26 Feb 07 jari 1127     fireValueChanged(event);  
2 26 Feb 07 jari 1128     if (isSignificant(i)) {
2 26 Feb 07 jari 1129         sigGenes.add(new Integer(i));
2 26 Feb 07 jari 1130                     fValuesVector.add(new Float(currentF));
2 26 Feb 07 jari 1131                     //pValuesVector.add(new Float(currentP));
2 26 Feb 07 jari 1132     } else {
2 26 Feb 07 jari 1133         nonSigGenes.add(new Integer(i));
2 26 Feb 07 jari 1134                     fValuesVector.add(new Float(currentF));
2 26 Feb 07 jari 1135                     //pValuesVector.add(new Float(currentP));                    
2 26 Feb 07 jari 1136     }                
2 26 Feb 07 jari 1137             }
2 26 Feb 07 jari 1138         }
2 26 Feb 07 jari 1139         
2 26 Feb 07 jari 1140   Vector sortedGenes = new Vector();
2 26 Feb 07 jari 1141   sortedGenes.add(sigGenes);
2 26 Feb 07 jari 1142   sortedGenes.add(nonSigGenes);   
2 26 Feb 07 jari 1143   
2 26 Feb 07 jari 1144         return sortedGenes;        
2 26 Feb 07 jari 1145     }
2 26 Feb 07 jari 1146      */
2 26 Feb 07 jari 1147     
2 26 Feb 07 jari 1148     private float[] getGene(int gene) {
2 26 Feb 07 jari 1149         float[] currentGene = new float[expMatrix.getColumnDimension()];
2 26 Feb 07 jari 1150         for (int i = 0; i < currentGene.length; i++) {
2 26 Feb 07 jari 1151             currentGene[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1152         }
2 26 Feb 07 jari 1153         return currentGene;
2 26 Feb 07 jari 1154     }
2 26 Feb 07 jari 1155     
2 26 Feb 07 jari 1156     private int getDfNum(int gene) {
2 26 Feb 07 jari 1157   //float[] geneValues = new float[numExps];
2 26 Feb 07 jari 1158         int n = 0;
2 26 Feb 07 jari 1159   for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1160       //geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1161             if ((!Float.isNaN(expMatrix.A[gene][i])) && (groupAssignments[i] != 0)) {
2 26 Feb 07 jari 1162                 n++;
2 26 Feb 07 jari 1163             }
2 26 Feb 07 jari 1164   } 
2 26 Feb 07 jari 1165          if (n == 0) {
2 26 Feb 07 jari 1166              return (-1); // will be exported as Float.NaN to OWAGUI
2 26 Feb 07 jari 1167          }
2 26 Feb 07 jari 1168         
2 26 Feb 07 jari 1169          return (numGroups - 1);
2 26 Feb 07 jari 1170     }
2 26 Feb 07 jari 1171     
2 26 Feb 07 jari 1172     private int getDfDenom(int gene) {
2 26 Feb 07 jari 1173         int n = 0;
2 26 Feb 07 jari 1174   for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1175       //geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1176             if ((!Float.isNaN(expMatrix.A[gene][i])) && (groupAssignments[i] != 0)) {
2 26 Feb 07 jari 1177                 n++;
2 26 Feb 07 jari 1178             }
2 26 Feb 07 jari 1179   } 
2 26 Feb 07 jari 1180          if (n == 0) {
2 26 Feb 07 jari 1181              return (-1); // will be exported as Float.NaN to OWAGUI
2 26 Feb 07 jari 1182          } 
2 26 Feb 07 jari 1183         return (n - numGroups);
2 26 Feb 07 jari 1184     }
2 26 Feb 07 jari 1185     
2 26 Feb 07 jari 1186     private int getDfNum(int gene, FloatMatrix permMatrix) {
2 26 Feb 07 jari 1187   //float[] geneValues = new float[numExps];
2 26 Feb 07 jari 1188         int n = 0;
2 26 Feb 07 jari 1189   for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1190       //geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1191             if ((!Float.isNaN(permMatrix.A[gene][i])) && (groupAssignments[i] != 0)) {
2 26 Feb 07 jari 1192                 n++;
2 26 Feb 07 jari 1193             }
2 26 Feb 07 jari 1194   } 
2 26 Feb 07 jari 1195          if (n == 0) {
2 26 Feb 07 jari 1196              return (-1); // will be exported as Float.NaN to OWAGUI
2 26 Feb 07 jari 1197          }
2 26 Feb 07 jari 1198         
2 26 Feb 07 jari 1199          return (numGroups - 1);
2 26 Feb 07 jari 1200     }
2 26 Feb 07 jari 1201     
2 26 Feb 07 jari 1202     private int getDfDenom(int gene, FloatMatrix permMatrix) {
2 26 Feb 07 jari 1203         int n = 0;
2 26 Feb 07 jari 1204   for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1205       //geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1206             if ((!Float.isNaN(permMatrix.A[gene][i])) && (groupAssignments[i] != 0)) {
2 26 Feb 07 jari 1207                 n++;
2 26 Feb 07 jari 1208             }
2 26 Feb 07 jari 1209   } 
2 26 Feb 07 jari 1210          if (n == 0) {
2 26 Feb 07 jari 1211              return (-1); // will be exported as Float.NaN to OWAGUI
2 26 Feb 07 jari 1212          } 
2 26 Feb 07 jari 1213         return (n - numGroups);
2 26 Feb 07 jari 1214     }    
2 26 Feb 07 jari 1215     
2 26 Feb 07 jari 1216     /*
2 26 Feb 07 jari 1217     private boolean isSignificant(int gene) {
2 26 Feb 07 jari 1218         boolean sig = false;
2 26 Feb 07 jari 1219   float[] geneValues = new float[numExps];
2 26 Feb 07 jari 1220         int n = 0;
2 26 Feb 07 jari 1221   for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1222       geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1223             if (!Float.isNaN(geneValues[i])) {
2 26 Feb 07 jari 1224                 n++;
2 26 Feb 07 jari 1225             }
2 26 Feb 07 jari 1226   }  
2 26 Feb 07 jari 1227         
2 26 Feb 07 jari 1228         if (n == 0) {
2 26 Feb 07 jari 1229             currentF = Float.NaN;
2 26 Feb 07 jari 1230             currentP = Float.NaN;
2 26 Feb 07 jari 1231             return false;
2 26 Feb 07 jari 1232         }
2 26 Feb 07 jari 1233         
2 26 Feb 07 jari 1234         constant = getConstant(geneValues);
2 26 Feb 07 jari 1235         
2 26 Feb 07 jari 1236         double totalSS = getTotalSS(geneValues);
2 26 Feb 07 jari 1237         double groupsSS = getGroupsSS(geneValues);
2 26 Feb 07 jari 1238         double errorSS = totalSS - groupsSS;
2 26 Feb 07 jari 1239         
2 26 Feb 07 jari 1240         if ((Double.isNaN(totalSS))||(Double.isNaN(groupsSS))||(Double.isNaN(errorSS))) {
2 26 Feb 07 jari 1241             currentF = Float.NaN;
2 26 Feb 07 jari 1242             currentP = Float.NaN;
2 26 Feb 07 jari 1243             return false;
2 26 Feb 07 jari 1244         }
2 26 Feb 07 jari 1245         
2 26 Feb 07 jari 1246         //int totalDF = validN - 1;
2 26 Feb 07 jari 1247         int groupsDF = getDfNum(gene);
2 26 Feb 07 jari 1248         int errorDF = getDfDenom(gene);
2 26 Feb 07 jari 1249         
2 26 Feb 07 jari 1250         double groupsMS = groupsSS / groupsDF;
2 26 Feb 07 jari 1251         double errorMS = errorSS / errorDF;
2 26 Feb 07 jari 1252         
2 26 Feb 07 jari 1253         double fValue = groupsMS/errorMS;
2 26 Feb 07 jari 1254         currentF = (float)(fValue);
2 26 Feb 07 jari 1255         FDistribution fDist = new FDistribution(groupsDF, errorDF);
2 26 Feb 07 jari 1256         //System.out.println("Gene " + gene + ": fValue = " + fValue + ", dfNum = " + groupsDF + ", dfDenom = " + errorDF);        
2 26 Feb 07 jari 1257         double cumulProb = fDist.cumulative(fValue);
2 26 Feb 07 jari 1258         double pValue = 1 - cumulProb; //
2 26 Feb 07 jari 1259         
2 26 Feb 07 jari 1260         if (pValue > 1) {
2 26 Feb 07 jari 1261             pValue = 1.0d;
2 26 Feb 07 jari 1262         }
2 26 Feb 07 jari 1263         
2 26 Feb 07 jari 1264         currentP = (float)pValue;
2 26 Feb 07 jari 1265         
2 26 Feb 07 jari 1266         double criticalPValue = 0;
2 26 Feb 07 jari 1267         
2 26 Feb 07 jari 1268         if (correctionMethod == OneWayANOVAInitBox.JUST_ALPHA) {
2 26 Feb 07 jari 1269             criticalPValue = (double)alpha;
2 26 Feb 07 jari 1270         } else if (correctionMethod == OneWayANOVAInitBox.STD_BONFERRONI) {
2 26 Feb 07 jari 1271             criticalPValue = (double)alpha/numGenes;
2 26 Feb 07 jari 1272         }
2 26 Feb 07 jari 1273         
2 26 Feb 07 jari 1274         //double criticalFValue = fDist.inverse(criticalPValue);
2 26 Feb 07 jari 1275         
2 26 Feb 07 jari 1276         //System.out.println("critical P = " + criticalPValue + ", dfNum = " + groupsDF + ", dfDenom = " + errorDF + ", critical F = " + criticalFValue);
2 26 Feb 07 jari 1277         
2 26 Feb 07 jari 1278         if (pValue <= criticalPValue) {
2 26 Feb 07 jari 1279             sig = true;
2 26 Feb 07 jari 1280         } else {
2 26 Feb 07 jari 1281             sig = false;
2 26 Feb 07 jari 1282         }
2 26 Feb 07 jari 1283         
2 26 Feb 07 jari 1284         /*
2 26 Feb 07 jari 1285         if (currentF >= criticalFValue) {
2 26 Feb 07 jari 1286             sig = true;
2 26 Feb 07 jari 1287         } else {
2 26 Feb 07 jari 1288             sig = false;
2 26 Feb 07 jari 1289         }
2 26 Feb 07 jari 1290          */
2 26 Feb 07 jari 1291         /*
2 26 Feb 07 jari 1292         FDistribution testF = new FDistribution(11, 15);
2 26 Feb 07 jari 1293         double testP = testF.cumulative(2.51);
2 26 Feb 07 jari 1294         
2 26 Feb 07 jari 1295         System.out.println("F(11, 15)  = 2.51, cum. probability = " + (1 - testP));
2 26 Feb 07 jari 1296          /
2 26 Feb 07 jari 1297         //System.out.println("Gene " + gene + ": GroupsSS = " + groupsSS + ", errorSS = " + errorSS + ", groupsDF = " + groupsDF + ", errorDF = " + errorDF + ", F = " + fValue + ", p = " + pValue);
2 26 Feb 07 jari 1298         
2 26 Feb 07 jari 1299         return sig;
2 26 Feb 07 jari 1300     }
2 26 Feb 07 jari 1301          */
2 26 Feb 07 jari 1302
2 26 Feb 07 jari 1303     
2 26 Feb 07 jari 1304     private double getConstant(float[] geneValues) {
2 26 Feb 07 jari 1305         double sum = 0.0d;
2 26 Feb 07 jari 1306         double cons;
2 26 Feb 07 jari 1307         int n = 0;
2 26 Feb 07 jari 1308         for (int i = 0; i < geneValues.length; i++) {
2 26 Feb 07 jari 1309             if ((!Float.isNaN(geneValues[i])) && (groupAssignments[i] != 0)) {
2 26 Feb 07 jari 1310                 sum = sum + geneValues[i];
2 26 Feb 07 jari 1311                 n++;
2 26 Feb 07 jari 1312             }
2 26 Feb 07 jari 1313         }
2 26 Feb 07 jari 1314         
2 26 Feb 07 jari 1315         if (n == 0) {
2 26 Feb 07 jari 1316             return Double.NaN;
2 26 Feb 07 jari 1317         } else {
2 26 Feb 07 jari 1318             cons = (Math.pow((double)sum, 2d))/n;
2 26 Feb 07 jari 1319         }
2 26 Feb 07 jari 1320         return cons;
2 26 Feb 07 jari 1321        
2 26 Feb 07 jari 1322     }
2 26 Feb 07 jari 1323     
2 26 Feb 07 jari 1324     private double getTotalSS(float[] geneValues) {
2 26 Feb 07 jari 1325         double ss = 0;
2 26 Feb 07 jari 1326         int n = 0;
2 26 Feb 07 jari 1327         for (int i = 0; i < geneValues.length; i++) {
2 26 Feb 07 jari 1328             if ((!Float.isNaN(geneValues[i])) && (groupAssignments[i] != 0)) {
2 26 Feb 07 jari 1329                 ss = ss + Math.pow(geneValues[i], 2);
2 26 Feb 07 jari 1330                 n++;
2 26 Feb 07 jari 1331             }
2 26 Feb 07 jari 1332         }  
2 26 Feb 07 jari 1333         
2 26 Feb 07 jari 1334         if (n == 0) {
2 26 Feb 07 jari 1335             return Double.NaN;
2 26 Feb 07 jari 1336         } else {
2 26 Feb 07 jari 1337             ss = ss - constant;
2 26 Feb 07 jari 1338         }
2 26 Feb 07 jari 1339                        
2 26 Feb 07 jari 1340         return ss;
2 26 Feb 07 jari 1341     }
2 26 Feb 07 jari 1342     
2 26 Feb 07 jari 1343     private double getGroupsSS(float[] geneValues) {
2 26 Feb 07 jari 1344         float[][] geneValuesByGroups = new float[numGroups][];
2 26 Feb 07 jari 1345         
2 26 Feb 07 jari 1346         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1347             geneValuesByGroups[i] = getGeneValuesForGroup(geneValues, i+1);
2 26 Feb 07 jari 1348         }
2 26 Feb 07 jari 1349         
2 26 Feb 07 jari 1350         double[] avSquareArray = new double[numGroups];
2 26 Feb 07 jari 1351         
2 26 Feb 07 jari 1352         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1353             avSquareArray[i] = getAvSquare(geneValuesByGroups[i]);
2 26 Feb 07 jari 1354         }
2 26 Feb 07 jari 1355         
2 26 Feb 07 jari 1356         double ss = 0;
2 26 Feb 07 jari 1357         
2 26 Feb 07 jari 1358         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1359             ss = ss + avSquareArray[i];
2 26 Feb 07 jari 1360         }
2 26 Feb 07 jari 1361         
2 26 Feb 07 jari 1362         return (ss - constant);
2 26 Feb 07 jari 1363         
2 26 Feb 07 jari 1364         //double ss = 0;
2 26 Feb 07 jari 1365         //return ss;
2 26 Feb 07 jari 1366     }
2 26 Feb 07 jari 1367     
2 26 Feb 07 jari 1368     private float[] getGeneGroupMeans(int gene) {
2 26 Feb 07 jari 1369   float[] geneValues = new float[numExps];
2 26 Feb 07 jari 1370         for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1371       geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1372   } 
2 26 Feb 07 jari 1373         
2 26 Feb 07 jari 1374         float[][] geneValuesByGroups = new float[numGroups][];
2 26 Feb 07 jari 1375         
2 26 Feb 07 jari 1376         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1377             geneValuesByGroups[i] = getGeneValuesForGroup(geneValues, i+1);
2 26 Feb 07 jari 1378         } 
2 26 Feb 07 jari 1379         
2 26 Feb 07 jari 1380         float[] geneGroupMeans = new float[numGroups];
2 26 Feb 07 jari 1381         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1382             geneGroupMeans[i] = getMean(geneValuesByGroups[i]);
2 26 Feb 07 jari 1383         }
2 26 Feb 07 jari 1384         
2 26 Feb 07 jari 1385         return geneGroupMeans;
2 26 Feb 07 jari 1386     }
2 26 Feb 07 jari 1387     
2 26 Feb 07 jari 1388     private float[] getGeneGroupSDs(int gene) {
2 26 Feb 07 jari 1389   float[] geneValues = new float[numExps];        
2 26 Feb 07 jari 1390         for (int i = 0; i < numExps; i++) {
2 26 Feb 07 jari 1391       geneValues[i] = expMatrix.A[gene][i];
2 26 Feb 07 jari 1392   }   
2 26 Feb 07 jari 1393         
2 26 Feb 07 jari 1394         float[][] geneValuesByGroups = new float[numGroups][];
2 26 Feb 07 jari 1395         
2 26 Feb 07 jari 1396         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1397             geneValuesByGroups[i] = getGeneValuesForGroup(geneValues, i+1);
2 26 Feb 07 jari 1398         }  
2 26 Feb 07 jari 1399         
2 26 Feb 07 jari 1400         float[] geneGroupSDs = new float[numGroups];
2 26 Feb 07 jari 1401         for (int i = 0; i < numGroups; i++) {
2 26 Feb 07 jari 1402             geneGroupSDs[i] = getStdDev(geneValuesByGroups[i]);
2 26 Feb 07 jari 1403         }
2 26 Feb 07 jari 1404         
2 26 Feb 07 jari 1405         return geneGroupSDs;        
2 26 Feb 07 jari 1406     }
2 26 Feb 07 jari 1407     
2 26 Feb 07 jari 1408     private float[] getGeneValuesForGroup(float[] geneValues, int group) {
2 26 Feb 07 jari 1409         Vector groupValuesVector = new Vector();
2 26 Feb 07 jari 1410         
2 26 Feb 07 jari 1411         for (int i = 0; i < groupAssignments.length; i++) {
2 26 Feb 07 jari 1412             if (groupAssignments[i] == group) {
2 26 Feb 07 jari 1413                 groupValuesVector.add(new Float(geneValues[i]));
2 26 Feb 07 jari 1414             }
2 26 Feb 07 jari 1415         }
2 26 Feb 07 jari 1416         
2 26 Feb 07 jari 1417         float[] groupGeneValues = new float[groupValuesVector.size()];
2 26 Feb 07 jari 1418         
2 26 Feb 07 jari 1419         for (int i = 0; i < groupValuesVector.size(); i++) {
2 26 Feb 07 jari 1420             groupGeneValues[i] = ((Float)(groupValuesVector.get(i))).floatValue();
2 26 Feb 07 jari 1421         }
2 26 Feb 07 jari 1422         
2 26 Feb 07 jari 1423         return groupGeneValues;
2 26 Feb 07 jari 1424     }
2 26 Feb 07 jari 1425     
2 26 Feb 07 jari 1426     private FloatMatrix getAllGeneGroupMeans() {
2 26 Feb 07 jari 1427         FloatMatrix means = new FloatMatrix(numGenes, numGroups);
2 26 Feb 07 jari 1428         for (int i = 0; i < means.getRowDimension(); i++) {
2 26 Feb 07 jari 1429             means.A[i] = getGeneGroupMeans(i);
2 26 Feb 07 jari 1430         }
2 26 Feb 07 jari 1431         return means;
2 26 Feb 07 jari 1432     }
2 26 Feb 07 jari 1433     
2 26 Feb 07 jari 1434     private FloatMatrix getAllGeneGroupSDs() {
2 26 Feb 07 jari 1435         FloatMatrix sds = new FloatMatrix(numGenes, numGroups);
2 26 Feb 07 jari 1436         for (int i = 0; i < sds.getRowDimension(); i++) {
2 26 Feb 07 jari 1437             sds.A[i] = getGeneGroupSDs(i);
2 26 Feb 07 jari 1438         }
2 26 Feb 07 jari 1439         return sds;        
2 26 Feb 07 jari 1440     }
2 26 Feb 07 jari 1441     
2 26 Feb 07 jari 1442     private double getAvSquare(float[] values) {
2 26 Feb 07 jari 1443         double ss = 0;
2 26 Feb 07 jari 1444         double sum = 0;
2 26 Feb 07 jari 1445         int n = 0;
2 26 Feb 07 jari 1446         for (int i = 0; i < values.length; i++) {
2 26 Feb 07 jari 1447             if (!Float.isNaN(values[i])) {
2 26 Feb 07 jari 1448                 sum  = sum + values[i];
2 26 Feb 07 jari 1449                 n++;
2 26 Feb 07 jari 1450             }
2 26 Feb 07 jari 1451         }
2 26 Feb 07 jari 1452         
2 26 Feb 07 jari 1453         if (n == 0) {
2 26 Feb 07 jari 1454             return Double.NaN;
2 26 Feb 07 jari 1455         } else {
2 26 Feb 07 jari 1456             ss = (Math.pow(sum, 2)) / n;
2 26 Feb 07 jari 1457         }
2 26 Feb 07 jari 1458         
2 26 Feb 07 jari 1459         return ss;
2 26 Feb 07 jari 1460     }
2 26 Feb 07 jari 1461     
2 26 Feb 07 jari 1462     private float getMean(float[] group) {
2 26 Feb 07 jari 1463   float sum = 0;
2 26 Feb 07 jari 1464   int n = 0;
2 26 Feb 07 jari 1465   
2 26 Feb 07 jari 1466   for (int i = 0; i < group.length; i++) {
2 26 Feb 07 jari 1467       //System.out.println("getMean(): group[" + i + "] = " + group[i]);
2 26 Feb 07 jari 1468       if (!Float.isNaN(group[i])) {
2 26 Feb 07 jari 1469     sum = sum + group[i];
2 26 Feb 07 jari 1470     n++;
2 26 Feb 07 jari 1471       }
2 26 Feb 07 jari 1472   }
2 26 Feb 07 jari 1473   
2 26 Feb 07 jari 1474   //System.out.println("getMean(): sum = " +sum);
2 26 Feb 07 jari 1475   if (n == 0) {
2 26 Feb 07 jari 1476             return Float.NaN;
2 26 Feb 07 jari 1477         }
2 26 Feb 07 jari 1478   float mean =  sum / (float)n;
2 26 Feb 07 jari 1479         
2 26 Feb 07 jari 1480         if (Float.isInfinite(mean)) {
2 26 Feb 07 jari 1481             return Float.NaN;
2 26 Feb 07 jari 1482         }
2 26 Feb 07 jari 1483         
2 26 Feb 07 jari 1484   return mean;
2 26 Feb 07 jari 1485     }
2 26 Feb 07 jari 1486     
2 26 Feb 07 jari 1487     private float getStdDev(float[] group) {
2 26 Feb 07 jari 1488   float mean = getMean(group);
2 26 Feb 07 jari 1489   int n = 0;
2 26 Feb 07 jari 1490   
2 26 Feb 07 jari 1491   float sumSquares = 0;
2 26 Feb 07 jari 1492   
2 26 Feb 07 jari 1493   for (int i = 0; i < group.length; i++) {
2 26 Feb 07 jari 1494       if (!Float.isNaN(group[i])) {
2 26 Feb 07 jari 1495     sumSquares = (float)(sumSquares + Math.pow((group[i] - mean), 2));
2 26 Feb 07 jari 1496     n++;
2 26 Feb 07 jari 1497       }
2 26 Feb 07 jari 1498   }
2 26 Feb 07 jari 1499         
2 26 Feb 07 jari 1500         if (n < 2) {
2 26 Feb 07 jari 1501             return Float.NaN;
2 26 Feb 07 jari 1502         }
2 26 Feb 07 jari 1503   
2 26 Feb 07 jari 1504   float var = sumSquares / (float)(n - 1);
2 26 Feb 07 jari 1505   if (Float.isInfinite(var)) {
2 26 Feb 07 jari 1506             return Float.NaN;
2 26 Feb 07 jari 1507         }
2 26 Feb 07 jari 1508   return (float)(Math.sqrt((double)var));
2 26 Feb 07 jari 1509     }    
2 26 Feb 07 jari 1510     
2 26 Feb 07 jari 1511 }

Generated on Thu Jul 20 09:25:34 2023 (UTC) by svndigest 0.11.2pre (dev build 1672)