[R] extend summary.lm for hccm?
John Fox
jfox at mcmaster.ca
Sun Dec 24 17:30:57 CET 2006
Dear Frank,
If I remember Freedman's recent paper correctly, he argues that sandwich
variance estimator, though problematic in general, is not problematic in the
case that White described -- an otherwise correctly specified linear model
with heteroscedasticity estimated by least-squares.
Regards,
John
--------------------------------
John Fox
Department of Sociology
McMaster University
Hamilton, Ontario
Canada L8S 4M4
905-525-9140x23604
http://socserv.mcmaster.ca/jfox
--------------------------------
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch
> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Frank
> E Harrell Jr
> Sent: Sunday, December 24, 2006 10:29 AM
> To: John Fox
> Cc: r-help at stat.math.ethz.ch; 'ivo welch'; 'Dirk Eddelbuettel'
> Subject: Re: [R] extend summary.lm for hccm?
>
> John Fox wrote:
> > Dear Dirk and Ivo,
> >
> > It's true that the sandwich package has more extensive
> facilities for
> > sandwich variance estimation than the hccm function in the car
> > package, but I think that the thrust of Ivo's question is to get a
> > model summary that automatically uses the adjusted standard errors
> > rather than having to compute them and use them "manually." Here's
> > one approach, which could be modified slightly if Ivo wants
> "hc0" as
> > the default; it could also be adapted to use the sandwich package.
> >
> > I hope this helps,
> > John
>
> Another approach:
>
> library(Design) # also requires library(Hmisc) f <- ols(y ~
> x1 + x2, x=TRUE, y=TRUE)
> f <- robcov(f) # sandwich; also allows clustering. Also see bootcov
> anova(f) # all later steps use sandwich variance matrix
> summmary(f)
> contrast(f, list(x1=.5), list(x1=.2))
>
> BUT note that sandwich covariance matrix estimators can have
> poor mean squared error (a paper by Bill Gould in Stata
> Technical Bulletin years ago related to logistic regression
> showed an example with a 100-fold increase in the variance of
> a variance estimate) and can give you the right estimate of
> the wrong quantity (reference below).
>
> Frank Harrell
>
> @Article{free06so,
> author = {Freedman, David A.},
> title = {On the so-called ``{Huber} sandwich
> estimator'' and
> ``robust standard errors''},
> journal = The American Statistician,
> year = 2006,
> volume = 60,
> pages = {299-302},
> annote = {nice summary of derivation of sandwich
> estimators;questions why we should be interested in getting
> the right variance of the wrong parameters when the model
> doesn't fit} }
>
>
> >
> > ----------- snip --------------
> >
> > summaryHCCM <- function(model, ...) UseMethod("summaryHCCM")
> >
> > summaryHCCM.lm <- function(model, type=c("hc3", "hc0",
> "hc1", "hc2",
> > "hc4"),
> >
> > ...){
> > if (!require(car)) stop("Required car package is missing.")
> > type <- match.arg(type)
> > V <- hccm(model, type=type)
> > sumry <- summary(model)
> > table <- coef(sumry)
> > table[,2] <- sqrt(diag(V))
> > table[,3] <- table[,1]/table[,2]
> > table[,4] <- 2*pt(abs(table[,3]), df.residual(model),
> lower.tail=FALSE)
> > sumry$coefficients <- table
> > print(sumry)
> > cat("Note: Heteroscedasticity-consistant standard errors using
> > adjustment",
> > type, "\n")
> > }
> >
> > --------------------------------
> > John Fox
> > Department of Sociology
> > McMaster University
> > Hamilton, Ontario
> > Canada L8S 4M4
> > 905-525-9140x23604
> > http://socserv.mcmaster.ca/jfox
> > --------------------------------
> >
> >> -----Original Message-----
> >> From: r-help-bounces at stat.math.ethz.ch
> >> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Dirk
> >> Eddelbuettel
> >> Sent: Saturday, December 23, 2006 11:16 PM
> >> To: ivo welch
> >> Cc: r-help at stat.math.ethz.ch
> >> Subject: Re: [R] extend summary.lm for hccm?
> >>
> >>
> >> On 23 December 2006 at 20:46, ivo welch wrote:
> >> | I wonder whether it is possible to extend the summary
> >> method for the
> >> | lm function, so that it uses an option "hccm" (well, model
> >> "hc0"). In
> >> | my line of work, it is pretty much required in reporting of
> >> almost all
> >> | linear regressions these days, which means that it would be
> >> very nice
> >> | not to have to manually library car, then sqrt the diagonal, and
> >> | recompute T-stats; instead, I would love to get everything
> >> in the same
> >> | format as the current output---except errors heteroskedasticity
> >> | adjusted.
> >> |
> >> | easy or hard?
> >>
> >> Did you consider the 'sandwich' package? A simple
> >>
> >> > install.packages("sandwich")
> >> > library(sandwich)
> >> > ?vcovHC
> >> > ?vcovHAC
> >>
> >> should get you there.
> >>
> >> Dirk
> >>
> >> --
>
>
>
> --
> Frank E Harrell Jr Professor and Chair School of Medicine
> Department of Biostatistics
> Vanderbilt University
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
More information about the R-help
mailing list