Examples with DoubleMatrix1D cern.colt.matrix.DoubleMatrix1D used on opensource projects

Example 21 with DoubleMatrix1D

use of cern.colt.matrix.DoubleMatrix1D in project Gemma by PavlidisLab.

the class ExpressionDataSVD method winnow.

/**
 * Implements method described in Skillicorn et al., "Strategies for winnowing microarray data" (also section 3.5.5
 * of his book)
 *
 * @param thresholdQuantile Enter 0.5 for median. Value must be > 0 and < 1.
 * @return a filtered matrix
 */
public ExpressionDataDoubleMatrix winnow(double thresholdQuantile) {
    if (thresholdQuantile <= 0 || thresholdQuantile >= 1) {
        throw new IllegalArgumentException("Threshold quantile should be a value between 0 and 1 exclusive");
    }
    class NormCmp implements Comparable<NormCmp> {

        private Double norm;

        private int rowIndex;

        private NormCmp(int rowIndex, Double norm) {
            super();
            this.rowIndex = rowIndex;
            this.norm = norm;
        }

        @Override
        public int compareTo(NormCmp o) {
            return this.norm.compareTo(o.norm);
        }

        public int getRowIndex() {
            return rowIndex;
        }

        @Override
        public int hashCode() {
            final int prime = 31;
            int result = 1;
            result = prime * result + ((norm == null) ? 0 : norm.hashCode());
            return result;
        }

        @Override
        public boolean equals(Object obj) {
            if (this == obj)
                return true;
            if (obj == null)
                return false;
            if (this.getClass() != obj.getClass())
                return false;
            NormCmp other = (NormCmp) obj;
            if (norm == null) {
                return other.norm == null;
            } else
                return norm.equals(other.norm);
        }
    }
    // order rows by distance from the origin. This is proportional to the 1-norm.
    Algebra a = new Algebra();
    List<NormCmp> os = new ArrayList<>();
    for (int i = 0; i < this.expressionData.rows(); i++) {
        double[] row = this.getU().getRow(i);
        DoubleMatrix1D rom = new DenseDoubleMatrix1D(row);
        double norm1 = a.norm1(rom);
        os.add(new NormCmp(i, norm1));
    }
    Collections.sort(os);
    int quantileLimit = (int) Math.floor(this.expressionData.rows() * thresholdQuantile);
    quantileLimit = Math.max(0, quantileLimit);
    List<CompositeSequence> keepers = new ArrayList<>();
    for (int i = 0; i < quantileLimit; i++) {
        NormCmp x = os.get(i);
        CompositeSequence d = this.expressionData.getDesignElementForRow(x.getRowIndex());
        keepers.add(d);
    }
    // remove genes which are near the origin in SVD space. FIXME: make sure the missing values are still masked.
    return new ExpressionDataDoubleMatrix(this.expressionData, keepers);
}

Example 22 with DoubleMatrix1D

use of cern.colt.matrix.DoubleMatrix1D in project Gemma by PavlidisLab.

the class ComBat method getBatchData.

/**
 * @param sdata   data to be sliced
 * @param batchId which batch
 */
private DoubleMatrix2D getBatchData(DoubleMatrix2D sdata, String batchId) {
    Collection<C> sampleNames = batches.get(batchId);
    DoubleMatrix2D result = new DenseDoubleMatrix2D(sdata.rows(), sampleNames.size());
    int i = 0;
    for (C sname : sampleNames) {
        DoubleMatrix1D colInBatch = sdata.viewColumn(data.getColIndexByName(sname));
        for (int k = 0; k < colInBatch.size(); k++) {
            result.set(k, i, colInBatch.get(k));
        }
        i++;
    }
    // log.info( result );
    return result;
}

Example 23 with DoubleMatrix1D

use of cern.colt.matrix.DoubleMatrix1D in project Gemma by PavlidisLab.

the class ComBat method itSol.

