Why logLik works for PLM models in R?

Question

I'm estimating linear panel models with fixed effects in R using plm of the plm-package. As I was looking for ways to compare different models, I came across AIC (Akaike's ‘An Information Criterion’). AIC of the stats-package works according to the documentation for model objects for which a log-likelihood value can be obtained. According to the comments given the OLS-estimator of a linear model is the maximum likelihood estimator. But is this also true for fixed-effects panel models (within-estimator)? Here is a minimal example:

data("Produc", package = "plm")
zz <- plm::plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp, model="pooling",
          data = Produc, index = c("state","year"))
yy <- plm::plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp, model="within",
          data = Produc, index = c("state","year"))
stats::logLik(zz)
stats::logLik(yy)
stats::AIC(zz)
stats::AIC(yy) 

Thanks a lot for your help.

Answer 1

If anyone is looking to actually calculate the logLik for a plm object, consider this, which comes from a Stack Overflow question on the same topic here.

object here is a plm object, like zz above.

logLik.plm <- function(object){
  out <- -plm::nobs(object) * log(2 * var(object$residuals) * pi)/2 - deviance(object)/(2 * var(object$residuals))
  
  attr(out,"df") <- nobs(object) - object$df.residual
  attr(out,"nobs") <- plm::nobs(object)
  return(out)
}

Once this is defined, this should then also work with the AIC command:

data("Produc", package = "plm")
zz <- plm::plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp, model="pooling",
               data = Produc, index = c("state","year"))
yy <- plm::plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp, model="within",
               data = Produc, index = c("state","year"))


logLik.plm <- function(object){
  out <- -plm::nobs(object) * log(2 * var(object$residuals) * pi)/2 - deviance(object)/(2 * var(object$residuals))
  
  attr(out,"df") <- nobs(object) - object$df.residual
  attr(out,"nobs") <- plm::nobs(summary(object))
  return(out)
}

stats::logLik(zz)
[1] 826.9814
attr(,"df")
[1] 5
attr(,"nobs")
[1] 816

stats::logLik(yy)
[1] 1534.531
attr(,"df")
[1] 52
attr(,"nobs")
[1] 816

stats::AIC(zz)
[1] -1643.963

stats::AIC(yy) 
[1] -2965.063

I have created a public gist here.