How to predict random and fixed effect models?

Question

I just recently made a change from STATA to R and have some troubles implementing the R equivalent of the STATA commands xtlogit,fe or reand predict. May I ask for some assistance to adjust the following scenario:

  data <- read.table("http://people.stern.nyu.edu/wgreene/Econometrics/healthcare.csv",header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE)

   require(caret) # for confusionMatrix

   #### subset into test & train according to the panel nature (split  individuals rather then observations)
   nID <- length(unique(data$id))
   p = 0.50# partition

   inTrain <- sample(unique(data$id), round(nID * p), replace=FALSE)

   training <- data[data$id %in% inTrain, ] 

   testing <- data[!data$id %in% inTrain, ] 


   pooled <- glm(WORKING~WHITEC+FEMALE+BLUEC+HHNINC+AGE+AGESQ+EDUC+DOCVIS,data=training, family=binomial(link="logit"))

   prediction.working= round(predict(pooled,newdata=testing,type="response"))

   confusionMatrix(prediction.working,testing$WORKING) # Accuracy between both

Additionally, I would like to do these procedure for random effects and fixed effects. So I tried random effects first unsuccessfully:

   library(glmmML)
   RE <- glmmML(WORKING~WHITEC+FEMALE+BLUEC+HHNINC+AGE+AGESQ+EDUC+DOCVIS, family=binomial(link="logit"), data=training, cluster=id, method="ghq", n.points=12)



    prediction.working= round(predict(RE,newdata=testing,type="response"))

But that doesnt seem to work. May I ask for how to adjust the glmmodel regarding random effects and fixed effects in order to use the predictfunction.

Answer 1

Welcome to R. I am also a STATA convert.

This is a tricky problem, but its answer is essential to understand. To understand why the predict function does not work for glmmML you need to understand S3 methods (see http://adv-r.had.co.nz/OO-essentials.html).

Let me explain a bit. R is an object oriented language. This means that everything in R is an object (i.e. vector, function, data.frame). Each object contain attributes. An attribute is essentially meta-data about the object itself. For example the names of variables in a data.frame are attributes. One attribute all objects have is a class. To see the class of any object just call the class() function.

A lot, but not all, functions use S3 methods. The predict function is one of these function. This means when you call predict, the predict function looks at the class of the object. Then depending on the class picks which other predict function should be used. For example, if your object is class lm, the predict function will call the predict.lm function. Other predict functions include: predict.glm for objects of glm class, predict.loess for objects of loess class, predict.nls for objects of nls class, etc. (to see the full list read the predict help). Unfortunately, no predict.glmmML function exists. Consequently, when you call the predict function on a object of class glmmML you get an error.

id <- factor(rep(1:20, rep(5, 20)))
y <- rbinom(100, prob = rep(runif(20), rep(5, 20)), size = 1)
x <- rnorm(100)
dat <- data.frame(y = y, x = x, id = id)
fit.2 <- glmmML(y ~ x, data = dat, cluster = id)
predict(fit.2)
Error in UseMethod("predict") : 
  no applicable method for 'predict' applied to an object of class "glmmML"
class(fit.2)
 [1] "glmmML"

The error is very informative. It basically says that R tried to use S3 methods, however, there was no 'predict.glmmML'

What about the mclogit function which user227710 suggests. Let's see

data(Transport)

fit <- mclogit(
  cbind(resp,suburb)~distance+cost,
  data=Transport
)

class(fit)
[1] "mclogit" "lm"

The class of fit is mclogit and lm. Will predict work? Yes! When you call predict(fit) the predict function will first look for a predict.mclogit, which does not exist. It will next look for predict.lm. Which does exist.