Olivier Blaser <olivier.blaser <at> unil.ch> writes:
> I am having trouble with lmer. I am looking at recombinant versus non
> recombinant individuals. In the response variable recombinant
> individuals are coded as 1's and non-recombinant as 0's. I built a model
> with 2 fixed factors and 1 random effect. Sex (males/females) is the
> first fixed effect and sexual genotype (XY, YY, WX and WY) the second
> one. Sexual Genotype is nested in sex. XY and YY individuals are males
> and WX and WY females. I crossed 8 XY males with 8 WY females and 8 YY
> males with 8 WX females. Each cross corresponds to a family (i.e. family
> 1 to 8 for XY x WY and family 9 to 16 for YY WX). For each family I have
> 20 offspring. Family is nested in sexual genotype as a random factor.
>
> My data looks like:
>
reps<-factor(sample(c(0,1),640,replace=T,prob=c(0.95,0.05)))
sex<-factor(rep(c("M","F"),each=320))
geno<-factor(rep(c("XY","YY","WY","WX"),each=160))
fam<-factor(rep(c(1:8,9:16,1:8,9:16),each=20))
dat<-data.frame(reps,sex,geno,fam)
with(dat,table(sex,geno))
with(dat,table(sex,geno,fam))
Follow-ups to this question should go to r-sig-mixed-models
The slightly obscure error message is caused by a glitch
in lme4's printing code, where it tries to print a model
that hasn't been successfully fitted in the first place.
I actually get a **slightly** more informative error message
(maybe with a slightly later version of the package):
g1 <- glmer(reps~sex+geno+(1|fam),data=dat,family=binomial)
Error in mer_finalize(ans) : Downdated X'X is not positive definite, 1.
[although lmer() automatically calls glmer() for GLMMs, I
prefer to use glmer() explicitly]
To diagnose/show what the problem is, let's try fitting in glm.
It works, but the answer indicates that
your design matrix is rank-deficient/overfitted:
(g0 <- glm(reps~sex+geno,data=dat,family=binomial))
i.e., the 'genoYY' estimate is NA, because once 'sexM',
'genoWY', and 'genoXY' parameters are fitted the model is
completely determined.
Unlike [g]lm, [g]lmer can't handle rank-deficient/overfitted
fixed-effect design matrices (as I understand it, because the
somewhat fancier/more efficient linear algebra primitives used
in lme4 don't handle this case). Therefore, you have to figure
out a way to code your design matrix that will work. At the moment
I don't know of an elegant way to do that, but here's a hack:
## extract model matrix and drop intercept and genoYY dummy variables
mm <- model.matrix(g0)[,2:4]
dat2 <-data.frame(reps,mm,fam)
g2 <- glmer(reps~sexM+genoWY+genoXY+(1|fam),data=dat2,family=binomial)
I don't know if this gives you the answer you want or not ...
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
