Examples with Algebra cern.colt.matrix.linalg.Algebra used on opensource projects

Example 6 with Algebra

use of cern.colt.matrix.linalg.Algebra in project tdq-studio-se by Talend.

the class QRTest method main.

public static void main(String[] args) {
    // For COLT
    DoubleMatrix2D xmatrix, ymatrix, zmatrix;
    DoubleFactory2D myfactory;
    myfactory = DoubleFactory2D.dense;
    xmatrix = myfactory.make(8, 2);
    ymatrix = myfactory.make(8, 1);
    xmatrix.set(0, 0, 1);
    xmatrix.set(1, 0, 1);
    xmatrix.set(2, 0, 1);
    xmatrix.set(3, 0, 1);
    xmatrix.set(4, 0, 1);
    xmatrix.set(5, 0, 1);
    xmatrix.set(6, 0, 1);
    xmatrix.set(7, 0, 1);
    xmatrix.set(0, 1, 80);
    xmatrix.set(1, 1, 220);
    xmatrix.set(2, 1, 140);
    xmatrix.set(3, 1, 120);
    xmatrix.set(4, 1, 180);
    xmatrix.set(5, 1, 100);
    xmatrix.set(6, 1, 200);
    xmatrix.set(7, 1, 160);
    ymatrix.set(0, 0, 0.6);
    ymatrix.set(1, 0, 6.70);
    ymatrix.set(2, 0, 5.30);
    ymatrix.set(3, 0, 4.00);
    ymatrix.set(4, 0, 6.55);
    ymatrix.set(5, 0, 2.15);
    ymatrix.set(6, 0, 6.60);
    ymatrix.set(7, 0, 5.75);
    Algebra myAlgebra = new Algebra();
    zmatrix = myAlgebra.solve(xmatrix, ymatrix);
    System.err.println(xmatrix);
    System.err.println(ymatrix);
    System.err.println(zmatrix);
/*
       // For JAMA
       Matrix amatrix,bmatrix,cmatrix;
       amatrix = new Matrix(8,2);
       bmatrix = new Matrix(8,1);

       amatrix.set(0,0,1);
       amatrix.set(1,0,1);
       amatrix.set(2,0,1);
       amatrix.set(3,0,1);
       amatrix.set(4,0,1);
       amatrix.set(5,0,1);
       amatrix.set(6,0,1);
       amatrix.set(7,0,1);

       amatrix.set(0,1,80);
       amatrix.set(1,1,220);
       amatrix.set(2,1,140);
       amatrix.set(3,1,120);
       amatrix.set(4,1,180);
       amatrix.set(5,1,100);
       amatrix.set(6,1,200);
       amatrix.set(7,1,160);

       bmatrix.set(0,0,0.6);
       bmatrix.set(1,0,6.70);
       bmatrix.set(2,0,5.30);
       bmatrix.set(3,0,4.00);
       bmatrix.set(4,0,6.55);
       bmatrix.set(5,0,2.15);
       bmatrix.set(6,0,6.60);
       bmatrix.set(7,0,5.75);

       cmatrix = amatrix.solve(bmatrix);
       amatrix.print(8,5);
       bmatrix.print(8,5);
       cmatrix.print(8,5);
		*/
}

Example 7 with Algebra

use of cern.colt.matrix.linalg.Algebra in project Gemma by PavlidisLab.

the class ExpressionDataSVD method equalize.

/**
 * Implements the method described in the SPELL paper, alternative interpretation as related by Q. Morris. Set all
 * components to have equal weight (set all singular values to 1)
 *
 * @return the reconstructed matrix; values that were missing before are re-masked.
 */
public ExpressionDataDoubleMatrix equalize() {
    DoubleMatrix<Integer, Integer> copy = svd.getS().copy();
    for (int i = 0; i < copy.columns(); i++) {
        copy.set(i, i, 1.0);
    }
    double[][] rawU = svd.getU().getRawMatrix();
    double[][] rawS = copy.getRawMatrix();
    double[][] rawV = svd.getV().getRawMatrix();
    DoubleMatrix2D u = new DenseDoubleMatrix2D(rawU);
    DoubleMatrix2D s = new DenseDoubleMatrix2D(rawS);
    DoubleMatrix2D v = new DenseDoubleMatrix2D(rawV);
    Algebra a = new Algebra();
    DoubleMatrix<CompositeSequence, BioMaterial> reconstructed = new DenseDoubleMatrix<>(a.mult(a.mult(u, s), a.transpose(v)).toArray());
    reconstructed.setRowNames(this.expressionData.getMatrix().getRowNames());
    reconstructed.setColumnNames(this.expressionData.getMatrix().getColNames());
    // re-mask the missing values.
    for (int i = 0; i < reconstructed.rows(); i++) {
        for (int j = 0; j < reconstructed.columns(); j++) {
            if (Double.isNaN(this.missingValueInfo.get(i, j))) {
                reconstructed.set(i, j, Double.NaN);
            }
        }
    }
    return new ExpressionDataDoubleMatrix(this.expressionData, reconstructed);
}

