I'm trying to analyze some binary data that have been consolidated into probability data over time. The probability data is the dependent variable, pertaining to P(event), and I'd like analyze this as a mixed model logit or probit regression. To that end I've done some research and found incredibly helpful answers like this and this.
However, before I used this advice, I wanted to make sure I understood it. So I created a dummy dataset available here and tested how four different models turned out. To my surprise, there were some substantial differences between the models. Shouldn't these models be comparable?
Here are the results:
Coefficient p-value family link dependent_var data_format
4.2028 0.0074 Binomial logit Binary raw
*0.8101* 0.4126 Binomial logit P(Binary) consolidated
2.3332 0.0044 Binomial probit Binary raw
2.3458 0.0004 Binomial probit P(Binary) consolidated
Although I wouldn't expect logit and probit models to be identical, I would expect both logit models to be similar. It seems to be the case in the probit models!
Note that I've asterisked the coefficient for consolidated data analyzed with a logit link. It should be much closer to the raw data analyzed with the same logit model.
Here's my code:
### Read in raw data (available via link above) ###
test1 = read.csv("User/pathway/dropbox_data.csv")
### Consolidate the raw data ###
require(plyr)
tm_consolidate <- ddply(test1, .(subj, group), summarise, y=mean(y), count=length(subj))
### Run 4 models that should all be pretty similar ###
require(MASS)
logit_raw <- glmmPQL(y~group, random=~1|subj, family=binomial(link='logit'),
data=test1)
logit_cons <- glmmPQL(y~group, random=~1|subj, family=binomial(link='logit'),
data=tm_consolidate, weights=count)
probit_raw <- glmmPQL(y~group, random=~1|subj, family=binomial(link='probit'),
data=test1)
probit_cons <- glmmPQL(y~group, random=~1|subj, family=binomial(link='probit'),
data=tm_consolidate, weights=count)
Note that I've also tried running with a quasibinomial (for consolidated data) but this does not substantially change the results.
- Does anyone know why I'm finding this discrepancy?
- Shouldn't the binary data yield the same results as the P(event) data?
- If not, then why does this only apply to logit models?
- What am I missing?