[R] heteroskedasticity-robust standard errors

Ott Toomet siim at obs.ee
Fri Mar 22 10:02:38 CET 2002


Hi,

I attach a version of small program for White stderrors.  Basically, it
calculates a new covariance matrix and prints the results with new stderrors
(without any tests and additional information).  It should be rewritten to
use the object-oriented architecture, but I have not had need so long.

you should call it as:
summaryw(lm(y~x))

Regards,
Ott Toomet

P.S.  If you get/write a better version, let me know.

--------------------------------------------------------------------
### ols with White' heteroscedasticity consistent stderrors
summaryw <- function( model) {
  s <- summary( model)
  X <- model.matrix( model)
  u2 <- residuals( model)^2
  XDX <- 0
  ## here one needs essentially to calculate X'DX.  But due to the fact that D
  ## is huge (NxN), it is better to do it with a cycle.
  for( i in 1:nrow( X)) {
    XDX <- XDX + u2[i]*X[i,]%*%t( X[i,])
  }
  XX1 <- solve( t( X)%*%X)
  varcovar <- XX1 %*% XDX %*% XX1
  stdh <- sqrt( diag( varcovar))
  t <- model$coefficients/stdh
  p <- 2*pnorm( -abs( t))
  results <- cbind( model$coefficients, stdh, t, p)
  dimnames(results) <- dimnames( s$coefficients)
  results
}




On Thu, 21 Mar 2002, Grant Farnsworth wrote:

  |I am trying to compute the white heteroskedasticity-robust standard errors
  |(also called the Huber standard errors) in a linear model, but I can't seem
  |to find a function to do it.  I know that the design library in S+ has
  |something like this (robcov?), but I have not yet seen this library ported
  |to R.
  |
  |Anyone know if there is already a function built into R to do this
  |relatively simple job?
  |
  |Thanks,
  |Grant

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._



More information about the R-help mailing list