-
Notifications
You must be signed in to change notification settings - Fork 47
Watanabe Akaike information criterion (WAIC)
Waic example from Vehtari & Gelman 2014.
''' y = [44.6 57.76 49.91 61.34 49.60 61.79 48.95 44.70 59.17 53.94... 46.55 54.74 50.27 51.24 46.32]'; growth = [2.4 2.89 .85 4.21 3.02 3.62 1.08 -.39 3.86 2.27 .38 1.04 2.36 1.72 .1]'; X = [ones(size(y)) , growth];
model = { 'data {' ' int N;' ' int J;' ' vector[N] y;' ' matrix[N,J] X;' '}' 'parameters {' ' vector[J] b;' ' real<lower=0> sigma;' '}' 'model {' ' y ~ normal(X*b, sigma); // data model' ' increment_log_prob(-log(sigma)); // log prior for p(sigma) propto 1/sigma' '}' 'generated quantities {' ' vector[N] log_lik;' ' for (n in 1:N){' ' log_lik[n] <- normal_log(y[n], X[n]*b, sigma);' ' }' '}' };
data = struct('N',numel(y),'J',size(X,2),'X',X,'y',y);
hibbs_fit = stan('model_code',model,'data',data,'verbose',true);
log_lik = hibbs_fit.extract.log_lik;
mstan.waic(log_lik) '''