3

I would like to plot a model with ggplot2. I have estimated a robust variance-covariance matrix which I would like to use when estimating the confidence interval.

Can I tell ggplot2 to use my VCOV, or, alternatively, can I somehow force predict.lm to use my VCOV matrix? A dummy example:

source("http://people.su.se/~ma/clmclx.R")
df <- data.frame(x1 = rnorm(100), x2 = rnorm(100), y = rnorm(100), group = as.factor(sample(1:10, 100, replace=T))) 
lm1 <- lm(y ~ x1 + x2, data = df)
coeftest(lm1)
## outputs coef.test, but can be modified to output VCOV
clx(lm1, 1, df$group)

It would be relatively easy to add to a ggplot, if I could get 'correct' predictions given my augmented VCOV-matrix.

1 Answer 1

4

Only the standard errors, not the predictions, should change -- right?

getvcov <- function(fm,dfcw,cluster) {
  library(sandwich);library(lmtest)
  M <- length(unique(cluster))   
  N <- length(cluster)           
  K <- fm$rank                        
  dfc <- (M/(M-1))*((N-1)/(N-K))  
  uj  <- apply(estfun(fm),2, function(x) tapply(x, cluster, sum));
  dfc*sandwich(fm, meat=crossprod(uj)/N)*dfcw
}

V <- getvcov(lm1,1,df$group)
X <- as.matrix(model.frame(lm1))
se <- predict(lm1,se=TRUE)$se.fit
se_robust <- sqrt(diag(X %*% V %*% t(X)))
Sign up to request clarification or add additional context in comments.

1 Comment

This is awesome. Thanks Ben. Core R functions are awesome and powerful, but I often find myself not knowing them well enough to utilize the full potential.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.