Solved – Computing the Standard Error of the Estimate from the ANOVA table

rregression

My question is quite straightforward, but I did not find a clear answer anywhere.

I'm computing the Standard Error of the Estimate (SEE) by doing the square root of the Residuals Mean Square output of the anova table:

anovatable<-anova(lm(carb~hp,data=mtcars))

anovatable

Analysis of Variance Table

Response: carb

          Df Sum Sq Mean Sq F value    Pr(>F)    
hp         1 45.469  45.469  38.527 7.828e-07 ***
Residuals 30 35.406   1.180                      

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

SEE<-sqrt(anovatable$`Mean Sq`[2])
SEE
[1] 1.086363    

Is it the correct way of doing it?

Is there any already implemented way in R to obtain the SEE? If so, It will be better than accessing the Mean Sq term for Residuals, since its position depends upon the number of predictors.

Best Answer

The standard error of the estimate (SEE) is the following where $SSE$ is the sum of squares of the ordinary residuals (this sum of squares is also called the deviance) and $n$ is the number of observations and $k$ is the number of coefficients in the model. (The intercept counts as a coefficient so $k=2$ in the case of the example shown in the question.)

$ \sqrt{SSE / (n - k)} $

In R, it can also be calculated from a model object using the sigma function so any of these work (assuming no NA's):

fm <- lm(carb ~ hp, data = mtcars)

sigma(fm)
## [1] 1.086363

sqrt(sum(resid(fm)^2) / (nrow(mtcars) - 2))
## [1] 1.086363

sqrt(deviance(fm) / (nobs(fm) - length(coef(fm))))
## [1] 1.086363

summary(fm)$sigma
## [1] 1.086363

sqrt(anova(fm)["Residuals", "Mean Sq"])
## [1] 1.086363

If what you meant was the standard errors of the coefficient estimates then there would be one for each coefficient and those standard errors would be any of the following where the last one makes use of an estimate of $var(\hat{\beta})$ being $\hat{\sigma}^2 (X'X)^{-1}$

sqrt(diag(vcov(fm)))
## (Intercept)          hp 
## 0.459500176 0.002845806 

coef(summary(fm))[, "Std. Error"]
## (Intercept)          hp 
## 0.459500176 0.002845806 

sigma(fm) * sqrt(diag(solve(crossprod(model.matrix(fm)))))
## (Intercept)          hp 
## 0.459500176 0.002845806