Java Code Examples for org.apache.commons.math3.linear.LUDecomposition#getDeterminant()

/**
 * Returns the log posterior PDF of a particular point xi, to belong to this
 * cluster.
 * 
 * @param xi    The point for which we want to estimate the PDF.
 * @return      The log posterior PDF
 */
@Override
public double posteriorLogPdf(Point xi) {
    RealVector x_mu = xi.data.subtract(mean);

    if(cache_covariance_determinant==null || cache_covariance_inverse==null) {
        LUDecomposition lud = new LUDecomposition(covariance);
        cache_covariance_determinant = lud.getDeterminant();
        cache_covariance_inverse = lud.getSolver().getInverse();
        lud =null;
    }
    Double determinant=cache_covariance_determinant;
    RealMatrix invCovariance=cache_covariance_inverse;

    double x_muInvSx_muT = (invCovariance.preMultiply(x_mu)).dotProduct(x_mu);

    double normConst = 1.0/( Math.pow(2*Math.PI, dimensionality/2.0) * Math.pow(determinant, 0.5) );


    //double pdf = Math.exp(-0.5 * x_muInvSx_muT)*normConst;
    double logPdf = -0.5 * x_muInvSx_muT + Math.log(normConst);
    return logPdf;
}
 

/**
 * pdf(x, x_hat) = exp(-0.5 * (x-x_hat) * inv(Σ) * (x-x_hat)T) / ( 2π^0.5d * det(Σ)^0.5)
 * 
 * @return value of probabilistic density function
 * @link https://en.wikipedia.org/wiki/Multivariate_normal_distribution#Density_function
 */
public static double pdf(@Nonnull final RealVector x, @Nonnull final RealVector x_hat,
        @Nonnull final RealMatrix sigma) {
    final int dim = x.getDimension();
    Preconditions.checkArgument(x_hat.getDimension() == dim,
        "|x| != |x_hat|, |x|=" + dim + ", |x_hat|=" + x_hat.getDimension());
    Preconditions.checkArgument(sigma.getRowDimension() == dim,
        "|x| != |sigma|, |x|=" + dim + ", |sigma|=" + sigma.getRowDimension());
    Preconditions.checkArgument(sigma.isSquare(), "Sigma is not square matrix");

    LUDecomposition LU = new LUDecomposition(sigma);
    final double detSigma = LU.getDeterminant();
    double denominator = Math.pow(2.d * Math.PI, 0.5d * dim) * Math.pow(detSigma, 0.5d);
    if (denominator == 0.d) { // avoid divide by zero
        return 0.d;
    }

    final RealMatrix invSigma;
    DecompositionSolver solver = LU.getSolver();
    if (solver.isNonSingular() == false) {
        SingularValueDecomposition svd = new SingularValueDecomposition(sigma);
        invSigma = svd.getSolver().getInverse(); // least square solution
    } else {
        invSigma = solver.getInverse();
    }
    //EigenDecomposition eigen = new EigenDecomposition(sigma);
    //double detSigma = eigen.getDeterminant();
    //RealMatrix invSigma = eigen.getSolver().getInverse();

    RealVector diff = x.subtract(x_hat);
    RealVector premultiplied = invSigma.preMultiply(diff);
    double sum = premultiplied.dotProduct(diff);
    double numerator = Math.exp(-0.5d * sum);

    return numerator / denominator;
}
 

public GaussianProcess fit(double[][] X, double[][] Y){
	int ns = X.length;
	this.X = X;
	this.Y = Y;

	// Compute Gram Matrix (Symmetric matrix)
	K = computeGramMatrix(X);
	RealMatrix KMat = MatrixUtils.createRealMatrix(K);
	RealMatrix sigma_n2I = MatrixUtils.createRealIdentityMatrix(ns).scalarMultiply(sigmaN*sigmaN);
	RealMatrix Ky = KMat.add(sigma_n2I);
	this.Ky = Ky.getData();
	LUDecomposition LUDecomp = new  LUDecomposition(Ky);
	detKy = LUDecomp.getDeterminant();

	RealMatrix invKyMat = LUDecomp.getSolver().getInverse();
	invKy = invKyMat.getData();

	// Precompute dY = Y - mX
	int ny=Y[0].length;
	this.mX = new double[ns][ny];
	this.dY = new double[ns][ny];
	for(int i=0; i<ns; i++){
		for(int j=0; j<ny; j++){
			mX[i][j] = meanFunc(X[i],j);
			dY[i][j] = Y[i][j]-mX[i][j];
		}
	}
	
	if(optConstVar==1){
		invKyDY = computeInvKyDY(invKy, dY);
	}
	
	return this;
}
 

@Override
public double getDeterminant(Matrix m) {
  ArgChecker.notNull(m, "m");
  if (m instanceof DoubleMatrix) {
    RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
    LUDecomposition lud = new LUDecomposition(temp);
    return lud.getDeterminant();
  }
  throw new IllegalArgumentException("Can only find determinant of DoubleMatrix; have " + m.getClass());
}
 

/**
 * Creates an instance.
 * 
 * @param lu The result of the LU decomposition, not null. $\mathbf{L}$ cannot be singular.
 */
public LUDecompositionCommonsResult(LUDecomposition lu) {
  ArgChecker.notNull(lu, "LU decomposition");
  ArgChecker.notNull(lu.getL(), "Matrix is singular; could not perform LU decomposition");
  _determinant = lu.getDeterminant();
  _l = CommonsMathWrapper.unwrap(lu.getL());
  _p = CommonsMathWrapper.unwrap(lu.getP());
  _pivot = lu.getPivot();
  _solver = lu.getSolver();
  _u = CommonsMathWrapper.unwrap(lu.getU());
}
 

public static double det(@Nonnull final RealMatrix m) {
    LUDecomposition LU = new LUDecomposition(m);
    return LU.getDeterminant();
}
 

@Override
public Double getDeterminant(final IScope scope) throws GamaRuntimeException {
	final RealMatrix rm = toApacheMatrix(scope);
	final LUDecomposition ld = new LUDecomposition(rm);
	return ld.getDeterminant();
}
 

@Override
public Double getDeterminant(final IScope scope) throws GamaRuntimeException {
	final RealMatrix rm = toApacheMatrix(scope);
	final LUDecomposition ld = new LUDecomposition(rm);
	return ld.getDeterminant();
}