src/jalview/math/SparseMatrix.java

   1 /*
   2  * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$)
   3  * Copyright (C) $$Year-Rel$$ The Jalview Authors
   4  *
   5  * This file is part of Jalview.
   6  *
   7  * Jalview is free software: you can redistribute it and/or
   8  * modify it under the terms of the GNU General Public License
   9  * as published by the Free Software Foundation, either version 3
  10  * of the License, or (at your option) any later version.
  11  *
  12  * Jalview is distributed in the hope that it will be useful, but
  13  * WITHOUT ANY WARRANTY; without even the implied warranty
  14  * of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  15  * PURPOSE.  See the GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
  19  * The Jalview Authors are detailed in the 'AUTHORS' file.
  20  */
  21 package jalview.math;
  22
  23 import jalview.ext.android.SparseDoubleArray;
  24
  25 /**
  26  * A variant of Matrix intended for use for sparse (mostly zero) matrices. This
  27  * class uses a SparseDoubleArray to hold each row of the matrix. The sparse
  28  * array only stores non-zero values. This gives a smaller memory footprint, and
  29  * fewer matrix calculation operations, for mostly zero matrices.
  30  *
  31  * @author gmcarstairs
  32  */
  33 public class SparseMatrix extends Matrix
  34 {
  35   /*
  36    * we choose columns for the sparse arrays as this allows
  37    * optimisation of the preMultiply() method used in PCA.run()
  38    */
  39   SparseDoubleArray[] sparseColumns;
  40
  41   /**
  42    * Constructor given data in [row][column] order
  43    *
  44    * @param v
  45    */
  46   public SparseMatrix(double[][] v)
  47   {
  48     rows = v.length;
  49     if (rows > 0)
  50     {
  51       cols = v[0].length;
  52     }
  53     sparseColumns = new SparseDoubleArray[cols];
  54
  55     /*
  56      * transpose v[row][col] into [col][row] order
  57      */
  58     for (int col = 0; col < cols; col++)
  59     {
  60       SparseDoubleArray sparseColumn = new SparseDoubleArray();
  61       sparseColumns[col] = sparseColumn;
  62       for (int row = 0; row < rows; row++)
  63       {
  64         double value = v[row][col];
  65         if (value != 0d)
  66         {
  67           sparseColumn.put(row, value);
  68         }
  69       }
  70     }
  71   }
  72
  73   /**
  74    * Answers the value at row i, column j
  75    */
  76   @Override
  77   public double getValue(int i, int j)
  78   {
  79     return sparseColumns[j].get(i);
  80   }
  81
  82   /**
  83    * Sets the value at row i, column j to val
  84    */
  85   @Override
  86   public void setValue(int i, int j, double val)
  87   {
  88     if (val == 0d)
  89     {
  90       sparseColumns[j].delete(i);
  91     }
  92     else
  93     {
  94       sparseColumns[j].put(i, val);
  95     }
  96   }
  97
  98   @Override
  99   public double[] getColumn(int i)
 100   {
 101     double[] col = new double[height()];
 102
 103     SparseDoubleArray vals = sparseColumns[i];
 104     for (int nonZero = 0; nonZero < vals.size(); nonZero++)
 105     {
 106       col[vals.keyAt(nonZero)] = vals.valueAt(nonZero);
 107     }
 108     return col;
 109   }
 110
 111   @Override
 112   public MatrixI copy()
 113   {
 114     double[][] vals = new double[height()][width()];
 115     for (int i = 0; i < height(); i++)
 116     {
 117       vals[i] = getRow(i);
 118     }
 119     return new SparseMatrix(vals);
 120   }
 121
 122   @Override
 123   public MatrixI transpose()
 124   {
 125     double[][] out = new double[cols][rows];
 126
 127     /*
 128      * for each column...
 129      */
 130     for (int i = 0; i < cols; i++)
 131     {
 132       /*
 133        * put non-zero values into the corresponding row
 134        * of the transposed matrix
 135        */
 136       SparseDoubleArray vals = sparseColumns[i];
 137       for (int nonZero = 0; nonZero < vals.size(); nonZero++)
 138       {
 139         out[i][vals.keyAt(nonZero)] = vals.valueAt(nonZero);
 140       }
 141     }
 142
 143     return new SparseMatrix(out);
 144   }
 145
 146   /**
 147    * Answers a new matrix which is the product in.this. If the product contains
 148    * less than 20% non-zero values, it is returned as a SparseMatrix, else as a
 149    * Matrix.
 150    * <p>
 151    * This method is optimised for the sparse arrays which store column values
 152    * for a SparseMatrix. Note that postMultiply is not so optimised. That would
 153    * require redundantly also storing sparse arrays for the rows, which has not
 154    * been done. Currently only preMultiply is used in Jalview.
 155    */
 156   @Override
 157   public MatrixI preMultiply(MatrixI in)
 158   {
 159     if (in.width() != rows)
 160     {
 161       throw new IllegalArgumentException("Can't pre-multiply " + this.rows
 162               + " rows by " + in.width() + " columns");
 163     }
 164     double[][] tmp = new double[in.height()][this.cols];
 165
 166     long count = 0L;
 167     for (int i = 0; i < in.height(); i++)
 168     {
 169       for (int j = 0; j < this.cols; j++)
 170       {
 171         /*
 172          * result[i][j] is the vector product of
 173          * in.row[i] and this.column[j]
 174          * we only need to use non-zero values from the column
 175          */
 176         SparseDoubleArray vals = sparseColumns[j];
 177         boolean added = false;
 178         for (int nonZero = 0; nonZero < vals.size(); nonZero++)
 179         {
 180           int myRow = vals.keyAt(nonZero);
 181           double myValue = vals.valueAt(nonZero);
 182           tmp[i][j] += (in.getValue(i, myRow) * myValue);
 183           added = true;
 184         }
 185         if (added && tmp[i][j] != 0d)
 186         {
 187           count++; // non-zero entry in product
 188         }
 189       }
 190     }
 191
 192     /*
 193      * heuristic rule - if product is more than 80% zero
 194      * then construct a SparseMatrix, else a Matrix
 195      */
 196     if (count * 5 < in.height() * cols)
 197     {
 198       return new SparseMatrix(tmp);
 199     }
 200     else
 201     {
 202       return new Matrix(tmp);
 203     }
 204   }
 205
 206   @Override
 207   protected double divideValue(int i, int j, double divisor)
 208   {
 209     if (divisor == 0d)
 210     {
 211       return getValue(i, j);
 212     }
 213     double v = sparseColumns[j].divide(i, divisor);
 214     return v;
 215   }
 216
 217   @Override
 218   protected double addValue(int i, int j, double addend)
 219   {
 220     double v = sparseColumns[j].add(i, addend);
 221     return v;
 222   }
 223
 224   /**
 225    * Returns the fraction of the whole matrix size that is actually modelled in
 226    * sparse arrays (normally, the non-zero values)
 227    *
 228    * @return
 229    */
 230   public float getFillRatio()
 231   {
 232     long count = 0L;
 233     for (SparseDoubleArray col : sparseColumns)
 234     {
 235       count += col.size();
 236     }
 237     return count / (float) (height() * width());
 238   }
 239 }