Example usage for org.apache.commons.math3.util FastMath signum

List of usage examples for org.apache.commons.math3.util FastMath signum

Introduction

In this page you can find the example usage for org.apache.commons.math3.util FastMath signum.

Prototype

public static float signum(final float a) 

Source Link

Document

Compute the signum of a number.

Usage

From source file:gamlss.distributions.PE.java

/**  First derivative dldm = dl/dmu, where l - log-likelihood function.
 * @param y - vector of values of response variable
 * @return  a vector of first derivative dldm = dl/dmu
 *///from ww w.j  ava2s . c om
public final ArrayRealVector dldm(final ArrayRealVector y) {

    dldm = new double[size];
    for (int i = 0; i < size; i++) {

        //dldm <- (sign(z)*nu)/(2*sigma*abs(z));
        dldm[i] = (FastMath.signum(z[i]) * nuV.getEntry(i)) / (2 * sigmaV.getEntry(i) * FastMath.abs(z[i]));

        //dldm <- dldm*((abs(z/c))^nu) 
        dldm[i] = dldm[i] * (FastMath.pow(FastMath.abs(z[i] / c[i]), nuV.getEntry(i)));
    }
    logC = null;
    c = null;
    z = null;
    return new ArrayRealVector(dldm, false);
}

From source file:gamlss.distributions.GT.java

/**  First derivative dldm = dl/dmu, where l - log-likelihood function.
 * @param y - vector of values of response variable
 * @return  a vector of first derivative dldm = dl/dmu
 *//*from   w ww  .  j a va 2 s .  co  m*/
public final ArrayRealVector dldm(final ArrayRealVector y) {

    dldm = new double[size];
    for (int i = 0; i < size; i++) {

        //dldm <- w*((abs(z))^(tau-1))*sign(z)/sigma
        dldm[i] = w[i] * (FastMath.pow(FastMath.abs(z[i]), tauV.getEntry(i) - 1)) * (FastMath.signum(z[i]))
                / sigmaV.getEntry(i);
    }
    w = null;
    z = null;
    zt = null;
    return new ArrayRealVector(dldm, false);
}

From source file:gamlss.smoothing.PB.java

/**
 * The main fitting method, initiate appropriate smoothing 
 * function according to incoming parameters.
 * @param additiveFit -object of AdditiveFoit class
 * @return reiduals/*from   ww  w  .j a  v a2s  .  c  o  m*/
 */
