bbolker / glmmadmb

generalized linear mixed models with AD Model Builder
Other
11 stars 2 forks source link

Getting: Error in glmmadmb The function maximizer failed (couldn't find parameter file) when trying to fit a zero-inflated poisson model with random effects and an offset #9

Closed ndrubins closed 7 years ago

ndrubins commented 7 years ago

Hi,

I have count data of gene expression (measured by microscopy: number of gene molecules per cell, collected from many cells), from three groups (e.g., conditions). In each group the gene counts come from several different animals, and in all I have about 10% of the cells having counts of zero, where the total number of cells per group is balanced at ~800. Two of the groups have two animals and the third has three. Within the first group the #cells per animal is not balanced (~250 and ~550), but within the other two groups the numbers are pretty balanced.

Here's a histogram of the data: gene counts hist

And here's a boxplot: gene counts box

What I want is to estimate the effect of group on the gene counts. For this I thought that a zero-inflated poisson mixed-effects model would be appropriate. I followed the owls example and tried out the glmmadmb model: glmmadmb(gene.count ~ group + offset(log.n.cells) + (1|mouse.id), data=df, zeroInflation=TRUE, family="poisson") where log.n.cells is the log of the number of cells collected per each group. This model crashes with this error:

Parameters were estimated, but standard errors were not: the most likely problem is that the curvature at MLE was zero or negative

Error in glmmadmb(gene.count ~ group + offset(log.n.cells) + (1 | mouse.id),  : 
  The function maximizer failed (couldn't find parameter file) Troubleshooting steps include (1) run with 'save.dir' set and inspect output files; (2) change run parameters: see '?admbControl';(3) re-run with debug=TRUE for more information on failure mode
In addition: Warning message:
running command './glmmadmb -maxfn 500 -maxph 5 -noinit -shess' had status 1 

When I remove the offset term from the model: glmmadmb(gene.count ~ group + (1|mouse.id), data=df, zeroInflation=TRUE, family="poisson") it doesn't crash, however looking at the summary:

Call:
glmmadmb(formula = gene.count ~ group + (1 | mouse.id), data = df, 
    family = "poisson", zeroInflation = TRUE, save.dir = paste0(out.dir, 
        "/glmmadmb"))

AIC: 18691.6 

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)    1.682      0.192    8.76   <2e-16 ***
groupB         0.079      0.271    0.29     0.77    
groupC         0.231      0.248    0.93     0.35    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Number of observations: total=2497, mouse.id=7 
Random effect variance(s):
Group=mouse.id
            Variance StdDev
(Intercept)  0.07325 0.2706

Zero-inflation: 0.097608  (std. err.:  0.0061656 )

Log-likelihood: -9340.8 

the effects of groups B and C, relative to A, are weak and not significant. But looking at the histogram, at least I, got the feeling they are significant.

I also tried fitting the zipme model (from here):

zipme(cformula = gene.count ~ group + offset(log.n.cells) + (1|mouse.id),zformula=z ~ 1,data=df,maxitr=20,tol=1e-6)

It converges in 7 iterations and the cfit is:

Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: poisson  ( log )
Formula: gene.count ~ group + offset(log.n.cells) + (1 | mouse.id)
   Data: bydataw
Weights: (1 - z)
      AIC       BIC    logLik  deviance  df.resid 
17625.251 17648.542 -8808.626 17617.251      2493 
Random effects:
 Groups   Name        Std.Dev.
 mouse.id (Intercept) 0.2863  
Number of obs: 2497, groups:  mouse.id, 7
Fixed Effects:
(Intercept)   groupB   groupC  
   -4.96953      0.02189      0.24924  

While the effects have no errors or p-values, their estimates seem pretty comparable to that of the glmmadmb fit without the offset term.

So my questions are:

  1. Am I fitting the right models?
  2. What's the error I'm getting and how to fix it?
  3. Do the results of the last two models make sense (relative to the figures) and hence the conclusion is that there is no significant group effect on gene.count ?

Thanks a lot

bbolker commented 7 years ago

These seem reasonable.

Follow-ups to r-sig-mixed-models@r-project.org please ....

== @Article{ murtaugh_simplicity_2007, title = {Simplicity and Complexity in Ecological Data Analysis}, volume = {88}, url = {http://www.esajournals.org/doi/abs/10.1890/0012-9658%282007%2988%5B56%3ASACIED%5D2.0.CO%3B2} , number = {1}, journal = {Ecology}, author = {Murtaugh, Paul A}, year = {2007}, pages = {56--62} }

ndrubins commented 7 years ago

Thanks a lot

On Sat, Oct 21, 2017 at 5:07 PM, Ben Bolker notifications@github.com wrote:

Closed #9 https://github.com/bbolker/glmmadmb/issues/9.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/bbolker/glmmadmb/issues/9#event-1304350547, or mute the thread https://github.com/notifications/unsubscribe-auth/ABGsHilK7aceXaztAjWoFPXBJxTj7XPsks5suodbgaJpZM4QBv7n .