Bayesian Model Selection for Determining the Biokinetics of Zirconium in the Human Body after Ingestion

Multi-compartment models based on ordinary differential equations (ODEs) are well established tools for simulating complex systems evolving over time [5]. However, although they are very straightforward and easily interpretable, determining the compartment structure and interaction mechanisms can be a very daunting task. In this context, Bayesian model selection methods are an extremely useful tool for evaluating different kinds of models. In contrast to much of the frequentists' methodology, which is generally based on asymptotic approximations for large sample sizes and single best parameter value evaluations [1,10], the Bayesian framework provides a fully probabilistic approach that easily incorporates both model and parameter uncertainty [3]. Bayesian analysis is well-grounded on the so-called posterior distribution, i.e. the probability distribution of a problem specific parameter space conditioned on given data. This specifies a measure of belief for all possible parameter values. The posterior distribution is proportional to the product of the data likelihood and the parameter prior distributions, allowing to easily include a priori information into the modeling process. We here present a method for computing Bayes Factors between pairs of models for performing the task of model selection, i.e. we calculate the marginal likelihood for each of the models. This is particularly useful when there are just a few models to choose from. Due to the evaluation of each model on the whole parameter space, the Bayes factor also naturally penalizes model complexity and thus prevents overfitting issues [16,17]. The selection process is based on the inference of the model's parameters that represent the reaction rates governing the mass transport between the compartments. As it is computationally very difficult to compute the marginal likelihood of the model based on the model's parameters, we split the computation into several intermediate steps where not the full marginal likelihood has to be evaluated, but we only have to draw independent samples from the so-called power posterior [2,10]. The power posterior is similar to the general posterior density except for the introduction of an auxiliary “temperature” parameter that governs the influence of the parameter likelihood. This approach is known as thermodynamic integration [4,11]. It stabilizes the numerical evaluation of the marginal likelihood for each model by computing the expectation with respect to the temperature based power posterior. This expectation was calculated via Monte Carlo sampling using a copula based Metropolis-Hastings algorithm [18].