Dear Maria,
On 08-04-2009 23:06, [email protected] wrote:
The DV of my expt is a target response (poresp) which can be either
PO (coded as 1) or DO (coded as 0). So PO is a success and DO is a
failure. The predictors are prime and nounrep (nounrepetition), each
with 2 levels. I have effect coded the predictors (-.5 and +.5
respectively for each of the two levels of the predictors, following
Barr, 2008): with this coding, which I believe is a way of centering,
the intercept should correspond to the grand mean and the coefficients
to the differences between the means etc., like in an ANOVA. Please
correct me if I am wrong. I tried this coding in the past where the
outcome was a continuous variable and indeed when I checked the means
and differences in my data , I found a correspondence between the ANOVA
measures and the regression intercept and coefficients; however, now I
am dealing with binomially distributed data, where the model parameters
are in log odds space, so it is a bit more complicated, because I have
to do transformations).
The problem is that when I check the output of the mixed logit
regression (using lmer and family = binomial ) and try to make sense of
the coefficients by checking, for example, whether the intercept
corresponds to the grand mean, in log odds, of course, things do not
match. I checked and re-checked, and I do not know where I am going
wrong.
==========
head(verbdiff)
X subject cond item myscore poresp nounrep prime primec nounrepc
1 1069 1 2 19 p 1 0 1 0.5 -0.5
2 1070 1 3 6 p 1 1 0 -0.5 0.5
3 1071 1 3 5 p 1 1 0 -0.5 0.5
4 1072 1 2 12 p 1 0 1 0.5 -0.5
5 1073 1 4 16 p 1 1 1 0.5 0.5
6 1074 1 3 14 p 1 1 0 -0.5 0.5
#primec and nounrepc are effect coded
COUNTS OF 0s and 1s:
xtabs (~ poresp +prime +nounrep, data=verbdiff)
nounrep = 0
prime
poresp 0 1
0 89 58
1 207 235
, , nounrep = 1
prime
poresp 0 1
0 91 64
1 202 228
MODEL:
verbdiff.lmer = lmer(poresp ~ primec *nounrepc + (1|subject) +
(1|item), data= verbdiff, family="binomial")
print (verbdiff.lmer)
Generalized linear mixed model fit using Laplace
Formula: poresp ~ primec * nounrepc + (1 | subject) + (1 | item)
Data: verbdiff
Family: binomial(logit link)
AIC BIC logLik deviance
1092 1122 -540 1080
Random effects:
Groups Name Variance Std.Dev.
item (Intercept) 1.2143 1.1020
subject (Intercept) 1.8470 1.3590
number of obs: 1174, groups: item, 40; subject, 32
Estimated scale (compare to 1 ) 0.917418
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.6605 0.3121 5.320 1.04e-07 ***
primec 0.7854 0.1645 4.773 1.81e-06 ***
nounrepc -0.1054 0.1612 -0.653 0.513
primec:nounrepc -0.2138 0.3224 -0.663 0.507
---
Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.?
0.1 ? ? 1
Correlation of Fixed Effects:
(Intr) primec nonrpc
primec 0.053
nounrepc -0.009 -0.023
primc:nnrpc -0.012 -0.028 0.128
=============
QUESTIONS:
1. THE ESTIMATED INTERCEPT IN THE LOGIT MODEL IS 1.66 IN LOG ODDS,
WHICH should correspond to grand mean, i.e. likelihood of successes (=PO
responses coded as "1") independent of the predictors.
THERE WERE OVERALL 872 SUCCESSES AND 302 FAILURES IN THE EXPT, SO ODDS
SHOULD BE 872/302=2.88 or (in probability space) .74/.26 = 2.85;
THIS SHOULD GIVE A LOG OF ODDS OF APPROX 1.05, BUT THE INTERCEPT
PREDICTED BY THE MODEL IS MUCH HIGHER (1.66)
Indeed:
> logit(872/(872+302))
[1] 1.060362
This may be because with your coding, you are estimating the
intercept if primec==0 & nounrepc==0. This means that if prime and
DV are related (e.g. because there are more hits if prime==1), or
likewise if noun and DV are related, then the estimated intercept is
no longer the grand mean, but rather the estimated value of your DV
if all predictors would be zero.
SAME THING WHEN I USE THE CORRESPONDING DUMMY CODED (0-1) V
ARIABLES
-INTERCEPT IS MUCH HIGHER THAN WHAT MY DATA SUGGEST
In my experience, models are far easier to interpret if you do NOT
centralize the binary dummy predictors. The estimated intercept then
equals the mean logit for the particular cell where all dummies are
zero. So, my advice is to use your dummy coded variables and not the
centered dummies.
The main point I guess is that the intercept does not reflect the
grand mean, in your model, but it reflects the estimated DV if all
predictors would be zero, which is different from the grand mean IF
predictors and DV are related.
Did you try this model...
>> verbdiff.m0.lmer = lmer(poresp ~ 1+(1|subject) +
> (1|item), data= verbdiff, family="binomial")
It should give you ONLY the intercept which should then be about 1.06.
More background and details are available at
http://www.hugoquene.nl/mixedeffects/
and references in the articles mentioned there.
Hope this helps, with kind regards, Hugo Quene
--
Dr Hugo Quené | assoc prof Phonetics | Utrecht inst of Linguistics
OTS | Utrecht University | Trans 10 | 3512 JK Utrecht | The
Netherlands | T +31 30 253 6070 | F +31 30 253 6000 | [email protected]
| www.hugoquene.nl | www.hum.uu.nl
_______________________________________________
R-lang mailing list
[email protected]
http://pidgin.ucsd.edu/mailman/listinfo/r-lang
