[R-sig-ME] zero-inflated models in MCMCglmm
Jarrod Hadfield
j.hadfield at ed.ac.uk
Sat Nov 26 08:59:39 CET 2016
Hi,
If X1 was has two levels then it is odd that you have an intercept AND 2
X1 contrasts. Have you used singular.ok=TRUE? If not then are you sure
there is not a third (rare) level to X1? Try unique(your_data$X1):
trailing whitespace from excel spreadsheets is a common cause of such a
problem.
Regarding the predictions, the issue may be with the above. However,
note that there was a bug with posterior="mean" in the call to
predict.MCMCglmm (the first posterior iteration of the fixed effects was
used rather than there mean). This bug is fixed in the current version,
and there was no bug in the default posterior="all".
Cheers,
Jarrod
On 25/11/2016 13:16, Gustaf Granath wrote:
> Hi,
> A few questions regarding zero-inflated models using MCMCglmm.
>
> 1. Odd contrasts are created when using, e.g. :
> y ~ trait-1 + at.level(trait, 1):(X1+X2), family = "zipoisson" #X1 and
> X2 are factors with 2 levels
> COEF TABLE
> traity
> traitzi_y
> at.level(trait, 1):X1lev1
> at.level(trait, 1):X1lev2
> at.level(trait, 1):X2lev2
>
> It doesnt look like this in the course notes (what I can see). Is
> at.level(trait, 1):X1lev1 the reference level for everything below? I
> also get very high (>1000) estimates for traity, at.level(trait,
> 1):X1lev1 and at.level(trait, 1):X1lev2.
>
> 2. I made four models
> a) y ~ X1*X2, family = "poisson"
> Large overdispersion (units ~ 5) but everything looks fine
> (traceplots, random effects (including units) ). Predictive checks
> show that the models predict too few zeros though. And alarming is
> that predictions are really bad, e.g. a treatment mean of 40 is
> predicted to have 300 counts.
>
> b) y ~ trait-1 + at.level(trait, 1):(X1+X2), family = "zipoisson"
> Model looks OK, but again, predictions are way too high for the higher
> tretament means.
>
> c:d) y ~ trait-1 + at.level(trait, 1):(X1+X2), family = "zapoisson" or
> "hupoisson"
> Predicted values are much better. But is it the same way get
> predictions here? Possible to get predictions from the separate models
> (poisson/binomial) in the hurdle case? (code available?)
> #get e.g. treatment predictions?
> predict.MCMCglmm(model, marginal = ~ random_factor, type="response",
> posterior = "mean")
>
> The "zapoisson" model performs best but I dont understand why the
> "zipoisson" is so bad. Any typical things to look for when it looks
> like this?
> Interpreting "zapoisson" isnt easy, any good literature/tutorials on
> this model?
>
>
> Cheers
>
> Gustaf
>
>
--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.
More information about the R-sig-mixed-models
mailing list