Removed dependency on "GaussianRandomGenerator" (unit test).

This commit is contained in:
Gilles Sadowski 2021-05-30 01:10:18 +02:00
parent 665ae9cbad
commit 8a756d763d
3 changed files with 52 additions and 18 deletions

View File

@ -36,8 +36,8 @@ import org.apache.commons.math4.legacy.linear.RectangularCholeskyDecomposition;
* Multivariate Normal Distribution</a>. The approach using a Cholesky
* decomposition is quite usual in this case. However, it can be extended
* to other cases as long as the underlying random generator provides
* {@link NormalizedRandomGenerator normalized values} like {@link
* GaussianRandomGenerator} or {@link UniformRandomGenerator}.</p>
* {@link NormalizedRandomGenerator normalized values} like
* {@link UniformRandomGenerator}.</p>
* <p>Sometimes, the covariance matrix for a given simulation is not
* strictly positive definite. This means that the correlations are
* not all independent from each other. In this case, however, the non

View File

@ -23,7 +23,9 @@ import org.apache.commons.math4.legacy.TestUtils;
import org.apache.commons.math4.legacy.linear.Array2DRowRealMatrix;
import org.apache.commons.math4.legacy.linear.MatrixUtils;
import org.apache.commons.math4.legacy.linear.RealMatrix;
import org.apache.commons.rng.UniformRandomProvider;
import org.apache.commons.rng.simple.RandomSource;
import org.apache.commons.rng.sampling.distribution.ZigguratNormalizedGaussianSampler;
import org.apache.commons.math4.legacy.stat.correlation.StorelessCovariance;
import org.apache.commons.math4.legacy.stat.descriptive.moment.VectorialCovariance;
import org.apache.commons.math4.legacy.stat.descriptive.moment.VectorialMean;
@ -57,11 +59,10 @@ public class CorrelatedRandomVectorGeneratorTest {
}
}
GaussianRandomGenerator rawGenerator = new GaussianRandomGenerator(RandomSource.create(RandomSource.WELL_1024_A, 17399225432l));
generator = new CorrelatedRandomVectorGenerator(mean,
covariance,
1.0e-12 * covariance.getNorm(),
rawGenerator);
1e-12 * covariance.getNorm(),
gaussianRandom(RandomSource.create(RandomSource.WELL_1024_A)));
}
@Test
@ -79,9 +80,9 @@ public class CorrelatedRandomVectorGeneratorTest {
{ 6, 2, -1, 197 }
};
RealMatrix covRM = MatrixUtils.createRealMatrix(cov);
NormalizedRandomGenerator rg = new GaussianRandomGenerator(RandomSource.create(RandomSource.WELL_1024_A, 5322145245211l));
CorrelatedRandomVectorGenerator sg =
new CorrelatedRandomVectorGenerator(mean, covRM, 0.00001, rg);
new CorrelatedRandomVectorGenerator(mean, covRM, 0.00001,
gaussianRandom(RandomSource.create(RandomSource.WELL_1024_A)));
double[] min = new double[mean.length];
Arrays.fill(min, Double.POSITIVE_INFINITY);
@ -168,12 +169,11 @@ public class CorrelatedRandomVectorGeneratorTest {
private CorrelatedRandomVectorGenerator createSampler(double[][] cov) {
RealMatrix matrix = new Array2DRowRealMatrix(cov);
double small = 10e-12 * matrix.getNorm();
return new CorrelatedRandomVectorGenerator(
new double[cov.length],
matrix,
small,
new GaussianRandomGenerator(RandomSource.create(RandomSource.WELL_1024_A, 0x366a26b94e520f41l)));
double small = 1e-12 * matrix.getNorm();
return new CorrelatedRandomVectorGenerator(new double[cov.length],
matrix,
small,
gaussianRandom(RandomSource.create(RandomSource.WELL_1024_A)));
}
private void testSampler(final double[][] covMatrix, int samples, double epsilon) {
@ -189,4 +189,20 @@ public class CorrelatedRandomVectorGeneratorTest {
TestUtils.assertEquals(covMatrix[r], sampleCov[r], epsilon);
}
}
/**
* @param rng RNG.
* @return a N(0,1) sampler.
*/
private NormalizedRandomGenerator gaussianRandom(final UniformRandomProvider rng) {
final ZigguratNormalizedGaussianSampler n = new ZigguratNormalizedGaussianSampler(rng);
return new NormalizedRandomGenerator() {
/** {@inheritDoc} */
@Override
public double nextNormalizedDouble() {
return n.sample();
}
};
}
}

