Java Code Examples for org.apache.commons.math.linear.RealMatrix#multiply()

public void testUnderdetermined() {
    final Random r    = new Random(42185006424567123l);
    int          p    = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
    int          q    = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
    RealMatrix   a    = createTestMatrix(r, p, q);
    RealMatrix   xRef = createTestMatrix(r, q, BlockRealMatrix.BLOCK_SIZE + 3);
    RealMatrix   b    = a.multiply(xRef);
    RealMatrix   x = new QRDecompositionImpl(a).getSolver().solve(b);

    // too many equations, the system cannot be solved at all
    assertTrue(x.subtract(xRef).getNorm() / (p * q) > 0.01);

    // the last unknown should have been set to 0
    assertEquals(0.0, x.getSubMatrix(p, q - 1, 0, x.getColumnDimension() - 1).getNorm());

}
 

public void testUnderdetermined() {
    final Random r    = new Random(42185006424567123l);
    int          p    = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
    int          q    = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
    RealMatrix   a    = createTestMatrix(r, p, q);
    RealMatrix   xRef = createTestMatrix(r, q, BlockRealMatrix.BLOCK_SIZE + 3);
    RealMatrix   b    = a.multiply(xRef);
    RealMatrix   x = new QRDecompositionImpl(a).getSolver().solve(b);

    // too many equations, the system cannot be solved at all
    assertTrue(x.subtract(xRef).getNorm() / (p * q) > 0.01);

    // the last unknown should have been set to 0
    assertEquals(0.0, x.getSubMatrix(p, q - 1, 0, x.getColumnDimension() - 1).getNorm());

}
 

public void testOverdetermined() {
    final Random r    = new Random(5559252868205245l);
    int          p    = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
    int          q    = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
    RealMatrix   a    = createTestMatrix(r, p, q);
    RealMatrix   xRef = createTestMatrix(r, q, BlockRealMatrix.BLOCK_SIZE + 3);

    // build a perturbed system: A.X + noise = B
    RealMatrix b = a.multiply(xRef);
    final double noise = 0.001;
    b.walkInOptimizedOrder(new DefaultRealMatrixChangingVisitor() {
        @Override
        public double visit(int row, int column, double value) {
            return value * (1.0 + noise * (2 * r.nextDouble() - 1));
        }
    });

    // despite perturbation, the least square solution should be pretty good
    RealMatrix x = new QRDecompositionImpl(a).getSolver().solve(b);
    assertEquals(0, x.subtract(xRef).getNorm(), 0.01 * noise * p * q);

}
 

public void testOverdetermined() {
    final Random r    = new Random(5559252868205245l);
    int          p    = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
    int          q    = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
    RealMatrix   a    = createTestMatrix(r, p, q);
    RealMatrix   xRef = createTestMatrix(r, q, BlockRealMatrix.BLOCK_SIZE + 3);

    // build a perturbed system: A.X + noise = B
    RealMatrix b = a.multiply(xRef);
    final double noise = 0.001;
    b.walkInOptimizedOrder(new DefaultRealMatrixChangingVisitor() {
        @Override
        public double visit(int row, int column, double value) {
            return value * (1.0 + noise * (2 * r.nextDouble() - 1));
        }
    });

    // despite perturbation, the least square solution should be pretty good
    RealMatrix x = new QRDecompositionImpl(a).getSolver().solve(b);
    assertEquals(0, x.subtract(xRef).getNorm(), 0.01 * noise * p * q);

}
 

public void testOverdetermined() {
    final Random r    = new Random(5559252868205245l);
    int          p    = (7 * BlockRealMatrix.BLOCK_SIZE) / 4;
    int          q    = (5 * BlockRealMatrix.BLOCK_SIZE) / 4;
    RealMatrix   a    = createTestMatrix(r, p, q);
    RealMatrix   xRef = createTestMatrix(r, q, BlockRealMatrix.BLOCK_SIZE + 3);

    // build a perturbed system: A.X + noise = B
    RealMatrix b = a.multiply(xRef);
    final double noise = 0.001;
    b.walkInOptimizedOrder(new DefaultRealMatrixChangingVisitor() {
        @Override
        public double visit(int row, int column, double value) {
            return value * (1.0 + noise * (2 * r.nextDouble() - 1));
        }
    });

    // despite perturbation, the least square solution should be pretty good
    RealMatrix x = new QRDecompositionImpl(a).getSolver().solve(b);
    assertEquals(0, x.subtract(xRef).getNorm(), 0.01 * noise * p * q);

}
 

@Test
public void testRootMatrix() {
    RealMatrix b = generator.getRootMatrix();
    RealMatrix bbt = b.multiply(b.transpose());
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            Assert.assertEquals(covariance.getEntry(i, j), bbt.getEntry(i, j), 1.0e-12);
        }
    }
}
 

public void testRootMatrix() {
    RealMatrix b = generator.getRootMatrix();
    RealMatrix bbt = b.multiply(b.transpose());
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            assertEquals(covariance.getEntry(i, j), bbt.getEntry(i, j), 1.0e-12);
        }
    }
}
 

