Generating log_lik for Regime Switching Model help

Modeling
1

Hello,
I have stan code for Regime Switching Model from Jim Savage site. Here this site:
https://khakieconomics.github.io/2018/02/24/Regime-switching-models.html.
I have an issue in creating the log-likelihood in the generated quantities. I am not sure how to define the target syntax in the generated quantities block.
Here this code.

data {
    int T;
    vector[T] y;
}
parameters {
    vector<lower = 0, upper = 1>[2] p;
    real<lower = 0> phi;
    vector[2] delta;
    vector<lower = 0>[2] sigma;
    real<lower = 0, upper = 1> xi1_init; 
    real y_tm1_init;
}
transformed parameters {
    matrix[T, 2] eta;
    matrix[T, 2] xi;
    vector[T] f;

// fill in etas
for(t in 1:T) {
    eta[t,1] = exp(normal_lpdf(y[t]| delta[1], sigma[1]));
if(t==1) {
    eta[t,2] = exp(normal_lpdf(y[t]| delta[2] + phi * y_tm1_init, sigma[2]));
} else {
    eta[t,2] = exp(normal_lpdf(y[t]| delta[2] + phi * y[t-1], sigma[2]));
}
}

// work out likelihood contributions

for(t in 1:T) {
// for the first observation
if(t==1) {
    f[t] = p[1]*xi1_init*eta[t,1] + // stay in state 1
            (1 - p[1])*xi1_init*eta[t,2] + // transition from 1 to 2
            p[2]*(1 - xi1_init)*eta[t,2] + // stay in state 2 
            (1 - p[2])*(1 - xi1_init)*eta[t,1]; // transition from 2 to 1

xi[t,1] = (p[1]*xi1_init*eta[t,1] +(1 - p[2])*(1 - xi1_init)*eta[t,1])/f[t];
xi[t,2] = 1.0 - xi[t,1];

} else {
// and for the rest

f[t] = p[1]*xi[t-1,1]*eta[t,1] + // stay in state 1
        (1 - p[1])*xi[t-1,1]*eta[t,2] + // transition from 1 to 2
        p[2]*xi[t-1,2]*eta[t,2] + // stay in state 2 
        (1 - p[2])*xi[t-1,2]*eta[t,1]; // transition from 2 to 1

// work out xi

xi[t,1] = (p[1]*xi[t-1,1]*eta[t,1] +(1 - p[2])*xi[t-1,2]*eta[t,1])/f[t];

// there are only two states so the probability of the other state is 1 - prob of the first
xi[t,2] = 1.0 - xi[t,1];
}
}

}
model {
// priors
p ~ beta(10, 2);
phi ~ normal(1, .1);
delta ~ normal(0, .1);
sigma ~ cauchy(0, 1);
xi1_init ~ beta(2, 2);
y_tm1_init ~ normal(0, .1);


// likelihood is really easy here!
target += sum(log(f));
}
"

Thanks for any help you can provide.

1 Like
2

Sorry, can’t answer now, maybe @mcol is not busy and can help?

3

I’m sorry, I’m not familiar with this type of models and I’m very short on time. Perhaps the answer to this post will set you on the right path. Otherwise I suppose @James_Savage will have the answer ready.

4

If total log likelihood is sum(log(f)) then individual terms of the log likelihood are log(f)

generated quantities {
  vector [T] log_lik = log(f);
}
5

Thank you for the attention Mr. @avehtari.
I have read your article entitled “WAIC and cross-validation in Stan”. Is the log_lik that you proposed in that formula (log_lik = log(f)) has the same mean with log-predictive density (for computing waic)?

6

In the first message you mention log-likelihood and the code wrote computes log-likelihood.

The computed term is \log p(y_i|\theta^{(s)}) which both the log likelihood evaluated with given certain parameter values \theta^{(s)} and the log predictive density given certain parameter values \theta^{(s)} and certain observation y_i. See more in Practical Bayesian model evaluation using leave-one-out cross-validation and
WAIC. It’s better to use PSIS-LOO than WAIC, as PSIS-LOO is more robust and has better diagnostics as demonstrated in that paper.