mev-4.0.01/source/org/tigr/microarray/mev/cluster/algorithm/impl/RN.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: RN.java,v $
2 26 Feb 07 jari 7  * $Revision: 1.5 $
2 26 Feb 07 jari 8  * $Date: 2006/04/24 17:13:43 $
2 26 Feb 07 jari 9  * $Author: eleanorahowe $
2 26 Feb 07 jari 10  * $State: Exp $
2 26 Feb 07 jari 11  */
2 26 Feb 07 jari 12 package org.tigr.microarray.mev.cluster.algorithm.impl;
2 26 Feb 07 jari 13
2 26 Feb 07 jari 14 import java.util.Arrays;
2 26 Feb 07 jari 15
2 26 Feb 07 jari 16 import org.tigr.microarray.mev.cluster.Cluster;
2 26 Feb 07 jari 17 import org.tigr.microarray.mev.cluster.Node;
2 26 Feb 07 jari 18 import org.tigr.microarray.mev.cluster.NodeList;
2 26 Feb 07 jari 19 import org.tigr.microarray.mev.cluster.algorithm.AbortException;
2 26 Feb 07 jari 20 import org.tigr.microarray.mev.cluster.algorithm.AbstractAlgorithm;
2 26 Feb 07 jari 21 import org.tigr.microarray.mev.cluster.algorithm.Algorithm;
2 26 Feb 07 jari 22 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmData;
2 26 Feb 07 jari 23 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmEvent;
2 26 Feb 07 jari 24 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmException;
2 26 Feb 07 jari 25 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmListener;
2 26 Feb 07 jari 26 import org.tigr.microarray.mev.cluster.algorithm.AlgorithmParameters;
2 26 Feb 07 jari 27 import org.tigr.microarray.mev.cluster.algorithm.impl.util.FloatArray;
2 26 Feb 07 jari 28 import org.tigr.microarray.mev.cluster.algorithm.impl.util.IntArray;
2 26 Feb 07 jari 29 import org.tigr.microarray.mev.cluster.algorithm.impl.util.IntSorter;
2 26 Feb 07 jari 30 import org.tigr.util.FloatMatrix;
2 26 Feb 07 jari 31
2 26 Feb 07 jari 32
2 26 Feb 07 jari 33 public class RN extends AbstractAlgorithm {
2 26 Feb 07 jari 34
2 26 Feb 07 jari 35     private static final int c_DecileCount = 10;
2 26 Feb 07 jari 36     public static final double LOG2 = Math.log(2.0);
2 26 Feb 07 jari 37
2 26 Feb 07 jari 38     private Algorithm permAlgo;
2 26 Feb 07 jari 39     private boolean stop = false;
2 26 Feb 07 jari 40     private int number_of_samples;
2 26 Feb 07 jari 41     private FloatMatrix expMatrix;
2 26 Feb 07 jari 42     
2 26 Feb 07 jari 43     public AlgorithmData execute(AlgorithmData data) throws AlgorithmException {
2 26 Feb 07 jari 44
2 26 Feb 07 jari 45         expMatrix = data.getMatrix("experiment");
2 26 Feb 07 jari 46         if (expMatrix == null)
2 26 Feb 07 jari 47             return null;
2 26 Feb 07 jari 48
2 26 Feb 07 jari 49         AlgorithmParameters map = data.getParams();
2 26 Feb 07 jari 50         number_of_samples = expMatrix.getColumnDimension();
2 26 Feb 07 jari 51         boolean use_permutation = map.getBoolean("use-permutation");
2 26 Feb 07 jari 52         float perm_threshold = 0.8f;
2 26 Feb 07 jari 53         if (use_permutation) {
2 26 Feb 07 jari 54             this.permAlgo = new PermutationTest();
2 26 Feb 07 jari 55             permAlgo.addAlgorithmListener(new SubAlgoListener());
2 26 Feb 07 jari 56             AlgorithmData permResult = permAlgo.execute(data);
2 26 Feb 07 jari 57             perm_threshold = permResult.getParams().getFloat("threshold", perm_threshold);
2 26 Feb 07 jari 58         }
2 26 Feb 07 jari 59
2 26 Feb 07 jari 60         int function = map.getInt("distance-function", PEARSON);
2 26 Feb 07 jari 61         float factor = map.getFloat("distance-factor", 1.0f);
2 26 Feb 07 jari 62         boolean absolute = map.getBoolean("distance-absolute", true);
2 26 Feb 07 jari 63         float min_threshold = use_permutation ? perm_threshold : map.getFloat("min-threshold", 0.8f);
2 26 Feb 07 jari 64         float max_threshold = map.getFloat("max-threshold", 1.0f);
2 26 Feb 07 jari 65         boolean bFilterByEntropy = map.getBoolean("filter-by-entropy");
2 26 Feb 07 jari 66         float fltTopNPercent=map.getFloat("top-n-percent", 100f);
2 26 Feb 07 jari 67         if (fltTopNPercent < 0f || fltTopNPercent>100f) {
2 26 Feb 07 jari 68             throw new AlgorithmException("Filter value is out of range (0, 100)%");
2 26 Feb 07 jari 69         }
2 26 Feb 07 jari 70
2 26 Feb 07 jari 71         int nGenes = expMatrix.getRowDimension();
2 26 Feb 07 jari 72         int filteredSize = nGenes;
2 26 Feb 07 jari 73
2 26 Feb 07 jari 74         int[] entropyIndices = new int[nGenes];
2 26 Feb 07 jari 75         for (int i=0; i<entropyIndices.length; i++) {
2 26 Feb 07 jari 76             entropyIndices[i] = i;
2 26 Feb 07 jari 77         }
2 26 Feb 07 jari 78         if (bFilterByEntropy) {
2 26 Feb 07 jari 79             double[] entropyValues = new double[nGenes];
2 26 Feb 07 jari 80             for (int i=0; i<entropyValues.length; i++) {
2 26 Feb 07 jari 81                 entropyValues[i] = getEntropy(expMatrix.A[i]);
2 26 Feb 07 jari 82             }
2 26 Feb 07 jari 83             IntSorter.sort(entropyIndices, new RelNetComparator(entropyValues));
2 26 Feb 07 jari 84             filteredSize = (int)((float)nGenes*fltTopNPercent/100f);
2 26 Feb 07 jari 85         }
2 26 Feb 07 jari 86
2 26 Feb 07 jari 87         AlgorithmEvent event = new AlgorithmEvent(this, AlgorithmEvent.SET_UNITS, 100, "Calculating Network");
2 26 Feb 07 jari 88         fireValueChanged(event);
2 26 Feb 07 jari 89
2 26 Feb 07 jari 90         event.setId(AlgorithmEvent.PROGRESS_VALUE);
2 26 Feb 07 jari 91
2 26 Feb 07 jari 92         int progress = 0;
2 26 Feb 07 jari 93         int links = 0;
2 26 Feb 07 jari 94         int sum = filteredSize*(filteredSize+1)/2;
2 26 Feb 07 jari 95         int step = sum/100+1;
2 26 Feb 07 jari 96         float value;
2 26 Feb 07 jari 97
2 26 Feb 07 jari 98         IntArray[] indices = new IntArray[nGenes];
2 26 Feb 07 jari 99         FloatArray[] weights = new FloatArray[nGenes];
2 26 Feb 07 jari 100         for (int i=0; i<indices.length; i++) {
2 26 Feb 07 jari 101             indices[i] = new IntArray();
2 26 Feb 07 jari 102             indices[i].add(i);
2 26 Feb 07 jari 103             weights[i] = new FloatArray();
2 26 Feb 07 jari 104             weights[i].add(0f);
2 26 Feb 07 jari 105         }
2 26 Feb 07 jari 106         for (int nCurIndOuter=0; nCurIndOuter<filteredSize; nCurIndOuter++) {
2 26 Feb 07 jari 107             if (this.stop) {
2 26 Feb 07 jari 108                 throw new AbortException();
2 26 Feb 07 jari 109             }
2 26 Feb 07 jari 110             progress++;
2 26 Feb 07 jari 111             // calculation
2 26 Feb 07 jari 112             for (int nCurIndInner=nCurIndOuter+1; nCurIndInner<filteredSize; nCurIndInner++) {
2 26 Feb 07 jari 113                 value = ExperimentUtil.geneDistance(expMatrix, null, entropyIndices[nCurIndOuter], entropyIndices[nCurIndInner], function, factor, absolute);
2 26 Feb 07 jari 114                 value = value*value; // = abs(r^2)
2 26 Feb 07 jari 115                 if (value >= min_threshold && value <= max_threshold) {
2 26 Feb 07 jari 116                     links++;
2 26 Feb 07 jari 117                     indices[entropyIndices[nCurIndOuter]].add(entropyIndices[nCurIndInner]);
2 26 Feb 07 jari 118                     indices[entropyIndices[nCurIndInner]].add(entropyIndices[nCurIndOuter]);
2 26 Feb 07 jari 119                     weights[entropyIndices[nCurIndOuter]].add(value);
2 26 Feb 07 jari 120                     weights[entropyIndices[nCurIndInner]].add(value);
2 26 Feb 07 jari 121                 }
2 26 Feb 07 jari 122                 // progress events handling
2 26 Feb 07 jari 123                 progress++;
2 26 Feb 07 jari 124                 if (progress%step == 0) {
2 26 Feb 07 jari 125                     event.setIntValue(progress/step);
2 26 Feb 07 jari 126                     event.setDescription("Calculating Network ("+String.valueOf(links)+" links found)");
2 26 Feb 07 jari 127                     fireValueChanged(event);
2 26 Feb 07 jari 128                 }
2 26 Feb 07 jari 129             }
2 26 Feb 07 jari 130         }
2 26 Feb 07 jari 131         // return the result
2 26 Feb 07 jari 132         AlgorithmData result = new AlgorithmData();
2 26 Feb 07 jari 133
2 26 Feb 07 jari 134         Cluster cluster;
2 26 Feb 07 jari 135         NodeList clusterNodeList;
2 26 Feb 07 jari 136         // copy nodes to the cluster structure
2 26 Feb 07 jari 137         cluster = new Cluster();
2 26 Feb 07 jari 138         clusterNodeList = cluster.getNodeList();
2 26 Feb 07 jari 139         clusterNodeList.ensureCapacity(nGenes);
2 26 Feb 07 jari 140         
2 26 Feb 07 jari 141         FloatMatrix means = getMeans(indices);        
2 26 Feb 07 jari 142         result.addMatrix("means", means);
2 26 Feb 07 jari 143         result.addMatrix("variances", this.getVariances(indices, means));
2 26 Feb 07 jari 144         
2 26 Feb 07 jari 145         for (int i=0; i<nGenes; i++) {
2 26 Feb 07 jari 146             clusterNodeList.addNode(new Node(indices[i].toArray()));
2 26 Feb 07 jari 147             indices[i] = null; // gc
2 26 Feb 07 jari 148         }
2 26 Feb 07 jari 149         result.addCluster("cluster", cluster);
2 26 Feb 07 jari 150         // copy weights to the cluster structure
2 26 Feb 07 jari 151         cluster = new Cluster();
2 26 Feb 07 jari 152         clusterNodeList = cluster.getNodeList();
2 26 Feb 07 jari 153         clusterNodeList.ensureCapacity(nGenes);
2 26 Feb 07 jari 154         for (int i=0; i<nGenes; i++) {
2 26 Feb 07 jari 155             clusterNodeList.addNode(new Node(float2int(weights[i].toArray())));
2 26 Feb 07 jari 156             weights[i] = null; // gc
2 26 Feb 07 jari 157         }
2 26 Feb 07 jari 158         result.addCluster("weights", cluster);
2 26 Feb 07 jari 159
2 26 Feb 07 jari 160         result.addParam("links", String.valueOf(links));
2 26 Feb 07 jari 161         result.addParam("min_threshold", String.valueOf(min_threshold));
2 26 Feb 07 jari 162         return result;
2 26 Feb 07 jari 163     }
2 26 Feb 07 jari 164
2 26 Feb 07 jari 165     public void abort() {
2 26 Feb 07 jari 166         this.stop = true;
2 26 Feb 07 jari 167         if (this.permAlgo != null)
2 26 Feb 07 jari 168             this.permAlgo.abort();
2 26 Feb 07 jari 169     }
2 26 Feb 07 jari 170
2 26 Feb 07 jari 171     private static int[] float2int(float[] floats) {
2 26 Feb 07 jari 172         if (floats == null)
2 26 Feb 07 jari 173             return null;
2 26 Feb 07 jari 174         int[] ints = new int[floats.length];
2 26 Feb 07 jari 175         for (int i=0; i<ints.length; i++)
2 26 Feb 07 jari 176             ints[i] = Float.floatToRawIntBits(floats[i]);
2 26 Feb 07 jari 177         return ints;
2 26 Feb 07 jari 178     }
2 26 Feb 07 jari 179
2 26 Feb 07 jari 180     private double getEntropy(float[] pVector) {
2 26 Feb 07 jari 181         double fltMin = Double.MAX_VALUE;
2 26 Feb 07 jari 182         double fltMax = -Double.MAX_VALUE;
2 26 Feb 07 jari 183         int i=0;
2 26 Feb 07 jari 184         int[] arrDeciles = new int[c_DecileCount]; 
2 26 Feb 07 jari 185
2 26 Feb 07 jari 186         final int iSize = pVector.length;
2 26 Feb 07 jari 187         int iValCount = 0;
2 26 Feb 07 jari 188         for (i=0; i<iSize; i++) {
2 26 Feb 07 jari 189             if (Double.isNaN(pVector[i]))
2 26 Feb 07 jari 190                 continue;
2 26 Feb 07 jari 191             fltMin = Math.min(fltMin, pVector[i]);
2 26 Feb 07 jari 192             fltMax = Math.max(fltMax, pVector[i]);
2 26 Feb 07 jari 193             iValCount++;
2 26 Feb 07 jari 194         }
2 26 Feb 07 jari 195
2 26 Feb 07 jari 196         double fltStep = (fltMax-fltMin)/(c_DecileCount);
2 26 Feb 07 jari 197         if (fltStep == 0d) {
2 26 Feb 07 jari 198             return -1.0*Math.log(1.0)/LOG2;
2 26 Feb 07 jari 199         }
2 26 Feb 07 jari 200
2 26 Feb 07 jari 201         if (fltMin == Double.MAX_VALUE)
2 26 Feb 07 jari 202             return 0d;
2 26 Feb 07 jari 203
2 26 Feb 07 jari 204         Arrays.fill(arrDeciles, 0);
2 26 Feb 07 jari 205         for (i=0; i<iSize; i++) {
2 26 Feb 07 jari 206             if (Double.isNaN(pVector[i]))
2 26 Feb 07 jari 207                 continue;
2 26 Feb 07 jari 208             int iDecileInd = (int)Math.ceil((pVector[i]-fltMin)/fltStep)-1;
2 26 Feb 07 jari 209             if (iDecileInd < 0) {
2 26 Feb 07 jari 210                 iDecileInd = 0;
2 26 Feb 07 jari 211             }
2 26 Feb 07 jari 212             arrDeciles[iDecileInd]++;
2 26 Feb 07 jari 213         }
2 26 Feb 07 jari 214         if (iValCount == 0)
2 26 Feb 07 jari 215             return 0d;
2 26 Feb 07 jari 216
2 26 Feb 07 jari 217         // finally, calculate entropy
2 26 Feb 07 jari 218         double dblEntropy=0;
2 26 Feb 07 jari 219
2 26 Feb 07 jari 220         for (i=0; i<c_DecileCount; i++) {
2 26 Feb 07 jari 221             if (arrDeciles[i] == 0) {
2 26 Feb 07 jari 222                 continue;
2 26 Feb 07 jari 223             }
2 26 Feb 07 jari 224             double dblPx=((double)arrDeciles[i])/iValCount;
2 26 Feb 07 jari 225             dblEntropy += dblPx*Math.log(dblPx)/LOG2; // log2(x)==log(x)/log(2)
2 26 Feb 07 jari 226         }
2 26 Feb 07 jari 227         return -dblEntropy;
2 26 Feb 07 jari 228     }
2 26 Feb 07 jari 229
2 26 Feb 07 jari 230     
2 26 Feb 07 jari 231
2 26 Feb 07 jari 232     
2 26 Feb 07 jari 233     private FloatMatrix getMeans(IntArray [] clusters) {
2 26 Feb 07 jari 234   FloatMatrix means = new FloatMatrix(clusters.length, number_of_samples);
2 26 Feb 07 jari 235   FloatMatrix mean; 
2 26 Feb 07 jari 236   for (int i=0; i<clusters.length; i++) {
2 26 Feb 07 jari 237       mean = getMean(clusters[i]);
2 26 Feb 07 jari 238       means.A[i] = mean.A[0];
2 26 Feb 07 jari 239   }
2 26 Feb 07 jari 240   return means;
2 26 Feb 07 jari 241     }
2 26 Feb 07 jari 242     
2 26 Feb 07 jari 243     private FloatMatrix getMean(IntArray cluster) {
2 26 Feb 07 jari 244   FloatMatrix mean = new FloatMatrix(1, number_of_samples);
2 26 Feb 07 jari 245   float currentMean;
2 26 Feb 07 jari 246   int n = cluster.getSize();
2 26 Feb 07 jari 247   int denom = 0;
2 26 Feb 07 jari 248   float value;
2 26 Feb 07 jari 249   for (int i=0; i<number_of_samples; i++) {
2 26 Feb 07 jari 250       currentMean = 0f;
2 26 Feb 07 jari 251       denom  = 0;
2 26 Feb 07 jari 252       for (int j=0; j<n; j++) {
2 26 Feb 07 jari 253     value = expMatrix.get(cluster.get(j), i);
2 26 Feb 07 jari 254     if (!Float.isNaN(value)) {
2 26 Feb 07 jari 255         currentMean += value;
2 26 Feb 07 jari 256         denom++;
2 26 Feb 07 jari 257     }
2 26 Feb 07 jari 258       }
2 26 Feb 07 jari 259       mean.set(0, i, currentMean/(float)denom);
2 26 Feb 07 jari 260   }
2 26 Feb 07 jari 261   
2 26 Feb 07 jari 262   return mean;
2 26 Feb 07 jari 263     }
2 26 Feb 07 jari 264     
2 26 Feb 07 jari 265         private FloatMatrix getVariances(IntArray [] clusters, FloatMatrix means) {
2 26 Feb 07 jari 266   final int rows = means.getRowDimension();
2 26 Feb 07 jari 267   final int columns = means.getColumnDimension();
2 26 Feb 07 jari 268   FloatMatrix variances = new FloatMatrix(rows, columns);
2 26 Feb 07 jari 269   for (int row=0; row<rows; row++) {
2 26 Feb 07 jari 270       for (int column=0; column<columns; column++) {
2 26 Feb 07 jari 271     variances.set(row, column, getSampleVariance(clusters[row], column, means.get(row, column)));
2 26 Feb 07 jari 272       }
2 26 Feb 07 jari 273   }
2 26 Feb 07 jari 274   return variances;
2 26 Feb 07 jari 275     }
2 26 Feb 07 jari 276     
2 26 Feb 07 jari 277     int validN;
2 26 Feb 07 jari 278     
2 26 Feb 07 jari 279     private float getSampleNormalizedSum(IntArray cluster, int column, float mean) {
2 26 Feb 07 jari 280   final int size = cluster.getSize();
2 26 Feb 07 jari 281   float sum = 0f;
2 26 Feb 07 jari 282   validN = 0;
2 26 Feb 07 jari 283   float value;
2 26 Feb 07 jari 284   for (int i=0; i<size; i++) {
2 26 Feb 07 jari 285       value = expMatrix.get(cluster.get(i), column);
2 26 Feb 07 jari 286       if (!Float.isNaN(value)) {
2 26 Feb 07 jari 287     sum += Math.pow(value-mean, 2);
2 26 Feb 07 jari 288     validN++;
2 26 Feb 07 jari 289       }
2 26 Feb 07 jari 290   }
2 26 Feb 07 jari 291   return sum;
2 26 Feb 07 jari 292     }
2 26 Feb 07 jari 293     
2 26 Feb 07 jari 294     private float getSampleVariance(IntArray cluster, int column, float mean) {
2 26 Feb 07 jari 295   return(float)Math.sqrt(getSampleNormalizedSum(cluster, column, mean)/(float)(validN-1));
2 26 Feb 07 jari 296     }
2 26 Feb 07 jari 297     
2 26 Feb 07 jari 298     
2 26 Feb 07 jari 299     class SubAlgoListener implements AlgorithmListener {
2 26 Feb 07 jari 300         public void valueChanged(AlgorithmEvent event) {
2 26 Feb 07 jari 301             fireValueChanged(event);
2 26 Feb 07 jari 302         }
2 26 Feb 07 jari 303     }
2 26 Feb 07 jari 304 }