/22 R MCMC R R MCMC? 3. Gibbs sampler : kubo/

Size: px

Start display at page:

Download "/22 R MCMC R R MCMC? 3. Gibbs sampler : kubo/"

こおがふじした
4 years ago
Views:

1 /22 R MCMC R R MCMC? 3. Gibbs sampler : kubo@ees.hokudai.ac.jp kubo/

2 /22 : ( ) : : ( ) : (?) community ( )

3 /22 :? 1. ( ) 2. ( ) 3. or?? MCMC? 1.? 2. MCMC? 3.? ( ) R?!

4 /22 ( ) [ ] Bayes Bayes [ ] fixed effects Bayes (GLMM) [ (GLM)] +, fixed effects [ ] +

5 /22? random effects ( ) nest random effects Plot A Plot B : MCMC : Gaussian Random Field

6 /22 MCMC?! Markov Chain Monte Carlo? Gibbs (sampling) : (sampling!) sample random sample set log likelihood (non constant par leaf area index step (burn-in) (K MCMC step) sampling point (K MCMC step)

7 /22 ( ) ( ) ( ) ( ) p(β, α y) p(y β)p(β α)p(α) : Markov Chain Monte Carlo (MCMC) MCMC 1: Metropolis-Hastings MCMC 2: Gibbs sampler ( ) Likelihood(α y) p(y β)p(β α)dβ : α ( ) : (GLMM)? ( : )

8 /22 R MCMC / Gibbs sampler MCMC? ( ) R package: library(mcmcpack) ( ) Gibbs sampler (R ) WinBUGS OpenBUGS JAGS WinBUGS OpenBUGS WinBUGS, 2004 OpenBUGS WinBUGS project, GPL

9 /22 : WinBUGS Gibbs sampler BUGS adaptive rejection sampler ( OpenBUGS ) Windows Linux WINE MacOS X Darwine GUI (Linux ) R R2WinBUGS ( )

10 /22 BUGS Spiegelhalter et al BUGS: Bayesian Using Gibbs Sampling version model { mu dnorm(0, 1.0E-2) tau dgamma(1.0e-3, 1.0E-3) for (i in 1:n.samples) { re[i] dnorm(0.0, tau) p[i] <- 1.0 / (1.0 + exp(-(mu + re[i]))) n.seeds[i] dbin(p[i], n.ovules[i]) } } JAGS ;

11 /22 BUGS (?) : zero-inflated Poisson (ZIP) model BUGS :?

12 /22 GPL WinBUGS : OpenBUGS Thomas Andrew WinBUGS OpenBUGS is still in development and suffers frequent crashes. Component Pascal BlackBox Component Builder Windows Linux R library(brugs) Windows R!!

13 /22 R (?) Gibbs sampler: JAGS 0.97 R core team Martyn Plummer Just Another Gibbs Sampler C++ R Vine Linux RPM package Windows JAGS MacOS X R : library(rjags)

14 /22 JAGS JAGS directed cycle ( ) model { x dnorm(y, tau) y dnorm(x, tau) }! WinBUGS? R

15 /22 : JAGS WinBUGS

16 /22 : [ ] i {1,, 100} n.ovules[i] = 8 n.seeds[i] p[i]: x = 5 : : logit(p[i]) = 0 + re[i] random effects: re[i] N(0, 1/0.5) = 0.5, = 4

17 /22 : ( ) [ ] = p = overdispersion ( )! random effects x = 5 ( )

18 /22 JAGS 1. BUGS model 2. R 3. JAGS (foo.cmd) 4. jags foo.cmd 5. JAGS R library(coda) read.coda() (mcmc ) 6. mcmc plot(), summary(),

19 /22 library(coda): MCMC Convergence Diagnosis and Output Analysis for MCMC R package ( S-plus ) Martyn Plummer MCMC mcmc, mcmc.list MCMC

20 /22 mcmc summary > summary(r.mcmc) 1. Empirical mean and standard deviation for each variable, plus standard error of the mean: Mean SD Naive SE Time-series SE mu tau Quantiles for each variable: 2.5% 25% 50% 75% 97.5% mu tau

21 /22 WinBUGS : R2WinBUGS 1. BUGS model 2. R2WinBUGS package R 3. R bugs 5. plot() summary() 6. mcmc / mcmc.list

22 /22 : R MCMC? ( MCMC ) WinBUGS R2WinBUGS? WinBUGS R! ( ) JAGS R ( Plummer?) MCMC geor : Gaussian Random Field GRF MCMC : Langevin-Hastings? : OpenBUGS?... or... Umacs + rv B?

/ 1 1 c Mike Gonzalez, October 14, Wikimedia Commons.

/ *1 *1 c Mike Gonzalez, October 14, Wikimedia Commons. 2010 05 22 1/ 35 2010 2010 05 22 *1 kubo@ees.hokudai.ac.jp *1 c Mike Gonzalez, October 14, 2007. Wikimedia Commons. 2010 05 22 2/ 35 1. 2. 3. 2010 05 22 3/ 35 : 1.? 2. 2010 05 22 4/ 35 1. 2010 05 22 5/