Calculating R-Squared Using the 'Residuals' Band from Google Earth Engine's linearRegression() Reducer

Question

I am trying to calculate an R-Squared for a harmonic linear regression I did with two dependent variables across an image collection, using the linearRegression reducer on Google Earth Engine. The reduction produces an image with two bands, 'coefficients' and 'residuals'. I can get the coefficients as bands by projecting it as an array, and then splitting the array by the number of coefficients. However, I am really struggling to calculate the R-Squared by the same method.  I've consulted a stack exchange message that says that I can do it, but does not explain how: 

Calculating residual sum of squares and standard error from linearFit model

And there is this Google Earth Engine help file, but that also does not how to do it with an image collection: 

https://developers.google.com/earth-engine/guides/reducers_regression

I would post the file, but GEE is not allowing me to do so, and I do not know how to do that either. You can see the code for it here:

https://code.earthengine.google.com/accept_repo=users/elijahtangenberg/Shareable

Answer 1

Thought I would add the answer to my question here for the next person who has this problem. GEE developer Nick Clinton presented his solution to it in the Google Earth Developer Forum (https://groups.google.com/g/google-earth-engine-developers) back in 2017. You can see his post here:

https://groups.google.com/g/google-earth-engine-developers/c/TrKSNicCT_0/m/qLMRjt0UGwAJ

and the example here:

https://code.earthengine.google.com/00ba3bb9daea70f0b403b230050c6aa0

In essence, you extract the 'residuals' band from the reduced image, which has the root mean squared residuals, the same way you would extract the coefficents, through an array, which you then flatten. After it is flattened, you can apply the R-Squared equation. You can see the code for the procedure here:

// The output of the regression reduction is a 4x1 array image.
var linearTrend = imageCollection
  .select(Independents.add(dependent))
  .reduce(ee.Reducer.linearRegression(Independents.length(), 1));

// Compute the number of observations in each pixel.
var n = imageCollection.select(dependent).count();

// There are n-p degrees of freedom
var dof = n.subtract(Independents.length());

// Get the root mean square residuals in a pixel.
var rmsr = linearTrend.select('residuals')
    .arrayProject([0])
    .arrayFlatten([['rmsr']]);

// Residual sum of squares.
var rss = rmsr.pow(2).multiply(n);

// Estimated variance = RSS/(n-p).
var sSquared = rss.divide(dof);

// Sample variance of the dependent variable.
var yVariance = imageCollection.select(dependent)
    .reduce(ee.Reducer.sampleVariance());

// Adjusted, squared, multiple regression coefficient.
var rSquareAdj = ee.Image(1).subtract(sSquared.divide(yVariance));

Map.addLayer(rSquareAdj, {min: 0, max: 1}, 'R squared');}