public void testRootMatrix() {
    RealMatrix b = generator.getRootMatrix();
    RealMatrix bbt = b.multiply(b.transpose());
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            assertEquals(covariance.getEntry(i, j), bbt.getEntry(i, j), 1.0e-12);
        }
    }
}
 

public void testRootMatrix() {
    RealMatrix b = generator.getRootMatrix();
    RealMatrix bbt = b.multiply(b.transpose());
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            assertEquals(covariance.getEntry(i, j), bbt.getEntry(i, j), 1.0e-12);
        }
    }
}
 

public void testRootMatrix() {
    RealMatrix b = generator.getRootMatrix();
    RealMatrix bbt = b.multiply(b.transpose());
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            assertEquals(covariance.getEntry(i, j), bbt.getEntry(i, j), 1.0e-12);
        }
    }
}
 

public void testRootMatrix() {
    RealMatrix b = generator.getRootMatrix();
    RealMatrix bbt = b.multiply(b.transpose());
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            assertEquals(covariance.getEntry(i, j), bbt.getEntry(i, j), 1.0e-12);
        }
    }
}
 

public CorrelatedRandomVectorGeneratorTest() {
    mean = new double[] { 0.0, 1.0, -3.0, 2.3 };

    RealMatrix b = MatrixUtils.createRealMatrix(4, 3);
    int counter = 0;
    for (int i = 0; i < b.getRowDimension(); ++i) {
        for (int j = 0; j < b.getColumnDimension(); ++j) {
            b.setEntry(i, j, 1.0 + 0.1 * ++counter);
        }
    }
    RealMatrix bbt = b.multiply(b.transpose());
    covariance = MatrixUtils.createRealMatrix(mean.length, mean.length);
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        covariance.setEntry(i, i, bbt.getEntry(i, i));
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            double s = bbt.getEntry(i, j);
            covariance.setEntry(i, j, s);
            covariance.setEntry(j, i, s);
        }
    }

    RandomGenerator rg = new JDKRandomGenerator();
    rg.setSeed(17399225432l);
    GaussianRandomGenerator rawGenerator = new GaussianRandomGenerator(rg);
    generator = new CorrelatedRandomVectorGenerator(mean,
                                                    covariance,
                                                    1.0e-12 * covariance.getNorm(),
                                                    rawGenerator);
}
 

public CorrelatedRandomVectorGeneratorTest() {
    mean = new double[] { 0.0, 1.0, -3.0, 2.3 };

    RealMatrix b = MatrixUtils.createRealMatrix(4, 3);
    int counter = 0;
    for (int i = 0; i < b.getRowDimension(); ++i) {
        for (int j = 0; j < b.getColumnDimension(); ++j) {
            b.setEntry(i, j, 1.0 + 0.1 * ++counter);
        }
    }
    RealMatrix bbt = b.multiply(b.transpose());
    covariance = MatrixUtils.createRealMatrix(mean.length, mean.length);
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        covariance.setEntry(i, i, bbt.getEntry(i, i));
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            double s = bbt.getEntry(i, j);
            covariance.setEntry(i, j, s);
            covariance.setEntry(j, i, s);
        }
    }

    RandomGenerator rg = new JDKRandomGenerator();
    rg.setSeed(17399225432l);
    GaussianRandomGenerator rawGenerator = new GaussianRandomGenerator(rg);
    generator = new CorrelatedRandomVectorGenerator(mean,
                                                    covariance,
                                                    1.0e-12 * covariance.getNorm(),
                                                    rawGenerator);
}
 

public CorrelatedRandomVectorGeneratorTest() {
    mean = new double[] { 0.0, 1.0, -3.0, 2.3 };

    RealMatrix b = MatrixUtils.createRealMatrix(4, 3);
    int counter = 0;
    for (int i = 0; i < b.getRowDimension(); ++i) {
        for (int j = 0; j < b.getColumnDimension(); ++j) {
            b.setEntry(i, j, 1.0 + 0.1 * ++counter);
        }
    }
    RealMatrix bbt = b.multiply(b.transpose());
    covariance = MatrixUtils.createRealMatrix(mean.length, mean.length);
    for (int i = 0; i < covariance.getRowDimension(); ++i) {
        covariance.setEntry(i, i, bbt.getEntry(i, i));
        for (int j = 0; j < covariance.getColumnDimension(); ++j) {
            double s = bbt.getEntry(i, j);
            covariance.setEntry(i, j, s);
            covariance.setEntry(j, i, s);
        }
    }

    RandomGenerator rg = new JDKRandomGenerator();
    rg.setSeed(17399225432l);
    GaussianRandomGenerator rawGenerator = new GaussianRandomGenerator(rg);
    generator = new CorrelatedRandomVectorGenerator(mean,
                                                    covariance,
                                                    1.0e-12 * covariance.getNorm(),
                                                    rawGenerator);
}
 

