[R-sig-ME] Poor mixing and autocorrelation of ZIP data in MCMCglmm
Pierre de Villemereuil
p|erre@dev|||emereu|| @end|ng |rom ephe@p@|@eu
Tue Jun 9 18:30:37 CEST 2020
Hi,
I'm very not sure about the prior you used: why fix all the random effects variances to 1 in the ZI part? I suspect at least a part of your problem might stem from there.
Here's the kind of prior I would suggest for such a model:
list(R = list(V = diag(2), nu = 2, fix = 2),
G = list(G1 = list(V = diag(2), nu = 2, alpha.mu = c(0,0), alpha.V = diag(2)),
G2 = list(V = diag(2), nu = 2, alpha.mu = c(0,0), alpha.V = diag(2)))
It's informative toward small variances, but you a priori know variances will be small in such models, because of the link functions involved.
Hope this helps,
Pierre.
Le mardi 9 juin 2020, 12:49:09 CEST Vital Heim a écrit :
> Good morning,
>
> Currently, I am working on a Zero-Inflated Poisson Generalized Linear Mixed
> Model using the MCMCglmm package in R. I am looking at the food uptake of
> animal at a tourism provisioning site. I try to understand what variables
> impact the amount of food/bait the animals consume.
>
> The response variable (bait uptake) is zero-inflated. I have 6 fixed
> effects (nr. tourists, nr. operators at the site, nr. conspecifics,
> temperature, presence time of the animal at the site, and gender of the
> animal). Presence time of the animal and temperature are only entered
> at.level(trait, 1), all the other fixed effects are entered into the ZI
> part and the count part of the model. I have two random effects for ID of
> the animal and the provisioning event.
>
> Here is the model specification:
>
>
> modelZI.9 <- MCMCglmm(bait ~ trait - 1
>
> + trait:conspecifics + trait:tourists +
> trait:operators + at.level(trait,1):presence +
> at.level(trait,1):temperature
>
> + trait:gender
>
> , random = ~idh(trait):id + idh(trait):event
>
> , rcov =~ idh(trait):units
>
> , data=df, family="zipoisson", prior=prior6
>
> , nitt=NITT*MULT, burnin=BURN*MULT, thin = THIN*MULT,
> verbose=FALSE)
>
>
>
>
> I started with the chain parameters: nitt = 120000, burn = 20000, thin =
> 100 to get an effective sample size of 1000. My first prior was:
>
>
>
> prior.base<-list(R=list(V=diag(c(1,1)), nu=0.002, fix=2)
>
> ,G=list(G1=list(V=diag(c(1,1)), nu=0.002, fix=2)
>
> ,G2=list(V=diag(c(1,1)), nu=0.002, fix=2)))
>
>
>
> The chosen parameters and prior resulted in a very poor mixing of the ZI
> effects as well as the both random effects. The effective sample sizes for
> the ZI effects were very low (less than 100). I adjusted the chain
> parameters (see below) I also adjusted the prior to account for
> overdispersion in the fixed effects and making the prior stronger for the
> random intercepts. There was autocorrelation in my ZI fixed effects so I
> increased the thinning while increasing nitt and burn accordingly to get an
> effective size of 1100. These measures improved the model in some parts.
> However, I tried to validate my prior choice and parameters settings by
> reading the literature to see if I am not doing anything “illegal” with my
> model. By now the parameters and prior gives me slightly better mixing. I
> was able to get the random effect ID mixing nicely but only when I set the
> nu = 100.002, which seems very high to me.
>
> These are my current chain parameters and my prior:
>
> NITT <- 185000
>
> BURN <- 15000
>
> THIN <- 150
>
> MULT <- 10
>
>
>
> prior6<-list(R=list(V=diag(c(1e-6,1)), nu=10.002, fix=2)
>
> ,G=list(G1=list(V=diag(c(1,1)), nu=100.002, fix=2)
>
> ,G2=list(V=diag(c(1,1)), nu=0.002, fix=2)))
>
>
>
> The autocorrelation diagnostics are:
>
> 0.075907049 NaN
>
> > autocorr.diag(modelZI.9$Sol)
> bait zi_bait bait:conspecifics
> zi_bait:conspecifics bait:tourists zi_bait:tourists bait:operators
> Lag 0 1.00000000 1.00000000 1.0000000000 1.00000000
> 1.00000000 1.000000000 1.000000000
> Lag 1500 -0.03158767 0.45992662 -0.0216380764 0.52372246
> -0.02650593 0.546651545 -0.027195064
> Lag 7500 0.02284612 0.15705391 0.0114730973 0.12436415
> 0.03518200 0.122140444 -0.010875426
> Lag 15000 0.04558767 0.06901158 -0.0131587465 0.05700078
> 0.05371373 0.054479522 0.024634961
> Lag 75000 0.03567369 0.04099224 -0.0005695797 0.02014238
> -0.02546435 -0.006299296 -0.005569035
> zi_bait:operators presence temperature
> bait:gendermale zi_bait:gendermale
> Lag 0 1.00000000 1.000000000 1.00000000
> 1.000000000 1.00000000
> Lag 1500 0.43279808 -0.004789598 -0.03151510
> -0.002708028 0.33387297
> Lag 7500 0.05697190 0.012426457 0.02221331
> -0.025070273 0.11236772
> Lag 15000 0.05472752 -0.043356418 0.04499879
> -0.003589775 0.05618167
> Lag 75000 0.03342174 -0.007475513 0.03201374
> 0.015790374 -0.03167048
>
> > autocorr.diag(modelZI.9$VCV)
> bait.id zi_bait.id bait.dive
> zi_bait.dive bait.units zi_bait.units
> Lag 0 1.000000000 NaN 1.00000000 NaN
> 1.000000000 NaN
> Lag 1500 -0.026735489 NaN 0.03018431 NaN
> 0.003349623 NaN
> Lag 7500 -0.051666715 NaN -0.02011165 NaN
> -0.048462744 NaN
> Lag 15000 0.003918555 NaN -0.02777388 NaN
> -0.001752476 NaN
> Lag 75000 0.028123317 NaN -0.01640397 NaN
> 0.075907049 NaN
>
>
>
> I therefore wanted to ask if there are suggestions on how I could improve
> the mixing of my model without choosing these drastic changes or are these
> changes in the range of what I am allowed to do to my model?
>
>
> Thank you,
>
> Vital
> --
> *Vital Heim *
>
> *PhD student*
> University of Basel, Switzerland
>
> Bimini Biological Field Station Foundation
> Bimini, Bahamas
>
> +41 (0)79 732 05 57
> vital.heim using gmail.com
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models using 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