commons-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From er...@apache.org
Subject svn commit: r1370782 - in /commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution: AbstractMultivariateRealDistribution.java MultivariateNormalDistribution.java MultivariateRealDistribution.java
Date Wed, 08 Aug 2012 14:18:16 GMT
Author: erans
Date: Wed Aug  8 14:18:16 2012
New Revision: 1370782

URL: http://svn.apache.org/viewvc?rev=1370782&view=rev
Log:
MATH-815
Initial commit. Units test are yet to be added. This version also
contains a bug. Code contributed by Jared Becksfort, included with
modifications. [Cf. comments on JIRA.]

Added:
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/AbstractMultivariateRealDistribution.java
  (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateNormalDistribution.java
  (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateRealDistribution.java
  (with props)

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/AbstractMultivariateRealDistribution.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/AbstractMultivariateRealDistribution.java?rev=1370782&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/AbstractMultivariateRealDistribution.java
(added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/AbstractMultivariateRealDistribution.java
Wed Aug  8 14:18:16 2012
@@ -0,0 +1,104 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.math3.distribution;
+
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.exception.util.LocalizedFormats;
+import org.apache.commons.math3.random.RandomGenerator;
+
+/**
+ * Base class for probability distributions on the multivariate reals.
+ * Default implementations are provided for some of the methods that do
+ * not vary from distribution to distribution.
+ * 
+ */
+public abstract class AbstractMultivariateRealDistribution
+    implements MultivariateRealDistribution {
+    /** The number of dimensions or columns in the multivariate distribution. */
+    private final int numDimensions;
+    /** RNG instance used to generate samples from the distribution. */
+    protected final RandomGenerator random;
+
+    /**
+     * @param rng Random number generator.
+     * @param n Number of dimensions.
+     */
+    protected AbstractMultivariateRealDistribution(RandomGenerator rng,
+                                                   int n) {
+        random = rng;
+        numDimensions = n;
+    }
+
+    /** {@inheritDoc} */
+    public void reseedRandomGenerator(long seed) {
+        random.setSeed(seed);
+    }
+
+    /**
+     * 
+     * @return the number of dimensions in the multivariate distribution .
+     */
+    public int getDimensions() {
+        return numDimensions;
+    }
+
+    /** {@inheritDoc} */
+    public abstract double[] sample();
+
+    /** {@inheritDoc} */
+    public double[][] sample(final int sampleSize) {
+        if (sampleSize <= 0) {
+            throw new NotStrictlyPositiveException(LocalizedFormats.NUMBER_OF_SAMPLES,
+                                                   sampleSize);
+        }
+        final double[][] out = new double[sampleSize][numDimensions];
+        for (int i = 0; i < sampleSize; i++) {
+            out[i] = sample();
+        }
+        return out;
+    }
+
+    /** {@inheritDoc} */
+    public double probability(double[] x) {
+        return 0;
+    }
+
+    /** {@inheritDoc} */
+    public double getSupportLowerBound() {
+        return Double.NEGATIVE_INFINITY;
+    }
+
+    /** {@inheritDoc} */
+    public double getSupportUpperBound() {
+        return Double.POSITIVE_INFINITY;
+    }
+
+    /** {@inheritDoc} */
+    public boolean isSupportLowerBoundInclusive() {
+        return false;
+    }
+
+    /** {@inheritDoc} */
+    public boolean isSupportUpperBoundInclusive() {
+        return false;
+    }
+
+    /** {@inheritDoc} */
+    public boolean isSupportConnected() {
+        return false;
+    }
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/AbstractMultivariateRealDistribution.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateNormalDistribution.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateNormalDistribution.java?rev=1370782&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateNormalDistribution.java
(added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateNormalDistribution.java
Wed Aug  8 14:18:16 2012
@@ -0,0 +1,211 @@
+package org.apache.commons.math3.distribution;
+
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.linear.Array2DRowRealMatrix;
+import org.apache.commons.math3.linear.EigenDecomposition;
+import org.apache.commons.math3.linear.NonPositiveDefiniteMatrixException;
+import org.apache.commons.math3.linear.RealMatrix;
+import org.apache.commons.math3.linear.SingularMatrixException;
+import org.apache.commons.math3.random.RandomGenerator;
+import org.apache.commons.math3.random.Well19937c;
+import org.apache.commons.math3.stat.correlation.Covariance;
+import org.apache.commons.math3.util.FastMath;
+import org.apache.commons.math3.util.MathArrays;
+
+/**
+ * Implementation of the multivariate normal (Gaussian) distribution.
+ * 
+ * @see <a href="http://en.wikipedia.org/wiki/Multivariate_normal_distribution">
+ * Multivariate normal distribution (Wikipedia)</a>
+ * @see <a href="http://mathworld.wolfram.com/MultivariateNormalDistribution.html">
+ * Multivariate normal distribution (MathWorld)</a>
+ */
+public class MultivariateNormalDistribution
+    extends AbstractMultivariateRealDistribution {
+    /** Vector of means. */
+    private final double[] means;
+    /** Covariance matrix. */
+    private final RealMatrix covarianceMatrix;
+    /** The matrix inverse of the covariance matrix. */
+    private final RealMatrix covarianceMatrixInverse;
+    /** The determinant of the covariance matrix. */
+    private final double covarianceMatrixDeterminant;
+    /** Matrix used in computation of samples. */
+    private final RealMatrix samplingMatrix;
+
+    /**
+     * Creates a multivariate normal distribution with the given mean vector and
+     * covariance matrix.
+     * <br/>
+     * The number of dimensions is equal to the length of the mean vector
+     * and to the number of rows and columns of the covariance matrix.
+     * It is frequently written as "p" in formulae.
+     * 
+     * @param means Vector of means.
+     * @param covariances Covariance matrix.
+     */
+    public MultivariateNormalDistribution(final double[] means,
+                                          final double[][] covariances)
+        throws SingularMatrixException,
+               DimensionMismatchException,
+               NonPositiveDefiniteMatrixException {
+        this(new Well19937c(), means, covariances);
+    }
+
+    /**
+     * Creates a multivariate normal distribution with the given mean vector and
+     * covariance matrix.
+     * <br/>
+     * The number of dimensions is equal to the length of the mean vector
+     * and to the number of rows and columns of the covariance matrix.
+     * It is frequently written as "p" in formulae.
+     * 
+     * @param rng Random Number Generator.
+     * @param means Vector of means.
+     * @param covariances Covariance matrix.
+     */
+    public MultivariateNormalDistribution(RandomGenerator rng,
+                                          final double[] means,
+                                          final double[][] covariances)
+            throws SingularMatrixException,
+                   DimensionMismatchException,
+                   NonPositiveDefiniteMatrixException {
+        super(rng, means.length);
+
+        final int dim = means.length;
+
+        if (covariances.length != dim) {
+            throw new DimensionMismatchException(covariances.length, dim);
+        }
+
+        for (int i = 0; i < dim; i++) {
+            if (dim != covariances[i].length) {
+                throw new DimensionMismatchException(covariances[i].length, dim);
+            }
+        }
+
+        this.means = MathArrays.copyOf(means);
+
+        covarianceMatrix = new Array2DRowRealMatrix(covariances);
+
+        // Covariance matrix eigen decomposition.
+        final EigenDecomposition covMatDec = new EigenDecomposition(covarianceMatrix);
+
+        // Compute and store the inverse.
+        covarianceMatrixInverse = covMatDec.getSolver().getInverse();
+        // Compute and store the determinant.
+        covarianceMatrixDeterminant = covMatDec.getDeterminant();
+
+        // Eigenvalues of the covariance matrix.
+        final double[] covMatEigenvalues = covMatDec.getRealEigenvalues();
+
+        for (int i = 0; i < covMatEigenvalues.length; i++) {
+            if (covMatEigenvalues[i] < 0) {
+                throw new NonPositiveDefiniteMatrixException(covMatEigenvalues[i], i, 0);
+            }
+        }
+
+        // Matrix where each column is an eigenvector of the covariance matrix.
+        final Array2DRowRealMatrix covMatEigenvectors = new Array2DRowRealMatrix(dim, dim);
+        for (int v = 0; v < dim; v++) {
+            final double[] evec = covMatDec.getEigenvector(v).toArray();
+            covMatEigenvectors.setColumn(v, evec);
+        }
+
+        final RealMatrix tmpMatrix = covMatEigenvectors.transpose();
+
+        // Scale each eigenvector by the square root of its eigenvalue.
+        for (int row = 0; row < dim; row++) {
+            final double factor = FastMath.sqrt(covMatEigenvalues[row]);
+            for (int col = 0; col < dim; col++) {
+                tmpMatrix.multiplyEntry(row, col, factor);
+            }
+        }
+
+        samplingMatrix = covMatEigenvectors.multiply(tmpMatrix);
+    }
+
+    /**
+     * Gets the mean vector.
+     * 
+     * @return the mean vector.
+     */
+    public double[] getMeans() {
+        return MathArrays.copyOf(means);
+    }
+
+    /**
+     * Gets the covariance matrix.
+     * 
+     * @return the covariance matrix.
+     */
+    public RealMatrix getCovariances() {
+        return covarianceMatrix.copy();
+    }
+
+    /** {@inheritDoc} */
+    public double density(final double[] vals) throws DimensionMismatchException {
+        final int dim = getDimensions();
+        if (vals.length != dim) {
+            throw new DimensionMismatchException(vals.length, dim);
+        }
+
+        final double kernel = getKernel(vals);
+
+        return FastMath.pow(2 * FastMath.PI, -dim / 2) *
+            FastMath.pow(covarianceMatrixDeterminant, -0.5) *
+            FastMath.exp(kernel);
+    }
+
+    /**
+     * Gets the square root of each element on the diagonal of the covariance
+     * matrix.
+     * 
+     * @return the standard deviations.
+     */
+    public double[] getStandardDeviations() {
+        final int dim = getDimensions();
+        final double[] std = new double[dim];
+        final double[][] s = covarianceMatrix.getData();
+        for (int i = 0; i < dim; i++) {
+            std[i] = FastMath.sqrt(s[i][i]);
+        }
+        return std;
+    }
+
+    /** {@inheritDoc} */
+    public double[] sample() {
+        final int dim = getDimensions();
+        final double[] normalVals = new double[dim];
+
+        for (int i = 0; i < dim; i++) {
+            normalVals[i] = random.nextGaussian();
+        }
+
+        final double[] vals = samplingMatrix.operate(normalVals);
+
+        for (int i = 0; i < dim; i++) {
+            vals[i] += means[i];
+        }
+
+        return vals;
+    }
+
+    /**
+     * Precomputes some of the multiplications used for determining densities.
+     * 
+     * @param values Values at which to compute density.
+     * @return the multiplication factor of density calculations.
+     */
+    private double getKernel(final double[] values) {
+        double k = 0;
+        for (int col = 0; col < values.length; col++) {
+            for (int v = 0; v < values.length; v++) {
+                k += covarianceMatrixInverse.getEntry(v, col)
+                    * FastMath.pow(values[v] - means[v], 2);
+            }
+        }
+        return -0.5 * k;
+    }
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateNormalDistribution.java
------------------------------------------------------------------------------
    svn:eol-style = native

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateRealDistribution.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateRealDistribution.java?rev=1370782&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateRealDistribution.java
(added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateRealDistribution.java
Wed Aug  8 14:18:16 2012
@@ -0,0 +1,130 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.math3.distribution;
+
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.exception.NumberIsTooLargeException;
+import org.apache.commons.math3.exception.OutOfRangeException;
+
+/**
+ * Base interface for multivariate distributions on the reals.
+ * 
+ * This is based largely on the RealDistribution interface, but cumulative
+ * distribution functions are not required because they are often quite
+ * difficult to compute for multivariate distributions.
+ */
+public interface MultivariateRealDistribution {
+    /**
+     * For a random variable {@code X} whose values are distributed according to
+     * this distribution, this method returns {@code P(X = x)}. In other words,
+     * this method represents the probability mass function (PMF) for the
+     * distribution.
+     * 
+     * @param x Point at which the PMF is evaluated.
+     * @return the value of the probability mass function at point {@code x}.
+     */
+    double probability(double[] x);
+
+    /**
+     * Returns the probability density function (PDF) of this distribution
+     * evaluated at the specified point {@code x}. In general, the PDF is the
+     * derivative of the cumulative distribution function. If the derivative
+     * does not exist at {@code x}, then an appropriate replacement should be
+     * returned, e.g. {@code Double.POSITIVE_INFINITY}, {@code Double.NaN}, or
+     * the limit inferior or limit superior of the difference quotient.
+     * 
+     * @param x Point at which the PDF is evaluated.
+     * @return the value of the probability density function at point {@code x}.
+     */
+    double density(double[] x) throws DimensionMismatchException;
+
+    /**
+     * Access the lower bound of the support.
+     * This method must return the same value as {@code inverseCumulativeProbability(0)}.
+     * In other words, this method must return
+     * <p>
+     * <code>inf {x in R | P(X <= x) > 0}</code>.
+     * </p>
+     * 
+     * @return the lower bound of the support (might be
+     * {@code Double.NEGATIVE_INFINITY}).
+     */
+    double getSupportLowerBound();
+
+    /**
+     * Access the upper bound of the support.
+     * This method must return the same value as {@code inverseCumulativeProbability(1)}.
+     * In other words, this method must return
+     * <p>
+     * <code>inf {x in R | P(X <= x) = 1}</code>.
+     * </p>
+     * 
+     * @return the upper bound of the support (might be
+     * {@code Double.POSITIVE_INFINITY}).
+     */
+    double getSupportUpperBound();
+
+    /**
+     * Gets information about whether the lower bound of the support is
+     * inclusive or not.
+     * 
+     * @return whether the lower bound of the support is inclusive or not.
+     */
+    boolean isSupportLowerBoundInclusive();
+
+    /**
+     * gets information about whether the upper bound of the support is
+     * inclusive or not.
+     * 
+     * @return whether the upper bound of the support is inclusive or not.
+     */
+    boolean isSupportUpperBoundInclusive();
+
+    /**
+     * Gets information about whether the support is connected (i.e. all
+     * values between the lower and upper bound of the support are included
+     * in the support).
+     * 
+     * @return whether the support is connected or not.
+     */
+    boolean isSupportConnected();
+
+    /**
+     * Reseeds the random generator used to generate samples.
+     * 
+     * @param seed Seed with which to initialize the random number generator.
+     */
+    void reseedRandomGenerator(long seed);
+
+    /**
+     * Generates a random value vector sampled from this distribution.
+     * 
+     * @return a random value vector.
+     */
+    double[] sample();
+
+    /**
+     * Generates a list of a random value vectors from the distribution.
+     * 
+     * @param sampleSize the number of random vectors to generate.
+     * @return an array representing the random samples.
+     * @throws org.apache.commons.math3.exception.NotStrictlyPositiveException
+     * if {@code sampleSize} is not positive.
+     */
+    double[][] sample(int sampleSize) throws NotStrictlyPositiveException;
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/distribution/MultivariateRealDistribution.java
------------------------------------------------------------------------------
    svn:eol-style = native



Mime
View raw message