/**
 * <p>Calculates the variance on the beta by OLS.
 * </p>
 * <p>Var(b) = (X<sup>T</sup>X)<sup>-1</sup>
 * </p>
 * <p>Uses QR decomposition to reduce (X<sup>T</sup>X)<sup>-1</sup>
 * to (R<sup>T</sup>R)<sup>-1</sup>, with only the top p rows of
 * R included, where p = the length of the beta vector.</p> 
 * 
 * @return The beta variance
 */
@Override
protected RealMatrix calculateBetaVariance() {
    int p = X.getColumnDimension();
    RealMatrix Raug = qr.getR().getSubMatrix(0, p - 1 , 0, p - 1);
    RealMatrix Rinv = new LUDecompositionImpl(Raug).getSolver().getInverse();
    return Rinv.multiply(Rinv.transpose());
}
 

/**
 * <p>Calculates the variance-covariance matrix of the regression parameters.
 * </p>
 * <p>Var(b) = (X<sup>T</sup>X)<sup>-1</sup>
 * </p>
 * <p>Uses QR decomposition to reduce (X<sup>T</sup>X)<sup>-1</sup>
 * to (R<sup>T</sup>R)<sup>-1</sup>, with only the top p rows of
 * R included, where p = the length of the beta vector.</p>
 *
 * @return The beta variance-covariance matrix
 */
@Override
protected RealMatrix calculateBetaVariance() {
    int p = X.getColumnDimension();
    RealMatrix Raug = qr.getR().getSubMatrix(0, p - 1 , 0, p - 1);
    RealMatrix Rinv = new LUDecompositionImpl(Raug).getSolver().getInverse();
    return Rinv.multiply(Rinv.transpose());
}
 

/**
 * <p>Calculates the variance on the beta by OLS.
 * </p>
 * <p>Var(b) = (X<sup>T</sup>X)<sup>-1</sup>
 * </p>
 * <p>Uses QR decomposition to reduce (X<sup>T</sup>X)<sup>-1</sup>
 * to (R<sup>T</sup>R)<sup>-1</sup>, with only the top p rows of
 * R included, where p = the length of the beta vector.</p>
 *
 * @return The beta variance
 */
@Override
protected RealMatrix calculateBetaVariance() {
    int p = X.getColumnDimension();
    RealMatrix Raug = qr.getR().getSubMatrix(0, p - 1 , 0, p - 1);
    RealMatrix Rinv = new LUDecompositionImpl(Raug).getSolver().getInverse();
    return Rinv.multiply(Rinv.transpose());
}
 

/**
 * <p>Calculates the variance on the beta by OLS.
 * </p>
 * <p>Var(b) = (X<sup>T</sup>X)<sup>-1</sup>
 * </p>
 * <p>Uses QR decomposition to reduce (X<sup>T</sup>X)<sup>-1</sup>
 * to (R<sup>T</sup>R)<sup>-1</sup>, with only the top p rows of
 * R included, where p = the length of the beta vector.</p>
 *
 * @return The beta variance
 */
@Override
protected RealMatrix calculateBetaVariance() {
    int p = X.getColumnDimension();
    RealMatrix Raug = qr.getR().getSubMatrix(0, p - 1 , 0, p - 1);
    RealMatrix Rinv = new LUDecompositionImpl(Raug).getSolver().getInverse();
    return Rinv.multiply(Rinv.transpose());
}
 

/**
 * <p>Calculates the variance on the beta by OLS.
 * </p>
 * <p>Var(b) = (X<sup>T</sup>X)<sup>-1</sup>
 * </p>
 * <p>Uses QR decomposition to reduce (X<sup>T</sup>X)<sup>-1</sup>
 * to (R<sup>T</sup>R)<sup>-1</sup>, with only the top p rows of
 * R included, where p = the length of the beta vector.</p>
 *
 * @return The beta variance
 */
@Override
protected RealMatrix calculateBetaVariance() {
    int p = X.getColumnDimension();
    RealMatrix Raug = qr.getR().getSubMatrix(0, p - 1 , 0, p - 1);
    RealMatrix Rinv = new LUDecompositionImpl(Raug).getSolver().getInverse();
    return Rinv.multiply(Rinv.transpose());
}
 

/**
 * <p>Calculates the variance-covariance matrix of the regression parameters.
 * </p>
 * <p>Var(b) = (X<sup>T</sup>X)<sup>-1</sup>
 * </p>
 * <p>Uses QR decomposition to reduce (X<sup>T</sup>X)<sup>-1</sup>
 * to (R<sup>T</sup>R)<sup>-1</sup>, with only the top p rows of
 * R included, where p = the length of the beta vector.</p>
 *
 * @return The beta variance-covariance matrix
 */
@Override
protected RealMatrix calculateBetaVariance() {
    int p = X.getColumnDimension();
    RealMatrix Raug = qr.getR().getSubMatrix(0, p - 1 , 0, p - 1);
    RealMatrix Rinv = new LUDecompositionImpl(Raug).getSolver().getInverse();
    return Rinv.multiply(Rinv.transpose());
}