Search in sources :

Example 11 with TetradVector

use of edu.cmu.tetrad.util.TetradVector in project tetrad by cmu-phil.

the class LingamPattern2 method getScore1.

// =============================PRIVATE METHODS=========================//
private Score getScore1(Graph dag, List<TetradMatrix> data, List<Node> variables) {
    // System.out.println("Scoring DAG: " + dag);
    List<Regression> regressions = new ArrayList<>();
    for (TetradMatrix _data : data) {
        regressions.add(new RegressionDataset(_data, variables));
    }
    int totalSampleSize = 0;
    for (TetradMatrix _data : data) {
        totalSampleSize += _data.rows();
    }
    int numCols = data.get(0).columns();
    List<Node> nodes = dag.getNodes();
    double score = 0.0;
    double[] pValues = new double[nodes.size()];
    TetradMatrix absoluteStandardizedResiduals = new TetradMatrix(totalSampleSize, numCols);
    for (int i = 0; i < nodes.size(); i++) {
        List<Double> _absoluteStandardizedResiduals = new ArrayList<>();
        for (int j = 0; j < data.size(); j++) {
            Node _target = nodes.get(i);
            List<Node> _regressors = dag.getParents(_target);
            Node target = getVariable(variables, _target.getName());
            List<Node> regressors = new ArrayList<>();
            for (Node _regressor : _regressors) {
                Node variable = getVariable(variables, _regressor.getName());
                regressors.add(variable);
            }
            RegressionResult result = regressions.get(j).regress(target, regressors);
            TetradVector residualsColumn = result.getResiduals();
            DoubleArrayList _absoluteStandardizedResidualsColumn = new DoubleArrayList(residualsColumn.toArray());
            double mean = Descriptive.mean(_absoluteStandardizedResidualsColumn);
            double std = Descriptive.standardDeviation(Descriptive.variance(_absoluteStandardizedResidualsColumn.size(), Descriptive.sum(_absoluteStandardizedResidualsColumn), Descriptive.sumOfSquares(_absoluteStandardizedResidualsColumn)));
            for (int i2 = 0; i2 < _absoluteStandardizedResidualsColumn.size(); i2++) {
                _absoluteStandardizedResidualsColumn.set(i2, (_absoluteStandardizedResidualsColumn.get(i2) - mean) / std);
                _absoluteStandardizedResidualsColumn.set(i2, Math.abs(_absoluteStandardizedResidualsColumn.get(i2)));
            }
            for (int k = 0; k < _absoluteStandardizedResidualsColumn.size(); k++) {
                _absoluteStandardizedResiduals.add(_absoluteStandardizedResidualsColumn.get(k));
            }
        }
        DoubleArrayList absoluteStandardResidualsList = new DoubleArrayList(absoluteStandardizedResiduals.getColumn(i).toArray());
        for (int k = 0; k < _absoluteStandardizedResiduals.size(); k++) {
            absoluteStandardizedResiduals.set(k, i, _absoluteStandardizedResiduals.get(k));
        }
        double _mean = Descriptive.mean(absoluteStandardResidualsList);
        double diff = _mean - Math.sqrt(2.0 / Math.PI);
        score += diff * diff;
    }
    for (int j = 0; j < absoluteStandardizedResiduals.columns(); j++) {
        double[] x = absoluteStandardizedResiduals.getColumn(j).toArray();
        double p = new AndersonDarlingTest(x).getP();
        pValues[j] = p;
    }
    return new Score(score, pValues);
}
Also used : Regression(edu.cmu.tetrad.regression.Regression) DoubleArrayList(cern.colt.list.DoubleArrayList) ArrayList(java.util.ArrayList) TetradMatrix(edu.cmu.tetrad.util.TetradMatrix) DoubleArrayList(cern.colt.list.DoubleArrayList) RegressionDataset(edu.cmu.tetrad.regression.RegressionDataset) TetradVector(edu.cmu.tetrad.util.TetradVector) AndersonDarlingTest(edu.cmu.tetrad.data.AndersonDarlingTest) RegressionResult(edu.cmu.tetrad.regression.RegressionResult)

Example 12 with TetradVector

use of edu.cmu.tetrad.util.TetradVector in project tetrad by cmu-phil.

the class MVPLikelihood method multipleRegression.