//gamlss.pb <- function(x, y, w, xeval = NULL, ...)      
public ArrayRealVector solve(final AdditiveFit additiveFit) {

    Double lambda = Controls.LAMBDAS_USER_DEFINED;
    colNum = additiveFit.getColNum();
    y = additiveFit.getZ();
    w = additiveFit.getW();
    whichDistParameter = additiveFit.getWhichDistParameter();
    basisM = (BlockRealMatrix) getBasisM().get(colNum);
    penaltyM = (BlockRealMatrix) getPenaltyM().get(colNum);
    df = (Integer) getDfValues().get(colNum);

    //n <- nrow(X) # the no of observations
    n = basisM.getRowDimension();

    //lambdaS <-  get(startLambdaName, envir=gamlss.env)
    //geting the starting value
    lambdaS = getLambdas().get(whichDistParameter);

    //if (lambdaS>=1e+07) lambda <- 1e+07
    if (lambdaS >= 1e+07) {
        lambda = 1e+07;
    }
    //if (lambdaS<=1e-07) lambda <- 1e-07
    if (lambdaS <= 1e-07) {
        lambda = 1e-07;
    }

    //if (is.null(df)&&!is.null(lambda)||!is.null(df)&&!is.null(lambda))
    if (lambda != null) {

        //fit <- regpen(y, X, w, lambda, D)
        beta = regpen(lambda);

        //fv <- X %*% fit$beta
        fv = (ArrayRealVector) basisM.operate(beta);

        //else if (is.null(df)&&is.null(lambda))
    } else if (df == null) {

        //lambda <- lambdaS
        lambda = lambdaS;

        //switch(control$c," ML"={
        switch (Controls.SMOOTH_METHOD) {
        case DistributionSettings.GAIC:
            fv = functionGAIC(lambda);
            break;
        case DistributionSettings.ML:
            fv = functionML(lambda);
            break;
        case DistributionSettings.ML1:
            fv = functionML1(lambda);
            break;
        case DistributionSettings.EM:
            System.err.println("EM has not been implemented yet");
            break;
        case DistributionSettings.GCV:
            fv = functionGCV(lambda);
            break;
        default:
            System.err.println(
                    " Cannot recognise the " + "smoothing method or it has" + "not been implemented yet");
        }
    } else {

        //QR <- qr(sqrt(w)*X)
        //Rinv <- solve(qr.R(QR))
        rM = (BlockRealMatrix) new QRDecomposition(
                MatrixFunctions.multVectorMatrix(MatrixFunctions.sqrtVec(w), basisM)).getR();
        rM = rM.getSubMatrix(0, rM.getColumnDimension() - 1, 0, rM.getColumnDimension() - 1);
        rMinv = (BlockRealMatrix) new QRDecomposition(rM).getSolver().getInverse();

        //S   <- t(D)%*%D
        sM = penaltyM.transpose().multiply(penaltyM);

        //UDU <- eigen(t(Rinv)%*%S%*%Rinv)
        uduV = new EigenDecomposition(rMinv.transpose().multiply(sM).multiply(rMinv)).getRealEigenvalues();

        //lambda <- if (sign(edf1_df(0))==sign(edf1_df(100000))) 100000
        //in case they have the some sign
        edfTemp1 = edf1_df(0, df);
        edfTemp2 = edf1_df(100000.0, df);

        if (FastMath.signum(edfTemp1) == FastMath.signum(edfTemp2)) {
            lambda = 100000.0;
        } else {

            //else uniroot(edf1_df, c(0,100000))$root
            uniRootObj.setDf(df);
            final double relativeAccuracy = 1.0e-12;
            final double absoluteAccuracy = 1.0e-8;
            UnivariateSolver uniRootSolver = new BrentSolver(relativeAccuracy, absoluteAccuracy);
            lambda = uniRootSolver.solve(1000, uniRootObj, 0.0, 100000.0);
        }

        //fit <- regpen(y, X, w, lambda, D)
        beta = regpen(lambda);
        fv = (ArrayRealVector) basisM.operate(beta);
    }
    if (Controls.IF_NEED_THIS) {

        //waug <- as.vector(c(w, rep(1,nrow(D))))
        waug = w.append(MatrixFunctions.repV(1.0, penaltyM.getRowDimension()));

        //xaug <- as.matrix(rbind(X,sqrt(lambda)*D))
        xaug = MatrixFunctions.appendMatricesRows(basisM,
                (BlockRealMatrix) penaltyM.scalarMultiply(FastMath.sqrt(lambda)));

        //lev <- hat(sqrt(waug)*xaug,intercept=FALSE)[1:n]
        //get the hat matrix
        lev = (ArrayRealVector) MatrixFunctions.getMainDiagonal(new BlockRealMatrix(MatrixFunctions
                .calculateHat(MatrixFunctions.multVectorMatrix(MatrixFunctions.sqrtVec(waug), xaug)).getData()))
                .getSubVector(0, n);

        //lev <- (lev-.hat.WX(w,x))
        //subtract  the linear since is already fitted 
        lev = lev.subtract(
                hatWX((ArrayRealVector) getSmootherMatrices().get(whichDistParameter).getColumnVector(colNum)));

        // var <- lev/w
        //the variance of the smootherz
        var = lev.ebeDivide(w);
    }

    //if (is.null(xeval)) # if no prediction
    if (Controls.XEVAL_USER_DEFINED == null) {
        // Residuals
        return y.subtract(fv);

        //else # for prediction
    } else {

        //ll <- dim(as.matrix(attr(x,"X")))[1]                 
        //nx <- as.matrix(attr(x,"X"))[seq(length(y)+1,ll),]
        tempArr = ArithmeticSeries.getSeries(y.getDimension() + 1, basisM.getRowDimension(), 1);
        BlockRealMatrix nx = new BlockRealMatrix(tempArr.length, basisM.getColumnDimension());
        for (int i = 0; i < tempArr.length; i++) {
            nx.setRowVector(i, basisM.getRowVector((int) tempArr[i]));
        }

        //pred <- drop(nx %*% fit$beta)
        ArrayRealVector pred = (ArrayRealVector) nx.operate(beta);

        // Residuals
        return y.subtract(fv);
    }
}

From source file:gamlss.distributions.BCPE.java

/**
 * Supportive function to compute first derivatives.
 * @param t -  value//w  ww. j  a va 2s  .  c o m
 * @param tau -  value of tau distribution parameter
 * @return modified value of regularized gamma function 
 * @param i - loop count
 */
