Matrix common = X.transposeTimesTranspose(F).times(F);
Matrix b = common.times(X).inverse().times(common.times(Y));
// Estimate sigma_0 and sigma:
// sigma_sum_square = sigma_0*sigma_0 + sigma*sigma
Matrix sigmaMat = F.times(X.times(b).minus(F.times(Y)));
final double sigma_sum_square = sigmaMat.normF() / (relationx.size() - 6 - 1);
final double norm = 1 / Math.sqrt(sigma_sum_square);
// calculate the absolute values of standard residuals
Matrix E = F.times(Y.minus(X.times(b))).timesEquals(norm);