Example 8 with Algebra

use of cern.colt.matrix.linalg.Algebra 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 9 with Algebra

use of cern.colt.matrix.linalg.Algebra 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;
}

Example 10 with Algebra

use of cern.colt.matrix.linalg.Algebra in project midpoint by Evolveum.

the class AbstractOrgClosureTest method checkClosureMatrix.

protected boolean checkClosureMatrix() throws SchemaException {
    try (Session session = openSession()) {
        // we compute the closure table "by hand" as 1 + A + A^2 + A^3 + ... + A^n where n is the greatest expected path length
        int vertices = getVertices().size();
        long start = System.currentTimeMillis();
        // used to give indices to vertices
        List<String> vertexList = new ArrayList<>(getVertices());
        if (DUMP_TC_MATRIX_DETAILS) {
            logger.info("Vertex list = {}", vertexList);
        }
        DoubleMatrix2D a = new SparseDoubleMatrix2D(vertices, vertices);
        // }
        for (DefaultEdge edge : orgGraph.edgeSet()) {
            a.set(vertexList.indexOf(orgGraph.getEdgeSource(edge)), vertexList.indexOf(orgGraph.getEdgeTarget(edge)), 1.0);
        }
        DoubleMatrix2D result = new SparseDoubleMatrix2D(vertices, vertices);
        for (int i = 0; i < vertices; i++) {
            result.setQuick(i, i, 1.0);
        }
        DoubleMatrix2D power = result.copy();
        Algebra alg = new Algebra();
        for (int level = 1; level <= maxLevel; level++) {
            power = alg.mult(power, a);
            result.assign(power, Functions.plus);
        // System.out.println("a=" + a);
        // System.out.println("a^"+level+"="+power);
        }
        logger.info("TC matrix computed in {} ms", System.currentTimeMillis() - start);
        if (DUMP_TC_MATRIX_DETAILS) {
            logger.info("TC matrix expected = {}", result);
        }
        Query q = session.createNativeQuery("select descendant_oid, ancestor_oid, val from m_org_closure").addScalar("descendant_oid", StringType.INSTANCE).addScalar("ancestor_oid", StringType.INSTANCE).addScalar("val", LongType.INSTANCE);
        List<Object[]> list = q.list();
        logger.info("OrgClosure has {} rows", list.size());
        DoubleMatrix2D closureInDatabase = new SparseDoubleMatrix2D(vertices, vertices);
        for (Object[] item : list) {
            int val = Integer.parseInt(item[2].toString());
            if (val == 0) {
                throw new IllegalStateException("Row with val == 0 in closure table: " + list);
            }
            closureInDatabase.set(vertexList.indexOf(item[0]), vertexList.indexOf(item[1]), val);
        }
        if (DUMP_TC_MATRIX_DETAILS) {
            logger.info("TC matrix fetched from db = {}", closureInDatabase);
        }
        double zSumResultBefore = result.zSum();
        double zSumClosureInDb = closureInDatabase.zSum();
        result.assign(closureInDatabase, Functions.minus);
        double zSumResultAfter = result.zSum();
        logger.info("Summary of items in closure computed: {}, in DB-stored closure: {}, delta: {}", zSumResultBefore, zSumClosureInDb, zSumResultAfter);
        if (DUMP_TC_MATRIX_DETAILS) {
            logger.info("Difference matrix = {}", result);
        }
        boolean problem = false;
        for (int i = 0; i < vertices; i++) {
            for (int j = 0; j < vertices; j++) {
                double delta = result.get(i, j);
                if (Math.round(delta) != 0) {
                    System.err.println("delta(" + vertexList.get(i) + "," + vertexList.get(j) + ") = " + delta + " (closureInDB=" + closureInDatabase.get(i, j) + ", expected=" + (result.get(i, j) + closureInDatabase.get(i, j)) + ")");
                    logger.error("delta(" + vertexList.get(i) + "," + vertexList.get(j) + ") = " + delta);
                    problem = true;
                }
            }
        }
        if (problem) {
            checkOrgGraph();
        }
        return problem;
    }
}