3

I'm using plm package to estimate a random effects model on panel data. Reading this question about the prediction in plm package gives me some doubts. How exactly it works? I tried 3 alternatives methods and they give different solutions. Why ?

library(data.table); library(plm)

set.seed(100)
DT <- data.table(CJ(id=c(1,2,3,4), time=c(1:10)))
DT[, x1:=rnorm(40)]
DT[, x2:=rnorm(40)]
DT[, y:=x1 + 2*x2 + rnorm(40)/10 + id]
DT <- DT[!(id=="a" & time==4)] # just to make it an unbalanced panel
setkey(DT, id, time)    

summary(plmFEit <- plm(data=DT, id=c("id","time"), formula=y ~ x1 + x2, model="random"))    
    ###################
    #method 1
    ###################
    # Extract the fitted values from the plm object
    FV <- data.table(plmFEit$model, residuals=as.numeric(plmFEit$residuals))
    FV[, y := as.numeric(y)]
    FV[, x1 := as.numeric(x1)]
    FV[, x2 := as.numeric(x2)]

    DT <- merge(x=DT, y=FV, by=c("y","x1","x2"), all=TRUE)
    DT[, fitted.plm_1 := as.numeric(y) - as.numeric(residuals)]                
    ###################
    #method 2
    ###################        
    # calculate the fitted values 
    DT[, fitted.plm_2 := as.numeric(coef(plmFEit)[1]+coef(plmFEit)[2] * x1 + coef(plmFEit)[3]*x2)]                
    ###################
    #method 3
    ###################
    # using pmodel.response 
    DT$fitted.plm_3 <-pmodel.response(plmFEit,model='random') 
Community
  • 1
  • 1
dax90
  • 1,088
  • 14
  • 29

1 Answers1

2

Method 1: This can be done easier:

FV <- data.table(plmFEit$model, residuals=as.numeric(plmFEit$residuals))
FV[ , fitted := y - residuals]

However, it gives the same fitted values like in you approach with merge()

Method 2: You are fitting a random effects model (albeit you name it plmFEit, where FE usually implies fixed effects). Compared to method 1, you are missing the idiosyncratic error term.

Method 3: pmodel.response() gives you the response variable (in this case that would be y) but with the specified transformation (random effects transformation ("quasi demeaning")) applied to it (see pmodel.response()).

I think, what you want is given by method 1.

Helix123
  • 3,502
  • 2
  • 16
  • 36