View File

@ -25,10 +25,11 @@ import org.apache.commons.math4.legacy.exception.NullArgumentException;
import org.apache.commons.math4.legacy.linear.MatrixUtils;
import org.apache.commons.math4.legacy.linear.RealMatrix;
import org.apache.commons.math4.legacy.linear.RealVector;
import org.apache.commons.math4.legacy.random.NormalizedRandomGenerator;
import org.apache.commons.math4.legacy.random.CorrelatedRandomVectorGenerator;
import org.apache.commons.math4.legacy.random.GaussianRandomGenerator;
import org.apache.commons.rng.UniformRandomProvider;
import org.apache.commons.rng.simple.RandomSource;
import org.apache.commons.rng.sampling.distribution.ZigguratNormalizedGaussianSampler;
import org.apache.commons.statistics.distribution.ContinuousDistribution;
import org.apache.commons.statistics.distribution.NormalDistribution;
import org.apache.commons.math4.legacy.stat.correlation.Covariance;
@ -220,7 +221,7 @@ public class GLSMultipleLinearRegressionTest extends MultipleLinearRegressionAbs
*/
@Test
public void testGLSEfficiency() {
final UniformRandomProvider rg = RandomSource.create(RandomSource.MT, 123456789L);
final UniformRandomProvider rg = RandomSource.create(RandomSource.MT);
final ContinuousDistribution.Sampler gauss = new NormalDistribution(0, 1).createSampler(rg);
// Assume model has 16 observations (will use Longley data). Start by generating
@ -245,10 +246,11 @@ public class GLSMultipleLinearRegressionTest extends MultipleLinearRegressionAbs
RealMatrix cov = (new Covariance(errorSeeds)).getCovarianceMatrix();
// Create a CorrelatedRandomVectorGenerator to use to generate correlated errors
GaussianRandomGenerator rawGenerator = new GaussianRandomGenerator(rg);
double[] errorMeans = new double[nObs]; // Counting on init to 0 here
CorrelatedRandomVectorGenerator gen = new CorrelatedRandomVectorGenerator(errorMeans, cov,
1.0e-12 * cov.getNorm(), rawGenerator);
CorrelatedRandomVectorGenerator gen
= new CorrelatedRandomVectorGenerator(errorMeans, cov,
1e-12 * cov.getNorm(),
gaussianRandom(rg));
// Now start generating models. Use Longley X matrix on LHS
// and Longley OLS beta vector as "true" beta. Generate
@ -297,4 +299,20 @@ public class GLSMultipleLinearRegressionTest extends MultipleLinearRegressionAbs
assert(olsBetaStats.getMean() > 1.5 * glsBetaStats.getMean());
assert(olsBetaStats.getStandardDeviation() > glsBetaStats.getStandardDeviation());
}
/**
* @param rng RNG.
* @return a N(0,1) sampler.
*/
private NormalizedRandomGenerator gaussianRandom(final UniformRandomProvider rng) {
final ZigguratNormalizedGaussianSampler n = new ZigguratNormalizedGaussianSampler(rng);
return new NormalizedRandomGenerator() {
/** {@inheritDoc} */
@Override
public double nextNormalizedDouble() {
return n.sample();
}
};
}
}