Bayes Factor using NUTS sampler

I have created a hierarchical model that can be scaled in dimensions. Since there are multiple combinations of the right sets of parameters I tend to archive better results with NUTS sampler instead of SMC. (example below)
For Model Selection I want to use the Bayes Factor. I know that the marginal_likelihood is accessible for SMC sampled traces via trace.report.log_marginal_likelihood and that there is an old Post where it has been accessed via model.marginal_likelihood. Is it possible to adapt this while remaining with NUTS sampling? The last way seems not to work with 3.9.3

For completeness, my symplified model:

traces = []
time = np.linspace(0,2,2000) #any timegrid
for order in range(1:4)
    with pm.model as model:
            A = pm.Uniform("A",0,5,shape=order)
            B = pm.HalfNormal("T", 3,shape=order)
            C = pm.HalfNormal("C", 1)
            link = C
            for idx in range(order):
                link += A[idx] * np.exp(time / B[idx])
            like = pm.StudentT("like", nu=len(time), mu=link,observed=myData)  
            traces.append(pm.sample(10_000, target_accept=0.9, tune=2000))

The goal is now to compare models of different orders to find the most favored one.

I’m fairly new to pymc and Bayesian statistics at all, so every critique is appreciated! Also if there are better solutions to solve this “multimodal problem” since sampling jumps in between two valid but close parameters. (Nested sampling is mentioned in one of my papers for Evidence computations but not implemented in pymc3)

Is there a particular reason you have to use Bayes Factor? For model comparison method like pm.loo is much better Model comparison — PyMC3 3.10.0 documentation

Thx for that quick reply! I might be able to convince the supervisors from alternative methods (which I would offer anyways) - but to be able to compare with already published papers it would be helpful to provide a Bayes Factor.

You can try bridge sampling (i have discussed with @Olaf in Fix Bridge sampling for variables with shape>1?).

1 Like

THX again, it’s working very smooth with my data!

1 Like