private double multipleRegression(TetradVector Y, TetradMatrix X) {
    int n = X.rows();
    TetradVector r;
    if (X.columns() >= n) {
        TetradVector ones = new TetradVector(n);
        for (int i = 0; i < n; i++) ones.set(i, 1);
        r = ones.scalarMult(ones.dotProduct(Y) / (double) n).minus(Y);
    } else {
        try {
            TetradMatrix Xt = X.transpose();
            TetradMatrix XtX = Xt.times(X);
            r = X.times(XtX.inverse().times(Xt.times(Y))).minus(Y);
        } catch (Exception e) {
            TetradVector ones = new TetradVector(n);
            for (int i = 0; i < n; i++) ones.set(i, 1);
            r = ones.scalarMult(ones.dotProduct(Y) / (double) n).minus(Y);
        }
    }
    double sigma2 = r.dotProduct(r) / n;
    double lik;
    if (sigma2 < 0) {
        TetradVector ones = new TetradVector(n);
        for (int i = 0; i < n; i++) ones.set(i, 1);
        r = ones.scalarMult(ones.dotProduct(Y) / (double) Math.max(n, 2)).minus(Y);
        sigma2 = r.dotProduct(r) / n;
        lik = -(n / 2) * (Math.log(2 * Math.PI) + Math.log(sigma2) + 1);
    } else if (sigma2 == 0) {
        lik = 0;
    } else {
        lik = -(n / 2) * (Math.log(2 * Math.PI) + Math.log(sigma2) + 1);
    }
    if (Double.isInfinite(lik) || Double.isNaN(lik)) {
        System.out.println(lik);
    }
    return lik;
}
Also used : TetradVector(edu.cmu.tetrad.util.TetradVector) TetradMatrix(edu.cmu.tetrad.util.TetradMatrix)

Example 13 with TetradVector

use of edu.cmu.tetrad.util.TetradVector in project tetrad by cmu-phil.

the class MVPLikelihood method getLik.

public double getLik(int child_index, int[] parents) {
    double lik = 0;
    Node c = variables.get(child_index);
    List<ContinuousVariable> continuous_parents = new ArrayList<>();
    List<DiscreteVariable> discrete_parents = new ArrayList<>();
    if (c instanceof DiscreteVariable && discretize) {
        for (int p : parents) {
            Node parent = discreteVariables.get(p);
            discrete_parents.add((DiscreteVariable) parent);
        }
    } else {
        for (int p : parents) {
            Node parent = variables.get(p);
            if (parent instanceof ContinuousVariable) {
                continuous_parents.add((ContinuousVariable) parent);
            } else {
                discrete_parents.add((DiscreteVariable) parent);
            }
        }
    }
    int p = continuous_parents.size();
    List<List<Integer>> cells = adTree.getCellLeaves(discrete_parents);
    int[] continuousCols = new int[p];
    for (int j = 0; j < p; j++) continuousCols[j] = nodesHash.get(continuous_parents.get(j));
    for (List<Integer> cell : cells) {
        // for (int[] cell : cells) {
        int r = cell.size();
        // int r = cell.length;
        if (r > 1) {
            double[] mean = new double[p];
            double[] var = new double[p];
            for (int i = 0; i < p; i++) {
                for (int j = 0; j < r; j++) {
                    mean[i] += continuousData[continuousCols[i]][cell.get(j)];
                    var[i] += Math.pow(continuousData[continuousCols[i]][cell.get(j)], 2);
                }
                mean[i] /= r;
                var[i] /= r;
                var[i] -= Math.pow(mean[i], 2);
                var[i] = Math.sqrt(var[i]);
                if (Double.isNaN(var[i])) {
                    System.out.println(var[i]);
                }
            }
            int degree = fDegree;
            if (fDegree < 1) {
                degree = (int) Math.floor(Math.log(r));
            }
            TetradMatrix subset = new TetradMatrix(r, p * degree + 1);
            for (int i = 0; i < r; i++) {
                subset.set(i, p * degree, 1);
                for (int j = 0; j < p; j++) {
                    for (int d = 0; d < degree; d++) {
                        subset.set(i, p * d + j, Math.pow((continuousData[continuousCols[j]][cell.get(i)] - mean[j]) / var[j], d + 1));
                    }
                }
            }
            if (c instanceof ContinuousVariable) {
                TetradVector target = new TetradVector(r);
                for (int i = 0; i < r; i++) {
                    target.set(i, continuousData[child_index][cell.get(i)]);
                // target.set(i, continuousData[child_index][cell[i]]);
                }
                lik += multipleRegression(target, subset);
            } else {
                TetradMatrix target = new TetradMatrix(r, ((DiscreteVariable) c).getNumCategories());
                for (int i = 0; i < r; i++) {
                    target.set(i, discreteData[child_index][cell.get(i)], 1);
                }
                lik += approxMultinomialRegression(target, subset);
            }
        }
    }
    return lik;
}
Also used : Node(edu.cmu.tetrad.graph.Node) ArrayList(java.util.ArrayList) TetradMatrix(edu.cmu.tetrad.util.TetradMatrix) AtomicInteger(java.util.concurrent.atomic.AtomicInteger) TetradVector(edu.cmu.tetrad.util.TetradVector) List(java.util.List) ArrayList(java.util.ArrayList)