private double f2T(final double t, final double tau, final int i) {

    //s <- 0.5*((abs(t/c))^tau)
    final double s = 0.5 * ((FastMath.pow(FastMath.abs(t / c[i]), tau)));

    //F.s <- pgamma(s,shape = 1/tau, scale = 1)
    gammaDistr.setDistrParameters(1 / tau, 1.0);
    final double fS = gammaDistr.cumulativeProbability(s);

    //cdf <- 0.5*(1+F.s*sign(t))
    return 0.5 * (1 + fS * FastMath.signum(t));
}

From source file:gamlss.distributions.BCPE.java

/** First derivative dldn = dl/dnu,
 *  where l - log-likelihood function.//from   ww w . ja  v a2s.  c  om
 * @param y - vector of values of response variable
 * @return  a vector of First derivative dldn = dl/dnu
 */
public final ArrayRealVector dldn(final ArrayRealVector y) {

    dldn = new double[size];
    for (int i = 0; i < size; i++) {

        //h <- f.T(1/(sigma*abs(nu)))/F.T(1/(sigma*abs(nu)))
        final double h = f1T(1 / (sigmaV.getEntry(i) * FastMath.abs(nuV.getEntry(i))), tauV.getEntry(i), false,
                i) / f2T(1 / (sigmaV.getEntry(i) * FastMath.abs(nuV.getEntry(i))), tauV.getEntry(i), i);

        //dldv <-  -(tau/(2*nu*c^2))*((abs(z/c))^(tau-2))
        //*z*((nu*z+1/sigma) *log(y/mu)-z)
        dldn[i] = -(tauV.getEntry(i) / (2 * nuV.getEntry(i) * c[i] * c[i]))
                * (FastMath.pow((FastMath.abs(z[i] / c[i])), (tauV.getEntry(i) - 2))) * z[i]
                * ((nuV.getEntry(i) * z[i] + 1 / sigmaV.getEntry(i))
                        * FastMath.log(y.getEntry(i) / muV.getEntry(i)) - z[i]);

        //dldv <- dldv+log(y/mu)+sign(nu)*h/(sigma*nu^2)
        dldn[i] = dldn[i] + FastMath.log(y.getEntry(i) / muV.getEntry(i)) + FastMath.signum(nuV.getEntry(i)) * h
                / (sigmaV.getEntry(i) * nuV.getEntry(i) * nuV.getEntry(i));
    }
    c = null;
    logC = null;
    z = null;
    return new ArrayRealVector(dldn, false);
}

From source file:de.tuberlin.uebb.jbop.example.DerivativeStructure.java

/**
 * Compute the signum of the instance.//ww w. ja v  a2  s.com
 * The signum is -1 for negative numbers, +1 for positive numbers and 0 otherwise
 * 
 * @return -1.0, -0.0, +0.0, +1.0 or NaN depending on sign of a
 */
public DerivativeStructure signum() {
    return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), FastMath.signum(data[0]));
}

From source file:de.tuberlin.uebb.jbop.example.DerivativeStructureOnlyCompose.java

/**
 * Compute the signum of the instance./*from w w w.  j a va 2  s .c  o m*/
 * The signum is -1 for negative numbers, +1 for positive numbers and 0 otherwise
 * 
 * @return -1.0, -0.0, +0.0, +1.0 or NaN depending on sign of a
 */
public DerivativeStructureOnlyCompose signum() {
    return new DerivativeStructureOnlyCompose(compiler.getFreeParameters(), compiler.getOrder(),
            FastMath.signum(data[0]));
}

From source file:gamlss.distributions.PE.java

/** Computes the cumulative distribution 
 * function P(X <= q) for a random variable X .
 * whose values are distributed according to this distribution
 * @param q - value of quantile//  w w w .  ja  va  2  s.  c o  m
 * @param mu - value of mu distribution parameter
 * @param sigma - value of sigma distribution parameter
 * @param nu - value of nu distribution parameter 
 * @param lowerTail - logical, if TRUE (default), probabilities
 *  are P[X <= x] otherwise, P[X > x].
 * @param isLog - logical, if TRUE, probabilities p are given as log(p)
 * @return value of cumulative probability function values P(X <= q)
 */
