<p>我相信在这种情况下计算对数可能性的正确方法是:</p>
<pre><code>generated quantities {
vector[J] log_lik;
for (i in 1:J)
log_lik[i] = normal_lpdf(y[i] | theta, sigma);
}
</code></pre>
<p>之后,您可以运行以下PSI:</p>
^{pr2}$
<p>完整代码将变成:</p>
<pre><code>import pystan
from psis import psisloo
schools_code = """
data {
int<lower=0> J; // number of schools
real y[J]; // estimated treatment effects
real<lower=0> sigma[J]; // s.e. of effect estimates
}
parameters {
real mu;
real<lower=0> tau;
real eta[J];
}
transformed parameters {
real theta[J];
for (j in 1:J)
theta[j] = mu + tau * eta[j];
}
model {
eta ~ normal(0, 1);
y ~ normal(theta, sigma);
}
generated quantities {
vector[J] log_lik;
for (i in 1:J)
log_lik[i] = normal_lpdf(y[i] | theta, sigma);
}
"""
schools_dat = {'J': 8,
'y': [28, 8, -3, 7, -1, 1, 18, 12],
'sigma': [15, 10, 16, 11, 9, 11, 10, 18]}
sm = pystan.StanModel(model_code=schools_code)
fit = sm.sampling(data=schools_dat, iter=1000, chains=4)
loo, loos, ks = psisloo(fit['log_lik'])
print('PSIS-LOO value: {:.2f}'.format(loo))
</code></pre>