private DoubleMatrix1D[] itSol(DoubleMatrix2D matrix, DoubleMatrix1D gHat, DoubleMatrix1D dHat, double gbar, double t2b, double a, double b) throws ComBatException {
    DoubleMatrix1D n = this.rowNonMissingCounts(matrix);
    DoubleMatrix1D gold = gHat;
    DoubleMatrix1D dold = dHat;
    final double conv = 0.0001;
    double change = 1.0;
    int count = 0;
    int MAXITERS = 500;
    while (change > conv) {
        DoubleMatrix1D gnew = this.postMean(gHat, gbar, n, dold, t2b);
        DoubleMatrix1D sum2 = this.stepSum(matrix, gnew);
        DoubleMatrix1D dnew = this.postVar(sum2, n, a, b);
        DoubleMatrix1D gnewtmp = gnew.copy().assign(gold, Functions.minus).assign(Functions.abs).assign(gold, Functions.div);
        DoubleMatrix1D dnewtmp = dnew.copy().assign(dold, Functions.minus).assign(Functions.abs).assign(dold, Functions.div);
        double gnewmax;
        double dnewmax;
        if (hasMissing) {
            gnewmax = DescriptiveWithMissing.max(new DoubleArrayList(gnewtmp.toArray()));
            dnewmax = DescriptiveWithMissing.max(new DoubleArrayList(dnewtmp.toArray()));
        } else {
            gnewmax = gnewtmp.aggregate(Functions.max, Functions.identity);
            dnewmax = dnewtmp.aggregate(Functions.max, Functions.identity);
        }
        change = Math.max(gnewmax, dnewmax);
        gold = gnew;
        dold = dnew;
        if (count++ > MAXITERS) {
            /*
                 * For certain data sets, we just flail around; for example if there are only two samples. This is a
                 * bailout for exceptional circumstances.
                 */
            throw new ComBatException("Failed to converge within " + MAXITERS + " iterations, last delta was " + String.format("%.2g", change));
        }
    }
    return new DoubleMatrix1D[] { gold, dold };
}

Example 24 with DoubleMatrix1D

use of cern.colt.matrix.DoubleMatrix1D in project Gemma by PavlidisLab.

the class ComBat method runParametric.

private void runParametric(final DoubleMatrix2D sdata, DoubleMatrix2D gammastar, DoubleMatrix2D deltastar) {
    int batchIndex = 0;
    for (String batchId : batches.keySet()) {
        DoubleMatrix2D batchData = this.getBatchData(sdata, batchId);
        DoubleMatrix1D[] batchResults;
        batchResults = this.itSol(batchData, gammaHat.viewRow(batchIndex), deltaHat.viewRow(batchIndex), gammaBar.get(batchIndex), t2.get(batchIndex), aPrior.get(batchIndex), bPrior.get(batchIndex));
        for (int j = 0; j < batchResults[0].size(); j++) {
            gammastar.set(batchIndex, j, batchResults[0].get(j));
        }
        for (int j = 0; j < batchResults[1].size(); j++) {
            deltastar.set(batchIndex, j, batchResults[1].get(j));
        }
        batchIndex++;
    }
}

Example 25 with DoubleMatrix1D

use of cern.colt.matrix.DoubleMatrix1D in project Gemma by PavlidisLab.

the class ComBat method stepSum.

private DoubleMatrix1D stepSum(DoubleMatrix2D matrix, DoubleMatrix1D gnew) {
    Algebra s = new Algebra();
    DoubleMatrix2D g = new DenseDoubleMatrix2D(1, gnew.size());
    for (int i = 0; i < gnew.size(); i++) {
        g.set(0, i, gnew.get(i));
    }
    DoubleMatrix2D a = new DenseDoubleMatrix2D(1, matrix.columns());
    a.assign(1.0);
    /*
         * subtract column gnew from each column of data; square; then sum over each row.
         */
    DoubleMatrix2D deltas = matrix.copy().assign((s.mult(s.transpose(g), a)), Functions.minus).assign(Functions.square);
    DoubleMatrix1D sumsq = new DenseDoubleMatrix1D(deltas.rows());
    sumsq.assign(0.0);
    for (int i = 0; i < deltas.rows(); i++) {
        sumsq.set(i, DescriptiveWithMissing.sum(new DoubleArrayList(deltas.viewRow(i).toArray())));
    }
    return sumsq;
}