0
$\begingroup$

I try to calculate the marginal likelihood of the example in the article " a tutorial on bridge sampling", which is estimating the marginal likelihood for a binomial model assuming a uniform prior on the rate parameter θ (i.e., the beta-binomial model). Consider a single participant who answered k = 2 out of n = 10 true/false questions correctly. The analytical result is 1/(n+1), that is 1/11.

And Then I get a stanfit with rstan, and use bridge_sampler(stanfit) to get the log marginal likelihood. After transforming the log marginal likelihood to the marginal likelihood, I find it is totally different from 1/11. I don't know the reason for this. besides, there is a warning message: “effective sample size cannot be calculated, has been replaced by number of samples.”

stanmodelcode <- "
data {
    int<lower=0> N; 
    int y;
}

parameters {
    real<lower=0,upper=1> theta; 
}

model {
    theta ~ beta(1,1);
    y ~ binomial(N,theta) ;
}"

y <- 2
dat <- list(N = 15, y = y);

stanfit <- stan(model_code = stanmodelcode,data=dat,chains=1 , iter=50000 , warmup=500)

library(bridgesampling)
(Lik_marginal <- bridge_sampler(stanfit))
exp(logml(Lik_marginal))
$\endgroup$

1 Answer 1

0
$\begingroup$

The model specifications should all be in the form "target += X" instead of "y ~ X" because the latter drops constants.

So you should have "target += binomial_lpmf(y | N, theta);" instead of "y~binomial(N,theta);". And the same with the theta distribution.

see: https://rdrr.io/cran/bridgesampling/f/vignettes/bridgesampling_stan_ttest.Rmd

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.