public final double pPE(final double q, final double mu, final double sigma, final double nu,
        final boolean lowerTail, final boolean isLog) {

    // {  if (any(sigma < 0))stop(paste("sigma must be positive",))
    if (sigma < 0) {
        System.err.println("sigma must be positive");
        return -1.0;
    }

    //if (any(nu < 0))  stop(paste("nu must be positive", "\n", ""))
    if (nu < 0) {
        System.err.println("nu must be positive");
        return -1.0;
    }

    double out = 0;

    //ifelse(nu>10000, (q-(mu-sqrt(3)*sigma))/(sqrt(12)*sigma),cdf)
    if (nu > 10000) {
        out = (q - (mu - FastMath.sqrt(3) * sigma)) / (FastMath.sqrt(12) * sigma);
    } else {

        //log.c <- 0.5*(-(2/nu)*log(2)+lgamma(1/nu)-lgamma(3/nu))
        final double logC = 0.5
                * (-(2 / nu) * FastMath.log(2) + Gamma.logGamma(1 / nu) - Gamma.logGamma(3 / nu));

        //c <- exp(log.c)
        final double c = FastMath.exp(logC);

        //z <- (y-mu)/sigma
        final double z = (q - mu) / sigma;

        //s <- 0.5*((abs(z/c))^nu)
        final double s = 0.5 * (FastMath.pow(FastMath.abs(z / c), nu));

        //cdf <- 0.5*(1+pgamma(s,shape=1/nu,scale=1)*sign(z))
        out = 0.5 * (1 + Gamma.regularizedGammaP(1 / nu, s) * FastMath.signum(z));
    }

    //if(lower.tail==TRUE) cdf  <- cdf else  cdf <- 1-cdf 
    //if(log.p==FALSE) cdf  <- cdf else  cdf <- log(cdf) 
    if (!lowerTail) {
        if (isLog) {
            out = FastMath.log(1 - out);
        } else {
            out = 1 - out;
        }
    } else if (isLog) {
        //if(log.p==FALSE) cdf  <- cdf else  cdf <- log(cdf)
        out = FastMath.log(out);
    }
    return out;
}

From source file:gamlss.distributions.PE.java

/** Computes the quantile (inverse cumulative probability)
 *  function  of this distribution.//  w ww .j av a 2  s .c o m
* @param p - value of cumulative probability
* @param mu -  value of mu distribution parameters
* @param sigma -  value of sigma distribution parameters
* @param nu -  value of nu distribution parameters 
* @param lowerTail - logical; if TRUE (default), probabilities 
* are P[X <= x] otherwise, P[X > x]
* @param isLog - logical; if TRUE, probabilities p are given as log(p).
* @return value of quantile function
*/
public final double qPE(final double p, final double mu, final double sigma, final double nu,
        final boolean lowerTail, final boolean isLog) {

    // {  if (any(sigma < 0))stop(paste("sigma must be positive",))
    if (sigma < 0) {
        System.err.println("sigma must be positive");
        return -1.0;
    }

    //if (any(nu < 0))  stop(paste("nu must be positive", "\n", ""))
    if (nu < 0) {
        System.err.println("nu must be positive");
        return -1.0;
    }

    //if (log.p==TRUE) p <- exp(p) else p <- p
    double out = p;
    if (isLog) {

        //if (lower.tail==TRUE) p <- p else p <- 1-p
        if (!lowerTail) {

            out = 1 - FastMath.exp(out);
            //out = FastMath.exp(1 -out);
        } else {

            out = FastMath.exp(out);
        }
    } else if (!lowerTail) {

        out = 1 - out;
    }

    //if (any(p < 0)|any(p > 1))  
    //stop(paste("p must be between 0 and 1", "\n", "")
    if (out < 0 || out > 1) {
        System.err.println("p must be between 0 and 1");
        return -1.0;
    }

    //log.c <- 0.5*(-(2/nu)*log(2)+lgamma(1/nu)-lgamma(3/nu))
    final double logC = 0.5 * (-(2 / nu) * FastMath.log(2) + Gamma.logGamma(1 / nu) - Gamma.logGamma(3 / nu));

    //c <- exp(log.c)
    final double c = FastMath.exp(logC);

    //suppressWarnings(
    //s <- qgamma((2*p-1)*sign(p-0.5),shape=(1/nu),scale=1))
    gammaDistr.setDistrParameters(1 / nu, 1.0);
    final double s = gammaDistr.inverseCumulativeProbability((2 * out - 1) * FastMath.signum(out - 0.5));

    //z <- sign(p-0.5)*((2*s)^(1/nu))*c
    final double z = FastMath.signum(out - 0.5) * (FastMath.pow((2 * s), (1 / nu))) * c;

    //ya <- mu + sigma*z 
    out = mu + sigma * z;

    return out;
}

From source file:com.musicplayer.AudioDecoderThread.java

private double sign(double x) {
    // TODO Auto-generated method stub
    return (FastMath.signum(x) / 2.0) + 0.5;
}