JAGS(Just Another Gibbs Sampler) 사용법

R의 JAGS(Just Another Gibbs Sampler)의 사용법에 대해 알아볼 것이다. JAGS를 통한 데이터 생성 과정은 4단계로 나눌 수 있다.

 

1. Specify the model

2. Set up the model

3. Run the MCMC(Markov Chain Monte Carlo) sampler

4. Post processing

 

다음의 모델에 이를 단계별로 적용하여 사후분포로부터 데이터를 생성하여 보자.

 

 

1. Specify the model

 

In:

library(rjags)

mod_string = " model {
  for (i in 1:n) {
    y[i] ~ dnorm(mu, 1.0/sig2)
  }
  
  mu ~ dt(0.0, 1.0/1.0, 1)
  sig2 = 1.0
} "

 

▷ 위의 코드와 같이 가능도 함수와 사전분포에 대해 정의한다.

 

▷ R과 JAGS의 가장 주요한 차이는 분포 함수의 사용 방법에 있다. R에서는 평균과 표준편차를 rnorm()의 모수로 받지만, JAGS에서는 표준편차 대신 분산의 역수인 Precision을 모수로 받는다. 이외의 JAGS의 함수와 관련된 내용은 메뉴얼(https://web.sgh.waw.pl/~atoroj/ekonometria_bayesowska/jags_user_manual.pdf)을 통해 확인할 수 있다.

 

2. Set up the model

 

In:

set.seed(777)

y = c(1.2, 1.4, -0.5, 0.3, 0.9, 
      2.3, 1.0, 0.1, 1.3, 1.9)
n = length(y)

data = list(y = y, n = n)
params = c('mu')

inits = function() {
  inits = list('mu' = 0.0)
}

mod = jags.model(textConnection(mod_string), 
                 data = data, 
                 inits = inits)

 

Out:

Compiling model graph
   Resolving undeclared variables
   Allocating nodes
Graph information:
   Observed stochastic nodes: 10
   Unobserved stochastic nodes: 1
   Total graph size: 15

Initializing model

  |++++++++++++++++++++++++++++++++++++++++++++++++++| 100%

 

▷ 위의 코드와 같이 데이터와 초기값을 설정한 뒤, jags.model()을 이용하여 모델의 코드를 컴파일 한다.

 

3. Run the MCMC sampler

 

In:

update(mod, 500)

mod_sim = coda.samples(model = mod, 
                       variable.names = params, 
                       n.iter = 1000)

 

Out:

|**************************************************| 100%
|**************************************************| 100%

 

▷ update()를 이용하여 모델을 학습시킨다. 즉, 초기값부터 시작하여 주어진 Iteration 동안 데이터 생성과정을 거친다. 이때, 생성된 데이터는 저장되지 않는다. 이 과정은 마르코프 체인이 정상 상태(Stationary state)에 도달하도록 도움을 준다.

 

▷ coda.samples()를 이용하여 학습된 모델로부터 n.iter만큼 데이터를 생성한다.

 

4. Post processing

 

In:

summary(mod_sim)

plot(mod_sim)

 

Out:

Iterations = 1501:2500
Thinning interval = 1 
Number of chains = 1 
Sample size per chain = 1000 

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

          Mean             SD       Naive SE Time-series SE 
       0.88664        0.31683        0.01002        0.01315 

2. Quantiles for each variable:

  2.5%    25%    50%    75%  97.5% 
0.2810 0.6746 0.8952 1.0971 1.5277 

 

 

▷ summary()와 plot()을 이용하여 생성된 데이터에 대한 분석 결과를 얻을 수 있다.

 

---

Reference:

"Bayesian Statistics: From Concept to Data AnalysisTechniques and Models," Coursera, https://www.coursera.org/learn/bayesian-statistics/.