The statistical tools in Excel have always been black boxes. There's nothing for it but to do some forensic reverse-engineering. By performing a simple regression in Excel 2013, involving the data $x=(1,2,3,4,5,6,7,8,9)$ and $y=(2,1,4,3,6,5,9,8,7)$, and requesting "standardized residuals" in the dialog, I obtained output that states
The "Standard Error" is $1.3723\ldots$.
There are $9$ observations.
The residuals $r_i$ are $(0.5333\ldots, -1.35, \ldots, 0.35, -1.533\ldots)$.
The corresponding "Standard Residuals" are $(0.4154\ldots, -1.0516\ldots, \ldots, 0.2726\ldots, -1.1944\ldots)$.
Since "standardized" values are typically numbers divided by some estimate of their standard error, I compared these "Standard Residuals" to the residuals and to the "Standard Error." Knowing that various formulas for variances are sums of squares of residuals $r_i$ divided variously by $n$ (the number of data) or $n-p$ (the number of data reduced by the number of variables, in this case two: one for the intercept and a second for the slope), I squared everything in sight. It became immediately obvious that Excel is computing the "Standard Residual" as
$$\frac{r_i}{\sqrt{\frac{1}{n-1}\sum_{i=1}^n r_i^2}}.$$
This formula reproduced Excel's output exactly--not even a trace of floating point roundoff error.
The denominator is what would be computed by Excel's STDEV
function. For residuals from a mean, it is an unbiased estimate of their variance. For residuals in a regression, however, it has no standard meaning or value. It's garbage! But now you know how to compute it... .
Syre, you say about the linear regression
A linear regression residual close to zero means that the model is a
good fit for the observed value. A negative residual means that the
model overestimates the effect of the independent variables in that
particular case.
and I think this is where the misunderstanding starts - a linear regression where you have all residuals close to zero (close by units of the standard deviation of the regression) is actually NOT a good fit. In a perfectly fitting linear regression, you assume that residuals scatter around the mean predicted value with a normal distribution. Hence, you completely expect that some values are higher and some are lower. This is not an overestimation of the effect, but a requirement of the model.
The goal of the residual checks for the linear regression is thus not to see if residuals are close to zero, but if they scatter normally distributed around zero!
The same is true for DHARMa residuals. The only difference is that the expected distribution is uniform, not normal. I quote from the vignette:
As discussed above, for a correctly specified model we would expect
So, interpretation of the residuals is really like in a linear regression, only that the distribution is uniform, and that the mean expectation is at 0.5.
Addition in response to the question below:
Yes, you could look at patterns in the DHARMa residuals and attempt an interpretation of why they occur, in the same way as you might do this in a linear regression.
Note that the quote in the paper assumes the most simple linear regression, where a point that is further away from the regression line is also less likely. If you include the possibility in the model that the variance of the residuals changes (e.g. in a gls), such an interpretation of raw residuals doesn't make sense any more to define outliers or especially interesting points. The most basic solution is to divide residuals by expected variance (= Pearson residuals). The quantile residuals in DHARMa generalize this idea.
A special property of the quantile residuals is that you compare against a simulated distribution. In DHARMa, I call 0 / 1 outliers, because they are outside the simulation range. What's different compared to normal outliers is that we know they are outside, but you don't know HOW FAR they are outside (you get a value of zero, if the observed value is smaller than all simulations, regardless of how much smaller). That's why this type of outliers are extra highlighted in DHARMa.
Best Answer
One way to try to understand this is to try to make a simple example. But first, let us clear up your doubt in
One does exactly as is done in parametric bootstrapping, assuming that the fitted model is the true one, replacing the unknown parameters with the fitted values. So, with a normal (Gaussian) family, we assume normal errors, with a Poisson family we simulate response values from a Poisson distribution, and so on.
A simple example, with code in R. Let us first simulate some data:
To repeat, this is really parametric bootstrapping!
Then we can program the tree steps of the algorithm you reference as a simple R function, which will work for linear models (
lm
) andglm
(generalized linear models) objects. TheDHARMa
package must be more complicated, as it works for many other model types. Making a minimal code can be a good way to understand an algorithm, as usual working code as that in theDHARMa
package might be difficult to read, as such code is often dominated by error checks, guards and strange border cases. No such distractions below:Hopefully, this code should clarify your doubts:
The loop in the code goes through the dataset, for each observation $y_i, i=1, 2, \dotsc, N$ it takes the $B$ simulated (bootstrapped) responses, stored in row $i$ in object
ystars
, calculates its ecdf (empirical cumulative distribution function), and then applies this function to the actual observation $y_i$.The cdf is created using the simulated responses.