5
$\begingroup$

After having visited stats stack exchange countless times, I'm finally asking a question!

For my research, I am try to run a model of the form: $$ Y = f(X,B)+ g(X) + \epsilon$$

Where $f(X,B)$ is a zero-inflated poisson/negative binomial density, $g(X)$ is a Gaussian process representing spatial correlation, and $epsilon$ is your usual error. Sorry if my notation is a little loose. Additionally, I want my zero-inflated portion to be able to have random effects. In summary, I need a program/package that can run:

  1. Zero-inflated GLMs
  2. Covariance structure (gaussian process)
  3. Random Effects

However, I can only find programs that are capable of two out of these three items. R is highly preferred, but I am also open to SAS as a last resort. I already looked at looked at PROC GLMMIX in SAS, but that does not have all three of these items. I'm already using packages like pscl (which does not have all three items in the list above). Stan is too slow to handle the size of the dataset. I have around 100,000 observations, and it would appear that the processing time is on the order of days or even weeks (100 observations took 18 minutes, and this was with assigning each of the two chains I was running to its own core (i5). And this was before even using an random or spatial effects in the ZIP model! In my literature review, it would appear that R-INLA could potentially solve this issue if I stick with bayesian, but it would appear people still need computing clusters to make this work.

Can anyone suggest ways to resolve this ?

Background: I am modeling the abundance of a species of bird in Wisconsin. As you can imagine, most observations even for commons species have "zero" as the count. Hence the need for zero inflation. There is also observer data with repeated measures, hence the mixed model. And there is spatial correlation to deal with, hence the Gaussian process.

Improve this question
$\endgroup$

1 Answer 1

Reset to default
2
$\begingroup$

If you adopt a Bayesian approach then you should be able to fit such a model. This can be estimated using Markov Chain Monte Carlo for example with JAGS (Just Another Gibbs Sampler) or Hamilton Monte Carlo with Stan. Both can be set up and run from within R. I appreciate that you say it may take a long time to run, but there are ways to speed things up such as running multiple chains on multiple compute nodes.

Improve this answer
$\endgroup$
4
  • $\begingroup$ Hi Robert, thank you thoughts. I outlined in more detail the time its taking to run these procedures in my post. Even if I were to cut the time in half, or even a quarter, this would not be fast enough. In contrast, the frequentist ZIP model in PSCL can run in about 4 seconds. It seems to me my only hope is R-INLA or a computing cluster (is that the same as "multiple compute nodes?"). We do have one at my university. Because I need to tamper a lot with the model, I'm probably only willing to wait a few minutes for it to run. Do you feel this can be accomplished with multiple computer nodes? $\endgroup$
    – user271536
    Commented Jan 22, 2020 at 16:05
  • $\begingroup$ @MarkMiller yes I was referring to a HPC (high performance computing) cluster. You say that you do have one so I would try that out. Also are you in the UK ? If so then you should be to get access to this. Also, the stan users forum is monitored by the developers who can often help you speed up a model considerably with their knowledge. Finally have you seen Ben Bolker's brain dump about modelling autocorrelation in mixed models with lme4 and similar packages ? $\endgroup$
    – Robert Long
    Commented Jan 22, 2020 at 16:18
  • $\begingroup$ If this is the link you mean (bbolker.github.io/mixedmodels-misc/notes/corr_braindump.html), then it would appear I had visited it prior to making this post. Why did I pass it by?! Ben's brain dump appears to be a comprehensive response to the status of my question, even if none of them work out. I will mark this as the answer for your node suggestion and Ben's resource; you can edit however you please. $\endgroup$
    – user271536
    Commented Jan 22, 2020 at 17:55
  • $\begingroup$ @MarkMiller thank you. Yes, that's the one. Glad I could help. I will also update the answer a little later with the link so that future visitors to this thread won't need to read these comments. $\endgroup$
    – Robert Long
    Commented Jan 22, 2020 at 18:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.