use of ubic.basecode.math.linearmodels.LeastSquaresFit in project Gemma by PavlidisLab.
the class ComBat method gammaHat.
private void gammaHat(DoubleMatrix2D sdata) {
DoubleMatrix2D Xb = x.viewPart(0, 0, x.rows(), numBatches);
gammaHat = new LeastSquaresFit(Xb, sdata).getCoefficients();
}
use of ubic.basecode.math.linearmodels.LeastSquaresFit in project Gemma by PavlidisLab.
the class ComBat method standardize.
/**
* Special standardization: partial regression of covariates
*
* @param b b
* @param A A
* @return double matrix 2d
*/
DoubleMatrix2D standardize(DoubleMatrix2D b, DoubleMatrix2D A) {
DoubleMatrix2D beta = new LeastSquaresFit(A, b).getCoefficients();
// assertEquals( 3.7805, beta.get( 0, 0 ), 0.001 );
// assertEquals( 0.0541, beta.get( 2, 18 ), 0.001 );
int batchIndex = 0;
DoubleMatrix2D bba = new DenseDoubleMatrix2D(1, numBatches);
for (String batchId : batches.keySet()) {
bba.set(0, batchIndex++, (double) batches.get(batchId).size() / numSamples);
}
/*
* Weight the non-batch coefficients by the batch sizes.
*/
DoubleMatrix2D grandMeanM = solver.mult(bba, beta.viewPart(0, 0, numBatches, beta.columns()));
if (hasMissing) {
varpooled = y.copy().assign(solver.transpose(solver.mult(x, beta)), Functions.minus);
DoubleMatrix2D var = new DenseDoubleMatrix2D(varpooled.rows(), 1);
for (int i = 0; i < varpooled.rows(); i++) {
DoubleMatrix1D row = varpooled.viewRow(i);
double m = DescriptiveWithMissing.mean(new DoubleArrayList(row.toArray()));
double v = DescriptiveWithMissing.sampleVariance(new DoubleArrayList(row.toArray()), m);
var.set(i, 0, v);
}
varpooled = var;
} else {
varpooled = y.copy().assign(solver.transpose(solver.mult(x, beta)), Functions.minus).assign(Functions.pow(2));
DoubleMatrix2D scale = new DenseDoubleMatrix2D(numSamples, 1);
scale.assign(1.0 / numSamples);
varpooled = solver.mult(varpooled, scale);
}
DoubleMatrix2D size = new DenseDoubleMatrix2D(numSamples, 1);
size.assign(1.0);
/*
* The coefficients repeated for each sample.
*/
standMean = solver.mult(solver.transpose(grandMeanM), solver.transpose(size));
/*
* Erase the batch factors from a copy of the design matrix
*/
DoubleMatrix2D tmpX = x.copy();
for (batchIndex = 0; batchIndex < numBatches; batchIndex++) {
for (int j = 0; j < x.rows(); j++) {
tmpX.set(j, batchIndex, 0.0);
}
}
/*
* row means, adjusted "per group", and ignoring batch effects.
*/
standMean = standMean.assign(solver.transpose(solver.mult(tmpX, beta)), Functions.plus);
DoubleMatrix2D varsq = solver.mult(varpooled.copy().assign(Functions.sqrt), solver.transpose(size));
/*
* Subtract the mean and divide by the standard deviations.
*/
DoubleMatrix2D meansubtracted = y.copy().assign(standMean, Functions.minus);
return meansubtracted.assign(varsq, Functions.div);
}
Aggregations