diff options
Diffstat (limited to 'src/de/lmu/ifi/dbs/elki/math/linearalgebra')
60 files changed, 2135 insertions, 783 deletions
diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/AffineTransformation.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/AffineTransformation.java index 6200164b..5d185ec5 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/AffineTransformation.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/AffineTransformation.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -150,7 +150,7 @@ public class AffineTransformation { * @param v translation vector */ public void addTranslation(Vector v) { - assert (v.getRowDimensionality() == dim); + assert (v.getDimensionality() == dim); // reset inverse transformation - needs recomputation. inv = null; @@ -295,7 +295,7 @@ public class AffineTransformation { * @return vector of dim+1, with new column having the value 1.0 */ public Vector homogeneVector(Vector v) { - assert (v.getRowDimensionality() == dim); + assert (v.getDimensionality() == dim); double[] dv = new double[dim + 1]; for(int i = 0; i < dim; i++) { dv[i] = v.get(i); @@ -311,7 +311,7 @@ public class AffineTransformation { * @return vector of dim+1, with new column having the value 0.0 */ public Vector homogeneRelativeVector(Vector v) { - assert (v.getRowDimensionality() == dim); + assert (v.getDimensionality() == dim); // TODO: this only works properly when trans[dim][dim] == 1.0, right? double[] dv = new double[dim + 1]; for(int i = 0; i < dim; i++) { @@ -328,7 +328,7 @@ public class AffineTransformation { * @return vector of dimension dim */ public Vector unhomogeneVector(Vector v) { - assert (v.getRowDimensionality() == dim + 1); + assert (v.getDimensionality() == dim + 1); // TODO: this only works properly when trans[dim][dim] == 1.0, right? double[] dv = new double[dim]; double scale = v.get(dim); @@ -346,7 +346,7 @@ public class AffineTransformation { * @return vector of dimension dim */ public Vector unhomogeneRelativeVector(Vector v) { - assert (v.getRowDimensionality() == dim + 1); + assert (v.getDimensionality() == dim + 1); double[] dv = new double[dim]; double scale = v.get(dim); assert (Math.abs(scale) == 0.0); diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Centroid.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Centroid.java index 5ba3d4a6..cf3bb25b 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Centroid.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Centroid.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -41,11 +41,6 @@ import de.lmu.ifi.dbs.elki.utilities.DatabaseUtil; */ public class Centroid extends Vector { /** - * Serial version - */ - private static final long serialVersionUID = 1L; - - /** * The current weight */ protected double wsum; @@ -147,7 +142,7 @@ public class Centroid extends Vector { * @return the data */ public <F extends NumberVector<? extends F, ?>> F toVector(Relation<? extends F> relation) { - return DatabaseUtil.assumeVectorField(relation).getFactory().newInstance(elements); + return DatabaseUtil.assumeVectorField(relation).getFactory().newNumberVector(elements); } /** @@ -160,7 +155,7 @@ public class Centroid extends Vector { int n = mat.getColumnDimensionality(); for(int i = 0; i < n; i++) { // TODO: avoid constructing the vector objects? - c.put(mat.getColumnVector(i)); + c.put(mat.getCol(i)); } return c; } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CholeskyDecomposition.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CholeskyDecomposition.java index ff3b8cfb..5209468f 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CholeskyDecomposition.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CholeskyDecomposition.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CovarianceMatrix.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CovarianceMatrix.java index 24ea740b..5d9d27a1 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CovarianceMatrix.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/CovarianceMatrix.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -245,7 +245,7 @@ public class CovarianceMatrix { * @return Mean vector */ public <F extends NumberVector<? extends F, ?>> F getMeanVector(Relation<? extends F> relation) { - return DatabaseUtil.assumeVectorField(relation).getFactory().newInstance(mean); + return DatabaseUtil.assumeVectorField(relation).getFactory().newNumberVector(mean); } /** @@ -332,7 +332,7 @@ public class CovarianceMatrix { int n = mat.getColumnDimensionality(); for(int i = 0; i < n; i++) { // TODO: avoid constructing the vector objects? - c.put(mat.getColumnVector(i)); + c.put(mat.getCol(i)); } return c; } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenPair.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenPair.java index 6b747dae..132393c6 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenPair.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenPair.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenvalueDecomposition.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenvalueDecomposition.java index 67d77089..48ce9c7b 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenvalueDecomposition.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/EigenvalueDecomposition.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -253,7 +253,7 @@ public class EigenvalueDecomposition implements java.io.Serializable { double g = d[l]; double p = (d[l + 1] - g) / (2.0 * e[l]); - double r = MathUtil.hypotenuse(p, 1.0); + double r = MathUtil.fastHypot(p, 1.0); if(p < 0) { r = -r; } @@ -281,7 +281,7 @@ public class EigenvalueDecomposition implements java.io.Serializable { s2 = s; g = c * e[i]; h = c * p; - r = MathUtil.hypotenuse(p, e[i]); + r = MathUtil.fastHypot(p, e[i]); e[i + 1] = s * r; s = e[i] / r; c = p / r; diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LUDecomposition.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LUDecomposition.java index 5adafe10..08634279 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LUDecomposition.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LUDecomposition.java @@ -1,10 +1,11 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; + /* This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -23,7 +24,6 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; along with this program. If not, see <http://www.gnu.org/licenses/>. */ - /** * LU Decomposition. * <P> @@ -78,11 +78,21 @@ public class LUDecomposition implements java.io.Serializable { * @param A Rectangular matrix */ public LUDecomposition(Matrix A) { - // Use a "left-looking", dot-product, Crout/Doolittle algorithm. + this(A.getArrayCopy(), A.getRowDimensionality(), A.getColumnDimensionality()); + } - LU = A.getArrayCopy(); - m = A.getRowDimensionality(); - n = A.getColumnDimensionality(); + /** + * LU Decomposition + * + * @param LU Rectangular matrix + * @param m row dimensionality + * @param n column dimensionality + */ + public LUDecomposition(double[][] LU, int m, int n) { + this.LU = LU; + this.m = m; + this.n = n; + // Use a "left-looking", dot-product, Crout/Doolittle algorithm. piv = new int[m]; for(int i = 0; i < m; i++) { piv[i] = i; @@ -273,25 +283,57 @@ public class LUDecomposition implements java.io.Serializable { Matrix Xmat = B.getMatrix(piv, 0, nx - 1); double[][] X = Xmat.getArrayRef(); + solveInplace(X, nx); + return Xmat; + } + + /** + * Solve A*X = B + * + * @param B A Matrix with as many rows as A and any number of columns. + * @return X so that L*U*X = B(piv,:) + * @exception IllegalArgumentException Matrix row dimensions must agree. + * @exception RuntimeException Matrix is singular. + */ + public double[][] solve(double[][] B) { + int mx = B.length; + int nx = B[0].length; + if(mx != m) { + throw new IllegalArgumentException("Matrix row dimensions must agree."); + } + if(!this.isNonsingular()) { + throw new RuntimeException("Matrix is singular."); + } + double[][] Xmat = new Matrix(B).getMatrix(piv, 0, nx - 1).getArrayRef(); + solveInplace(Xmat, nx); + return Xmat; + } + + /** + * Solve A*X = B + * + * @param B A Matrix with as many rows as A and any number of columns. + * @param nx Number of columns + */ + private void solveInplace(double[][] B, int nx) { // Solve L*Y = B(piv,:) for(int k = 0; k < n; k++) { for(int i = k + 1; i < n; i++) { for(int j = 0; j < nx; j++) { - X[i][j] -= X[k][j] * LU[i][k]; + B[i][j] -= B[k][j] * LU[i][k]; } } } // Solve U*X = Y; for(int k = n - 1; k >= 0; k--) { for(int j = 0; j < nx; j++) { - X[k][j] /= LU[k][k]; + B[k][j] /= LU[k][k]; } for(int i = 0; i < k; i++) { for(int j = 0; j < nx; j++) { - X[i][j] -= X[k][j] * LU[i][k]; + B[i][j] -= B[k][j] * LU[i][k]; } } } - return Xmat; } }
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LinearEquationSystem.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LinearEquationSystem.java index decec2b6..b6a53aa2 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LinearEquationSystem.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/LinearEquationSystem.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Matrix.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Matrix.java index fefe7e16..f64b1129 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Matrix.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Matrix.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -24,7 +24,6 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; */ import java.io.BufferedReader; -import java.io.Serializable; import java.io.StreamTokenizer; import java.util.Arrays; import java.util.logging.Logger; @@ -50,12 +49,7 @@ import de.lmu.ifi.dbs.elki.utilities.FormatUtil; * @apiviz.uses Vector * @apiviz.landmark */ -public class Matrix implements MatrixLike<Matrix>, Serializable { - /** - * Serial version - */ - private static final long serialVersionUID = 1L; - +public class Matrix { /** * A small number to handle numbers near 0 as 0. */ @@ -278,7 +272,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return a new matrix containing the same values as this matrix */ - @Override public final Matrix copy() { final Matrix X = new Matrix(elements.length, columndimension); for(int i = 0; i < elements.length; i++) { @@ -324,7 +317,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return m, the number of rows. */ - @Override public final int getRowDimensionality() { return elements.length; } @@ -334,7 +326,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return n, the number of columns. */ - @Override public final int getColumnDimensionality() { return columndimension; } @@ -347,7 +338,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @return A(i,j) * @throws ArrayIndexOutOfBoundsException on bounds error */ - @Override public final double get(final int i, final int j) { return elements[i][j]; } @@ -361,7 +351,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @return modified matrix * @throws ArrayIndexOutOfBoundsException on bounds error */ - @Override public final Matrix set(final int i, final int j, final double s) { elements[i][j] = s; return this; @@ -376,7 +365,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @return modified matrix * @throws ArrayIndexOutOfBoundsException on bounds error */ - @Override public final Matrix increment(final int i, final int j, final double s) { elements[i][j] += s; return this; @@ -597,22 +585,12 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { } /** - * Returns the <code>i</code>th row of this matrix. - * - * @param i the index of the row to be returned - * @return the <code>i</code>th row of this matrix - */ - public final Matrix getRow(final int i) { - return getMatrix(i, i, 0, columndimension - 1); - } - - /** * Returns the <code>i</code>th row of this matrix as vector. * * @param i the index of the row to be returned * @return the <code>i</code>th row of this matrix */ - public final Vector getRowVector(final int i) { + public final Vector getRow(final int i) { double[] row = elements[i].clone(); return new Vector(row); } @@ -620,26 +598,10 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { /** * Sets the <code>j</code>th row of this matrix to the specified vector. * - * @param j the index of the row to be set - * @param row the value of the row to be set - */ - public final void setRow(final int j, final Matrix row) { - if(row.columndimension != columndimension) { - throw new IllegalArgumentException("Matrix must consist of the same no of columns!"); - } - if(row.elements.length != 1) { - throw new IllegalArgumentException("Matrix must consist of one row!"); - } - setMatrix(elements.length - 1, 0, j, j, row); - } - - /** - * Sets the <code>j</code>th row of this matrix to the specified vector. - * * @param j the index of the column to be set * @param row the value of the column to be set */ - public final void setRowVector(final int j, final Vector row) { + public final void setRow(final int j, final Vector row) { if(row.elements.length != columndimension) { throw new IllegalArgumentException("Matrix must consist of the same no of columns!"); } @@ -649,23 +611,12 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { } /** - * Returns the <code>j</code>th column of this matrix. - * - * @param j the index of the column to be returned - * @return the <code>j</code>th column of this matrix - */ - public final Matrix getColumn(final int j) { - return getMatrix(0, elements.length - 1, j, j); - } - - /** * Returns the <code>j</code>th column of this matrix as vector. * * @param j the index of the column to be returned * @return the <code>j</code>th column of this matrix */ - @Override - public final Vector getColumnVector(final int j) { + public final Vector getCol(final int j) { final Vector v = new Vector(elements.length); for(int i = 0; i < elements.length; i++) { v.elements[i] = elements[i][j]; @@ -679,23 +630,7 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param j the index of the column to be set * @param column the value of the column to be set */ - public final void setColumn(final int j, final Matrix column) { - if(column.elements.length != elements.length) { - throw new IllegalArgumentException("Matrix must consist of the same no of rows!"); - } - if(column.columndimension != 1) { - throw new IllegalArgumentException("Matrix must consist of one column!"); - } - setMatrix(0, elements.length - 1, j, j, column); - } - - /** - * Sets the <code>j</code>th column of this matrix to the specified column. - * - * @param j the index of the column to be set - * @param column the value of the column to be set - */ - public final void setColumnVector(final int j, final Vector column) { + public final void setCol(final int j, final Vector column) { if(column.elements.length != elements.length) { throw new IllegalArgumentException("Matrix must consist of the same no of rows!"); } @@ -709,7 +644,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return A<sup>T</sup> */ - @Override public final Matrix transpose() { final Matrix X = new Matrix(columndimension, elements.length); for(int i = 0; i < elements.length; i++) { @@ -726,7 +660,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param B another matrix * @return A + B in a new Matrix */ - @Override public final Matrix plus(final Matrix B) { return copy().plusEquals(B); } @@ -738,7 +671,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param s scalar * @return A + s * B in a new Matrix */ - @Override public final Matrix plusTimes(final Matrix B, final double s) { return copy().plusTimesEquals(B, s); } @@ -749,7 +681,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param B another matrix * @return A + B in this Matrix */ - @Override public final Matrix plusEquals(final Matrix B) { checkMatrixDimensions(B); for(int i = 0; i < elements.length; i++) { @@ -767,7 +698,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param s Scalar * @return A + s * B in this Matrix */ - @Override public final Matrix plusTimesEquals(final Matrix B, final double s) { checkMatrixDimensions(B); for(int i = 0; i < elements.length; i++) { @@ -784,7 +714,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param B another matrix * @return A - B in a new Matrix */ - @Override public final Matrix minus(final Matrix B) { return copy().minusEquals(B); } @@ -796,7 +725,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param s Scalar * @return A - s * B in a new Matrix */ - @Override public final Matrix minusTimes(final Matrix B, final double s) { return copy().minusTimesEquals(B, s); } @@ -807,7 +735,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param B another matrix * @return A - B in this Matrix */ - @Override public final Matrix minusEquals(final Matrix B) { checkMatrixDimensions(B); for(int i = 0; i < elements.length; i++) { @@ -825,7 +752,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param s Scalar * @return A - s * B in this Matrix */ - @Override public final Matrix minusTimesEquals(final Matrix B, final double s) { checkMatrixDimensions(B); for(int i = 0; i < elements.length; i++) { @@ -842,7 +768,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param s scalar * @return s*A */ - @Override public final Matrix times(final double s) { return copy().timesEquals(s); } @@ -853,7 +778,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @param s scalar * @return replace A by s*A */ - @Override public final Matrix timesEquals(final double s) { for(int i = 0; i < elements.length; i++) { for(int j = 0; j < columndimension; j++) { @@ -1036,91 +960,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { } /** - * Returns the scalar product of the colA column of this and the colB column - * of B. - * - * @param colA The column of A to compute scalar product for - * @param B second Matrix - * @param colB The column of B to compute scalar product for - * @return double The scalar product of the first column of this and B - */ - public double scalarProduct(int colA, Matrix B, int colB) { - double scalarProduct = 0.0; - for(int row = 0; row < getRowDimensionality(); row++) { - double prod = elements[row][colA] * B.elements[row][colB]; - scalarProduct += prod; - } - return scalarProduct; - } - - /** - * Returns the scalar product of the colA column of this and the colB column - * of B. - * - * @param colA The column of A to compute scalar product for - * @param B Vector - * @return double The scalar product of the first column of this and B - */ - public double scalarProduct(int colA, Vector B) { - double scalarProduct = 0.0; - for(int row = 0; row < getRowDimensionality(); row++) { - double prod = elements[row][colA] * B.elements[row]; - scalarProduct += prod; - } - return scalarProduct; - } - - /** - * LU Decomposition - * - * @return LUDecomposition - * @see LUDecomposition - */ - public final LUDecomposition lu() { - return new LUDecomposition(this); - } - - /** - * QR Decomposition - * - * @return QRDecomposition - * @see QRDecomposition - */ - public final QRDecomposition qr() { - return new QRDecomposition(this); - } - - /** - * Cholesky Decomposition - * - * @return CholeskyDecomposition - * @see CholeskyDecomposition - */ - public final CholeskyDecomposition chol() { - return new CholeskyDecomposition(this); - } - - /** - * Singular Value Decomposition - * - * @return SingularValueDecomposition - * @see SingularValueDecomposition - */ - public final SingularValueDecomposition svd() { - return new SingularValueDecomposition(this); - } - - /** - * Eigenvalue Decomposition - * - * @return EigenvalueDecomposition - * @see EigenvalueDecomposition - */ - public final EigenvalueDecomposition eig() { - return new EigenvalueDecomposition(this); - } - - /** * Solve A*X = B * * @param B right hand side @@ -1131,16 +970,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { } /** - * Solve X*A = B, which is also A'*X' = B' - * - * @param B right hand side - * @return solution if A is square, least squares solution otherwise. - */ - public final Matrix solveTranspose(final Matrix B) { - return transpose().solve(B.transpose()); - } - - /** * Matrix inverse or pseudoinverse * * @return inverse(A) if A is square, pseudoinverse otherwise. @@ -1194,7 +1023,7 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return maximum column sum. */ - public double norm1() { + public final double norm1() { double f = 0; for(int j = 0; j < columndimension; j++) { double s = 0; @@ -1220,7 +1049,7 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return maximum row sum. */ - public double normInf() { + public final double normInf() { double f = 0; for(int i = 0; i < elements.length; i++) { double s = 0; @@ -1237,56 +1066,20 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * * @return sqrt of sum of squares of all elements. */ - public double normF() { + public final double normF() { double f = 0; for(int i = 0; i < elements.length; i++) { for(int j = 0; j < columndimension; j++) { - f = MathUtil.hypotenuse(f, elements[i][j]); + f = MathUtil.fastHypot(f, elements[i][j]); } } return f; } /** - * distanceCov returns distance of two Matrices A and B, i.e. the root of the - * sum of the squared distances A<sub>ij</sub>-B<sub>ij</sub>. - * - * @param B Matrix to compute distance from this (A) - * @return distance of Matrices - */ - // TODO: unused - remove / move into a MatrixDistance helper? - public final double distanceCov(final Matrix B) { - double distance = 0.0; - double distIJ; - int row; - for(int col = 0; col < columndimension; col++) { - for(row = 0; row < elements.length; row++) { - distIJ = elements[row][col] - B.elements[row][col]; - distance += (distIJ * distIJ); - } - } - distance = Math.sqrt(distance); - return distance; - } - - /** - * getDiagonal returns array of diagonal-elements. - * - * @return double[] the values on the diagonal of the Matrix - */ - public final double[] getDiagonal() { - int n = Math.min(columndimension, elements.length); - final double[] diagonal = new double[n]; - for(int i = 0; i < n; i++) { - diagonal[i] = elements[i][i]; - } - return diagonal; - } - - /** * Normalizes the columns of this matrix to length of 1.0. */ - public void normalizeColumns() { + public final void normalizeColumns() { for(int col = 0; col < columndimension; col++) { double norm = 0.0; for(int row = 0; row < elements.length; row++) { @@ -1592,10 +1385,10 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { for(int i = 0; i < result.columndimension; i++) { // FIXME: optimize - excess copying! if(i < columndimension) { - result.setColumn(i, getColumn(i)); + result.setCol(i, getCol(i)); } else { - result.setColumn(i, columns.getColumn(i - columndimension)); + result.setCol(i, columns.getCol(i - columndimension)); } } return result; @@ -1607,19 +1400,19 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { * @return the orthonormalized matrix */ public final Matrix orthonormalize() { - Matrix v = getColumn(0); + Matrix v = copy(); // FIXME: optimize - excess copying! for(int i = 1; i < columndimension; i++) { - final Matrix u_i = getColumn(i); - final Matrix sum = new Matrix(elements.length, 1); + final Vector u_i = getCol(i); + final Vector sum = new Vector(elements.length); for(int j = 0; j < i; j++) { - final Matrix v_j = v.getColumn(j); - double scalar = u_i.scalarProduct(0, v_j, 0) / v_j.scalarProduct(0, v_j, 0); - sum.plusEquals(v_j.times(scalar)); + final Vector v_j = v.getCol(j); + double scalar = u_i.transposeTimes(v_j) / v_j.transposeTimes(v_j); + sum.plusTimesEquals(v_j, scalar); } - final Matrix v_i = u_i.minus(sum); - v = v.appendColumns(v_i); + final Vector v_i = u_i.minus(sum); + v.setCol(i, v_i); } v.normalizeColumns(); @@ -1716,7 +1509,7 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { /** * Check if size(A) == size(B) */ - protected void checkMatrixDimensions(MatrixLike<?> B) { + protected void checkMatrixDimensions(Matrix B) { if(B.getRowDimensionality() != getRowDimensionality() || B.getColumnDimensionality() != getColumnDimensionality()) { throw new IllegalArgumentException("Matrix dimensions must agree."); } @@ -1807,15 +1600,6 @@ public class Matrix implements MatrixLike<Matrix>, Serializable { } /** - * Returns the dimensionality of this matrix as a string. - * - * @return the dimensionality of this matrix as a string - */ - public String dimensionInfo() { - return getRowDimensionality() + " x " + getColumnDimensionality(); - } - - /** * toString returns String-representation of Matrix. */ @Override diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/MatrixLike.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/MatrixLike.java deleted file mode 100644 index ff1ec5ba..00000000 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/MatrixLike.java +++ /dev/null @@ -1,194 +0,0 @@ -package de.lmu.ifi.dbs.elki.math.linearalgebra; - -/* - This file is part of ELKI: - Environment for Developing KDD-Applications Supported by Index-Structures - - Copyright (C) 2011 - Ludwig-Maximilians-Universität München - Lehr- und Forschungseinheit für Datenbanksysteme - ELKI Development Team - - This program is free software: you can redistribute it and/or modify - it under the terms of the GNU Affero General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - This program is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU Affero General Public License for more details. - - You should have received a copy of the GNU Affero General Public License - along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - - -/** - * Common Interface for Matrix and Vector objects, where M is the actual type. - * - * The type M guarantees type safety for many operations. - * - * @param M the actual type - * - * @apiviz.landmark - * - * @author Elke Achtert - * @author Erich Schubert - */ -public interface MatrixLike<M extends MatrixLike<M>> extends Cloneable { - /** - * Make a deep copy of a matrix. - * - * @return a new matrix containing the same values as this matrix - */ - public M copy(); - - /** - * Clone the Matrix object. - */ - public Object clone(); - - /** - * Returns the dimensionality of the rows of this matrix. - * - * @return m, the number of rows. - */ - public int getRowDimensionality(); - - /** - * Returns the dimensionality of the columns of this matrix. - * - * @return n, the number of columns. - */ - public int getColumnDimensionality(); - - /** - * Get a single element. - * - * @param i Row index. - * @param j Column index. - * @return A(i,j) - * @throws ArrayIndexOutOfBoundsException on bounds error - */ - public double get(int i, int j); - - /** - * Set a single element. - * - * @param i Row index. - * @param j Column index. - * @param s A(i,j). - * @throws ArrayIndexOutOfBoundsException on bounds error - */ - public M set(int i, int j, double s); - - /** - * Increments a single element. - * - * @param i the row index - * @param j the column index - * @param s the increment value: A(i,j) = A(i.j) + s. - * @throws ArrayIndexOutOfBoundsException on bounds error - */ - public M increment(int i, int j, double s); - - /** - * Returns the <code>i</code>th column of this matrix as vector. - * - * @param i the index of the column to be returned - * @return the <code>i</code>th column of this matrix - */ - public Vector getColumnVector(int i); - - /** - * Matrix transpose. - * - * @return A<sup>T</sup> - */ - public Matrix transpose(); - - /** - * C = A + B - * - * @param B another matrix - * @return A + B in a new Matrix - */ - public M plus(M B); - - /** - * C = A + s*B - * - * @param B another matrix - * @param s scalar - * @return A + s*B in a new Matrix - */ - public M plusTimes(M B, double s); - - /** - * A = A + B - * - * @param B another matrix - * @return A + B in this Matrix - */ - public M plusEquals(M B); - - /** - * C = A + s*B - * - * @param B another matrix - * @param s scalar - * @return A + s*B in this Matrix - */ - public M plusTimesEquals(M B, double s); - - /** - * C = A - B - * - * @param B another matrix - * @return A - B in a new Matrix - */ - public M minus(M B); - - /** - * C = A - s*B - * - * @param B another matrix - * @param s Scalar - * @return A - s*B in a new Matrix - */ - public M minusTimes(M B, double s); - - /** - * A = A - B - * - * @param B another matrix - * @return A - B in this Matrix - */ - public M minusEquals(M B); - - /** - * C = A - s*B - * - * @param B another matrix - * @param s Scalar - * @return A - s*B in a new Matrix - */ - public M minusTimesEquals(M B, double s); - - /** - * Multiply a matrix by a scalar, C = s*A - * - * @param s scalar - * @return s*A - */ - public M times(double s); - - /** - * Multiply a matrix by a scalar in place, A = s*A - * - * @param s scalar - * @return replace A by s*A - */ - public M timesEquals(double s); -}
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectedCentroid.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectedCentroid.java index e54f75fe..02b5b424 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectedCentroid.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectedCentroid.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -43,11 +43,6 @@ import de.lmu.ifi.dbs.elki.utilities.DatabaseUtil; */ public class ProjectedCentroid extends Centroid { /** - * Serial version - */ - private static final long serialVersionUID = 1L; - - /** * The selected dimensions. */ private BitSet dims; diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectionResult.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectionResult.java index 040bd4b1..d8858657 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectionResult.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/ProjectionResult.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/QRDecomposition.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/QRDecomposition.java index c053de0d..5b52d837 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/QRDecomposition.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/QRDecomposition.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -74,13 +74,22 @@ public class QRDecomposition implements java.io.Serializable { * QR Decomposition, computed by Householder reflections. * * @param A Rectangular matrix - * */ public QRDecomposition(Matrix A) { - // Initialize. - QR = A.getArrayCopy(); - m = A.getRowDimensionality(); - n = A.getColumnDimensionality(); + this(A.getArrayRef(), A.getRowDimensionality(), A.getColumnDimensionality()); + } + + /** + * QR Decomposition, computed by Householder reflections. + * + * @param A Rectangular matrix + * @param m row dimensionality + * @param n column dimensionality + */ + public QRDecomposition(double[][] A, int m, int n) { + this.QR = new Matrix(A).getArrayCopy(); + this.m = QR.length; + this.n = QR[0].length; Rdiag = new double[n]; // Main loop. @@ -88,7 +97,7 @@ public class QRDecomposition implements java.io.Serializable { // Compute 2-norm of k-th column without under/overflow. double nrm = 0; for(int i = k; i < m; i++) { - nrm = MathUtil.hypotenuse(nrm, QR[i][k]); + nrm = MathUtil.fastHypot(nrm, QR[i][k]); } if(nrm != 0.0) { @@ -227,8 +236,38 @@ public class QRDecomposition implements java.io.Serializable { // Copy right hand side int nx = B.getColumnDimensionality(); - double[][] X = B.getArrayCopy(); + Matrix X = B.copy(); + + solveInplace(X.getArrayRef(), nx); + return X.getMatrix(0, n - 1, 0, nx - 1); + } + + /** + * Least squares solution of A*X = B + * + * @param B A Matrix with as many rows as A and any number of columns. + * @return X that minimizes the two norm of Q*R*X-B. + * @exception IllegalArgumentException Matrix row dimensions must agree. + * @exception RuntimeException Matrix is rank deficient. + */ + public double[][] solve(double[][] B) { + int rows = B.length; + int cols = B[0].length; + if(rows != m) { + throw new IllegalArgumentException("Matrix row dimensions must agree."); + } + if(!this.isFullRank()) { + throw new RuntimeException("Matrix is rank deficient."); + } + + // Copy right hand side + Matrix X = new Matrix(B).copy(); + + solveInplace(X.getArrayRef(), cols); + return X.getMatrix(0, n - 1, 0, cols - 1).getArrayRef(); + } + private void solveInplace(double[][] X, int nx) { // Compute Y = transpose(Q)*B for(int k = 0; k < n; k++) { for(int j = 0; j < nx; j++) { @@ -253,6 +292,5 @@ public class QRDecomposition implements java.io.Serializable { } } } - return (new Matrix(X).getMatrix(0, n - 1, 0, nx - 1)); } }
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SingularValueDecomposition.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SingularValueDecomposition.java index 38d70443..183a8034 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SingularValueDecomposition.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SingularValueDecomposition.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -41,12 +41,7 @@ import de.lmu.ifi.dbs.elki.math.MathUtil; * * @apiviz.uses Matrix - - transforms */ -@SuppressWarnings("serial") -public class SingularValueDecomposition implements java.io.Serializable { - /* - * ------------------------ Class variables ------------------------ - */ - +public class SingularValueDecomposition { /** * Arrays for internal storage of U and V. * @@ -70,23 +65,26 @@ public class SingularValueDecomposition implements java.io.Serializable { */ private int m, n; - /* - * ------------------------ Constructor ------------------------ - */ - /** * Construct the singular value decomposition * * @param Arg Rectangular matrix */ - public SingularValueDecomposition(Matrix Arg) { + this(Arg.getArrayRef()); + } + /** + * Constructor. + * + * @param Arg Rectangular input matrix + */ + public SingularValueDecomposition(double[][] Arg) { + double[][] A = new Matrix(Arg).getArrayCopy(); + this.m = A.length; + this.n = A[0].length; // Derived from LINPACK code. // Initialize. - double[][] A = Arg.getArrayCopy(); - m = Arg.getRowDimensionality(); - n = Arg.getColumnDimensionality(); int nu = Math.min(m, n); s = new double[Math.min(m + 1, n)]; U = new double[m][nu]; @@ -103,13 +101,12 @@ public class SingularValueDecomposition implements java.io.Serializable { int nrt = Math.max(0, Math.min(n - 2, m)); for(int k = 0; k < Math.max(nct, nrt); k++) { if(k < nct) { - // Compute the transformation for the k-th column and // place the k-th diagonal in s[k]. // Compute 2-norm of k-th column without under/overflow. s[k] = 0; for(int i = k; i < m; i++) { - s[k] = MathUtil.hypotenuse(s[k], A[i][k]); + s[k] = MathUtil.fastHypot(s[k], A[i][k]); } if(s[k] != 0.0) { if(A[k][k] < 0.0) { @@ -124,9 +121,7 @@ public class SingularValueDecomposition implements java.io.Serializable { } for(int j = k + 1; j < n; j++) { if((k < nct) & (s[k] != 0.0)) { - // Apply the transformation. - double t = 0; for(int i = k; i < m; i++) { t += A[i][k] * A[i][j]; @@ -143,7 +138,6 @@ public class SingularValueDecomposition implements java.io.Serializable { e[j] = A[k][j]; } if(wantu & (k < nct)) { - // Place the transformation in U for subsequent back // multiplication. @@ -152,13 +146,12 @@ public class SingularValueDecomposition implements java.io.Serializable { } } if(k < nrt) { - // Compute the k-th row transformation and place the // k-th super-diagonal in e[k]. // Compute 2-norm without under/overflow. e[k] = 0; for(int i = k + 1; i < n; i++) { - e[k] = MathUtil.hypotenuse(e[k], e[i]); + e[k] = MathUtil.fastHypot(e[k], e[i]); } if(e[k] != 0.0) { if(e[k + 1] < 0.0) { @@ -171,9 +164,7 @@ public class SingularValueDecomposition implements java.io.Serializable { } e[k] = -e[k]; if((k + 1 < m) & (e[k] != 0.0)) { - // Apply the transformation. - for(int i = k + 1; i < m; i++) { work[i] = 0.0; } @@ -190,10 +181,8 @@ public class SingularValueDecomposition implements java.io.Serializable { } } if(wantv) { - // Place the transformation in V for subsequent // back multiplication. - for(int i = k + 1; i < n; i++) { V[i][k] = e[i]; } @@ -254,7 +243,6 @@ public class SingularValueDecomposition implements java.io.Serializable { } // If required, generate V. - if(wantv) { for(int k = n - 1; k >= 0; k--) { if((k < nrt) & (e[k] != 0.0)) { @@ -336,14 +324,13 @@ public class SingularValueDecomposition implements java.io.Serializable { // Perform the task indicated by kase. switch(kase){ - + // Deflate negligible s(p). - case 1: { double f = e[p - 2]; e[p - 2] = 0.0; for(int j = p - 2; j >= k; j--) { - double t = MathUtil.hypotenuse(s[j], f); + double t = MathUtil.fastHypot(s[j], f); double cs = s[j] / t; double sn = f / t; s[j] = t; @@ -363,12 +350,11 @@ public class SingularValueDecomposition implements java.io.Serializable { break; // Split at negligible s(k). - case 2: { double f = e[k - 1]; e[k - 1] = 0.0; for(int j = k; j < p; j++) { - double t = MathUtil.hypotenuse(s[j], f); + double t = MathUtil.fastHypot(s[j], f); double cs = s[j] / t; double sn = f / t; s[j] = t; @@ -386,11 +372,8 @@ public class SingularValueDecomposition implements java.io.Serializable { break; // Perform one qr step. - case 3: { - // Calculate the shift. - double scale = Math.max(Math.max(Math.max(Math.max(Math.abs(s[p - 1]), Math.abs(s[p - 2])), Math.abs(e[p - 2])), Math.abs(s[k])), Math.abs(e[k])); double sp = s[p - 1] / scale; double spm1 = s[p - 2] / scale; @@ -411,9 +394,8 @@ public class SingularValueDecomposition implements java.io.Serializable { double g = sk * ek; // Chase zeros. - for(int j = k; j < p - 1; j++) { - double t = MathUtil.hypotenuse(f, g); + double t = MathUtil.fastHypot(f, g); double cs = f / t; double sn = g / t; if(j != k) { @@ -430,7 +412,7 @@ public class SingularValueDecomposition implements java.io.Serializable { V[i][j] = t; } } - t = MathUtil.hypotenuse(f, g); + t = MathUtil.fastHypot(f, g); cs = f / t; sn = g / t; s[j] = t; @@ -452,11 +434,8 @@ public class SingularValueDecomposition implements java.io.Serializable { break; // Convergence. - case 4: { - // Make the singular values positive. - if(s[k] <= 0.0) { s[k] = (s[k] < 0.0 ? -s[k] : 0.0); if(wantv) { @@ -467,7 +446,6 @@ public class SingularValueDecomposition implements java.io.Serializable { } // Order the singular values. - while(k < pp) { if(s[k] >= s[k + 1]) { break; @@ -499,16 +477,11 @@ public class SingularValueDecomposition implements java.io.Serializable { } } - /* - * ------------------------ Public Methods ------------------------ - */ - /** * Return the left singular vectors * * @return U */ - public Matrix getU() { return new Matrix(U); } @@ -518,7 +491,6 @@ public class SingularValueDecomposition implements java.io.Serializable { * * @return V */ - public Matrix getV() { return new Matrix(V); } @@ -528,7 +500,6 @@ public class SingularValueDecomposition implements java.io.Serializable { * * @return diagonal of S. */ - public double[] getSingularValues() { return s; } @@ -538,7 +509,6 @@ public class SingularValueDecomposition implements java.io.Serializable { * * @return S */ - public Matrix getS() { Matrix X = new Matrix(n, n); double[][] S = X.getArrayRef(); @@ -556,7 +526,6 @@ public class SingularValueDecomposition implements java.io.Serializable { * * @return max(S) */ - public double norm2() { return s[0]; } @@ -566,7 +535,6 @@ public class SingularValueDecomposition implements java.io.Serializable { * * @return max(S)/min(S) */ - public double cond() { return s[0] / s[Math.min(m, n) - 1]; } @@ -576,7 +544,6 @@ public class SingularValueDecomposition implements java.io.Serializable { * * @return Number of nonnegligible singular values. */ - public int rank() { double eps = Math.pow(2.0, -52.0); double tol = Math.max(m, n) * s[0] * eps; @@ -588,4 +555,4 @@ public class SingularValueDecomposition implements java.io.Serializable { } return r; } -} +}
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SortedEigenPairs.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SortedEigenPairs.java index b0264c35..5fa023ca 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SortedEigenPairs.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SortedEigenPairs.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -65,7 +65,7 @@ public class SortedEigenPairs { this.eigenPairs = new EigenPair[eigenvalues.length]; for(int i = 0; i < eigenvalues.length; i++) { double e = java.lang.Math.abs(eigenvalues[i]); - Vector v = eigenvectors.getColumnVector(i); + Vector v = eigenvectors.getCol(i); eigenPairs[i] = new EigenPair(v, e); } @@ -124,7 +124,7 @@ public class SortedEigenPairs { Matrix eigenVectors = new Matrix(eigenPairs.length, eigenPairs.length); for(int i = 0; i < eigenPairs.length; i++) { EigenPair eigenPair = eigenPairs[i]; - eigenVectors.setColumnVector(i, eigenPair.getEigenvector()); + eigenVectors.setCol(i, eigenPair.getEigenvector()); } return eigenVectors; } @@ -139,7 +139,7 @@ public class SortedEigenPairs { Matrix eigenVectors = new Matrix(eigenPairs.length, n); for(int i = 0; i < n; i++) { EigenPair eigenPair = eigenPairs[i]; - eigenVectors.setColumnVector(i, eigenPair.getEigenvector()); + eigenVectors.setCol(i, eigenPair.getEigenvector()); } return eigenVectors; } @@ -154,7 +154,7 @@ public class SortedEigenPairs { Matrix eigenVectors = new Matrix(eigenPairs.length, n); for(int i = 0; i < n; i++) { EigenPair eigenPair = eigenPairs[eigenPairs.length - 1 - i]; - eigenVectors.setColumnVector(i, eigenPair.getEigenvector()); + eigenVectors.setCol(i, eigenPair.getEigenvector()); } return eigenVectors; } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SubspaceProjectionResult.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SubspaceProjectionResult.java index 3e8815e3..3cab3b51 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SubspaceProjectionResult.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/SubspaceProjectionResult.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/VMath.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/VMath.java new file mode 100644 index 00000000..97466b20 --- /dev/null +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/VMath.java @@ -0,0 +1,1512 @@ +package de.lmu.ifi.dbs.elki.math.linearalgebra; + +/* + This file is part of ELKI: + Environment for Developing KDD-Applications Supported by Index-Structures + + Copyright (C) 2012 + Ludwig-Maximilians-Universität München + Lehr- und Forschungseinheit für Datenbanksysteme + ELKI Development Team + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU Affero General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU Affero General Public License for more details. + + You should have received a copy of the GNU Affero General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. + */ + +import java.util.Arrays; + +/** + * Class providing basic vector mathematics, for low-level vectors stored as + * {@code double[]}. While this is less nice syntactically, it reduces memory + * usage and VM overhead. + * + * @author Erich Schubert + * + * @apiviz.landmark + */ +public final class VMath { + /** + * A small number to handle numbers near 0 as 0. + */ + public static final double DELTA = 1E-5; + + /** + * Fake constructor. Static class. + */ + private VMath() { + // Cannot be instantiated + } + + /** + * Returns a randomly created vector of length 1.0 + * + * @param dimensionality dimensionality + * @return Random vector of length 1.0 + */ + public final static double[] randomNormalizedVector(final int dimensionality) { + final double[] v = new double[dimensionality]; + for(int i = 0; i < dimensionality; i++) { + v[i] = Math.random(); + } + double norm = euclideanLength(v); + if(norm != 0) { + for(int row = 0; row < v.length; row++) { + v[row] /= norm; + } + return v; + } + else { + return randomNormalizedVector(dimensionality); + } + } + + /** + * Returns the ith unit vector of the specified dimensionality. + * + * @param dimensionality the dimensionality of the vector + * @param i the index + * @return the ith unit vector of the specified dimensionality + */ + public final static double[] unitVector(final int dimensionality, final int i) { + final double[] v = new double[dimensionality]; + v[i] = 1; + return v; + } + + /** + * Returns a copy of this vector. + * + * @param v original vector + * @return a copy of this vector + */ + public final static double[] copy(final double[] v) { + return Arrays.copyOf(v, v.length); + } + + /** + * Transpose vector to a matrix. + * + * @param v Vector + * @return Matrix + */ + public final static double[][] transpose(final double[] v) { + double[][] re = new double[v.length][1]; + for(int i = 0; i < v.length; i++) { + re[i][0] = v[i]; + } + return re; + } + + /** + * Computes v1 + v2 for vectors. + * + * @param v1 first vector + * @param v2 second vector + * @return the sum v1 + v2 + */ + public final static double[] plus(final double[] v1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + final double[] result = new double[v1.length]; + for(int i = 0; i < result.length; i++) { + result[i] = v1[i] + v2[i]; + } + return result; + } + + /** + * Computes v1 + v2 * s2 + * + * @param v1 first vector + * @param v2 second vector + * @param s2 the scalar + * @return the result of v1 + v2 * s2 + */ + public final static double[] plusTimes(final double[] v1, final double[] v2, final double s2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + final double[] result = new double[v1.length]; + for(int i = 0; i < result.length; i++) { + result[i] = v1[i] + v2[i] * s2; + } + return result; + } + + /** + * Computes v1 * s1 + v2 + * + * @param v1 first vector + * @param s1 the scalar for v1 + * @param v2 second vector + * @return the result of v1 * s1 + v2 + */ + public final static double[] timesPlus(final double[] v1, final double s1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + final double[] result = new double[v1.length]; + for(int i = 0; i < result.length; i++) { + result[i] = v1[i] * s1 + v2[i]; + } + return result; + } + + /** + * Computes v1 * s1 + v2 * s2 + * + * @param v1 first vector + * @param s1 the scalar for v1 + * @param v2 second vector + * @param s2 the scalar for v2 + * @return the result of v1 * s1 + v2 * s2 + */ + public final static double[] timesPlusTimes(final double[] v1, final double s1, final double[] v2, final double s2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + final double[] result = new double[v1.length]; + for(int i = 0; i < result.length; i++) { + result[i] = v1[i] * s1 + v2[i] * s2; + } + return result; + } + + /** + * Computes v1 = v1 + v2, overwriting v1 + * + * @param v1 first vector (overwritten) + * @param v2 second vector + * @return v1 = v1 + v2 + */ + public final static double[] plusEquals(final double[] v1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] += v2[i]; + } + return v1; + } + + /** + * Computes v1 = v1 + v2 * s2, overwriting v1 + * + * @param v1 first vector + * @param v2 another vector + * @param s2 scalar vor v2 + * @return v1 = v1 + v2 * s2 + */ + public final static double[] plusTimesEquals(final double[] v1, final double[] v2, final double s2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] += s2 * v2[i]; + } + return v1; + } + + /** + * Computes v1 = v1 * s1 + v2, overwriting v1 + * + * @param v1 first vector + * @param s1 scalar for v1 + * @param v2 another vector + * @return v1 = v1 * s1 + v2 + */ + public final static double[] timesPlusEquals(final double[] v1, final double s1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] = v1[i] * s1 + v2[i]; + } + return v1; + } + + /** + * Computes v1 = v1 * s1 + v2 * s2, overwriting v1 + * + * @param v1 first vector + * @param s1 scalar for v1 + * @param v2 another vector + * @param s2 scalar for v2 + * @return v1 = v1 * s1 + v2 * s2 + */ + public final static double[] timesPlusTimesEquals(final double[] v1, final double s1, final double[] v2, final double s2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] = v1[i] * s1 + v2[i] * s2; + } + return v1; + } + + /** + * Computes v1 + d + * + * @param v1 vector to add to + * @param d value to add + * @return v1 + d + */ + public final static double[] plus(final double[] v1, final double d) { + final double[] result = new double[v1.length]; + for(int i = 0; i < result.length; i++) { + result[i] = v1[i] + d; + } + return result; + } + + /** + * Computes v1 = v1 + d, overwriting v1 + * + * @param v1 vector to add to + * @param d value to add + * @return Modified vector + */ + public final static double[] plusEquals(final double[] v1, final double d) { + for(int i = 0; i < v1.length; i++) { + v1[i] += d; + } + return v1; + } + + /** + * Computes v1 - v2 + * + * @param v1 first vector + * @param v2 the vector to be subtracted from this vector + * @return v1 - v2 + */ + public final static double[] minus(final double[] v1, final double[] v2) { + final double[] sub = new double[v1.length]; + for(int i = 0; i < v1.length; i++) { + sub[i] = v1[i] - v2[i]; + } + return sub; + } + + /** + * Computes v1 - v2 * s2 + * + * @param v1 first vector + * @param v2 the vector to be subtracted from this vector + * @param s2 the scaling factor for v2 + * @return v1 - v2 * s2 + */ + public final static double[] minusTimes(final double[] v1, final double[] v2, final double s2) { + final double[] sub = new double[v1.length]; + for(int i = 0; i < v1.length; i++) { + sub[i] = v1[i] - v2[i] * s2; + } + return sub; + } + + /** + * Computes v1 * s1 - v2 + * + * @param v1 first vector + * @param s1 the scaling factor for v1 + * @param v2 the vector to be subtracted from this vector + * @return v1 * s1 - v2 + */ + public final static double[] timesMinus(final double[] v1, final double s1, final double[] v2) { + final double[] sub = new double[v1.length]; + for(int i = 0; i < v1.length; i++) { + sub[i] = v1[i] * s1 - v2[i]; + } + return sub; + } + + /** + * Computes v1 * s1 - v2 * s2 + * + * @param v1 first vector + * @param s1 the scaling factor for v1 + * @param v2 the vector to be subtracted from this vector + * @param s2 the scaling factor for v2 + * @return v1 * s1 - v2 * s2 + */ + public final static double[] timesMinusTimes(final double[] v1, final double s1, final double[] v2, final double s2) { + final double[] sub = new double[v1.length]; + for(int i = 0; i < v1.length; i++) { + sub[i] = v1[i] * s1 - v2[i] * s2; + } + return sub; + } + + /** + * Computes v1 = v1 - v2, overwriting v1 + * + * @param v1 vector + * @param v2 another vector + * @return v1 = v1 - v2 + */ + public final static double[] minusEquals(final double[] v1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] -= v2[i]; + } + return v1; + } + + /** + * Computes v1 = v1 - v2 * s2, overwriting v1 + * + * @param v1 vector + * @param v2 another vector + * @param s2 scalar for v2 + * @return v1 = v1 - v2 * s2 + */ + public final static double[] minusTimesEquals(final double[] v1, final double[] v2, final double s2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] -= v2[i] * s2; + } + return v1; + } + + /** + * Computes v1 = v1 * s1 - v2, overwriting v1 + * + * @param v1 vector + * @param s1 scalar for v1 + * @param v2 another vector + * @return v1 = v1 * s1 - v2 + */ + public final static double[] timesMinusEquals(final double[] v1, final double s1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] = v1[i] * s1 - v2[i]; + } + return v1; + } + + /** + * Computes v1 = v1 * s1 - v2 * s2, overwriting v1 + * + * @param v1 vector + * @param s1 scalar for v1 + * @param v2 another vector + * @param s2 Scalar + * @return v1 = v1 * s1 - v2 * s2 + */ + public final static double[] timesMinusTimesEquals(final double[] v1, final double s1, final double[] v2, final double s2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + for(int i = 0; i < v1.length; i++) { + v1[i] = v1[i] * s1 - v2[i] * s2; + } + return v1; + } + + /** + * Compute v1 - d + * + * @param v1 original vector + * @param d Value to subtract + * @return v1 - d + */ + public final static double[] minus(final double[] v1, final double d) { + final double[] result = new double[v1.length]; + for(int i = 0; i < v1.length; i++) { + result[i] = v1[i] - d; + } + return result; + } + + /** + * Computes v1 = v1 - d, overwriting v1 + * + * @param v1 original vector + * @param d Value to subtract + * @return v1 = v1 - d + */ + public final static double[] minusEquals(final double[] v1, final double d) { + for(int i = 0; i < v1.length; i++) { + v1[i] -= d; + } + return v1; + } + + /** + * Computes v1 * s1 + * + * @param v1 original vector + * @param s1 the scalar to be multiplied + * @return v1 * s1 + */ + public final static double[] times(final double[] v1, final double s1) { + final double[] v = new double[v1.length]; + for(int i = 0; i < v1.length; i++) { + v[i] = v1[i] * s1; + } + return v; + } + + /** + * Computes v1 = v1 * s1, overwritings v1 + * + * @param v1 original vector + * @param s scalar + * @return v1 = v1 * s1 + */ + public final static double[] timesEquals(final double[] v1, final double s) { + for(int i = 0; i < v1.length; i++) { + v1[i] *= s; + } + return v1; + } + + /** + * Matrix multiplication: v1 * m2 + * + * @param v1 vector + * @param m2 other matrix + * @return Matrix product, v1 * m2 + */ + public final static double[][] times(final double[] v1, final double[][] m2) { + assert (m2.length == 1) : "Matrix inner dimensions must agree."; + final int columndimension = m2[0].length; + final double[][] re = new double[v1.length][columndimension]; + for(int j = 0; j < columndimension; j++) { + for(int i = 0; i < v1.length; i++) { + re[i][j] = v1[i] * m2[0][j]; + } + } + return re; + } + + /** + * Linear algebraic matrix multiplication, v1<sup>T</sup> * m2 + * + * @param v1 vector + * @param m2 other matrix + * @return Matrix product, v1<sup>T</sup> * m2 + */ + public final static double[][] transposeTimes(final double[] v1, final double[][] m2) { + assert (m2.length == v1.length) : "Matrix inner dimensions must agree."; + final int columndimension = m2[0].length; + final double[][] re = new double[1][columndimension]; + for(int j = 0; j < columndimension; j++) { + double s = 0; + for(int k = 0; k < v1.length; k++) { + s += v1[k] * m2[k][j]; + } + re[0][j] = s; + } + return re; + } + + /** + * Linear algebraic matrix multiplication, v1<sup>T</sup> * v2 + * + * @param v1 vector + * @param v2 other vector + * @return Matrix product, v1<sup>T</sup> * v2 + */ + public final static double transposeTimes(final double[] v1, final double[] v2) { + assert (v2.length == v1.length) : "Matrix inner dimensions must agree."; + double s = 0; + for(int k = 0; k < v1.length; k++) { + s += v1[k] * v2[k]; + } + return s; + } + + /** + * Linear algebraic matrix multiplication, v1 * m2^T + * + * @param v1 vector + * @param m2 other matrix + * @return Matrix product, v1 * m2^T + */ + public final static double[][] timesTranspose(final double[] v1, final double[][] m2) { + assert (m2[0].length == 1) : "Matrix inner dimensions must agree."; + + final double[][] re = new double[v1.length][m2.length]; + for(int j = 0; j < m2.length; j++) { + for(int i = 0; i < v1.length; i++) { + re[i][j] = v1[i] * m2[j][0]; + } + } + return re; + } + + /** + * Linear algebraic matrix multiplication, v1 * v2^T + * + * @param v1 vector + * @param v2 other vector + * @return Matrix product, v1 * v2^T + */ + public final static double[][] timesTranspose(final double[] v1, final double[] v2) { + final double[][] re = new double[v1.length][v2.length]; + for(int j = 0; j < v2.length; j++) { + for(int i = 0; i < v1.length; i++) { + re[i][j] = v1[i] * v2[j]; + } + } + return re; + } + + /** + * Returns the scalar product of this vector and the specified vector v. + * + * This is the same as transposeTimes. + * + * @param v1 vector + * @param v2 other vector + * @return double the scalar product of vectors v1 and v2 + */ + public final static double scalarProduct(final double[] v1, final double[] v2) { + assert (v1.length == v2.length) : "Vector dimensions must agree."; + double scalarProduct = 0.0; + for(int row = 0; row < v1.length; row++) { + scalarProduct += v1[row] * v2[row]; + } + return scalarProduct; + } + + /** + * Euclidean length of the vector + * + * @param v1 vector + * @return euclidean length of this vector + */ + public final static double euclideanLength(final double[] v1) { + double acc = 0.0; + for(int row = 0; row < v1.length; row++) { + final double v = v1[row]; + acc += v * v; + } + return Math.sqrt(acc); + } + + /** + * Normalizes v1 to the length of 1.0. + * + * @param v1 vector + * @return normalized copy of v1 + */ + public final static double[] normalize(final double[] v1) { + double norm = euclideanLength(v1); + double[] re = new double[v1.length]; + if(norm != 0) { + for(int row = 0; row < v1.length; row++) { + re[row] = v1[row] / norm; + } + } + return re; + } + + /** + * Normalizes v1 to the length of 1.0. + * + * @param v1 vector + * @return normalized v1 + */ + public final static double[] normalizeEquals(final double[] v1) { + double norm = euclideanLength(v1); + if(norm != 0) { + for(int row = 0; row < v1.length; row++) { + v1[row] /= norm; + } + } + return v1; + } + + /** + * Projects this row vector into the subspace formed by the specified matrix + * v. + * + * @param m2 the subspace matrix + * @return the projection of p into the subspace formed by v + */ + public final static double[] project(final double[] v1, final double[][] m2) { + assert (v1.length == m2.length) : "v1 and m2 differ in dimensionality!"; + final int columndimension = m2[0].length; + + double[] sum = new double[v1.length]; + for(int i = 0; i < columndimension; i++) { + // TODO: optimize - copy less. + double[] v_i = getCol(m2, i); + plusTimesEquals(sum, v_i, scalarProduct(v1, v_i)); + } + return sum; + } + + /** + * Compute the hash code for the vector + * + * @param v1 elements + * @return hash code + */ + public final static int hashCode(final double[] v1) { + return Arrays.hashCode(v1); + } + + /** + * Compare for equality. + * + * @param v1 first vector + * @param v2 second vector + * @return comparison result + */ + public final static boolean equals(final double[] v1, final double[] v2) { + return Arrays.equals(v1, v2); + } + + /** + * Reset the Vector to 0. + * + * @param v1 vector + */ + public final static void clear(final double[] v1) { + Arrays.fill(v1, 0.0); + } + + /** + * Rotate vector by 90 degrees. + * + * @param v1 first vector + * @return modified v1, rotated by 90 degrees + */ + public final static double[] rotate90Equals(final double[] v1) { + assert (v1.length == 2) : "rotate90Equals is only valid for 2d vectors."; + double temp = v1[0]; + v1[0] = v1[1]; + v1[1] = -temp; + return v1; + } + + // *********** MATRIX operations + + /** + * Returns the unit matrix of the specified dimension. + * + * @param dim the dimensionality of the unit matrix + * @return the unit matrix of the specified dimension + */ + public final static double[][] unitMatrix(final int dim) { + final double[][] e = new double[dim][dim]; + for(int i = 0; i < dim; i++) { + e[i][i] = 1; + } + return e; + } + + /** + * Returns the zero matrix of the specified dimension. + * + * @param dim the dimensionality of the unit matrix + * @return the zero matrix of the specified dimension + */ + public final static double[][] zeroMatrix(final int dim) { + final double[][] z = new double[dim][dim]; + return z; + } + + /** + * Generate matrix with random elements + * + * @param m Number of rows. + * @param n Number of columns. + * @return An m-by-n matrix with uniformly distributed random elements. + */ + public final static double[][] random(final int m, final int n) { + final double[][] A = new double[m][n]; + for(int i = 0; i < m; i++) { + for(int j = 0; j < n; j++) { + A[i][j] = Math.random(); + } + } + return A; + } + + /** + * Generate identity matrix + * + * @param m Number of rows. + * @param n Number of columns. + * @return An m-by-n matrix with ones on the diagonal and zeros elsewhere. + */ + public final static double[][] identity(final int m, final int n) { + final double[][] A = new double[m][n]; + for(int i = 0; i < Math.min(m, n); i++) { + A[i][i] = 1.0; + } + return A; + } + + /** + * Returns a quadratic Matrix consisting of zeros and of the given values on + * the diagonal. + * + * @param v1 the values on the diagonal + * @return the resulting matrix + */ + public final static double[][] diagonal(final double[] v1) { + final double[][] result = new double[v1.length][v1.length]; + for(int i = 0; i < v1.length; i++) { + result[i][i] = v1[i]; + } + return result; + } + + /** + * Make a deep copy of a matrix. + * + * @param m1 Input matrix + * @return a new matrix containing the same values as this matrix + */ + public static final double[][] copy(final double[][] m1) { + final int columndimension = m1[0].length; + final double[][] X = new double[m1.length][columndimension]; + for(int i = 0; i < m1.length; i++) { + System.arraycopy(m1[i], 0, X[i], 0, columndimension); + } + return X; + } + + /** + * Make a one-dimensional row packed copy of the internal array. + * + * @param m1 Input matrix + * @return Matrix elements packed in a one-dimensional array by rows. + */ + public static final double[] rowPackedCopy(final double[][] m1) { + final int columndimension = m1[0].length; + double[] vals = new double[m1.length * columndimension]; + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + vals[i * columndimension + j] = m1[i][j]; + } + } + return vals; + } + + /** + * Make a one-dimensional column packed copy of the internal array. + * + * @param m1 Input matrix + * @return Matrix elements packed in a one-dimensional array by columns. + */ + public static final double[] columnPackedCopy(final double[][] m1) { + final int columndimension = m1[0].length; + final double[] vals = new double[m1.length * columndimension]; + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + vals[i + j * m1.length] = m1[i][j]; + } + } + return vals; + } + + /** + * Get a submatrix. + * + * @param m1 Input matrix + * @param r0 Initial row index + * @param r1 Final row index + * @param c0 Initial column index + * @param c1 Final column index + * @return m1(r0:r1,c0:c1) + */ + public static final double[][] getMatrix(final double[][] m1, final int r0, final int r1, final int c0, final int c1) { + final double[][] X = new double[r1 - r0 + 1][c1 - c0 + 1]; + for(int i = r0; i <= r1; i++) { + for(int j = c0; j <= c1; j++) { + X[i - r0][j - c0] = m1[i][j]; + } + } + return X; + } + + /** + * Get a submatrix. + * + * @param m1 Input matrix + * @param r Array of row indices. + * @param c Array of column indices. + * @return m1(r(:),c(:)) + */ + public static final double[][] getMatrix(final double[][] m1, final int[] r, final int[] c) { + final double[][] X = new double[r.length][c.length]; + for(int i = 0; i < r.length; i++) { + for(int j = 0; j < c.length; j++) { + X[i][j] = m1[r[i]][c[j]]; + } + } + return X; + } + + /** + * Get a submatrix. + * + * @param m1 Input matrix + * @param r Array of row indices. + * @param c0 Initial column index + * @param c1 Final column index + * @return m1(r(:),c0:c1) + */ + public static final double[][] getMatrix(final double[][] m1, final int[] r, final int c0, final int c1) { + final double[][] X = new double[r.length][c1 - c0 + 1]; + for(int i = 0; i < r.length; i++) { + for(int j = c0; j <= c1; j++) { + X[i][j - c0] = m1[r[i]][j]; + } + } + return X; + } + + /** + * Get a submatrix. + * + * @param m1 Input matrix + * @param r0 Initial row index + * @param r1 Final row index + * @param c Array of column indices. + * @return m1(r0:r1,c(:)) + */ + public static final double[][] getMatrix(final double[][] m1, final int r0, final int r1, final int[] c) { + final double[][] X = new double[r1 - r0 + 1][c.length]; + for(int i = r0; i <= r1; i++) { + for(int j = 0; j < c.length; j++) { + X[i - r0][j] = m1[i][c[j]]; + } + } + return X; + } + + /** + * Set a submatrix. + * + * @param m1 Original matrix + * @param r0 Initial row index + * @param r1 Final row index + * @param c0 Initial column index + * @param c1 Final column index + * @param m2 New values for m1(r0:r1,c0:c1) + */ + public static final void setMatrix(final double[][] m1, final int r0, final int r1, final int c0, final int c1, final double[][] m2) { + for(int i = r0; i <= r1; i++) { + for(int j = c0; j <= c1; j++) { + m1[i][j] = m2[i - r0][j - c0]; + } + } + } + + /** + * Set a submatrix. + * + * @param m1 Original matrix + * @param r Array of row indices. + * @param c Array of column indices. + * @param m2 New values for m1(r(:),c(:)) + */ + public static final void setMatrix(final double[][] m1, final int[] r, final int[] c, final double[][] m2) { + for(int i = 0; i < r.length; i++) { + for(int j = 0; j < c.length; j++) { + m1[r[i]][c[j]] = m2[i][j]; + } + } + } + + /** + * Set a submatrix. + * + * @param m1 Input matrix + * @param r Array of row indices. + * @param c0 Initial column index + * @param c1 Final column index + * @param m2 New values for m1(r(:),c0:c1) + */ + public static final void setMatrix(final double[][] m1, final int[] r, final int c0, final int c1, final double[][] m2) { + for(int i = 0; i < r.length; i++) { + for(int j = c0; j <= c1; j++) { + m1[r[i]][j] = m2[i][j - c0]; + } + } + } + + /** + * Set a submatrix. + * + * @param m1 Input matrix + * @param r0 Initial row index + * @param r1 Final row index + * @param c Array of column indices. + * @param m2 New values for m1(r0:r1,c(:)) + */ + public static final void setMatrix(final double[][] m1, final int r0, final int r1, final int[] c, final double[][] m2) { + for(int i = r0; i <= r1; i++) { + for(int j = 0; j < c.length; j++) { + m1[i][c[j]] = m2[i - r0][j]; + } + } + } + + /** + * Returns the <code>r</code>th row of this matrix as vector. + * + * @param m1 Input matrix + * @param r the index of the row to be returned + * @return the <code>r</code>th row of this matrix + */ + public static final double[] getRow(final double[][] m1, final int r) { + return m1[r].clone(); + } + + /** + * Sets the <code>r</code>th row of this matrix to the specified vector. + * + * @param m1 Original matrix + * @param r the index of the column to be set + * @param row the value of the column to be set + */ + public static final void setRow(final double[][] m1, final int r, final double[] row) { + final int columndimension = getColumnDimensionality(m1); + assert (row.length == columndimension) : "Matrix must consist of the same no of columns!"; + for(int i = 0; i < columndimension; i++) { + m1[r][i] = row[i]; + } + } + + /** + * Get a column from a matrix as vector. + * + * @param m1 Matrix to extract the column from + * @param col Column number + * @return Column + */ + public final static double[] getCol(double[][] m1, int col) { + double[] ret = new double[m1.length]; + for(int i = 0; i < ret.length; i++) { + ret[i] = m1[i][col]; + } + return ret; + } + + /** + * Sets the <code>c</code>th column of this matrix to the specified column. + * + * @param m1 Input matrix + * @param c the index of the column to be set + * @param column the value of the column to be set + */ + public static final void setCol(final double[][] m1, final int c, final double[] column) { + assert (column.length == m1.length) : "Matrix must consist of the same no of rows!"; + for(int i = 0; i < m1.length; i++) { + m1[i][c] = column[i]; + } + } + + /** + * Matrix transpose + * + * @param m1 Input matrix + * @return m1<sup>T</sup> as copy + */ + public static final double[][] transpose(final double[][] m1) { + final int columndimension = getColumnDimensionality(m1); + final double[][] re = new double[columndimension][m1.length]; + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + re[j][i] = m1[i][j]; + } + } + return re; + } + + /** + * m3 = m1 + m2 + * + * @param m1 Input matrix + * @param m2 another matrix + * @return m1 + m1 in a new Matrix + */ + public static final double[][] plus(final double[][] m1, final double[][] m2) { + return plusEquals(copy(m1), m2); + } + + /** + * m3 = m1 + s2 * m2 + * + * @param m1 Input matrix + * @param m2 another matrix + * @param s2 scalar + * @return m1 + s2 * m2 in a new Matrix + */ + public static final double[][] plusTimes(final double[][] m1, final double[][] m2, final double s2) { + return plusTimesEquals(copy(m1), m2, s2); + } + + /** + * m1 = m1 + m2, overwriting m1 + * + * @param m1 input matrix + * @param m2 another matrix + * @return m1 = m1 + m2 + */ + public static final double[][] plusEquals(final double[][] m1, final double[][] m2) { + final int columndimension = getColumnDimensionality(m1); + assert (getRowDimensionality(m1) == getRowDimensionality(m2) && columndimension == getColumnDimensionality(m2)) : "Matrix dimensions must agree."; + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + m1[i][j] += m2[i][j]; + } + } + return m1; + } + + /** + * m1 = m1 + s2 * m2, overwriting m1 + * + * @param m1 input matrix + * @param m2 another matrix + * @param s2 scalar for s2 + * @return m1 = m1 + s2 * m2, overwriting m1 + */ + public static final double[][] plusTimesEquals(final double[][] m1, final double[][] m2, final double s2) { + final int columndimension = getColumnDimensionality(m1); + assert (getRowDimensionality(m1) == getRowDimensionality(m2) && columndimension == getColumnDimensionality(m2)) : "Matrix dimensions must agree."; + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + m1[i][j] += s2 * m2[i][j]; + } + } + return m1; + } + + /** + * m3 = m1 - m2 + * + * @param m1 Input matrix + * @param m2 another matrix + * @return m1 - m2 in a new matrix + */ + public static final double[][] minus(final double[][] m1, final double[][] m2) { + return minusEquals(copy(m1), m2); + } + + /** + * m3 = m1 - s2 * m2 + * + * @param m1 Input matrix + * @param m2 another matrix + * @param s2 Scalar + * @return m1 - s2 * m2 in a new Matrix + */ + public static final double[][] minusTimes(final double[][] m1, final double[][] m2, final double s2) { + return minusTimesEquals(copy(m1), m2, s2); + } + + /** + * m1 = m1 - m2, overwriting m1 + * + * @param m1 Input matrix + * @param m2 another matrix + * @return m1 - m2, overwriting m1 + */ + public static final double[][] minusEquals(final double[][] m1, final double[][] m2) { + final int columndimension = getColumnDimensionality(m1); + assert (getRowDimensionality(m1) == getRowDimensionality(m2) && columndimension == getColumnDimensionality(m2)) : "Matrix dimensions must agree."; + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + m1[i][j] -= m2[i][j]; + } + } + return m1; + } + + /** + * m1 = m1 - s2 * m2, overwriting m1 + * + * @param m1 Input matrix + * @param m2 another matrix + * @param s2 Scalar + * @return m1 = m1 - s2 * m2, overwriting m1 + */ + public static final double[][] minusTimesEquals(final double[][] m1, final double[][] m2, final double s2) { + assert (getRowDimensionality(m1) == getRowDimensionality(m2) && getColumnDimensionality(m1) == getColumnDimensionality(m2)) : "Matrix dimensions must agree."; + for(int i = 0; i < m1.length; i++) { + final double[] row1 = m1[i]; + final double[] row2 = m2[i]; + for(int j = 0; j < row1.length; j++) { + row1[j] -= s2 * row2[j]; + } + } + return m1; + } + + /** + * Multiply a matrix by a scalar, m3 = s1*m1 + * + * @param m1 Input matrix + * @param s1 scalar + * @return s1*m1, in a new matrix + */ + public static final double[][] times(final double[][] m1, final double s1) { + return timesEquals(copy(m1), s1); + } + + /** + * Multiply a matrix by a scalar in place, m1 = s1 * m1 + * + * @param m1 Input matrix + * @param s1 scalar + * @return m1 = s1 * m1, overwriting m1 + */ + public static final double[][] timesEquals(final double[][] m1, final double s1) { + for(int i = 0; i < m1.length; i++) { + final double[] row = m1[i]; + for(int j = 0; j < row.length; j++) { + row[j] *= s1; + } + } + return m1; + } + + /** + * Linear algebraic matrix multiplication, m1 * m2 + * + * @param m1 Input matrix + * @param m2 another matrix + * @return Matrix product, m1 * m2 + */ + public static final double[][] times(final double[][] m1, final double[][] m2) { + final int columndimension = getColumnDimensionality(m1); + final int bcolumndimension = getColumnDimensionality(m2); + // Optimized implementation, exploiting the storage layout + assert (m2.length == columndimension) : "Matrix inner dimensions must agree: " + getRowDimensionality(m1) + "," + getColumnDimensionality(m1) + " * " + getRowDimensionality(m2) + "," + getColumnDimensionality(m2); + final double[][] r2 = new double[m1.length][bcolumndimension]; + // Optimized ala Jama. jik order. + final double[] Bcolj = new double[columndimension]; + for(int j = 0; j < bcolumndimension; j++) { + // Make a linear copy of column j from B + // TODO: use column getter from B? + for(int k = 0; k < columndimension; k++) { + Bcolj[k] = m2[k][j]; + } + // multiply it with each row from A + for(int i = 0; i < m1.length; i++) { + final double[] Arowi = m1[i]; + double s = 0; + for(int k = 0; k < columndimension; k++) { + s += Arowi[k] * Bcolj[k]; + } + r2[i][j] = s; + } + } + return r2; + } + + /** + * Linear algebraic matrix multiplication, m1 * v2 + * + * @param m1 Input matrix + * @param v2 a vector + * @return Matrix product, m1 * v2 + */ + public static final double[] times(final double[][] m1, final double[] v2) { + assert (v2.length == getColumnDimensionality(m1)) : "Matrix inner dimensions must agree."; + final double[] re = new double[m1.length]; + // multiply it with each row from A + for(int i = 0; i < m1.length; i++) { + final double[] Arowi = m1[i]; + double s = 0; + for(int k = 0; k < Arowi.length; k++) { + s += Arowi[k] * v2[k]; + } + re[i] = s; + } + return re; + } + + /** + * Linear algebraic matrix multiplication, m1<sup>T</sup> * v2 + * + * @param m1 Input matrix + * @param v2 another matrix + * @return Matrix product, m1<sup>T</sup> * v2 + */ + public static final double[] transposeTimes(final double[][] m1, final double[] v2) { + final int columndimension = getColumnDimensionality(m1); + assert (v2.length == m1.length) : "Matrix inner dimensions must agree."; + final double[] re = new double[columndimension]; + // multiply it with each row from A + for(int i = 0; i < columndimension; i++) { + double s = 0; + for(int k = 0; k < m1.length; k++) { + s += m1[k][i] * v2[k]; + } + re[i] = s; + } + return re; + } + + /** + * Linear algebraic matrix multiplication, m1<sup>T</sup> * m2 + * + * @param m1 Input matrix + * @param m2 another matrix + * @return Matrix product, m1<sup>T</sup> * m2 + */ + public static final double[][] transposeTimes(final double[][] m1, final double[][] m2) { + final int coldim1 = getColumnDimensionality(m1); + final int coldim2 = getColumnDimensionality(m2); + assert (m2.length == m1.length) : "Matrix inner dimensions must agree."; + final double[][] re = new double[coldim1][coldim2]; + final double[] Bcolj = new double[m1.length]; + for(int j = 0; j < coldim2; j++) { + // Make a linear copy of column j from B + for(int k = 0; k < m1.length; k++) { + Bcolj[k] = m2[k][j]; + } + // multiply it with each row from A + for(int i = 0; i < coldim1; i++) { + double s = 0; + for(int k = 0; k < m1.length; k++) { + s += m1[k][i] * Bcolj[k]; + } + re[i][j] = s; + } + } + return re; + } + + /** + * Linear algebraic matrix multiplication, m1 * m2^T + * + * @param m1 Input matrix + * @param m2 another matrix + * @return Matrix product, m1 * m2^T + */ + public static final double[][] timesTranspose(final double[][] m1, final double[][] m2) { + assert (getColumnDimensionality(m2) == getColumnDimensionality(m1)) : "Matrix inner dimensions must agree."; + final double[][] re = new double[m1.length][m2.length]; + for(int j = 0; j < re.length; j++) { + final double[] Browj = m2[j]; + // multiply it with each row from A + for(int i = 0; i < m1.length; i++) { + final double[] Arowi = m1[i]; + double s = 0; + for(int k = 0; k < Browj.length; k++) { + s += Arowi[k] * Browj[k]; + } + re[i][j] = s; + } + } + return re; + } + + /** + * Linear algebraic matrix multiplication, m1^T * m2^T. Computed as (m2*m1)^T + * + * @param m1 Input matrix + * @param m2 another matrix + * @return Matrix product, m1^T * m2^T + */ + public static final double[][] transposeTimesTranspose(final double[][] m1, final double[][] m2) { + // Optimized implementation, exploiting the storage layout + assert (m1.length == getColumnDimensionality(m2)) : "Matrix inner dimensions must agree: " + getRowDimensionality(m1) + "," + getColumnDimensionality(m1) + " * " + getRowDimensionality(m2) + "," + getColumnDimensionality(m2); + final double[][] re = new double[getColumnDimensionality(m1)][m2.length]; + // Optimized ala Jama. jik order. + final double[] Acolj = new double[m1.length]; + for(int j = 0; j < re.length; j++) { + // Make a linear copy of column j from B + for(int k = 0; k < m1.length; k++) { + Acolj[k] = m1[k][j]; + } + final double[] Xrow = re[j]; + // multiply it with each row from A + for(int i = 0; i < m2.length; i++) { + final double[] Browi = m2[i]; + double s = 0; + for(int k = 0; k < m1.length; k++) { + s += Browi[k] * Acolj[k]; + } + Xrow[i] = s; + } + } + return re; + } + + /** + * getDiagonal returns array of diagonal-elements. + * + * @param m1 Input matrix + * @return values on the diagonal of the Matrix + */ + public final static double[] getDiagonal(final double[][] m1) { + final int dim = Math.min(getColumnDimensionality(m1), m1.length); + final double[] diagonal = new double[dim]; + for(int i = 0; i < dim; i++) { + diagonal[i] = m1[i][i]; + } + return diagonal; + } + + /** + * Normalizes the columns of this matrix to length of 1.0. + * + * @param m1 Input matrix + */ + public final static void normalizeColumns(final double[][] m1) { + final int columndimension = getColumnDimensionality(m1); + for(int col = 0; col < columndimension; col++) { + double norm = 0.0; + for(int row = 0; row < m1.length; row++) { + norm = norm + (m1[row][col] * m1[row][col]); + } + norm = Math.sqrt(norm); + if(norm != 0) { + for(int row = 0; row < m1.length; row++) { + m1[row][col] /= norm; + } + } + else { + // TODO: else: throw an exception? + } + } + } + + /** + * Returns a matrix which consists of this matrix and the specified columns. + * + * @param m1 Input matrix + * @param m2 the columns to be appended + * @return the new matrix with the appended columns + */ + public static final double[][] appendColumns(final double[][] m1, final double[][] m2) { + final int columndimension = getColumnDimensionality(m1); + final int ccolumndimension = getColumnDimensionality(m2); + assert (m1.length == m2.length) : "m.getRowDimension() != column.getRowDimension()"; + + final int rcolumndimension = columndimension + ccolumndimension; + final double[][] result = new double[m1.length][rcolumndimension]; + for(int i = 0; i < rcolumndimension; i++) { + // FIXME: optimize - excess copying! + if(i < columndimension) { + setCol(result, i, getCol(m1, i)); + } + else { + setCol(result, i, getCol(m2, i - columndimension)); + } + } + return result; + } + + /** + * Returns an orthonormalization of this matrix. + * + * @param m1 Input matrix + * @return the orthonormalized matrix + */ + public static final double[][] orthonormalize(final double[][] m1) { + final int columndimension = getColumnDimensionality(m1); + double[][] v = copy(m1); + + // FIXME: optimize - excess copying! + for(int i = 1; i < columndimension; i++) { + final double[] u_i = getCol(m1, i); + final double[] sum = new double[m1.length]; + for(int j = 0; j < i; j++) { + final double[] v_j = getCol(v, j); + double scalar = scalarProduct(u_i, v_j) / scalarProduct(v_j, v_j); + plusEquals(sum, times(v_j, scalar)); + } + final double[] v_i = minus(u_i, sum); + setCol(v, i, v_i); + } + + normalizeColumns(v); + return v; + } + + /** + * Compute hash code + * + * @param m1 Input matrix + * @return Hash code + */ + public static final int hashCode(final double[][] m1) { + return Arrays.hashCode(m1); + } + + /** + * Test for equality + * + * @param m1 Input matrix + * @param m2 Other matrix + * @return Equality + */ + public static final boolean equals(final double[][] m1, final double[][] m2) { + return Arrays.equals(m1, m2); + } + + /** + * Compare two matrices with a delta parameter to take numerical errors into + * account. + * + * @param m1 Input matrix + * @param m2 other matrix to compare with + * @param maxdelta maximum delta allowed + * @return true if delta smaller than maximum + */ + public static final boolean almostEquals(final double[][] m1, final double[][] m2, final double maxdelta) { + if(m1 == m2) { + return true; + } + if(m2 == null) { + return false; + } + if(m1.getClass() != m2.getClass()) { + return false; + } + if(m1.length != m2.length) { + return false; + } + final int columndimension = getColumnDimensionality(m1); + if(columndimension != getColumnDimensionality(m2)) { + return false; + } + for(int i = 0; i < m1.length; i++) { + for(int j = 0; j < columndimension; j++) { + if(Math.abs(m1[i][j] - m2[i][j]) > maxdelta) { + return false; + } + } + } + return true; + } + + /** + * Compare two matrices with a delta parameter to take numerical errors into + * account. + * + * @param m1 Input matrix + * @param m2 other matrix to compare with + * @return almost equals with delta {@link #DELTA} + */ + public static final boolean almostEquals(final double[][] m1, final double[][] m2) { + return almostEquals(m1, m2, DELTA); + } + + /** + * Returns the dimensionality of the rows of this matrix. + * + * @param m1 Input matrix + * @return the number of rows. + */ + public static final int getRowDimensionality(final double[][] m1) { + return m1.length; + } + + /** + * Returns the dimensionality of the columns of this matrix. + * + * @param m1 Input matrix + * @return the number of columns. + */ + public static final int getColumnDimensionality(final double[][] m1) { + return m1[0].length; + } +}
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Vector.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Vector.java index 8cf3cf78..0a674d87 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Vector.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/Vector.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -23,11 +23,12 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra; along with this program. If not, see <http://www.gnu.org/licenses/>. */ -import java.io.Serializable; import java.util.Arrays; -import de.lmu.ifi.dbs.elki.math.MathUtil; +import de.lmu.ifi.dbs.elki.data.NumberVector; import de.lmu.ifi.dbs.elki.utilities.FormatUtil; +import de.lmu.ifi.dbs.elki.utilities.datastructures.arraylike.ArrayAdapter; +import de.lmu.ifi.dbs.elki.utilities.datastructures.arraylike.NumberArrayAdapter; /** * Provides a vector object that encapsulates an m x 1 - matrix object. @@ -36,12 +37,7 @@ import de.lmu.ifi.dbs.elki.utilities.FormatUtil; * * @apiviz.landmark */ -public class Vector implements MatrixLike<Vector>, Serializable { - /** - * Serial version - */ - private static final long serialVersionUID = 1L; - +public class Vector implements NumberVector<Vector, Double> { /** * Array for internal storage of elements. * @@ -73,13 +69,18 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param dimensionality dimensionality * @return the dimensionality of the vector */ - // FIXME: may also return null vector by chance. public static final Vector randomNormalizedVector(final int dimensionality) { final Vector v = new Vector(dimensionality); - for(int i = 0; i < dimensionality; i++) { - v.elements[i] = Math.random(); + double norm = 0; + while(norm <= 0) { + for(int i = 0; i < dimensionality; i++) { + v.elements[i] = Math.random(); + } + norm = v.euclideanLength(); + } + for(int row = 0; row < dimensionality; row++) { + v.elements[row] /= norm; } - v.normalize(); return v; } @@ -101,7 +102,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { * * @return a copy of this vector */ - @Override public final Vector copy() { return new Vector(elements.clone()); } @@ -137,20 +137,11 @@ public class Vector implements MatrixLike<Vector>, Serializable { * * @return the dimensionality of this vector */ - public final int getDimensionality() { - return elements.length; - } - @Override - public final int getRowDimensionality() { + public final int getDimensionality() { return elements.length; } - @Override - public final int getColumnDimensionality() { - return 1; - } - /** * Returns the value at the specified row. * @@ -161,14 +152,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { return elements[i]; } - @Override - public final double get(final int i, final int j) { - if(j != 0) { - throw new ArrayIndexOutOfBoundsException(); - } - return elements[i]; - } - /** * Sets the value at the specified row. * @@ -182,37 +165,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { return this; } - @Override - public final Vector set(final int i, final int j, final double s) { - if(j != 0) { - throw new ArrayIndexOutOfBoundsException(); - } - elements[i] = s; - return this; - } - - @Override - public final Vector increment(final int i, final int j, final double s) { - if(j != 0) { - throw new ArrayIndexOutOfBoundsException(); - } - elements[i] += s; - return this; - } - - @Override - public final Vector getColumnVector(final int i) { - if(i != 0) { - throw new ArrayIndexOutOfBoundsException(); - } - return this; - } - - @Override - public final Matrix transpose() { - return new Matrix(this.elements, 1); - } - /** * Returns a new vector which is the result of this vector plus the specified * vector. @@ -220,9 +172,8 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param v the vector to be added * @return the resulting vector */ - @Override public final Vector plus(final Vector v) { - checkDimensions(v); + assert (this.elements.length == v.elements.length) : "Vector dimensions must agree."; final Vector result = new Vector(elements.length); for(int i = 0; i < elements.length; i++) { result.elements[i] = elements[i] + v.elements[i]; @@ -238,9 +189,8 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param s the scalar * @return the resulting vector */ - @Override public final Vector plusTimes(final Vector v, final double s) { - checkDimensions(v); + assert (this.elements.length == v.elements.length) : "Vector dimensions must agree."; final Vector result = new Vector(elements.length); for(int i = 0; i < elements.length; i++) { result.elements[i] = elements[i] + v.elements[i] * s; @@ -254,11 +204,10 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param B another matrix * @return A + B in this Matrix */ - @Override public final Vector plusEquals(final Vector B) { - checkDimensions(B); + assert (this.elements.length == B.elements.length) : "Vector dimensions must agree."; for(int i = 0; i < elements.length; i++) { - elements[i] += B.get(i, 0); + elements[i] += B.elements[i]; } return this; } @@ -270,11 +219,10 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param s Scalar * @return A + s * B in this Matrix */ - @Override public final Vector plusTimesEquals(final Vector B, final double s) { - checkDimensions(B); + assert (this.elements.length == B.elements.length) : "Vector dimensions must agree."; for(int i = 0; i < elements.length; i++) { - elements[i] += s * B.get(i, 0); + elements[i] += s * B.elements[i]; } return this; } @@ -298,7 +246,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param v the vector to be subtracted from this vector * @return this vector minus the specified vector v */ - @Override public final Vector minus(final Vector v) { final Vector sub = new Vector(elements.length); for(int i = 0; i < elements.length; i++) { @@ -314,7 +261,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param s the scaling factor * @return this vector minus the specified vector v */ - @Override public final Vector minusTimes(final Vector v, final double s) { final Vector sub = new Vector(elements.length); for(int i = 0; i < elements.length; i++) { @@ -329,11 +275,10 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param B another matrix * @return A - B in this Matrix */ - @Override public final Vector minusEquals(final Vector B) { - checkDimensions(B); + assert (this.elements.length == B.elements.length) : "Vector dimensions must agree."; for(int i = 0; i < elements.length; i++) { - elements[i] -= B.get(i, 0); + elements[i] -= B.elements[i]; } return this; } @@ -345,11 +290,10 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param s Scalar * @return A - s * B in this Matrix */ - @Override public final Vector minusTimesEquals(final Vector B, final double s) { - checkDimensions(B); + assert (this.elements.length == B.elements.length) : "Vector dimensions must agree."; for(int i = 0; i < elements.length; i++) { - elements[i] -= s * B.get(i, 0); + elements[i] -= s * B.elements[i]; } return this; } @@ -374,7 +318,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param s the scalar to be multiplied * @return the resulting vector */ - @Override public final Vector times(final double s) { final Vector v = new Vector(elements.length); for(int i = 0; i < elements.length; i++) { @@ -389,7 +332,6 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @param s scalar * @return replace A by s*A */ - @Override public final Vector timesEquals(final double s) { for(int i = 0; i < elements.length; i++) { elements[i] *= s; @@ -402,12 +344,9 @@ public class Vector implements MatrixLike<Vector>, Serializable { * * @param B another matrix * @return Matrix product, A * B - * @throws IllegalArgumentException Matrix inner dimensions must agree. */ public final Matrix times(final Matrix B) { - if(B.elements.length != 1) { - throw new IllegalArgumentException("Matrix inner dimensions must agree."); - } + assert (B.elements.length == 1) : "Matrix inner dimensions must agree."; final Matrix X = new Matrix(this.elements.length, B.columndimension); for(int j = 0; j < B.columndimension; j++) { for(int i = 0; i < this.elements.length; i++) { @@ -422,12 +361,9 @@ public class Vector implements MatrixLike<Vector>, Serializable { * * @param B another matrix * @return Matrix product, A<sup>T</sup> * B - * @throws IllegalArgumentException Matrix inner dimensions must agree. */ public final Matrix transposeTimes(final Matrix B) { - if(B.elements.length != this.elements.length) { - throw new IllegalArgumentException("Matrix inner dimensions must agree."); - } + assert (B.elements.length == this.elements.length) : "Matrix inner dimensions must agree."; final Matrix X = new Matrix(1, B.columndimension); for(int j = 0; j < B.columndimension; j++) { // multiply it with each row from A @@ -441,16 +377,34 @@ public class Vector implements MatrixLike<Vector>, Serializable { } /** + * Linear algebraic matrix multiplication, a<sup>T</sup> * B * c + * + * @param B matrix + * @param c vector on the right + * @return Matrix product, a<sup>T</sup> * B + */ + public final double transposeTimesTimes(final Matrix B, final Vector c) { + assert (B.elements.length == this.elements.length) : "Matrix inner dimensions must agree."; + double sum = 0.0; + for(int j = 0; j < B.columndimension; j++) { + // multiply it with each row from A + double s = 0; + for(int k = 0; k < this.elements.length; k++) { + s += this.elements[k] * B.elements[k][j]; + } + sum += s * c.elements[j]; + } + return sum; + } + + /** * Linear algebraic matrix multiplication, A<sup>T</sup> * B * * @param B another vector * @return Matrix product, A<sup>T</sup> * B - * @throws IllegalArgumentException Matrix inner dimensions must agree. */ public final double transposeTimes(final Vector B) { - if(B.elements.length != this.elements.length) { - throw new IllegalArgumentException("Matrix inner dimensions must agree."); - } + assert (B.elements.length == this.elements.length) : "Matrix inner dimensions must agree."; double s = 0; for(int k = 0; k < this.elements.length; k++) { s += this.elements[k] * B.elements[k]; @@ -463,61 +417,32 @@ public class Vector implements MatrixLike<Vector>, Serializable { * * @param B another matrix * @return Matrix product, A * B^T - * @throws IllegalArgumentException Matrix inner dimensions must agree. */ public final Matrix timesTranspose(final Matrix B) { - if(B.columndimension != 1) { - throw new IllegalArgumentException("Matrix inner dimensions must agree."); - } + assert (B.columndimension == 1) : "Matrix inner dimensions must agree."; final Matrix X = new Matrix(this.elements.length, B.elements.length); for(int j = 0; j < B.elements.length; j++) { for(int i = 0; i < this.elements.length; i++) { - X.elements[i][j] = elements[i] * B.elements[0][j]; + X.elements[i][j] = elements[i] * B.elements[j][0]; } } return X; } /** - * Returns the scalar product of this vector and the specified vector v. - * - * @param v the vector - * @return double the scalar product of this vector and v - */ - public final double scalarProduct(final Vector v) { - checkDimensions(v); - double scalarProduct = 0.0; - for(int row = 0; row < elements.length; row++) { - final double prod = elements[row] * v.elements[row]; - scalarProduct += prod; - } - return scalarProduct; - } - - /** - * Inverts every element of the vector. - * - * @return the resulting vector - */ - public final Vector inverseVector() { - final Vector inv = new Vector(elements.length); - for(int i = 0; i < elements.length; i++) { - inv.elements[i] = 1.0 / elements[i]; - } - return inv; - } - - /** - * Square roots every element of the vector. + * Linear algebraic matrix multiplication, A * B^T * - * @return the resulting vector + * @param B another matrix + * @return Matrix product, A * B^T */ - public final Vector sqrtVector() { - final Vector sqrt = new Vector(elements.length); - for(int i = 0; i < elements.length; i++) { - sqrt.elements[i] = Math.sqrt(elements[i]); + public final Matrix timesTranspose(final Vector B) { + final Matrix X = new Matrix(this.elements.length, B.elements.length); + for(int j = 0; j < B.elements.length; j++) { + for(int i = 0; i < this.elements.length; i++) { + X.elements[i][j] = elements[i] * B.elements[j]; + } } - return sqrt; + return X; } /** @@ -526,27 +451,15 @@ public class Vector implements MatrixLike<Vector>, Serializable { * @return the length of this vector */ public final double euclideanLength() { - double sqlen = 0.0; + double acc = 0.0; for(int row = 0; row < elements.length; row++) { - sqlen += elements[row] * elements[row]; + final double v = elements[row]; + acc += v * v; } - return Math.sqrt(sqlen); + return Math.sqrt(acc); } /** - * Frobenius norm - * - * @return sqrt of sum of squares of all elements. - */ - public double normF() { - double f = 0; - for(int i = 0; i < elements.length; i++) { - f = MathUtil.hypotenuse(f, elements[i]); - } - return f; - } - - /** * Normalizes this vector to the length of 1.0. */ public final Vector normalize() { @@ -565,40 +478,21 @@ public class Vector implements MatrixLike<Vector>, Serializable { * * @param v the subspace matrix * @return the projection of p into the subspace formed by v - * @throws IllegalArgumentException if this matrix is no row vector, i.e. this - * matrix has more than one column or this matrix and v have different - * length of rows */ public final Vector projection(final Matrix v) { - if(elements.length != v.elements.length) { - throw new IllegalArgumentException("p and v differ in row dimensionality!"); - } + assert (elements.length == v.elements.length) : "p and v differ in row dimensionality!"; Vector sum = new Vector(elements.length); for(int i = 0; i < v.columndimension; i++) { - // TODO: optimize - copy less. - Vector v_i = v.getColumnVector(i); - sum.plusEquals(v_i.times(scalarProduct(v_i))); + // TODO: optimize - copy less? + Vector v_i = v.getCol(i); + sum.plusTimesEquals(v_i, this.transposeTimes(v_i)); } return sum; } - /** - * Check if this.getDimensionality() == v.getDimensionality(). - * - * @throws IllegalArgumentException if the dimensions do not agree - */ - private final void checkDimensions(final Vector v) { - if(this.elements.length != v.elements.length) { - throw new IllegalArgumentException("Vector dimensions must agree."); - } - } - @Override public int hashCode() { - final int PRIME = 31; - int result = 1; - result = PRIME * result + Arrays.hashCode(this.elements); - return result; + return Arrays.hashCode(this.elements); } @Override @@ -616,12 +510,7 @@ public class Vector implements MatrixLike<Vector>, Serializable { if(this.elements.length != other.elements.length) { return false; } - for(int i = 0; i < this.elements.length; i++) { - if(this.elements[i] != other.elements[i]) { - return false; - } - } - return true; + return Arrays.equals(this.elements, other.elements); } /** @@ -650,4 +539,95 @@ public class Vector implements MatrixLike<Vector>, Serializable { public void setZero() { Arrays.fill(elements, 0.0); } + + /** + * Rotate vector by 90 degrees. + * + * @return self, for operation chaining. + */ + public Vector rotate90Equals() { + assert (elements.length == 2); + double temp = elements[0]; + elements[0] = elements[1]; + elements[1] = -temp; + return this; + } + + // ////// NumberVector API. A bit hackish. :-( + + @Override + public double getMin(int dimension) { + return elements[dimension - 1]; + } + + @Override + public double getMax(int dimension) { + return elements[dimension - 1]; + } + + @Override + public Double getValue(int dimension) { + return elements[dimension - 1]; + } + + @Override + public double doubleValue(int dimension) { + return elements[dimension - 1]; + } + + @Override + public float floatValue(int dimension) { + return (float) elements[dimension - 1]; + } + + @Override + public int intValue(int dimension) { + return (int) elements[dimension - 1]; + } + + @Override + public long longValue(int dimension) { + return (long) elements[dimension - 1]; + } + + @Override + public short shortValue(int dimension) { + return (short) elements[dimension - 1]; + } + + @Override + public byte byteValue(int dimension) { + return (byte) elements[dimension - 1]; + } + + @Override + public Vector getColumnVector() { + return copy(); + } + + @Override + public Vector newNumberVector(double[] values) { + return new Vector(values); + } + + @Override + public <A> Vector newNumberVector(A array, NumberArrayAdapter<?, A> adapter) { + double[] raw = new double[adapter.size(array)]; + for(int i = 0; i < raw.length; i++) { + raw[i] = adapter.getDouble(array, i); + } + return new Vector(raw); + } + + @Override + public <A> Vector newFeatureVector(A array, ArrayAdapter<Double, A> adapter) { + if(adapter instanceof NumberArrayAdapter) { + return newNumberVector(array, (NumberArrayAdapter<?, A>) adapter); + } + double[] raw = new double[adapter.size(array)]; + for(int i = 0; i < raw.length; i++) { + raw[i] = adapter.get(array, i); + } + return new Vector(raw); + } }
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunction.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunction.java index 0cb67acc..b04c64d5 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunction.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunction.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.fitting; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunctionResult.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunctionResult.java index 46245349..76f90138 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunctionResult.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/FittingFunctionResult.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.fitting; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/GaussianFittingFunction.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/GaussianFittingFunction.java index 327b92b4..42badb2f 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/GaussianFittingFunction.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/GaussianFittingFunction.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.fitting; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/LevenbergMarquardtMethod.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/LevenbergMarquardtMethod.java index 8b8253dc..d87e1208 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/LevenbergMarquardtMethod.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/LevenbergMarquardtMethod.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.fitting; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/package-info.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/package-info.java index f2ad8c35..2c4143a8 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/package-info.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/fitting/package-info.java @@ -7,7 +7,7 @@ This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures -Copyright (C) 2011 +Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/package-info.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/package-info.java index d38b9a7b..5c55dfd1 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/package-info.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/package-info.java @@ -43,7 +43,7 @@ can be found at <a href="http://math.nist.gov/javanumerics/jama/">http://math.ni This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures -Copyright (C) 2011 +Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/AbstractCovarianceMatrixBuilder.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/AbstractCovarianceMatrixBuilder.java index 197ad939..c14986bd 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/AbstractCovarianceMatrixBuilder.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/AbstractCovarianceMatrixBuilder.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -54,17 +54,17 @@ public abstract class AbstractCovarianceMatrixBuilder<V extends NumberVector<? e public abstract Matrix processIds(DBIDs ids, Relation<? extends V> database); @Override - public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<DistanceResultPair<D>> results, Relation<? extends V> database, int k) { + public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database, int k) { ModifiableDBIDs ids = DBIDUtil.newArray(k); int have = 0; - for(Iterator<DistanceResultPair<D>> it = results.iterator(); it.hasNext() && have < k; have++) { + for(Iterator<? extends DistanceResultPair<D>> it = results.iterator(); it.hasNext() && have < k; have++) { ids.add(it.next().getDBID()); } return processIds(ids, database); } @Override - final public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<DistanceResultPair<D>> results, Relation<? extends V> database) { + final public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database) { return processQueryResults(results, database, results.size()); } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CompositeEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CompositeEigenPairFilter.java index 112a6df4..bc9486c5 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CompositeEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CompositeEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CovarianceMatrixBuilder.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CovarianceMatrixBuilder.java index 7416abce..5098ffac 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CovarianceMatrixBuilder.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/CovarianceMatrixBuilder.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -67,7 +67,7 @@ public interface CovarianceMatrixBuilder<V extends NumberVector<? extends V, ?>> * @param k the number of entries to process * @return Covariance Matrix */ - public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<DistanceResultPair<D>> results, Relation<? extends V> database, int k); + public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database, int k); /** * Compute Covariance Matrix for a QueryResult Collection @@ -78,5 +78,5 @@ public interface CovarianceMatrixBuilder<V extends NumberVector<? extends V, ?>> * @param database the database used * @return Covariance Matrix */ - public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<DistanceResultPair<D>> results, Relation<? extends V> database); + public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database); }
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/EigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/EigenPairFilter.java index 3c9d8603..553a111f 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/EigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/EigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FilteredEigenPairs.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FilteredEigenPairs.java index 9fdac000..a2c83249 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FilteredEigenPairs.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FilteredEigenPairs.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FirstNEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FirstNEigenPairFilter.java index c4839109..08482b71 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FirstNEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/FirstNEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/LimitEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/LimitEigenPairFilter.java index c1edb99a..e8d2b844 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/LimitEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/LimitEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/NormalizingEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/NormalizingEigenPairFilter.java index 47dbb01a..29be965c 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/NormalizingEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/NormalizingEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -86,7 +86,7 @@ public class NormalizingEigenPairFilter implements EigenPairFilter { */ private void normalizeEigenPair(final EigenPair eigenPair) { final Vector eigenvector = eigenPair.getEigenvector(); - final double scaling = 1.0 / Math.sqrt(eigenPair.getEigenvalue()) * eigenvector.normF(); + final double scaling = 1.0 / Math.sqrt(eigenPair.getEigenvalue()) * eigenvector.euclideanLength(); eigenvector.timesEquals(scaling); } }
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredAutotuningRunner.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredAutotuningRunner.java new file mode 100644 index 00000000..8b53dc43 --- /dev/null +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredAutotuningRunner.java @@ -0,0 +1,232 @@ +package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; + +/* + This file is part of ELKI: + Environment for Developing KDD-Applications Supported by Index-Structures + + Copyright (C) 2012 + Ludwig-Maximilians-Universität München + Lehr- und Forschungseinheit für Datenbanksysteme + ELKI Development Team + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU Affero General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU Affero General Public License for more details. + + You should have received a copy of the GNU Affero General Public License + along with this program. If not, see <http://www.gnu.org/licenses/>. + */ + +import java.util.ArrayList; +import java.util.Collection; +import java.util.Collections; +import java.util.Iterator; +import java.util.LinkedList; +import java.util.List; + +import de.lmu.ifi.dbs.elki.data.NumberVector; +import de.lmu.ifi.dbs.elki.database.ids.DBID; +import de.lmu.ifi.dbs.elki.database.ids.DBIDs; +import de.lmu.ifi.dbs.elki.database.query.DistanceResultPair; +import de.lmu.ifi.dbs.elki.database.query.DoubleDistanceResultPair; +import de.lmu.ifi.dbs.elki.database.relation.Relation; +import de.lmu.ifi.dbs.elki.distance.distancefunction.EuclideanDistanceFunction; +import de.lmu.ifi.dbs.elki.distance.distancevalue.NumberDistance; +import de.lmu.ifi.dbs.elki.math.linearalgebra.EigenPair; +import de.lmu.ifi.dbs.elki.math.linearalgebra.EigenvalueDecomposition; +import de.lmu.ifi.dbs.elki.math.linearalgebra.Matrix; +import de.lmu.ifi.dbs.elki.math.linearalgebra.SortedEigenPairs; +import de.lmu.ifi.dbs.elki.utilities.DatabaseUtil; +import de.lmu.ifi.dbs.elki.utilities.documentation.Reference; + +/** + * Performs a self-tuning local PCA based on the covariance matrices of given + * objects. At most the closest 'k' points are used in the calculation and a + * weight function is applied. + * + * The number of points actually used depends on when the strong eigenvectors + * exhibit the clearest correlation. + * + * @author Erich Schubert + * @param <V> vector type + */ +@Reference(authors = "H.-P. Kriegel, P. Kröger, E. Schubert, A. Zimek", title = "A General Framework for Increasing the Robustness of PCA-based Correlation Clustering Algorithms", booktitle = "Proceedings of the 20th International Conference on Scientific and Statistical Database Management (SSDBM), Hong Kong, China, 2008", url = "http://dx.doi.org/10.1007/978-3-540-69497-7_27") +public class PCAFilteredAutotuningRunner<V extends NumberVector<? extends V, ?>> extends PCAFilteredRunner<V> { + /** + * Constructor. + * + * @param covarianceMatrixBuilder + * @param eigenPairFilter + * @param big + * @param small + */ + public PCAFilteredAutotuningRunner(CovarianceMatrixBuilder<V> covarianceMatrixBuilder, EigenPairFilter eigenPairFilter, double big, double small) { + super(covarianceMatrixBuilder, eigenPairFilter, big, small); + } + + @Override + public PCAFilteredResult processIds(DBIDs ids, Relation<? extends V> database) { + // Assume Euclidean distance. In the context of PCA, the neighborhood should + // be L2-spherical to be unbiased. + V center = DatabaseUtil.centroid(database, ids); + List<DoubleDistanceResultPair> dres = new ArrayList<DoubleDistanceResultPair>(ids.size()); + for(DBID id : ids) { + final double dist = EuclideanDistanceFunction.STATIC.doubleDistance(center, database.get(id)); + dres.add(new DoubleDistanceResultPair(dist, id)); + } + Collections.sort(dres); + return processQueryResult(dres, database); + } + + @Override + public <D extends NumberDistance<?, ?>> PCAFilteredResult processQueryResult(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database) { + assertSortedByDistance(results); + final int dim = DatabaseUtil.dimensionality(database); + + List<Matrix> best = new LinkedList<Matrix>(); + for(int i = 0; i < dim; i++) { + best.add(null); + } + double[] beststrength = new double[dim]; + for(int i = 0; i < dim; i++) { + beststrength[i] = -1; + } + int[] bestk = new int[dim]; + // 'history' + LinkedList<Matrix> prevM = new LinkedList<Matrix>(); + LinkedList<Double> prevS = new LinkedList<Double>(); + LinkedList<Integer> prevD = new LinkedList<Integer>(); + // TODO: starting parameter shouldn't be hardcoded... + int smooth = 3; + int startk = 4; + if(startk > results.size() - 1) { + startk = results.size() - 1; + } + // TODO: add smoothing options, handle border cases better. + for(int k = startk; k < results.size(); k++) { + // sorted eigenpairs, eigenvectors, eigenvalues + Matrix covMat = covarianceMatrixBuilder.processQueryResults(results, database); + EigenvalueDecomposition evd = new EigenvalueDecomposition(covMat); + SortedEigenPairs eigenPairs = new SortedEigenPairs(evd, false); + FilteredEigenPairs filteredEigenPairs = getEigenPairFilter().filter(eigenPairs); + + // correlationDimension = #strong EV + int thisdim = filteredEigenPairs.countStrongEigenPairs(); + + // FIXME: handle the case of no strong EVs. + assert ((thisdim > 0) && (thisdim <= dim)); + double thisexplain = computeExplainedVariance(filteredEigenPairs); + + prevM.add(covMat); + prevS.add(thisexplain); + prevD.add(thisdim); + assert (prevS.size() == prevM.size()); + assert (prevS.size() == prevD.size()); + + if(prevS.size() >= 2 * smooth + 1) { + // all the same dimension? + boolean samedim = true; + for(Iterator<Integer> it = prevD.iterator(); it.hasNext();) { + if(it.next().intValue() != thisdim) { + samedim = false; + } + } + if(samedim) { + // average their explain values + double avgexplain = 0.0; + for(Iterator<Double> it = prevS.iterator(); it.hasNext();) { + avgexplain += it.next().doubleValue(); + } + avgexplain /= prevS.size(); + + if(avgexplain > beststrength[thisdim - 1]) { + beststrength[thisdim - 1] = avgexplain; + best.set(thisdim - 1, prevM.get(smooth)); + bestk[thisdim - 1] = k - smooth; + } + } + prevM.removeFirst(); + prevS.removeFirst(); + prevD.removeFirst(); + assert (prevS.size() == prevM.size()); + assert (prevS.size() == prevD.size()); + } + } + // Try all dimensions, lowest first. + for(int i = 0; i < dim; i++) { + if(beststrength[i] > 0.0) { + // If the best was the lowest or the biggest k, skip it! + if(bestk[i] == startk + smooth) { + continue; + } + if(bestk[i] == results.size() - smooth - 1) { + continue; + } + Matrix covMat = best.get(i); + + // We stop at the lowest dimension that did the job for us. + // System.err.println("Auto-k: "+bestk[i]+" dim: "+(i+1)); + return processCovarMatrix(covMat); + } + } + // NOTE: if we didn't get a 'maximum' anywhere, we end up with the data from + // the last run of the loop above. I.e. PCA on the full data set. That is + // intended. + return processCovarMatrix(covarianceMatrixBuilder.processQueryResults(results, database)); + } + + /** + * Compute the explained variance for a FilteredEigenPairs + * + * @param filteredEigenPairs + * @return explained variance by the strong eigenvectors. + */ + private double computeExplainedVariance(FilteredEigenPairs filteredEigenPairs) { + double strongsum = 0.0; + double weaksum = 0.0; + for(EigenPair ep : filteredEigenPairs.getStrongEigenPairs()) { + strongsum += ep.getEigenvalue(); + } + for(EigenPair ep : filteredEigenPairs.getWeakEigenPairs()) { + weaksum += ep.getEigenvalue(); + } + return strongsum / (strongsum / weaksum); + } + + /** + * Ensure that the results are sorted by distance. + * + * @param results + */ + private <D extends NumberDistance<?, ?>> void assertSortedByDistance(Collection<? extends DistanceResultPair<D>> results) { + // TODO: sort results instead? + double dist = -1.0; + for(Iterator<? extends DistanceResultPair<D>> it = results.iterator(); it.hasNext();) { + double qr = it.next().getDistance().doubleValue(); + if(qr < dist) { + System.err.println("WARNING: results not sorted by distance!"); + } + dist = qr; + } + } + + /** + * Parameterization class. + * + * @author Erich Schubert + * + * @apiviz.exclude + */ + public static class Parameterizer<V extends NumberVector<? extends V, ?>> extends PCAFilteredRunner.Parameterizer<V> { + @Override + protected PCAFilteredAutotuningRunner<V> makeInstance() { + return new PCAFilteredAutotuningRunner<V>(covarianceMatrixBuilder, eigenPairFilter, big, small); + } + } +}
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredResult.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredResult.java index 30453f0c..4aa626a9 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredResult.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredResult.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -118,7 +118,7 @@ public class PCAFilteredResult extends PCAResult implements ProjectionResult { for(Iterator<EigenPair> it = strongEigenPairs.iterator(); it.hasNext(); i++) { EigenPair eigenPair = it.next(); strongEigenvalues[i] = eigenPair.getEigenvalue(); - strongEigenvectors.setColumnVector(i, eigenPair.getEigenvector()); + strongEigenvectors.setCol(i, eigenPair.getEigenvector()); sumStrongEigenvalues += strongEigenvalues[i]; } } @@ -131,7 +131,7 @@ public class PCAFilteredResult extends PCAResult implements ProjectionResult { for(Iterator<EigenPair> it = weakEigenPairs.iterator(); it.hasNext(); i++) { EigenPair eigenPair = it.next(); weakEigenvalues[i] = eigenPair.getEigenvalue(); - weakEigenvectors.setColumnVector(i, eigenPair.getEigenvector()); + weakEigenvectors.setCol(i, eigenPair.getEigenvector()); sumWeakEigenvalues += weakEigenvalues[i]; } } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredRunner.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredRunner.java index 799502fa..2391446d 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredRunner.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAFilteredRunner.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -141,7 +141,7 @@ public class PCAFilteredRunner<V extends NumberVector<? extends V, ?>> extends P * @return PCA result */ @Override - public <D extends NumberDistance<?, ?>> PCAFilteredResult processQueryResult(Collection<DistanceResultPair<D>> results, Relation<? extends V> database) { + public <D extends NumberDistance<?, ?>> PCAFilteredResult processQueryResult(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database) { return processCovarMatrix(covarianceMatrixBuilder.processQueryResults(results, database)); } @@ -153,7 +153,7 @@ public class PCAFilteredRunner<V extends NumberVector<? extends V, ?>> extends P @Override public PCAFilteredResult processCovarMatrix(Matrix covarMatrix) { // TODO: add support for a different implementation to do EVD? - EigenvalueDecomposition evd = covarMatrix.eig(); + EigenvalueDecomposition evd = new EigenvalueDecomposition(covarMatrix); return processEVD(evd); } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAResult.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAResult.java index 26497645..6969a3a3 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAResult.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCAResult.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCARunner.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCARunner.java index 21b4304f..661fa5c5 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCARunner.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PCARunner.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -113,7 +113,7 @@ public class PCARunner<V extends NumberVector<? extends V, ?>> implements Parame * @param database the database used * @return PCA result */ - public <D extends NumberDistance<?, ?>> PCAResult processQueryResult(Collection<DistanceResultPair<D>> results, Relation<? extends V> database) { + public <D extends NumberDistance<?, ?>> PCAResult processQueryResult(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database) { return processCovarMatrix(covarianceMatrixBuilder.processQueryResults(results, database)); } @@ -125,7 +125,7 @@ public class PCARunner<V extends NumberVector<? extends V, ?>> implements Parame */ public PCAResult processCovarMatrix(Matrix covarMatrix) { // TODO: add support for a different implementation to do EVD? - EigenvalueDecomposition evd = covarMatrix.eig(); + EigenvalueDecomposition evd = new EigenvalueDecomposition(covarMatrix); return processEVD(evd); } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PercentageEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PercentageEigenPairFilter.java index b0065eae..321c12cc 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PercentageEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/PercentageEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/ProgressiveEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/ProgressiveEigenPairFilter.java index f43ac424..f66a1e96 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/ProgressiveEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/ProgressiveEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/RelativeEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/RelativeEigenPairFilter.java index 5109b010..59b2b750 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/RelativeEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/RelativeEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/SignificantEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/SignificantEigenPairFilter.java index d09b7312..94636553 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/SignificantEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/SignificantEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/StandardCovarianceMatrixBuilder.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/StandardCovarianceMatrixBuilder.java index 393b26ec..2c88d490 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/StandardCovarianceMatrixBuilder.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/StandardCovarianceMatrixBuilder.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeakEigenPairFilter.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeakEigenPairFilter.java index 57d21ec1..fbca039d 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeakEigenPairFilter.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeakEigenPairFilter.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeightedCovarianceMatrixBuilder.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeightedCovarianceMatrixBuilder.java index 547737fe..db5e8702 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeightedCovarianceMatrixBuilder.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/WeightedCovarianceMatrixBuilder.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -158,7 +158,7 @@ public class WeightedCovarianceMatrixBuilder<V extends NumberVector<? extends V, * @return Covariance Matrix */ @Override - public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<DistanceResultPair<D>> results, Relation<? extends V> database, int k) { + public <D extends NumberDistance<?, ?>> Matrix processQueryResults(Collection<? extends DistanceResultPair<D>> results, Relation<? extends V> database, int k) { final int dim = DatabaseUtil.dimensionality(database); final CovarianceMatrix cmat = new CovarianceMatrix(dim); @@ -172,7 +172,7 @@ public class WeightedCovarianceMatrixBuilder<V extends NumberVector<? extends V, double stddev = 0.0; { int i = 0; - for(Iterator<DistanceResultPair<D>> it = results.iterator(); it.hasNext() && i < k; i++) { + for(Iterator<? extends DistanceResultPair<D>> it = results.iterator(); it.hasNext() && i < k; i++) { DistanceResultPair<D> res = it.next(); final double dist; if(res instanceof DoubleDistanceResultPair) { @@ -194,7 +194,7 @@ public class WeightedCovarianceMatrixBuilder<V extends NumberVector<? extends V, // calculate weighted PCA int i = 0; - for(Iterator<DistanceResultPair<D>> it = results.iterator(); it.hasNext() && i < k; i++) { + for(Iterator<? extends DistanceResultPair<D>> it = results.iterator(); it.hasNext() && i < k; i++) { DistanceResultPair<? extends NumberDistance<?, ?>> res = it.next(); final double dist; if(res instanceof DoubleDistanceResultPair) { diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/package-info.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/package-info.java index 94b0f5a5..f08a016c 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/package-info.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/package-info.java @@ -5,7 +5,7 @@ This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures -Copyright (C) 2011 +Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ConstantWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ConstantWeight.java index b31d15ea..ef089980 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ConstantWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ConstantWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcStddevWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcStddevWeight.java index 71096f05..656406c3 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcStddevWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcStddevWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -24,6 +24,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; */ import de.lmu.ifi.dbs.elki.math.MathUtil; +import de.lmu.ifi.dbs.elki.math.statistics.distribution.NormalDistribution; /** * Gaussian Error Function Weight function, scaled using stddev. This probably @@ -42,6 +43,6 @@ public final class ErfcStddevWeight implements WeightFunction { if(stddev <= 0) { return 1; } - return MathUtil.erfc(MathUtil.SQRTHALF * distance / stddev); + return NormalDistribution.erfc(MathUtil.SQRTHALF * distance / stddev); } -} +}
\ No newline at end of file diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcWeight.java index 1a6cb5b6..6a8783ed 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ErfcWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team @@ -23,7 +23,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; along with this program. If not, see <http://www.gnu.org/licenses/>. */ -import de.lmu.ifi.dbs.elki.math.MathUtil; +import de.lmu.ifi.dbs.elki.math.statistics.distribution.NormalDistribution; /** * Gaussian Error Function Weight function, scaled such that the result it 0.1 @@ -48,6 +48,6 @@ public final class ErfcWeight implements WeightFunction { double relativedistance = distance / max; // the scaling was picked such that getWeight(a,a,0) is 0.1 // since erfc(1.1630871536766736) == 1.0 - return MathUtil.erfc(1.1630871536766736 * relativedistance); + return NormalDistribution.erfc(1.1630871536766736 * relativedistance); } } diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialStddevWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialStddevWeight.java index 979ccaf5..e6cb9951 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialStddevWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialStddevWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialWeight.java index 96077581..0792e764 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/ExponentialWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussStddevWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussStddevWeight.java index 33ff90c5..f73d4a7d 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussStddevWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussStddevWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussWeight.java index 00eeb8b1..05c3389f 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/GaussWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseLinearWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseLinearWeight.java index f618c08b..d5843cf1 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseLinearWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseLinearWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalStddevWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalStddevWeight.java index da1c0a61..740ec8ca 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalStddevWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalStddevWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalWeight.java index ee17b25b..c4ae9a97 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/InverseProportionalWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/LinearWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/LinearWeight.java index 241c2673..38318133 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/LinearWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/LinearWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticStddevWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticStddevWeight.java index 7d9b41ad..372dc016 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticStddevWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticStddevWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticWeight.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticWeight.java index 9755c963..49cb1aca 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticWeight.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/QuadraticWeight.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/WeightFunction.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/WeightFunction.java index c3102350..79c9f8a5 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/WeightFunction.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/WeightFunction.java @@ -4,7 +4,7 @@ package de.lmu.ifi.dbs.elki.math.linearalgebra.pca.weightfunctions; This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures - Copyright (C) 2011 + Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team diff --git a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/package-info.java b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/package-info.java index 268a3e89..649a8fb7 100644 --- a/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/package-info.java +++ b/src/de/lmu/ifi/dbs/elki/math/linearalgebra/pca/weightfunctions/package-info.java @@ -5,7 +5,7 @@ This file is part of ELKI: Environment for Developing KDD-Applications Supported by Index-Structures -Copyright (C) 2011 +Copyright (C) 2012 Ludwig-Maximilians-Universität München Lehr- und Forschungseinheit für Datenbanksysteme ELKI Development Team |