use of cern.colt.matrix.linalg.Algebra in project tetrad by cmu-phil.
the class Ricf method lik.
private double lik(DoubleMatrix2D K, DoubleMatrix2D S, int n, int k) {
Algebra algebra = new Algebra();
DoubleMatrix2D SK = algebra.mult(S, K);
return (algebra.trace(SK) - Math.log(algebra.det(SK)) - k) * n;
}
use of cern.colt.matrix.linalg.Algebra in project tetrad by cmu-phil.
the class Glasso method fatmul.
// s = vv * x, or s12 = Theta.1 * Theta.12
private void fatmul(int it, int n, DoubleMatrix2D vv, DoubleMatrix1D x, DoubleMatrix1D s, DoubleMatrix1D z, int[] m) {
double fac = 0.2;
// z consists of the nonzero entries of x. m indexes these. If there are enough zeroes in x,
// use the simpler multiplication method.
int l = 0;
for (int j = 0; j < n; j++) {
if (x.get(j) == 0.0)
continue;
l = l + 1;
m[l] = j;
}
if (l < (int) (fac * n)) {
int[] indices = new int[l];
for (int i = 0; i < indices.length; i++) indices[i] = i;
if (it == 1) {
for (int j = 0; j < n; j++) {
if (Thread.currentThread().isInterrupted()) {
break;
}
double dotProduct = 0.0;
for (int i = 0; i < l; i++) {
// dotProduct += vv.get(m[i], j) * z.get(j);
dotProduct += vv.get(m[i], j) * x.get(m[j]);
}
s.set(j, dotProduct);
}
} else {
for (int j = 0; j < n; j++) {
if (Thread.currentThread().isInterrupted()) {
break;
}
double dotProduct = 0.0;
for (int i = 0; i < l; i++) {
dotProduct += vv.get(m[i], j) * x.get(m[j]);
}
s.set(j, s.get(j) - dotProduct);
}
}
} else if (it == 1) {
s.assign(new Algebra().mult(vv, x));
} else {
s.assign(new Algebra().mult(vv, x), PlusMult.plusMult(-1));
}
return;
}
use of cern.colt.matrix.linalg.Algebra in project tetrad by cmu-phil.
the class IndTestFisherZGeneralizedInverse method isIndependent.
/**
* Determines whether variable x is independent of variable y given a list of conditioning variables z.
*
* @param xVar the one variable being compared.
* @param yVar the second variable being compared.
* @param z the list of conditioning variables.
* @return true iff x _||_ y | z.
* @throws RuntimeException if a matrix singularity is encountered.
*/
public boolean isIndependent(Node xVar, Node yVar, List<Node> z) {
if (z == null) {
throw new NullPointerException();
}
for (Node node : z) {
if (node == null) {
throw new NullPointerException();
}
}
int size = z.size();
int[] zCols = new int[size];
int xIndex = getVariables().indexOf(xVar);
int yIndex = getVariables().indexOf(yVar);
for (int i = 0; i < z.size(); i++) {
zCols[i] = getVariables().indexOf(z.get(i));
}
int[] zRows = new int[data.rows()];
for (int i = 0; i < data.rows(); i++) {
zRows[i] = i;
}
DoubleMatrix2D Z = data.viewSelection(zRows, zCols);
DoubleMatrix1D x = data.viewColumn(xIndex);
DoubleMatrix1D y = data.viewColumn(yIndex);
DoubleMatrix2D Zt = new Algebra().transpose(Z);
DoubleMatrix2D ZtZ = new Algebra().mult(Zt, Z);
TetradMatrix _ZtZ = new TetradMatrix(ZtZ.toArray());
TetradMatrix ginverse = _ZtZ.inverse();
DoubleMatrix2D G = new DenseDoubleMatrix2D(ginverse.toArray());
DoubleMatrix2D Zt2 = Zt.like();
Zt2.assign(Zt);
DoubleMatrix2D GZt = new Algebra().mult(G, Zt2);
DoubleMatrix1D b_x = new Algebra().mult(GZt, x);
DoubleMatrix1D b_y = new Algebra().mult(GZt, y);
DoubleMatrix1D xPred = new Algebra().mult(Z, b_x);
DoubleMatrix1D yPred = new Algebra().mult(Z, b_y);
DoubleMatrix1D xRes = xPred.copy().assign(x, Functions.minus);
DoubleMatrix1D yRes = yPred.copy().assign(y, Functions.minus);
// Note that r will be NaN if either xRes or yRes is constant.
double r = StatUtils.correlation(xRes.toArray(), yRes.toArray());
if (Double.isNaN(thresh)) {
this.thresh = cutoffGaussian();
}
if (Double.isNaN(r)) {
if (verbose) {
TetradLogger.getInstance().log("independencies", SearchLogUtils.independenceFactMsg(xVar, yVar, z, getPValue()));
}
return true;
}
if (r > 1)
r = 1;
if (r < -1)
r = -1;
this.fishersZ = Math.sqrt(sampleSize() - z.size() - 3.0) * 0.5 * (Math.log(1.0 + r) - Math.log(1.0 - r));
if (Double.isNaN(this.fishersZ)) {
throw new IllegalArgumentException("The Fisher's Z " + "score for independence fact " + xVar + " _||_ " + yVar + " | " + z + " is undefined.");
}
boolean indFisher = true;
// if(Math.abs(fishersZ) > 1.96) indFisher = false; //Two sided with alpha = 0.05
if (Math.abs(fishersZ) > thresh) {
// Two sided
indFisher = false;
}
if (verbose) {
if (indFisher) {
TetradLogger.getInstance().log("independencies", SearchLogUtils.independenceFactMsg(xVar, yVar, z, getPValue()));
} else {
TetradLogger.getInstance().log("independencies", SearchLogUtils.independenceFactMsg(xVar, yVar, z, getPValue()));
}
}
return indFisher;
}
Aggregations