Example 14 with TetradVector

use of edu.cmu.tetrad.util.TetradVector in project tetrad by cmu-phil.

the class Ling method simulateCyclic.

// used to produce dataset if one is not provided as the input to the constructor
private static TetradMatrix simulateCyclic(GraphWithParameters dwp, TetradVector errorCoefficients, int n, Distribution distribution) {
    TetradMatrix reducedForm = reducedForm(dwp);
    TetradMatrix vectors = new TetradMatrix(dwp.getGraph().getNumNodes(), n);
    for (int j = 0; j < n; j++) {
        TetradVector vector = simulateReducedForm(reducedForm, errorCoefficients, distribution);
        vectors.assignColumn(j, vector);
    }
    return vectors;
}
Also used : TetradVector(edu.cmu.tetrad.util.TetradVector) TetradMatrix(edu.cmu.tetrad.util.TetradMatrix)

Example 15 with TetradVector

use of edu.cmu.tetrad.util.TetradVector in project tetrad by cmu-phil.

the class Ling method pruneEdgesByResampling.

// ==============================PRIVATE METHODS====================//
/**
 * This is the method used in Patrik's code.
 */
public TetradMatrix pruneEdgesByResampling(TetradMatrix data) {
    TetradMatrix X = new TetradMatrix(data.transpose().toArray());
    int npieces = 10;
    int cols = X.columns();
    int rows = X.rows();
    int piecesize = (int) Math.floor(cols / npieces);
    List<TetradMatrix> bpieces = new ArrayList<>();
    List<TetradVector> diststdpieces = new ArrayList<>();
    List<TetradVector> cpieces = new ArrayList<>();
    for (int p = 0; p < npieces; p++) {
        // % Select subset of data, and permute the variables to the causal order
        // Xp = X(k,((p-1)*piecesize+1):(p*piecesize));
        int p0 = (p) * piecesize;
        int p1 = (p + 1) * piecesize - 1;
        int[] range = range(p0, p1);
        TetradMatrix Xp = X;
        // % Remember to subract out the mean
        // Xpm = mean(Xp,2);
        // Xp = Xp - Xpm*ones(1,size(Xp,2));
        // 
        // % Calculate covariance matrix
        // cov = (Xp*Xp')/size(Xp,2);
        TetradVector Xpm = new TetradVector(rows);
        for (int i = 0; i < rows; i++) {
            double sum = 0.0;
            for (int j = 0; j < Xp.columns(); j++) {
                sum += Xp.get(i, j);
            }
            Xpm.set(i, sum / Xp.columns());
        }
        for (int i = 0; i < rows; i++) {
            for (int j = 0; j < Xp.columns(); j++) {
                Xp.set(i, j, Xp.get(i, j) - Xpm.get(i));
            }
        }
        TetradMatrix Xpt = Xp.transpose();
        TetradMatrix cov = Xp.times(Xpt);
        for (int i = 0; i < cov.rows(); i++) {
            for (int j = 0; j < cov.columns(); j++) {
                cov.set(i, j, cov.get(i, j) / Xp.columns());
            }
        }
        // % Do QL decomposition on the inverse square root of cov
        // [Q,L] = tridecomp(cov^(-0.5),'ql');
        boolean posDef = LingUtils.isPositiveDefinite(cov);
        if (!posDef) {
            System.out.println("Covariance matrix is not positive definite.");
        }
        TetradMatrix sqrt = cov.sqrt();
        ;
        TetradMatrix I = TetradMatrix.identity(rows);
        TetradMatrix AI = I.copy();
        TetradMatrix invSqrt = sqrt.inverse();
        QRDecomposition qr = new QRDecomposition(invSqrt.getRealMatrix());
        RealMatrix r = qr.getR();
        // % The estimated disturbance-stds are one over the abs of the diag of L
        // newestdisturbancestd = 1./diag(abs(L));
        TetradVector newestdisturbancestd = new TetradVector(rows);
        for (int t = 0; t < rows; t++) {
            newestdisturbancestd.set(t, 1.0 / Math.abs(r.getEntry(t, t)));
        }
        // 
        for (int s = 0; s < rows; s++) {
            for (int t = 0; t < min(s, cols); t++) {
                r.setEntry(s, t, r.getEntry(s, t) / r.getEntry(s, s));
            }
        }
        // % Calculate corresponding B
        // bnewest = eye(dims)-L;
        TetradMatrix bnewest = TetradMatrix.identity(rows);
        bnewest = bnewest.minus(new TetradMatrix(r));
        TetradVector cnewest = new TetradMatrix(r).times(Xpm);
        bpieces.add(bnewest);
        diststdpieces.add(newestdisturbancestd);
        cpieces.add(cnewest);
    }
    // 
    // for i=1:dims,
    // for j=1:dims,
    // 
    // themean = mean(Bpieces(i,j,:));
    // thestd = std(Bpieces(i,j,:));
    // if abs(themean)<prunefactor*thestd,
    // Bfinal(i,j) = 0;
    // else
    // Bfinal(i,j) = themean;
    // end
    // 
    // end
    // end
    TetradMatrix means = new TetradMatrix(rows, rows);
    TetradMatrix stds = new TetradMatrix(rows, rows);
    TetradMatrix BFinal = new TetradMatrix(rows, rows);
    for (int i = 0; i < rows; i++) {
        for (int j = 0; j < rows; j++) {
            double sum = 0.0;
            for (int y = 0; y < npieces; y++) {
                sum += bpieces.get(y).get(i, j);
            }
            double themean = sum / (npieces);
            double sumVar = 0.0;
            for (int y = 0; y < npieces; y++) {
                sumVar += Math.pow((bpieces.get(y).get(i, j)) - themean, 2);
            }
            double thestd = Math.sqrt(sumVar / (npieces));
            means.set(i, j, themean);
            stds.set(i, j, thestd);
            if (Math.abs(themean) < threshold * thestd) {
                // getPruneFactor() * thestd) {
                BFinal.set(i, j, 0);
            } else {
                BFinal.set(i, j, themean);
            }
        }
    }
    return BFinal;
}
Also used : TetradVector(edu.cmu.tetrad.util.TetradVector) QRDecomposition(org.apache.commons.math3.linear.QRDecomposition) RealMatrix(org.apache.commons.math3.linear.RealMatrix) ArrayList(java.util.ArrayList) TetradMatrix(edu.cmu.tetrad.util.TetradMatrix)

