View Javadoc
1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      https://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  
18  /*
19   * This is not the original file distributed by the Apache Software Foundation
20   * It has been modified by the Hipparchus project
21   */
22  
23  package org.hipparchus.distribution.multivariate;
24  
25  import java.util.Random;
26  
27  import org.hipparchus.UnitTestUtils;
28  import org.hipparchus.distribution.continuous.NormalDistribution;
29  import org.hipparchus.exception.LocalizedCoreFormats;
30  import org.hipparchus.exception.MathIllegalArgumentException;
31  import org.hipparchus.linear.Array2DRowRealMatrix;
32  import org.hipparchus.linear.RealMatrix;
33  import org.hipparchus.random.Well19937c;
34  import org.hipparchus.util.Precision;
35  import org.junit.Assert;
36  import org.junit.Test;
37  
38  /**
39   * Test cases for {@link MultivariateNormalDistribution}.
40   */
41  public class MultivariateNormalDistributionTest {
42      /**
43       * Test the ability of the distribution to report its mean value parameter.
44       */
45      @Test
46      public void testGetMean() {
47          final double[] mu = { -1.5, 2 };
48          final double[][] sigma = { { 2, -1.1 },
49                                     { -1.1, 2 } };
50          final MultivariateNormalDistribution d = new MultivariateNormalDistribution(mu, sigma);
51  
52          final double[] m = d.getMeans();
53          for (int i = 0; i < m.length; i++) {
54              Assert.assertEquals(mu[i], m[i], 0);
55          }
56      }
57  
58      /**
59       * Test the ability of the distribution to report its covariance matrix parameter.
60       */
61      @Test
62      public void testGetCovarianceMatrix() {
63          final double[] mu = { -1.5, 2 };
64          final double[][] sigma = { { 2, -1.1 },
65                                     { -1.1, 2 } };
66          final MultivariateNormalDistribution d = new MultivariateNormalDistribution(mu, sigma);
67  
68          final RealMatrix s = d.getCovariances();
69          final int dim = d.getDimension();
70          for (int i = 0; i < dim; i++) {
71              for (int j = 0; j < dim; j++) {
72                  Assert.assertEquals(sigma[i][j], s.getEntry(i, j), 0);
73              }
74          }
75      }
76  
77      /**
78       * Test the accuracy of sampling from the distribution.
79       */
80      @Test
81      public void testSampling() {
82          final double[] mu = { -1.5, 2 };
83          final double[][] sigma = { { 2, -1.1 },
84                                     { -1.1, 2 } };
85          final MultivariateNormalDistribution d = new MultivariateNormalDistribution(mu, sigma);
86          d.reseedRandomGenerator(50);
87  
88          final int n = 500000;
89  
90          final double[][] samples = d.sample(n);
91          final int dim = d.getDimension();
92          final double[] sampleMeans = new double[dim];
93  
94          for (int i = 0; i < samples.length; i++) {
95              for (int j = 0; j < dim; j++) {
96                  sampleMeans[j] += samples[i][j];
97              }
98          }
99  
100         final double sampledValueTolerance = 1e-2;
101         for (int j = 0; j < dim; j++) {
102             sampleMeans[j] /= samples.length;
103             Assert.assertEquals(mu[j], sampleMeans[j], sampledValueTolerance);
104         }
105 
106         //final double[][] sampleSigma = new Covariance(samples).getCovarianceMatrix().getData();
107         final RealMatrix sampleSigma = UnitTestUtils.covarianceMatrix(new Array2DRowRealMatrix(samples));
108         for (int i = 0; i < dim; i++) {
109             for (int j = 0; j < dim; j++) {
110                 Assert.assertEquals(sigma[i][j], sampleSigma.getEntry(i, j), sampledValueTolerance);
111             }
112         }
113     }
114 
115     /**
116      * Test the accuracy of the distribution when calculating densities.
117      */
118     @Test
119     public void testDensities() {
120         final double[] mu = { -1.5, 2 };
121         final double[][] sigma = { { 2, -1.1 },
122                                    { -1.1, 2 } };
123         final MultivariateNormalDistribution d = new MultivariateNormalDistribution(mu, sigma);
124 
125         final double[][] testValues = { { -1.5, 2 },
126                                         { 4, 4 },
127                                         { 1.5, -2 },
128                                         { 0, 0 } };
129         final double[] densities = new double[testValues.length];
130         for (int i = 0; i < densities.length; i++) {
131             densities[i] = d.density(testValues[i]);
132         }
133 
134         // From dmvnorm function in R 2.15 CRAN package Mixtools v0.4.5
135         final double[] correctDensities = { 0.09528357207691344,
136                                             5.80932710124009e-09,
137                                             0.001387448895173267,
138                                             0.03309922090210541 };
139 
140         for (int i = 0; i < testValues.length; i++) {
141             Assert.assertEquals(correctDensities[i], densities[i], 1e-16);
142         }
143     }
144 
145     /**
146      * Test the accuracy of the distribution when calculating densities.
147      */
148     @Test
149     public void testUnivariateDistribution() {
150         final double[] mu = { -1.5 };
151         final double[][] sigma = { { 1 } };
152 
153         final MultivariateNormalDistribution multi = new MultivariateNormalDistribution(mu, sigma);
154 
155         final NormalDistribution uni = new NormalDistribution(mu[0], sigma[0][0]);
156         final Random rng = new Random();
157         final int numCases = 100;
158         final double tol = Math.ulp(1d);
159         for (int i = 0; i < numCases; i++) {
160             final double v = rng.nextDouble() * 10 - 5;
161             Assert.assertEquals(uni.density(v), multi.density(new double[] { v }), tol);
162         }
163     }
164 
165     /**
166      * Test getting/setting custom singularMatrixTolerance
167      */
168     @Test
169     public void testGetSingularMatrixTolerance() {
170         final double[] mu = { -1.5 };
171         final double[][] sigma = { { 1 } };
172 
173         final double tolerance1 = 1e-2;
174         final MultivariateNormalDistribution mvd1 = new MultivariateNormalDistribution(mu, sigma, tolerance1);
175         Assert.assertEquals(tolerance1, mvd1.getSingularMatrixCheckTolerance(), Precision.EPSILON);
176 
177         final double tolerance2 = 1e-3;
178         final MultivariateNormalDistribution mvd2 = new MultivariateNormalDistribution(mu, sigma, tolerance2);
179         Assert.assertEquals(tolerance2, mvd2.getSingularMatrixCheckTolerance(), Precision.EPSILON);
180     }
181 
182     @Test
183     public void testNotPositiveDefinite() {
184         try {
185             new MultivariateNormalDistribution(new Well19937c(0x543l), new double[2],
186                                                new double[][] { { -1.0, 0.0 }, { 0.0, -2.0 } });
187             Assert.fail("an exception should have been thrown");
188         } catch (MathIllegalArgumentException miae) {
189             Assert.assertEquals(LocalizedCoreFormats.NOT_POSITIVE_DEFINITE_MATRIX, miae.getSpecifier());
190         }
191     }
192 
193     @Test
194     public void testStd() {
195         MultivariateNormalDistribution d = new MultivariateNormalDistribution(new Well19937c(0x543l), new double[2],
196                                                                               new double[][] { { 4.0, 0.0 }, { 0.0, 9.0 } });
197         double[] s = d.getStandardDeviations();
198         Assert.assertEquals(2, s.length);
199         Assert.assertEquals(2.0, s[0], 1.0e-15);
200         Assert.assertEquals(3.0, s[1], 1.0e-15);
201     }
202 
203     @Test
204     public void testWrongDensity() {
205         try {
206             MultivariateNormalDistribution d = new MultivariateNormalDistribution(new Well19937c(0x543l), new double[2],
207                                                                                   new double[][] { { 4.0, 0.0 }, { 0.0, 4.0 } });
208             d.density(new double[3]);
209             Assert.fail("an exception should have been thrown");
210         } catch (MathIllegalArgumentException miae) {
211             Assert.assertEquals(LocalizedCoreFormats.DIMENSIONS_MISMATCH, miae.getSpecifier());
212         }
213     }
214 
215     @Test
216     public void testWrongArguments() {
217         checkWrongArguments(new double[3], new double[6][6]);
218         checkWrongArguments(new double[3], new double[3][6]);
219     }
220 
221     private void checkWrongArguments(double[] means, double[][] covariances) {
222         try {
223             new MultivariateNormalDistribution(new Well19937c(0x543l), means, covariances);
224             Assert.fail("an exception should have been thrown");
225         } catch (MathIllegalArgumentException miae) {
226             Assert.assertEquals(LocalizedCoreFormats.DIMENSIONS_MISMATCH, miae.getSpecifier());
227         }
228     }
229 }