Thursday, May 13, 2010

A function that reports White Standard Errors

Update: Since writing this post, I found a natural extension written by John Fox that does what the function in this post does and I have extended it in my CRAN-contributed library tonymisc. In my package, this modified and extended function is called summaryR().

One of the frustrating things about using R for econometrics is that R has some cumbersome functionality for computing the standard (in econometrics, anyway) White standard errors. As one of my favorite references on R learning points out, there are two ways to retrieve these precious standard errors.
  1. You could install and load the sandwich library. Then, apply the function vcovHC to your lmHC0. That will give you the White covariance matrix. Next, take the square root of the diagonals, and you have retrieved your standard errors.
  2. Install and load the car library. Then, apply the function hccm to your lm object, specifying the option hc0. Again, this only gives the covariance matrix, which you have to transform to obtain standard errors.
Both options produce the exact same answers (I tried it myself), but both are clunky if you want to use the White standard errors as a routine check for whether it is important. To mitigate this issue (at least for myself), I wrote a function that does the work for you.

The function is called summarize.r(), and it works just like the summary command. In fact, my function calls summary.lm, and then modifies the contents that are printed using method 1 from above.

To run my function, (1) download and install the sandwich library, (2) run the code that defines the function (see below), and (3) after fitting your regression using lm(), you can use the function by calling summarize.r(yourregression.lm).

Without further ado, here's the code you should run.

summarize.r = function(reg.lm, ...){
z = summary.lm(reg.lm)
W.VCOV = vcovHC(reg.lm,type="HC0")
sig.sq = ((summary(reg.lm)$sigma)^2)
coef = z$coefficients[,1]

White.SE <- sqrt(diag(W.VCOV)) t.robust <- coef/White.SE ## Element-by-element division df <- reg.lm$df.residual p.val.robust <- 2*(1-pt(abs(t.robust),df)) z$cov.unscaled = W.VCOV/sig.sq z$coefficients[,2] = White.SE z$coefficients[,3] = t.robust z$coefficients[,4] = p.val.robust z }

And, here's an example of putting the code to use. First, note what the summary output gives (cutting the call, residuals, and significance codes). Next, see the output from my summarize.r command.

>  summary(mathcomp.lm)

Estimate Std. Error t value Pr(>|t|)
(Intercept) 643.59 2.06 312.388 <> summarize.r(mathcomp.lm)

Estimate Std. Error t value Pr(>|t|)
(Intercept) 643.587 2.053 313.495 <>


  1. Cool. I didn't realize that R had a package for this. I'm hoping to pick it up soon.

    I thought you might enjoy the following link: It's a wiki for SAS users. Maybe one exists for R-users? I also don't know enough sas to evaluate how quality the answers are. It seems like wikis with an attached forum would be a great medium for exchanging learning about statistical software. Forums by themselves don't have the same ability to create categories that contain threads that share a common thread.

    Wikis also have the following advantage being able to write "pre" tags, or some markups, that indicate code (as shown below). Plus, the ability to have screenshots of output, etc.

    I'd be happy to point you in the right direction if you ever decide you'd like to set up a wiki for this.

  2. Check out my package "apsrtable" for handling this situation (and so much more) in output to latex.