[R-sig-ME] glmmadmb- problems with explanatory variable with a lot of zeros?
Mabille, Geraldine
Geraldine.Mabille at nina.no
Tue Oct 15 14:33:32 CEST 2013
Thanks heaps to Ben Bolker and Dave Fournier for their answers and suggestions!
I have tried several of the suggestions made by Ben Bolker and here is what I get:
1) Using the zipme function, I get the following error message
zipme(cformula=Killed_Lamb~Dens_Lynx+(1|KOM),zformula=z~1,data=DATA,maxitr=20,tol=1e-6,verbose=FALSE,cfamily=poisson)
Error in lmerFrames(mc, formula, contrasts) : negative weights or weights of zero are not allowed
2) Trying to add a binary variable Lynx_Bin into the model:
DATA$Lynx_Bin<-ifelse(DATA$Dens_Lynx>0,1,0)
Mod1a <-glmmadmb(Killed_Lamb~I(scale(Dens_Lynx))+Lynx_Bin, random=~ 1|KOM, data=DATA,zeroInflation=TRUE,family="nbinom")
I get the following error messages:
Parameters were estimated, but not standard errors were not: the most likely problem is that the curvature at MLE was zero or negative
Error in glmmadmb(Killed_Lamb ~ I(scale(Dens_Lynx)) + Lynx_Bin, random = ~1 | :
The function maximizer failed (couldn't find STD file) Troubleshooting steps include (1) run with 'save.dir' set and inspect output files; (2) change run parameters: see '?admbControl'
In addition: Warning message:
running command 'C:\Windows\system32\cmd.exe /c "C:/Myname/R-2.15.2/library/glmmADMB/bin/windows64/glmmadmb.exe" -maxfn 500 -maxph 5 -noinit -shess' had status 1
3) Trying to fit Dens_Lynx as a spline:
library(splines)
Mod1b <-glmmadmb(Killed_Lamb~ns(Dens_Lynx,3), random=~ 1|KOM, data=DATA,zeroInflation=TRUE,family="nbinom")
Warning message:
In glmmadmb(Killed_Lamb ~ ns(Dens_Lynx, 3), random = ~1 | KOM, data = DATA, :
Convergence failed:log-likelihood of gradient= -510.054
summary(Mod1b)
Call:
glmmadmb(formula = Killed_Lamb ~ ns(Dens_Lynx, 3), data = DATA,
family = "nbinom", random = ~1 | KOM, zeroInflation = TRUE)
AIC: 29691.4
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 7.13e+12 1.09e+14 0.07 0.95
ns(Dens_Lynx, 3)1 1.31e+13 1.82e+14 0.07 0.94
ns(Dens_Lynx, 3)2 -2.83e+13 2.17e+14 -0.13 0.90
ns(Dens_Lynx, 3)3 2.79e+11 6.93e+13 0.00 1.00
Number of observations: total=2151, KOM=257
Random effect variance(s):
Group=KOM
Variance StdDev
(Intercept) 9.115 3.019
Negative binomial dispersion parameter: 0.90085 (std. err.: 0.076799)
Zero-inflation: 0.034505 (std. err.: 0.01553 )
Log-likelihood: -14838.7
4)Trying to specify my starting values (I'm not sure I understood Ben's suggestion properly here: specify starting values for the fixed effects only? And not sure either if I implemented that properly). Here is what I tried
Mod1c<-glmmadmb(Killed_Lamb~I(scale(Dens_Lynx)), random=~ 1|KOM, data=DATA,zeroInflation=TRUE,family="nbinom",start=list(fixed=c(Mod0$b,0)))
summary(Mod1c) ####### IS THE SAME AS summary(Mod1)
Call:
glmmadmb(formula = Killed_Lamb ~ I(scale(Dens_Lynx)), data = DATA,
family = "nbinom", start = list(fixed = c(Mod0$b, 0)), random = ~1 |
KOM, zeroInflation = TRUE)
AIC: 19722.8
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.9544 0.1374 21.5 < 2e-16 ***
I(scale(Dens_Lynx)) 0.1122 0.0229 4.9 9.5e-07 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Number of observations: total=2151, KOM=257
Random effect variance(s):
Group=KOM
Variance StdDev
(Intercept) 4.053 2.013
Negative binomial dispersion parameter: 1.7665 (std. err.: 0.075737)
Zero-inflation: 0.090062 (std. err.: 0.0083659 )
Log-likelihood: -9856.38
Mod0c<-glmmadmb(Killed_Lamb~1, random=~ 1|KOM, data=DATA,zeroInflation=TRUE,family="nbinom",start=list(fixed=Mod1$b[1]))
Parameters were estimated, but not standard errors were not: the most likely problem is that the curvature at MLE was zero or negative
Error in glmmadmb(Killed_Lamb ~ 1, random = ~1 | KOM, data = DATA, zeroInflation = TRUE, :
The function maximizer failed (couldn't find STD file) Troubleshooting steps include (1) run with 'save.dir' set and inspect output files; (2) change run parameters: see '?admbControl'
In addition: Warning message: running command
'C:\Windows\system32\cmd.exe /c "C:/Myname/R-2.15.2/library/glmmADMB/bin/windows64/glmmadmb.exe" -maxfn 500 -maxph 5 -noinit -shess' had status 1
5) Finally, I also tried fitting the same models as Mod0 and Mod1, but without zero inflation. This time the log likelihood behaves normally and AIC decreases when I add the Dens_Lynx variable in the model (Mod1d compared to Mod0d). However, if I test for the existence of zero-inflation by comparing AIC of Mod1 with AIC of Mod1d (or Mod0 with Mod0d) , I see that there seems to be strong evidence for existence of zero-inflation. Should I just forget about zero-inflation and fit the model with glmmadmb and no zero-inflation? I don't know how to interpret values for zero-inflation parameter we get from Mod0 and Mod1. Is there also strong evidence for zero-inflation according to those parameters?
Mod0d <-glmmadmb(Killed_Lamb~1, random=~ 1|KOM, data=DATA,family="nbinom")
Call: glmmadmb(formula = Killed_Lamb ~ 1, data = DATA, family = "nbinom", random = ~1 | KOM, save.dir = "M:/Geraldine/Sheep/Results/Folder")
AIC: 19888.8
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.900 0.135 21.4 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Number of observations: total=2151, KOM=257
Random effect variance(s):
Group=KOM
Variance StdDev
(Intercept) 4.187 2.046
Negative binomial dispersion parameter: 1.0361 (std. err.: 0.038589)
Log-likelihood: -9941.4
Mod1d <-glmmadmb(Killed_Lamb~I(scale(Dens_Lynx)), random=~ 1|KOM, data=DATA,family="nbinom")
summary(Mod1d)
Call: glmmadmb(formula = Killed_Lamb ~ I(scale(Dens_Lynx)), data = DATA, family = "nbinom", random = ~1 | KOM)
AIC: 19871.8
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.9085 0.1326 21.93 < 2e-16 ***
I(scale(Dens_Lynx)) 0.1259 0.0292 4.32 1.6e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Number of observations: total=2151, KOM=257
Random effect variance(s):
Group=KOM
Variance StdDev
(Intercept) 4.009 2.002
Negative binomial dispersion parameter: 1.0432 (std. err.: 0.038955)
Log-likelihood: -9931.89
Thanks again with any help on those issues,
Cheers,
Geraldine
-----Original Message-----
From: r-sig-mixed-models-bounces at r-project.org [mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of Ben Bolker
Sent: 11. oktober 2013 22:30
To: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] glmmadmb- problems with explanatory variable with a lot of zeros?
Mabille, Geraldine <Geraldine.Mabille at ...> writes:
> Edit...I have just understood the "problem" with the residuals that
> were different even though they had the same number of Lamb Killed and
> the same Dens_Lynx. It is only that they come from different
> municipalities! The rest of the problem, with log likelihood behaving
> weirdly is still their though...
[snip]
> I'm just starting to use glmmadmb to try to model the number of lambs
> killed by lynx (Killed_Lamb) in each Norwegian municipality. We have
> repeated data over 11 years so I use municipality as a random effect
> in the analysis. The data contains a lot of zeros (403 over a total of
> 2151 lines) and evidence for overdispersion (tested using the
> overdisp_fun() after fitting of a poisson model in lme4, ratio
> obtained=23.8). I therefore decided to try using a negative binomial
> distribution with zero inflation.
> I first modeled the number of killed lambs as a function of the
> density of lynx (Dens_Lynx) in the municipality. We also have a lot of
> zero in Dens_Lynx (1518 over a total of 2151 lines).
> Mod1 <-glmmadmb(Killed_Lamb ~I(scale(Dens_Lynx)), random=~ 1|KOM,
> data=DATA,zeroInflation=TRUE,family="nbinom")
> I obtain the following summary for Mod 1:
> Call:
> glmmadmb(formula = Killed_Lamb ~ I(scale(Dens_Lynx)), data = DATA,
> family = "nbinom", random = ~1 | KOM, zeroInflation = TRUE)
[snip]
> Number of observations: total=2151, KOM=257
> Random effect variance(s):
> Group=KOM
> Variance StdDev
> (Intercept) 4.053 2.013
>
> Negative binomial dispersion parameter: 1.7665 (std. err.: 0.075737)
> Zero-inflation: 0.090062 (std. err.: 0.0083659 )
>
> Log-likelihood: -9856.38
> I then tried to compare AIC for Mod1, with AIC for a base model
> containing only an intercept as fixed effect:
> Mod0 <-glmmadmb(Killed_Lamb ~1, random=~ 1|KOM,
> data=DATA,zeroInflation=TRUE,family="nbinom")
> anova(Mod0,Mod1).
> I get the following result and a warning because the more complex model
> (Mod1) has a lower log-likelihood
> than the base model (Mod0)
>
> Analysis of Deviance Table
>
> Model 1: Killed_Lamb ~ 1
>
> Model 2: Killed_Lamb ~ I(scale(Dens_Lynx))
>
> NoPar LogLik Df Deviance Pr(>Chi)
>
> 1 4 -9849.6
>
> 2 5 -9856.4 1 -13.56 1
>
I don't know exactly what's going on here, but I agree it's weird.
Some suggestions:
* try fitting the more and less complex models with starting
values set to the parameters of the other model (i.e. fit the
more complex model starting at the reduced model, and the reduced
model starting at the estimate of the full model [forcing the
slope variable to zero])
* see whether adding a binary lynx present/absent variable helps
(this would add flexibility in the lynx density variable -- you
could also add a spline based on the lynx density; e.g. library(splines),
response ~ ns(Dens_Lynx,3); when
you plot the data, do the fitted parameters seem reasonable?)
* you can try a zero-inflated model via expectation-maximization
with lme4:
https://groups.nceas.ucsb.edu/non-linear-modeling/projects/owls/
WRITEUP/owls.Rnw
https://groups.nceas.ucsb.edu/non-linear-modeling/projects/owls/R/
owls_R_funs.R
_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
More information about the R-sig-mixed-models
mailing list