獨斷論

R의 rjags를 이용한 베이시언Bayesian 입문 3 본문

과학과 기술/R 통계

R의 rjags를 이용한 베이시언Bayesian 입문 3

부르칸 2021. 2. 7. 14:46

앞서 1과 2에서는 JAGS모델을 별도의 파일에 저장하여 R script에서 불러들였다.

이는 모델 전체를 이해하는데 도움이 되지 않는다.

여기서는 별도의 *.jags 파일을 작성하지 않고 R script에 직접 입력하는 방법을 알아보자.

 

입문1을 아래와같이 R script만으로 JAGS를 이용할수 있다.

 

1
2
3
4
5

6

7
8
9
10
11
12
13
14
15
16

18
19
20
21
22
23
24
25

setwd("d:/tmp/rcode/")
rm(list = ls())
x <- read.csv("example1.csv", header=FALSE)
x <- as.vector(x[,1])
N <- length(x)

library(rjags)

modelstring = "
  model {
    mu ~ dnorm(0, 1/100^2)
    tau <- 1 / sigma^2
    sigma ~ dunif(0, 100)
    for(i in 1:N) {
      x[i] ~ dnorm(mu, tau)
    }
  }
"

jags <- jags.model(
  textConnection(modelstring),
  data = list('x' = x, 'N' = N),
  n.chains = 4,
  n.adapt = 100
)
update(jags, 1000)
jags.results <- jags.samples( jags, c('mu', 'tau'), 1000 ) 

example1.jags를 만들어 파일을 저장하지 않고

example1.jags의 코딩을 modelstring이라는 변수안에 쌍따옴표를 이용하여 넣는다.

그리고 jags.model( )의 argument로 textConnection(modelstring)이라고 넘겨주면 된다(이전에는 "example1.jags" 파일이름을 직접 사용하였었다.

 

입문2도 아래와 같다.

같은 방법으로 example2.jags의 내용을 그대로 modelstring에 넣으면 된다.

1
2
3
4
5

6

7
8
9
10
11
12
13
14
15
16
17

18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34

rm(list = ls())
setwd("d:/tmp/rcode/")
regdat <- read.csv('example2.csv', header=TRUE, strip.white=TRUE)
head(regdat)
dim(regdat)

library(rjags)

modelstring = "
model {
  a ~ dnorm(0, 1/100^2)
  b ~ dnorm(0, 1/100^2)
  sigma ~ dunif(0, 100)
  tau <- 1/sigma^2
  for(i in 1:N) {
    mu[i] <- a + b*x[i]
    y[i] ~ dnorm(mu[i], tau)
  }
}"

x <- regdat$X
y <- regdat$Y
N <- length(y)
jags <- jags.model(
  textConnection(modelstring),
  data = list(
    'x' = x,
    'y' = y,
    'N' = N
  ),
  n.chains = 4,
  n.adapt = 100
)
update(jags, 1000)
jags.results <- jags.samples(jags, c('a', 'b'), 1000)
summary(jags.results)
jags.results

 

Comments