How to do a random effect regression with clustered error in R?

Question

I want to calculate the clustered error of regression residuals. I know how to run the random effect model in R.

set.seed(50)
data <- data.frame(id=rep(1:10,10),income=rnorm(100),education=rnorm(100))
data$id <- factor(data$id)

library(nlme)
reg <- lme(income ~ education -1,data = data, random = ~ 1 | id)

I am not sure how to make the standard errors clustered at the id level as well.

Your request suggests that you expect there to be subsequent inferences about the "random" variables. That's not my understanding of the goals of mixed models. So I wouldn't expect to find such estimates in either the reg structure or in the output of `summary` applied to it. — IRTFM, Jan 29 '18 at 22:43
In addition, why do you want to both cluster SEs and have individual-level random effects? If you believe the random effects are capturing the heterogeneity in the data (which presumably you do, or you would use another model), what are you hoping to capture with the clustered errors? — Sean Norton, Jan 30 '18 at 01:38
@SeanNorton Yes, you are right. I should have thought about that. — Leonhardt Guass, Jan 30 '18 at 01:53

score 0 · Accepted Answer · answered Jan 30 '18 at 01:31

You can use plm package. See details here. But here is how to specify a random effects model:

set.seed(50)
data <- data.frame(id=rep(1:10,10),income=rnorm(100),education=rnorm(100))
data$id <- factor(data$id)

library(plm)
reg <- plm(income ~ education,
           data= data, 
           index = c('id'), 
           effect = "individual", 
           model ='random')

summary(reg)

Oneway (individual) effect Random Effect Model 
   (Swamy-Arora's transformation)

Call:
plm(formula = income ~ education, data = data, effect = "individual", 
    model = "random", index = c("id"))

Balanced Panel: n = 10, T = 10, N = 100

Effects:
                  var std.dev share
idiosyncratic 0.98161 0.99076 0.973
individual    0.02766 0.16630 0.027
theta: 0.1167

Residuals:
     Min.   1st Qu.    Median   3rd Qu.      Max. 
-2.777437 -0.626498  0.082894  0.602743  2.732686 

Coefficients:
             Estimate Std. Error t-value Pr(>|t|)
(Intercept) -0.116959   0.112068 -1.0436   0.2992
education   -0.033914   0.101471 -0.3342   0.7389

Total Sum of Squares:    95.329
Residual Sum of Squares: 95.221
R-Squared:      0.0011386
Adj. R-Squared: -0.0090539
F-statistic: 0.111707 on 1 and 98 DF, p-value: 0.73892

And to obtain clustered standard errors:

coeftest(reg, vcov=vcovHC(reg,type="HC0",cluster="group"))

That doesn't deliver what was asked, although it is possible that the questioner didn't frame his question clearly. There are no stratum-specific "standard errors" for the "random effect" — IRTFM, Jan 30 '18 at 02:38

How to do a random effect regression with clustered error in R?

1 Answers1