DIC(Deviance Information Criterion)

베이지안 모델에서는 모델 선택을 위한 정보의 기준으로써 DIC(Deviance Information Criterion)을 제시하고 있다. DIC의 공식은 다음과 같다.

 

 

▷ theta hat은 각 모수의 사후평균이고, 사후분포로부터 얻은 theta hat의 로그 가능도와 실질적인 모수의 갯수(Effective number of parameters)를 고려하여 DIC를 구할 수 있다.

 

▶ 실질적인 모수의 갯수는 모델의 추정치 사이에 상관(Correlation)을 고려하기 위한 것이다. 예를 들어, 모델의 추정치 사이에 0.99의 상관이 존재한다면 이를 독립적인 모수로 간주한다면 합리적이지 않을 것이다.

 

R의 car 패키지의 Leihardt 데이터를 이용한 두 모델을 DIC를 통해 비교하여보자.

 

In:

library(car)
library(rjags)

data('Leinhardt')

Leinhardt$log_infant = log(Leinhardt$infant)
Leinhardt$log_income = log(Leinhardt$income)

data = na.omit(Leinhardt)

###########
# model 1 #
###########

mod1_string = " model {
  for (i in 1:n) {
    y[i] ~ dnorm(mu[i], prec)
    mu[i] = b[1] + b[2]*log_income[i] 
  }
  
  for (i in 1:2) {
    b[i] ~ dnorm(0.0, 1.0/1.0e6)
  }
  
  prec ~ dgamma(5/2.0, 5*10.0/2.0)
  sig2 = 1.0 / prec
  sig = sqrt(sig2)
} "

data1_jags = list(y = data$log_infant, 
                  n = nrow(data), 
                  log_income = data$log_income)

params1 = c('b', 'sig')

inits1 = function() {
  inits = list('b' = rnorm(2, 0.0, 100.0), 'prec' = rgamma(1, 1.0, 1.0))
}

mod1 = jags.model(textConnection(mod1_string), 
                  data = data1_jags, 
                  inits = inits1, 
                  n.chains = 3)

###########
# model 2 #
###########

mod2_string = " model {
    for (i in 1:length(y)) {
        y[i] ~ dnorm(mu[i], prec)
        mu[i] = b[1] + b[2]*log_income[i] + b[3]*is_oil[i]
    }
    
    for (i in 1:3) {
        b[i] ~ dnorm(0.0, 1.0/1.0e6)
    }
    
    prec ~ dgamma(5/2.0, 5*10.0/2.0)
    sig = sqrt( 1.0 / prec )
} "

data2_jags = list(y = data$log_infant, 
                  log_income = data$log_income,
                  is_oil = as.numeric(data$oil == 'yes'))

params2 = c('b', 'sig')

inits2 = function() {
  inits = list('b' = rnorm(3, 0.0, 100.0), 
               'prec' = rgamma(1, 1.0, 1.0))
}

mod2 = jags.model(textConnection(mod2_string), 
                  data = data2_jags, 
                  inits = inits2, 
                  n.chains = 3)
                  
update(mod1, 1e3)
update(mod2, 1e3)

 

Out:

Compiling model graph
   Resolving undeclared variables
   Allocating nodes
Graph information:
   Observed stochastic nodes: 101
   Unobserved stochastic nodes: 3
   Total graph size: 404

Initializing model
Compiling model graph
   Resolving undeclared variables
   Allocating nodes
Graph information:
   Observed stochastic nodes: 101
   Unobserved stochastic nodes: 4
   Total graph size: 507

Initializing model
  |**************************************************| 100%
  |**************************************************| 100%

 

▷ 위의 두 모델은 연속형 변수인 infant를 예측하기 위한 선형회귀 모델이다. 첫 번째 모델은 연속형 변수인 income만을 이용하였고, 두 번째 모델은 income과 범주형 변수인 oil을 이용하여 모델링하였다.

 

▷ DIC를 구하기 전, 1,000회 데이터 생성을 통해 마르코프 체인이 수렴하도록 하였다. 수렴을 확인하는 과정은 아래의 포스팅을 통해서 확인할 수 있다.

 

MCMC(Markov Chain Monte-Carlo)의 수렴(Convergence)

 

In:

dic.samples(mod1, n.iter = 1e3)
dic.samples(mod2, n.iter = 1e3)

 

Out:

  |**************************************************| 100%
Mean deviance:  231.4 
penalty 2.945 
Penalized deviance: 234.3 
  |**************************************************| 100%
Mean deviance:  225.5 
penalty 3.999 
Penalized deviance: 229.5 

 

▷ dic.samples()는 사후분포로부터 n.iter만큼 데이터를 생성하고, 이 데이터를 이용하여 추정치를 구하고 DIC를 계산한다.

 

▷ Mean deviance는 추정치의 로그 가능도에 -2를 곱한 것으로, 낮을수록 모델이 데이터에 잘 적합한다고 볼 수 있다. 반면, penalty는 실질적인 모수의 갯수로, 첫 번째 모델의 결과는 2.945인 것으로 나타났다. 모델의 파라미터 b[1], b[2], prec때문에 약 3개로 나타난 것을 확인할 수 있다. 두 번째 모델은 oil에 대한 파라미터가 추가되었으므로, 약 4개인 것으로 나타난 것을 알 수 있다.

 

▶ Penalized deviance는 Mean deviance와 penalty의 합으로 DIC를 의미한다.

 

▷ 모델 선택시, DIC가 낮을수록 더 적합한 모델임을 의미한다. 따라서 DIC가 229.5로 더 낮은 두 번째 모델을 최종 모델로 선정할 수 있을 것이다.

 

---

Reference:

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

"DIC와 WAIC: 베이지안 모형선택을 위한 정보 기준들," 박준석, https://bayestour.github.io/blog/2019/07/12/DIC_WAIC.html.