Aggregations

TetradVector (edu.cmu.tetrad.util.TetradVector)54 TetradMatrix (edu.cmu.tetrad.util.TetradMatrix)44 ArrayList (java.util.ArrayList)19 Node (edu.cmu.tetrad.graph.Node)11 RegressionResult (edu.cmu.tetrad.regression.RegressionResult)9 Regression (edu.cmu.tetrad.regression.Regression)7 RegressionDataset (edu.cmu.tetrad.regression.RegressionDataset)7 DoubleArrayList (cern.colt.list.DoubleArrayList)5 AndersonDarlingTest (edu.cmu.tetrad.data.AndersonDarlingTest)5 SingularMatrixException (org.apache.commons.math3.linear.SingularMatrixException)5 List (java.util.List)4 Test (org.junit.Test)4 ContinuousVariable (edu.cmu.tetrad.data.ContinuousVariable)3 DepthChoiceGenerator (edu.cmu.tetrad.util.DepthChoiceGenerator)2 RandomUtil (edu.cmu.tetrad.util.RandomUtil)2 Vector (java.util.Vector)2 DoubleMatrix1D (cern.colt.matrix.DoubleMatrix1D)1 DoubleMatrix2D (cern.colt.matrix.DoubleMatrix2D)1 DenseDoubleMatrix1D (cern.colt.matrix.impl.DenseDoubleMatrix1D)1 DenseDoubleMatrix2D (cern.colt.matrix.impl.DenseDoubleMatrix2D)1