matrix.java

:<<数据挖掘--实用机器学习技术及java实现>>一书的配套源程序

  public final Matrix add(Matrix other) {

    int nr = m_Elements.length, nc = m_Elements[0].length;
    Matrix b;
    try {
      b = (Matrix)clone();
    } catch (CloneNotSupportedException ex) {
      b = new Matrix(nr, nc);
    }
    
    for(int i = 0;i < nc; i++) {
      for(int j = 0; j < nr; j++) {
        b.m_Elements[i][j] = m_Elements[j][i] + other.m_Elements[j][i];
      }
    }
    return b;
  }
  
  /**
   * Returns the transpose of a matrix.
   *
   * @return the transposition of this instance).
   */
  public final Matrix transpose() {

    int nr = m_Elements.length, nc = m_Elements[0].length;
    Matrix b = new Matrix(nc, nr);

    for(int i = 0;i < nc; i++) {
      for(int j = 0; j < nr; j++) {
	b.m_Elements[i][j] = m_Elements[j][i];
      }
    }

    return b;
  }
  
  /**
   * Reurns the multiplication of two matrices
   *
   * @param b the multiplication matrix 
   * @return the product matrix
   */
  public final Matrix multiply(Matrix b) {
   
    int nr = m_Elements.length, nc = m_Elements[0].length;
    int bnr = b.m_Elements.length, bnc = b.m_Elements[0].length;
    Matrix c = new Matrix(nr, bnc);

    for(int i = 0; i < nr; i++) {
      for(int j = 0; j< bnc; j++) {
	for(int k = 0; k < nc; k++) {
	  c.m_Elements[i][j] += m_Elements[i][k] * b.m_Elements[k][j];
	}
      }
    }

    return c;
  }

  /**
   * Performs a (ridged) linear regression.
   *
   * @param y the dependent variable vector
   * @return the coefficients 
   * @exception IllegalArgumentException if not successful
   */
  public final double[] regression(Matrix y) {

    if (y.numColumns() > 1) {
      throw new IllegalArgumentException("Only one dependent variable allowed");
    }
    int nc = m_Elements[0].length;
    double[] b = new double[nc];
    Matrix xt = this.transpose();

    boolean success = true;
    double ridge = 1e-8;

    do {
      Matrix ss = xt.multiply(this);
      
      // Set ridge regression adjustment
      for (int i = 0; i < nc; i++) {
	ss.setElement(i, i, ss.getElement(i, i) + ridge);
      }

      // Carry out the regression
      Matrix bb = xt.multiply(y);
      for(int i = 0; i < nc; i++) {
	b[i] = bb.m_Elements[i][0];
      }
      try {
	ss.lubksb(ss.ludcmp(), b);
	success = true;
      } catch (Exception ex) {
	ridge *= 10;
	success = false;
      }
    } while (!success);

    return b;
  }

  /**
   * Performs a weighted (ridged) linear regression. 
   *
   * @param y the dependent variable vector
   * @param w the array of data point weights
   * @return the coefficients 
   * @exception IllegalArgumentException if the wrong number of weights were
   * provided.
   */
  public final double[] regression(Matrix y, double [] w) {

    if (w.length != numRows()) {
      throw new IllegalArgumentException("Incorrect number of weights provided");
    }
    Matrix weightedThis = new Matrix(numRows(), numColumns());
    Matrix weightedDep = new Matrix(numRows(), 1);
    for (int i = 0; i < w.length; i++) {
      double weight = Math.sqrt(w[i]);
      for (int j = 0; j < numColumns(); j++) {
	weightedThis.setElement(i, j, getElement(i, j) * weight);
      }
      weightedDep.setElement(i, 0, y.getElement(i, 0) * weight);
    }
    return weightedThis.regression(weightedDep);
  }

  /**
   * Performs LU backward substitution.  Adapted from Numerical Recipes in C.
   *
   * @param indx the indices of the permutation
   * @param b the double vector, storing constant terms in the equation set; 
   * it later stores the computed coefficients' values
   */
  public final void lubksb(int[] indx, double b[]) {
    
    int nc = m_Elements[0].length;
    int ii = -1, ip;
    double sum;
    
    for (int i = 0; i < nc; i++) {
      ip = indx[i];
      sum = b[ip];
      b[ip] = b[i];
      if (ii != -1) {
	for (int j = ii; j < i; j++) {
	  sum -= m_Elements[i][j] * b[j];
	}
      } else if (sum != 0.0) {
	ii = i;
      }
      b[i] = sum;
    }
    for (int i = nc - 1; i >= 0; i--) {
      sum = b[i];
      for (int j = i + 1; j < nc; j++) {
	sum -= m_Elements[i][j] * b[j];
      }
      b[i] = sum / m_Elements[i][i];
    }
  }
  
  /**
   * Performs LU decomposition. Adapted from Numerical Recipes in C.
   *
   * @return the indices of the permutation
   * @exception Exception if the matrix is singular
   */
  public final int[] ludcmp() throws Exception {

    int nc = m_Elements[0].length;
    int[] indx = new int[m_Elements.length];
    int imax = -1;
    double big, dum, sum, temp;
    double vv[];
    
    vv = new double[nc];
    for (int i = 0; i < nc; i++) {
      big = 0.0;
      for (int j = 0; j < nc; j++) {
	if ((temp = Math.abs(m_Elements[i][j])) > big) {
	  big = temp;
	}
      }
      if (big < 0.000000001) {
	throw new Exception("Matrix is singular!");
      }
      vv[i] = 1.0 / big;
    }
    for (int j = 0; j < nc; j++) {
      for (int i = 0; i < j; i++) {
	sum = m_Elements[i][j];
	for (int k = 0; k < i; k++) {
	  sum -= m_Elements[i][k] * m_Elements[k][j];
	}
	m_Elements[i][j] = sum;
      }
      big = 0.0;
      for (int i = j; i < nc; i++) {
	sum = m_Elements[i][j];
	for (int k = 0; k < j; k++) {
	  sum -= m_Elements[i][k] * m_Elements[k][j];
	}
	m_Elements[i][j] = sum;
	if ((dum = vv[i] * Math.abs(sum)) >= big) {
	  big = dum;
	  imax = i;
	}
      } 
      if (j != imax) {
	for (int k = 0; k < nc; k++) {
	  dum = m_Elements[imax][k];
	  m_Elements[imax][k] = m_Elements[j][k];
	  m_Elements[j][k] = dum;
	}
	vv[imax] = vv[j];
      }
      indx[j] = imax;
      if (m_Elements[j][j] == 0.0) {
	throw new Exception("Matrix is singular");
      }
      if (j != nc - 1) {
	dum = 1.0 / (m_Elements[j][j]);
	for (int i = j + 1; i < nc; i++) {
	  m_Elements[i][j] *= dum;
	}
      }
    }
    return indx;
  }

  /**
   * Main method for testing this class.
   */
  public static void main(String[] ops) {

    double[] first = {2.3, 1.2, 5};
    double[] second = {5.2, 1.4, 9};
    double[] response = {4, 7, 8};
    double[] weights = {1, 2, 3};

    try {
      Matrix a = new Matrix(2, 3);
      Matrix b = new Matrix(3, 2);
      System.out.println("Number of columns for a: " + a.numColumns());
      System.out.println("Number of rows for a: " + a.numRows());
      a.setRow(0, first);
      a.setRow(1, second);
      b.setColumn(0, first);
      b.setColumn(1, second);
      System.out.println("a:\n " + a);
      System.out.println("b:\n " + b);
      System.out.println("a (0, 0): " + a.getElement(0, 0));
      System.out.println("a transposed:\n " + a.transpose());
      System.out.println("a * b:\n " + a.multiply(b));
      Matrix r = new Matrix(3, 1);
      r.setColumn(0, response);
      System.out.println("r:\n " + r);
      System.out.println("Coefficients of regression of b on r: ");
      double[] coefficients = b.regression(r);
      for (int i = 0; i < coefficients.length; i++) {
	System.out.print(coefficients[i] + " ");
      }
      System.out.println();
      System.out.println("Weights: ");
      for (int i = 0; i < weights.length; i++) {
	System.out.print(weights[i] + " ");
      }
      System.out.println();
      System.out.println("Coefficients of weighted regression of b on r: ");
      coefficients = b.regression(r, weights);
      for (int i = 0; i < coefficients.length; i++) {
	System.out.print(coefficients[i] + " ");
      }
      System.out.println();
      a.setElement(0, 0, 6);
      System.out.println("a with (0, 0) set to 6:\n " + a);
      a.write(new java.io.FileWriter("main.matrix"));
      System.out.println("wrote matrix to \"main.matrix\"\n" + a);
      a = new Matrix(new java.io.FileReader("main.matrix"));
      System.out.println("read matrix from \"main.matrix\"\n" + a);
    } catch (Exception e) {
      e.printStackTrace();
    }
  }
}