2. Deviance Information Criterion#

Deviance Information Criterion (DIC) is similar in spirit to Akaike Information Criteria (AIC)

\[ \text{AIC} = -2\log f(y|\hat{\theta}) + 2p \]

where \(\hat{\theta}\) is the MLE of \(\theta\), and \(p\) is the number of parameters. A model with lower AIC is favored.

The Bayesian version of DIC was introduced by Spiegelhalter, 2002.

\[\begin{split}\begin{align*} \bar{D} & = \text{E}^{\pi(\cdot|y)}(-2\log f(y|\theta)) \\ D(\bar{\theta}) & = -2 \log f(y|\bar{\theta}), \space \bar{\theta} \text{ Bayes Estimator of } \theta \\ p_D & = \bar{D} - D(\bar{\theta}) \equiv \text{ Effective number of parameters} \\ DIC & = \bar{D} + p_D = D(\bar{\theta}) + 2 p_D \end{align*}\end{split}\]

Widely Applicable Information Criteria (WAIC)#

The DIC described above was available as an automatic calculation in BUGS. In practice, we use a newer model fit metric called Widely Applicable Information Criteria (WAIC) available in Arviz. We won’t go into it’s formula, but you can read more about it at arviz.waic, with referenced papers.

Let’s see how it works with the same Weibull and exponential models. One extra thing we need to do is save the log-likehood in the trace object using the idata_kwargs=dict(log_likelihood=True) argument in pm.sample.

By default, az.waic provides WAIC on a log score scale, where the prefered model is the higher value. Alternatively, we can specifiy the deviance scale, where a lower score is prefered, like the AIC and DIC defined above.

import pymc as pm
import arviz as az
import numpy as np

y = np.array([1,1,2,2,3,4,4,5,5,8])

#Weibull Model
with pm.Model() as m1:
    #priors
    mu = pm.Normal("mu",0,tau=.001)
    gamma = pm.Gamma("gamma",0.001,0.001)

    #convert Weibull parameterization
    beta = gamma ** (-1 / mu)

    #likelihood
    like = pm.Weibull("like",mu,beta,observed=y)
    
    trace1 = pm.sample(idata_kwargs=dict(log_likelihood=True))

#Exponential Model
with pm.Model() as m2:
    #prior
    lambda_ = pm.Gamma("lambda_",0.001,0.001)

    #likelihood
    like = pm.Exponential("like",lambda_,observed=y)

    trace2 = pm.sample(idata_kwargs=dict(log_likelihood=True))

Initializing NUTS using jitter+adapt_diag...
Multiprocess sampling (4 chains in 4 jobs)
NUTS: [mu, gamma]



Sampling 4 chains for 1_000 tune and 1_000 draw iterations (4_000 + 4_000 draws total) took 0 seconds.
There were 13 divergences after tuning. Increase `target_accept` or reparameterize.
Initializing NUTS using jitter+adapt_diag...
Multiprocess sampling (4 chains in 4 jobs)
NUTS: [lambda_]







Sampling 4 chains for 1_000 tune and 1_000 draw iterations (4_000 + 4_000 draws total) took 0 seconds.













az.waic(trace1,scale='deviance')









/opt/anaconda3/envs/pymc_env/lib/python3.13/site-packages/arviz/stats/stats.py:1655: UserWarning: For one or more samples the posterior variance of the log predictive densities exceeds 0.4. This could be indication of WAIC starting to fail. 
See http://arxiv.org/abs/1507.04544 for details
  warnings.warn(





Computed from 4000 posterior samples and 10 observations log-likelihood matrix.

              Estimate       SE
deviance_waic    44.76     4.72
p_waic            1.73        -

There has been a warning during the calculation. Please check the results.













az.waic(trace2,scale='deviance')









Computed from 4000 posterior samples and 10 observations log-likelihood matrix.

              Estimate       SE
deviance_waic    46.52     3.82
p_waic            0.40        -













az.compare#


We see a lower deviance_waic for the Weibull model, as expected. Another way to compare 2 models is with az.compare






az.compare({"Weibull":trace1,"Exponential":trace2},
           ic='waic',
           scale='deviance')









/opt/anaconda3/envs/pymc_env/lib/python3.13/site-packages/arviz/stats/stats.py:1655: UserWarning: For one or more samples the posterior variance of the log predictive densities exceeds 0.4. This could be indication of WAIC starting to fail. 
See http://arxiv.org/abs/1507.04544 for details
  warnings.warn(









  
    

      
      rank
      elpd_waic
      p_waic
      elpd_diff
      weight
      se
      dse
      warning
      scale
    

  
  
    

      Weibull
      0
      44.759743
      1.728281
      0.000000
      0.828696
      4.715606
      0.000000
      True
      deviance
    

    

      Exponential
      1
      46.515467
      0.401046
      1.755724
      0.171304
      3.815093
      3.252716
      False
      deviance