Prev 15950 / 20628 Next

ZIP MCMCglmm - how to increase effective sample size?

dani

Tue, Oct 31, 2017 5:39 PM

Hello Pierre,


Thank you so much for your message! I truly appreciate your help with the prior, I have exhausted all the other options and I was going to give up using MCMCglmm. I will give it a try and hopefully it works out!


Cheers,

Dani

From: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> on behalf of Pierre de Villemereuil <pierre.de.villemereuil at mailoo.org>
Sent: Tuesday, October 31, 2017 2:08 PM
To: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] ZIP MCMCglmm - how to increase effective sample size?

Hi,

There are about three way to increase effective sample size:
- increase the number of iterations
- use a prior with better properties
- change your model somehow (you might not always want to use that one...)

In your case, using a slightly more informative prior and the extended parameters prior might help? Something like:

priori <- list(R=list(V=diag(2), nu=1, fix=2),
               G=list(G1=list(V=diag(2)/2, nu=2, alpha.mu=c(0,0), alpha.V=diag(2)*1000)),
                                                 G2=list(V=diag(2)/2, nu=2, alpha.mu=c(0,0), alpha.V=diag(2)*1000))))

Hope this helps,
Pierre.

On Wednesday, 25 October 2017 18:45:46 NZDT dani wrote:
> Dear list members,
>
> I need some advice regarding this ZIP MCMCglmm model:
>
> library(MCMCglmm)
>
> priori <- list(R=list(V=diag(2), nu=0.002,fix=2),
>                G=list(G1=list(V=diag(2), n=2),G2=list(V=diag(2), n=2)))
>
> mj <- MCMCglmm(y ~ trait - 1 + at.level(trait,1):(x1+x2+x3+x4+ x5 +x6+x7+ offset),
>                   random = ~idh(trait):group1 + idh(trait):group2,
>                   family = "zipoisson",
>                   prior = priori,
>                   rcov = ~idh(trait):units,
>                verbose=FALSE,
>                thin   = 100,
>                burnin = 3000,
>                nitt   = 103000,
>                saveX=TRUE, saveZ=TRUE, saveXL=TRUE, pr=TRUE, pl=FALSE,
>                data = s25h)
>
> summary(mj)
>
> # Iterations = 3001:102901
> # Thinning interval  = 100
> # Sample size  = 1000
> #
> # DIC: 4811.791
> #
> # G-structure:  ~idh(trait):group1
> #
> #                                post.mean l-95% CI u-95% CI eff.samp
> # traity.group1       0.4307   0.1351   0.9281    10.17
> # traitzi_y. group1  4.3196   2.1216   7.4310    31.26
> #
> # ~idh(trait):group2
> #
> #                              post.mean l-95% CI u-95% CI eff.samp
> # traity. group2       0.4233   0.2341   0.6781    30.81
> # traitzi_y. group2    3.5497   1.2365   6.1525    26.39
> #
> # R-structure:  ~idh(trait):units
> #
> #                            post.mean l-95% CI u-95% CI eff.samp
> # traity.units      0.02393 0.002833  0.06621    10.58
> # traitzi_y.units   1.00000 1.000000  1.00000     0.00
> #
> # Location effects: y ~ trait - 1 + at.level(trait, 1):(x1 + x2 + x3 + x4 + x5 + x6 + x7 + offset)
> #
> #                                                            post.mean   l-95% CI   u-95% CI eff.samp  pMCMC
> # traity                                               -4.3823820 -6.1496186 -2.6424402   23.592 <0.001 ***
> # traitzi_y                                           3.4696204  2.6430476  4.1392235    1.922 <0.001 ***
> # at.level(trait, 1):x1                         -0.0498043 -0.2192051  0.1097667   16.979  0.522
> # at.level(trait, 1):x2M                     -0.2088408 -0.4535085  0.0440055    8.727  0.088 .
> # at.level(trait, 1):x31                       0.1422342 -0.1473884  0.4199985   11.521  0.288
> # at.level(trait, 1):x4                         0.0007054 -0.0030953  0.0043456   24.299  0.680
> # at.level(trait, 1):x5                         0.1131704  0.0647184  0.1676469   26.621 <0.001 ***
> # at.level(trait, 1):x6                        -0.0128734 -0.0483344  0.0306350   13.599  0.588
> # at.level(trait, 1):x7                         0.0102356 -0.0276141  0.0540893   40.746  0.680
> # at.level(trait, 1):offset                   1.3511873  0.6963525  2.1299075   13.216 <0.001 ***
> #   ---
> #   Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
>
> I would like to increase my effective samples, but I am not sure which way to go. I tried increasing the NITT to 503000, but the effective samples actually got worse. Is there anything else I could do? I plan on dropping some variables from the model, but if I were to proceed with the model above, what could I have done better?
>
> Thanks in advance!
> DNM
> <http://aka.ms/weboutlook>
>
>        [[alternative HTML version deleted]]
>

_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
R-sig-mixed-models Info Page - stat.ethz.ch<https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>
stat.ethz.ch
Your email address: Your name (optional): You may enter a privacy password below. This provides only mild security, but should prevent others from messing ...

Thread (13 messages)

Pierre de Villemereuil ZIP MCMCglmm - how to increase effective sample size? Oct 31 dani ZIP MCMCglmm - how to increase effective sample size? Oct 31 dani ZIP MCMCglmm - how to increase effective sample size? Oct 31 Pierre de Villemereuil ZIP MCMCglmm - how to increase effective sample size? Oct 31 dani ZIP MCMCglmm - how to increase effective sample size? Nov 1 Pierre de Villemereuil ZIP MCMCglmm - how to increase effective sample size? Nov 1 dani ZIP MCMCglmm - how to increase effective sample size? Nov 1 Jarrod Hadfield ZIP MCMCglmm - how to increase effective sample size? Nov 1 dani ZIP MCMCglmm - how to increase effective sample size? Nov 2 Jarrod Hadfield ZIP MCMCglmm - how to increase effective sample size? Nov 2 dani ZIP MCMCglmm - how to increase effective sample size? Nov 3 dani ZIP MCMCglmm - how to increase effective sample size? Nov 3 dani ZIP MCMCglmm - how to increase effective sample size? Nov 3