# On integral priors for multiple comparison in Bayesian model selection

Diego Salmerón∗, Juan Antonio Cano∗111∗ Juan Antonio Cano (1956-2018) was
Professor in the Department of Statistics and Operations Research at the
University of Murcia. Dr. Cano was a dear mentor and friend who contributed
substantially to the theory of integral priors and was instrumental to the
developments of this article.
and Christian P. Robert†‡

∗Universidad de Murcia, Spain, †Université Paris-Dauphine, France
‡University of Warwick, United Kindom

###### Abstract

Noninformative priors constructed for estimation purposes are usually not
appropriate for model selection and testing. The methodology of integral
priors was developed to get prior distributions for Bayesian model selection
when comparing two models, modifying initial improper reference priors. We
propose a generalization of this methodology to more than two models. Our
approach adds an artificial copy of each model under comparison by
compactifying the parametric space and creating an ergodic Markov chain across
all models that returns the integral priors as marginals of the stationary
distribution. Besides the garantee of their existance and the lack of
paradoxes attached to estimation reference priors, an additional advantage of
this methodology is that the simulation of this Markov chain is
straightforward as it only requires simulations of imaginary training samples
for all models and from the corresponding posterior distributions. This
renders its implementation automatic and generic, both in the nested case and
in the nonnested case.

Keywords: Integral priors, Bayesian model selection, Objective Bayes factor,
Markov chains.

##  1 Introduction

There exists an extensive literature on the choice and variety of
noninformative prior distributions for estimation. The most famous of these
prior distributions are the Jeffreys prior (Jeffreys (1961)), and the
reference priors (Bernardo (1979), Berger and Bernardo (1992)). These default
or reference methods usually yield an improper prior distribution, which is
determined up to a positive multiplicative constant. Very often this is not a
serious issue for estimation because the constant cancels in the posterior
distribution. However, in model selection problems where we have two models
under consideration, $M_{i}:$ $\\{f_{i}(x\mid\theta_{i})$,
$\theta_{i}\in\Theta_{i}\\}$, $i=1,2$, and improper priors
$\pi_{i}^{N}(\theta_{i})=c_{i}h_{i}(\theta_{i})$, $i=1,2$, the Bayes factor

$B_{21}^{N}(x)=\frac{c_{2}\int
f_{2}(x\mid\theta_{2})h_{2}(\theta_{2})\mathrm{d}\theta_{2}}{c_{1}\int
f_{1}(x\mid\theta_{1})h_{1}(\theta_{1})\mathrm{d}\theta_{1}},$

depends on the arbitrary ratio $c_{2}/c_{1}$. An attempt towards solving this
issue consists of using intrinsic priors (Berger and Pericchi (1996)) that use
part of the data to turn the prior into a proper posterior. Intrinsic priors
are the solution $\\{\pi_{1}^{I}(\theta_{1}),\pi_{2}^{I}(\theta_{2})\\}$ of a
system of two functional equations. Under certain conditions, and assuming
that $M_{1}$ is nested in $M_{2}$, the system of functional equations reduces
to a single equation with two functional unknowns, in such a way that once an
arbitrary choice is made for $\pi_{1}^{I}(\theta_{1})$, the dual prior
$\pi_{2}^{I}(\theta_{2})$ is automatically determined. The solution
$\\{\pi_{1}^{I}(\theta_{1}),\pi_{2}^{I}(\theta_{2})\\}$ proposed by Moreno et
al. (1998) consists of choosing
$\pi_{1}^{I}(\theta_{1})=\pi_{1}^{N}(\theta_{1})$ and hence

$\pi_{2}^{I}(\theta_{2})=\pi_{2}^{N}(\theta_{2})\mathbb{E}_{f_{2}(z\mid\theta_{2})}(B_{12}^{N}(Z))=\int\pi_{2}^{N}(\theta_{2}\mid
z)m_{1}^{N}(z)\mathrm{d}z,$

where $z$ denotes an imaginary training sample,
$B_{12}^{N}(z)=m_{1}^{N}(z)/m_{2}^{N}(z)$, and $m_{i}^{N}(z)=\int
f_{i}(z\mid\theta_{i})\pi_{i}^{N}(\theta_{i})\mathrm{d}\theta_{i}$, $i=1,2$.
See Pérez and Berger (2002) for more details. However, intrinsic Bayes factors
are often not suitable in the non-nested case, see Berger and Mortera (1999)
and Moreno (2005). Another proposal consists of using fractional priors, which
are again defined as the solution to a system of functional equations (see
Moreno (1997)). However, this methodology has received less attention and also
faces difficulties (Robert, 2001).

Pérez and Berger (2002) have proposed a novel approach under the name of
expected posterior priors distributions. This method can be applied in the
general scenario where we have $q\geq 2$ models under consideration. The
expected posterior priors are defined as

$\pi_{i}^{*}(\theta_{i})=\int\pi_{i}^{N}(\theta_{i}\mid
z)m^{*}(z)\mathrm{d}z,\,\,\,i=1,\dots,q,$

where $m^{*}(z)$ is an arbitrary predictive density, which can be improper,
for the imaginary training sample $z$. When $M_{1}$, say, is nested in the
other models, Pérez and Berger (2002) propose $m^{*}(z)=m_{1}^{N}(z)$, and
therefore this method applied to two nested-encompassing models produces the
intrinsic priors introduced by Moreno et al. (1998).

Other approaches have been proposed for prior construction towards Bayesian
model selection, among which Bayarri and García-Donato (2008), Bayarri et al.
(2012), and Fouskakis et al. (2015). A review of reference priors for model
selection can be found in Consonni et al. (2018).

###  1.1 Integral priors

Integral priors have been introduced in Cano et al. (2008) for model selection
when two models, $M_{i}:$ $\\{f_{i}(x\mid\theta_{i})$,
$\theta_{i}\in\Theta_{i}\\}$, $i=1,2$, are under consideration. These prior
distributions, $\pi_{1}(\theta_{1})$ and $\pi_{2}(\theta_{2})$, are defined as
the solution of a system of integral equations, and can be considered as a
generalization of the expected posterior priors introduced in Pérez and Berger
(2002). Concretely, integral priors
$\\{\pi_{1}(\theta_{1}),\pi_{2}(\theta_{2})\\}$ satisfy the integral equations

$\pi_{1}(\theta_{1})=\int\pi_{1}^{N}(\theta_{1}\mid
z_{1})m_{2}(z_{1})\mathrm{d}z_{1},$
$\pi_{2}(\theta_{2})=\int\pi_{2}^{N}(\theta_{2}\mid
z_{2})m_{1}(z_{2})\mathrm{d}z_{2},$

where $\pi_{i}^{N}(\theta_{i})$ is a prior distribution for estimation
(typically a noninformative prior), $\pi_{i}^{N}(\theta_{i}\mid z)\propto
f_{i}(z\mid\theta_{i})\pi_{i}^{N}(\theta_{i})$, $m_{i}(z)=\int
f_{i}(z\mid\theta_{i})\pi_{i}(\theta_{i})\mathrm{d}\theta_{i}$, $i=1,2$, and
$z_{1}$ and $z_{2}$ are (minimal) imaginary training samples for $\theta_{1}$
and $\theta_{2}$, respectively. When $\theta_{i}$ or $z_{i}$ are discrete, the
corresponding integrals are replaced by sums with no loss of generality. These
integral equations balance each model with respect to the other one since the
prior $\pi_{i}(\theta_{i})$ is derived from the marginal $m_{j}(z_{i})$, and
therefore from $\pi_{j}(\theta_{j})$, $j\neq i$, as an unknown generalized
expected posterior prior.

Integral priors can be used to compare both nested and non-nested models
without the need to make this distinction, and the application of this
methodology does not need to specify a predictive distribution as with the
expected posterior priors approach.

Cano et al. (2008) showed that $\pi_{1}(\theta_{1})$ and $\pi_{2}(\theta_{2})$
are integral priors if and only if $\pi_{1}(\theta_{1})$ and
$\pi_{2}(\theta_{2})$ are the invariant $\sigma$-finite measures of two Markov
chains, $(\theta_{1}^{t})$ and $(\theta_{2}^{t})$, respectively. For instance,
the transition $\theta_{1}^{\prime}\rightarrow\theta_{1}$ of the Markov chain
$(\theta_{1}^{t})$ is made of the following steps:

  1. 1.

Generate a training sample $z_{2}\sim f_{1}(z_{2}\mid\theta_{1}^{\prime})$

  2. 2.

Generate $\theta_{2}\sim\pi_{2}^{N}(\theta_{2}\mid z_{2})$

  3. 3.

Generate a training sample $z_{1}\sim f_{2}(z_{1}\mid\theta_{2})$

  4. 4.

Generate $\theta_{1}\sim\pi_{1}^{N}(\theta_{1}\mid z_{1})$.

Therefore, the study of integral priors can be performed considering the
stochastic stability of these Markov chains, as in Cano et al. (2007a), Cano
et al. (2007b), Cano and Salmerón (2013), Salmerón et al. (2015), Cano and
Salmerón (2016), and Cano et al. (2018).

In some situations the Markov chains for integral priors are recurrent Markov
chains, even positive Harris recurrent, in which case the integral priors are
proper prior distributions. This implies that the Bayes factor $B_{21}(x)$ is
defined without ambiguity (Robert, 2001) and can be estimated by the Monte
Carlo estimator

$\hat{B}_{21}(x)=\frac{\sum_{t=1}^{T}f_{2}(x\mid\theta_{2}^{t})}{\sum_{t=1}^{T}f_{1}(x\mid\theta_{1}^{t})}$

associated with the above Markov chains.

In the general case Cano and Salmerón (2013) have proposed to consider
constrained imaginary training samples $z$ when the simulation of the chains
is performed, which guarantees the existence and the uniqueness of invariant
probability measures. This approach has been implemented successfully in
several situations, see Cano and Salmerón (2013), and Cano et al. (2018).
However, the simulation of constrained imaginary training samples can involve
high computing time, specially when the accept-reject method is implemented.
In this article we therefore propose a new approach.

The plan of this paper is as follows. In Section 2 we generalize the
definition of integral priors for more than two models. In Section 3 we
propose to add copies of the models under consideration but with compact
parametric spaces, and we show that integral priors are proper priors. Section
4 is dedicated to show how the methodology applies to the problem of testing
if the mean of a normal population is negative, zero, or positive. This
multiple non-nested hypothesis testing problem has previously been studied in
Berger and Mortera (1999) assuming known variance when using both intrinsic
and fractional Bayes factor. However, as pointed out by Berger and Mortera
(1999), these approaches do not correspond to genuine Bayes factors attached
to specific prior distributions. In Section 5 we apply our methodology to the
one way ANOVA model. Finally, in Section 6 we present some relevant
conclusions and outline incoming research.

##  2 The Markov chains associated with integral priors

When considering only two models, the integral priors are the invariant
$\sigma$-finite measures of two Markov chains, $(\theta_{1}^{t})$ and
$(\theta_{2}^{t})$. In this section we extend the definition of the integral
priors when comparing more than two models, through a generalization of the
above Markov chains.

As a starter, it proves convenient to represent the simulation of the dual
chains $(\theta_{1}^{t})$ and $(\theta_{2}^{t})$ as shown in Figure 1, where
each single step $\theta_{i}\rightarrow\theta_{j}$, $i\neq j$, is carried out
simulating $z_{j}$ (an imaginary training sample for $\theta_{j}$ given
$\theta_{i}$), and then simulating $\theta_{j}\sim\pi_{j}^{N}(\theta_{j}\mid
z_{j})$, that is, simulating from the transition

$\pi_{ij}(\theta_{i},\theta_{j})=\int
f_{i}(z_{j}\mid\theta_{i})\pi_{j}^{N}(\theta_{j}\mid z_{j})\mathrm{d}z_{j}.$
(1)

More specifically, given an initial point $\theta_{1}^{1}\in\Theta_{1}$, we
carry out the transitions
$\theta_{1}^{1}\rightarrow\theta_{2}^{1}\rightarrow\theta_{1}^{2}\rightarrow\theta_{2}^{2}\,\cdots$
obtaining two marginal Markov chains, $(\theta_{1}^{t})$ and
$(\theta_{2}^{t})$.

$\theta_{2}^{t}$$\theta_{1}^{t}$$\theta_{2}^{t+1}$$\theta_{1}^{t+1}$ Figure 1:
Transition of the Markov chains for two models under comparison.

When starting with an initial point $\theta_{2}^{1}\in\Theta_{2}$, the process
maintains the same stochastic behavior.

However, with more than two models under comparison, there exist several
possibilities for combining simulations from $\pi_{ij}(\theta_{i},\theta_{j})$
in a sequence of simulation steps, towards obtaining Markov chains on each
parametric space. For instance, consider the following three models involving
the exponential distribution. Under model $M_{i}$ the data
$x=(x_{1},\dots,x_{n})$ are drawn from the exponential distribution

$M_{i}:\,f_{i}(x_{j}\mid\theta_{i})=\theta_{i}^{-1}\exp(-x_{j}/\theta_{i}),\,\,\,j=1,\dots,n,\,\,\,\theta_{i}\in\Theta_{i},\,\,\,i=1,2,3,$
(2)

the difference being that $\Theta_{1}=(0,1)$, $\Theta_{2}=(1,+\infty)$, and
$\Theta_{3}=(0,+\infty)$, and $\pi_{i}^{N}(\theta_{i})\propto\theta_{i}^{-1}$,
$i=1,2,3$. Figure 2 shows two possibilities to combine single steps. Each
possibility produces a Markov chain on $\Theta_{i}$, $i=1,2,3$.

Possibility
1$\theta_{1}^{t}$$\theta_{2}^{t}$$\theta_{3}^{t}$$\theta_{1}^{t+1}$$\theta_{2}^{t+1}$$\theta_{3}^{t+1}$Possibility
2$\theta_{1}^{t}$$\theta_{3}^{t}$$\theta_{2}^{t}$$\theta_{1}^{t+1}$$\theta_{3}^{t+1}$$\theta_{2}^{t+1}$
Figure 2: Two possibilities to combine single steps when comparing three
models.

For $i=1,2,3$, we have simulated both the Markov chain $(\theta_{i,1}^{t})$
corresponding to the first possibility (left), and the Markov chain
$(\theta_{i,2}^{t})$ corresponding to the second possibility (right). If
$(\theta_{3,1}^{t})$ and $(\theta_{3,2}^{t})$ had the same stationary
distribution, then the kernel density estimates from the simulations
$(\log\theta_{3,1}^{t})$ and $(\log\theta_{3,2}^{t})$ should be quite similar,
but Figure 3 demonstrates a clear difference between both. Hence the specific
way in which the single steps are combined can produce different answers, with
the consequent problem of defining integral priors for more than two models in
a consistent manner.

Figure 3: Kernel density estimates from 100,000 simulations of the Markov
chains $(\log\theta_{3,1}^{t})$ and $(\log\theta_{3,2}^{t})$.

A straightforward proposal that makes this issue vanish is to randomly choose
the ordering in which the steps are performed: if we have carried out the
steps
$\theta_{i_{1}}^{t}\rightarrow\theta_{i_{2}}^{t}\rightarrow\theta_{i_{3}}^{t}$
for a given sequence $(i_{1},i_{2},i_{3})$, then we select at random a
sequence $(j_{1},j_{2},j_{3})$ with the condition that $j_{1}\neq i_{3}$, and
then we perform the steps
$\theta_{i_{3}}^{t}\rightarrow\theta_{j_{1}}^{t+1}\rightarrow\theta_{j_{2}}^{t+1}\rightarrow\theta_{j_{3}}^{t+1}$,
as shown in Figure (4).

$\theta_{i_{1}}^{t}$$\theta_{i_{2}}^{t}$$\theta_{i_{3}}^{t}$$\theta_{j_{1}}^{t+1}$$\theta_{j_{2}}^{t+1}$$\theta_{j_{3}}^{t+1}$
Figure 4: Definition for three models.

This produces three Markov chains, $(\theta_{1}^{t})$, $(\theta_{2}^{t})$ and
$(\theta_{3}^{t})$, and leads us to the following definition.

######  Definition.

The integral priors $\pi_{1}(\theta_{1})$, $\pi_{2}(\theta_{2})$ and
$\pi_{3}(\theta_{3})$, are the $\sigma$-finite invariant measures for the
Markov chains $(\theta_{1}^{t})$, $(\theta_{2}^{t})$ and $(\theta_{3}^{t})$
defined above, respectively.

This definition reproduces the original one (with two models) and naturally
extends to more than three models.

Based on this new definition we have simulated the Markov chains associated
with the corresponding integral priors for the models
$\\{M_{1},M_{2},M_{3}\\}$ in (2). The histograms corresponding to the chains
$(\log\theta_{i}^{t})$, $i=1,2,3$, are shown in Figure (5). It indicates that
the chains on the contrained parameter spaces accumulate some weight near the
boundary, while the unconstrained chain explores both positive and negative
parts of its support.

Figure 5: Histograms based on 100,000 simulations from the integral priors for
$\log\theta_{i}$, $i=1,2,3$, form left to right, using the associated Markov
chains.

##  3 Proper integral priors

As indicated above, in some model comparisons cases involving two models, the
Markov chains associated with integral priors are recurrent Markov chains. In
the general case Cano and Salmerón (2013) have proposed to consider
constrained imaginary training samples $z$ when the simulation of the chains
is performed, which guarantees the existence and the uniqueness of invariant
probability measures. Briefly, their practical recommendation consists in
keeping the imaginary training samples constrained within an interval of $\pm
5s$ about the sample mean, where $s$ is the sample standard deviation. This
approach has been implemented successfully in several situations comparing two
models, see Cano and Salmerón (2013) and Cano et al. (2018). However, since
the simulation of constrained imaginary training samples is carried out using
some accept-reject method (for instance, a truncated normal distribution, see
Robert (1995)), this induces a consequent computational overhead. The
computational alternative developed in this article consists in adding a copy
of each model with a corresponding compact parametric space.

Consider for instance the case of two models under comparison,
$M_{i}:f_{i}(x\mid\theta_{i})$, $\theta_{i}\in\Theta_{i}$, $i=1,2$, and
attached improper priors $\pi_{1}^{N}(\theta_{1})$ and
$\pi_{2}^{N}(\theta_{2})$. We introduce model $M_{3}$, a copy of model $M_{1}$
with compact parametric space, that is,

$M_{3}:f_{1}(x\mid\theta_{3}),\,\,\,\theta_{3}\in\Theta_{3},$

where $\Theta_{3}$ is a compact subset of $\Theta_{1}$ such that
$\pi_{3}^{N}(\theta_{3})\propto\pi_{1}^{N}(\theta_{3})\mathbb{I}_{\Theta_{3}}(\theta_{3})$
is a proper prior. Our goal is to obtain integral priors for both $\theta_{1}$
and $\theta_{2}$, but we can call the generalized methodology of constructing
integral priors for the three models $\\{M_{1},M_{2},M_{3}\\}$, which produces
the Markov chains $(\theta_{i}^{t})$ with transition densities
$Q_{i}(\theta_{i}\mid\theta_{i}^{\prime})$, and resulting integral priors
$\pi_{i}(\theta_{i})$, that is,

$\pi_{i}(\theta_{i})=\int
Q_{i}(\theta_{i}\mid\theta_{i}^{\prime})\pi_{i}(\theta_{i}^{\prime})\mathrm{d}\theta_{i}^{\prime},\,\,\,i=1,2,3.$

Note that we are using the conditional notation for the transition density
$Q_{i}(\theta_{i}\mid\theta_{i}^{\prime})$, while the alternative notation
$\pi_{ij}(\theta_{i},\theta_{j})$ corresponds for the single step transition
$\theta_{i}\rightarrow\theta_{j}$, $i\neq j$.

Since $\pi_{3}^{N}(\theta_{3})$ is a proper prior, we need no imaginary
training sample for $\theta_{3}$. In this case

$\pi_{13}(\theta_{1},\theta_{3})=\pi_{23}(\theta_{2},\theta_{3})=\pi_{3}^{N}(\theta_{3}).$

Thus, the Markov chain $(\theta_{3}^{t})$ is a sequence of $\mathrm{i.i.d.}$
random variables, $\theta_{3}^{t}\sim\pi_{3}^{N}(\theta_{3})$, and therefore
$\pi_{3}(\theta_{3})=\pi_{3}^{N}(\theta_{3})$ is the integral prior for
$\theta_{3}$.

$\theta_{i_{1}}^{t}$$\theta_{i_{2}}^{t}$$\theta_{i_{3}}^{t}$$\theta_{j_{1}}^{t+1}$$\theta_{j_{2}}^{t+1}$$\theta_{j_{3}}^{t+1}$$\theta_{k_{1}}^{t+2}$$\theta_{k_{2}}^{t+2}$$\theta_{k_{3}}^{t+2}$
Figure 6: 2-step transition for three models

The study of the stochastic stability of the correponding Markov chains
$(\theta_{1}^{t})$ and $(\theta_{2}^{t})$ entails a slightly greater level of
complexity. The simpler approach is based on the 2-step transition densities

$Q_{i}^{2}(\theta_{i}\mid\theta_{i}^{\prime})=\int
Q_{i}(\theta_{i}\mid\tilde{\theta_{i}})Q_{i}(\tilde{\theta}_{i}\mid\theta_{i}^{\prime})\mathrm{d}\tilde{\theta}_{i},\,\,\,i=1,2.$

For $i=1,2$, both the transition density
$Q_{i}(\theta_{i}\mid\theta_{i}^{\prime})$ and the 2-step transition density
$Q_{i}^{2}(\theta_{i}\mid\theta_{i}^{\prime})$, are mixtures of conditional
densities due to the randomness on the order of the sequences, see Figure 6.
Each component of the mixture $Q_{i}^{2}(\theta_{i}\mid\theta_{i}^{\prime})$
corresponds to a specific realisation of the sequences $(i_{1},i_{2},i_{3})$,
$(j_{1},j_{2},j_{3})$, and $(k_{1},k_{2},k_{3})$. However, all components of
this mixture share the very helpful property that a simulation from
$\pi_{3}^{N}(\theta_{3})$ occurs at some point, since $j=3$ is in the sequence
$(j_{1},j_{2},j_{3})$. For instance, with regards to the mixture
$Q_{1}^{2}(\theta_{1}\mid\theta_{1}^{\prime})$, the configuration 1 in Figure
7 corresponds to the component

$\int\pi_{3}^{N}(\theta_{3})\pi_{31}(\theta_{3},\theta_{1})\mathrm{d}\theta_{3},$

and the configuration 2 corresponds to the component

$\int\pi_{3}^{N}(\theta_{3})\pi_{32}(\theta_{3},\theta_{2})\pi_{21}(\theta_{2},\tilde{\theta}_{1})\pi_{12}(\tilde{\theta}_{1},\theta_{2}^{\prime})\pi_{21}(\theta_{2}^{\prime},\theta_{1})\mathrm{d}\theta_{2}^{\prime}\mathrm{d}\tilde{\theta}_{1}\mathrm{d}\theta_{2}\mathrm{d}\theta_{3}.$

Configuration 1Configuration
2$\theta_{3}^{t}$$\theta_{2}^{t}$$\theta_{1}^{t}$$\theta_{2}^{t+1}$$\theta_{3}^{t+1}$$\theta_{1}^{t+1}$$\theta_{3}^{t+2}$$\theta_{1}^{t+2}$$\theta_{2}^{t+2}$$\theta_{3}^{t}$$\theta_{2}^{t}$$\theta_{1}^{t}$$\theta_{3}^{t+1}$$\theta_{2}^{t+1}$$\theta_{1}^{t+1}$$\theta_{2}^{t+2}$$\theta_{1}^{t+2}$$\theta_{3}^{t+2}$
Figure 7: Examples of 2-step transitions.

More precisely, it is straightforward to derive that all components of these
mixtures are necessarily one of four cases that can be schematically
represented as follows

$1\rightarrow 3\rightarrow 1,\,\,\,1\rightarrow 3\rightarrow 2\rightarrow 1,$
$1\rightarrow 3\rightarrow 1\rightarrow 2\rightarrow
1,\,\,\,\mathrm{and}\,\,\,1\rightarrow 3\rightarrow 2\rightarrow 1\rightarrow
2\rightarrow 1.$

Then $Q^{2}_{1}(\theta_{1}\mid\theta_{1}^{\prime})$ does not depend on
$\theta_{1}^{\prime}$, and hence
$\pi_{1}(\theta_{1})=Q^{2}_{1}(\theta_{1}\mid\theta_{1}^{\prime})$ is the
integral prior for $\theta_{1}$:

$\int\pi_{i}(\theta_{i})Q_{i}(\theta_{i}^{*}\mid\theta_{i})\mathrm{d}\theta_{i}=\int
Q_{i}^{2}(\theta_{i}\mid\theta_{i}^{\prime})Q_{i}(\theta_{i}^{*}\mid\theta_{i})\mathrm{d}\theta_{i}$
$=\int
Q_{i}(\theta_{i}\mid\tilde{\theta}_{i})Q_{i}(\tilde{\theta}_{i}\mid\theta_{i}^{\prime})Q_{i}(\theta_{i}^{*}\mid\theta_{i})\mathrm{d}\theta_{i}\mathrm{d}\tilde{\theta}_{i}$
$=\int
Q_{i}(\tilde{\theta}_{i}\mid\theta_{i}^{\prime})Q_{i}^{2}(\theta_{i}^{*}\mid\tilde{\theta}_{i})\mathrm{d}\tilde{\theta}_{i}=Q_{i}^{2}(\theta_{i}^{*}\mid\theta_{i}^{\prime})=\pi_{i}(\theta_{i}^{*}).$

Similarly, the 2-step transition for $(\theta_{2}^{t})$ is the integral prior
for $\theta_{2}$, that is,
$\pi_{2}(\theta_{2})=Q^{2}_{2}(\theta_{2}\mid\theta_{2}^{\prime})$.

On the other hand, it is clear that the resulting priors depend on the model
that is copied in a compact version. In order to overcome the resulting
arbitrariness, we propose to duplicate each model into a compact version. For
instance, when facing two models, we add to the model $M_{3}$ the model

$M_{4}:f_{2}(x\mid\theta_{4}),\,\,\,\theta_{4}\in\Theta_{4},$

where $\Theta_{4}$ is a compact subset of $\Theta_{2}$ such that
$\pi_{4}^{N}(\theta_{4})\propto\pi_{2}^{N}(\theta_{4})\mathbb{I}_{\Theta_{4}}(\theta_{4})$
is a proper prior. We consider that the training sample size for $\theta_{4}$
is zero as well. Then we apply the general methodology of integral priors to
the collection of models $\\{M_{1},M_{2},M_{3},$ $M_{4}\\}$. Thus the integral
priors for the copied models are
$\pi_{i}(\theta_{i})=\pi_{i}^{N}(\theta_{i})$, $i=3,4$. The 2-step transition
$Q^{2}_{i}(\theta_{i}\mid\theta_{i}^{\prime})$ does not depend on
$\theta_{i}^{\prime}$ and hence
$\pi_{i}(\theta_{i})=Q^{2}_{i}(\theta_{i}\mid\theta_{i}^{\prime})$ is the
integral prior for $\theta_{i}$, $i=1,2$. Besides,
$Q^{2}_{1}(\theta_{1}\mid\theta_{1}^{\prime})$ is a mixture of densities and
the components of this mixture are reduced to eight cases that can be
schematically represented as follows:

$1\rightarrow 3\rightarrow 1,\,\,\,1\rightarrow 3\rightarrow 2\rightarrow 1,$
$1\rightarrow 3\rightarrow 1\rightarrow 2\rightarrow 1,\,\,\,1\rightarrow
3\rightarrow 2\rightarrow 1\rightarrow 2\rightarrow 1,$ $1\rightarrow
4\rightarrow 1,\,\,\,1\rightarrow 4\rightarrow 2\rightarrow 1,$ $1\rightarrow
4\rightarrow 1\rightarrow 2\rightarrow 1,\,\,\,\mathrm{and}\,\,\,1\rightarrow
4\rightarrow 2\rightarrow 1\rightarrow 2\rightarrow 1,$

which means that there exist functions $H_{1}^{1}(\theta_{1},\theta_{3})$ and
$H_{1}^{2}(\theta_{1},\theta_{4})$, defined for
$\theta_{1},\theta_{3}\in\Theta_{1}$ and $\theta_{4}\in\Theta_{2}$, such that
the integral prior $\pi_{1}(\theta_{1})$ is

$\pi_{1}(\theta_{1})=Q_{1}^{2}(\theta_{1}\mid\theta_{1}^{\prime})=\int\pi_{3}^{N}(\theta_{3})H_{1}^{1}(\theta_{1},\theta_{3})\mathrm{d}\theta_{3}+\int\pi_{4}^{N}(\theta_{4})H_{1}^{2}(\theta_{1},\theta_{4})\mathrm{d}\theta_{4}.$

Similarly, there exist functions $H_{2}^{1}(\theta_{2},\theta_{3})$ and
$H_{2}^{2}(\theta_{2},\theta_{4})$, defined for
$\theta_{2},\theta_{4}\in\Theta_{2}$ and $\theta_{3}\in\Theta_{1}$, such that

$\pi_{2}(\theta_{2})=Q_{2}^{2}(\theta_{2}\mid\theta_{2}^{\prime})=\int\pi_{3}^{N}(\theta_{3})H_{2}^{1}(\theta_{2},\theta_{3})\mathrm{d}\theta_{3}+\int\pi_{4}^{N}(\theta_{4})H_{2}^{2}(\theta_{2},\theta_{4})\mathrm{d}\theta_{4},$

is the integral prior for $\theta_{2}$.

Obviously, this approach of adding copycat models with a compact parametric
space applies to more than two models as follows. We star from $q\geq 2$
models

$M_{i}:f_{i}(x\mid\theta_{i}),\,\,\,\theta_{i}\in\Theta_{i},\,\,\,i=1,\dots,q,$

under comparison and associated improper prior distributions
$\pi_{i}^{N}(\theta_{i})$ $i=1,\dots,q$. Given a compact subset
$\Theta_{q+i}\subseteq\Theta_{i}$ such that
$\int_{\Theta_{q+i}}\pi_{i}^{N}(\theta_{i})\mathrm{d}\theta_{i}<+\infty$, we
consider a copy of $M_{i}$ but with parametric space $\Theta_{q+i}$, that is,

$M_{q+i}:f_{i}(x\mid\theta_{q+i}),\,\,\,\theta_{q+i}\in\Theta_{q+i},$

and
$\pi_{q+i}^{N}(\theta_{q+i})\propto\pi_{i}^{N}(\theta_{q+i})\mathbb{I}_{\Theta_{q+i}}(\theta_{q+i})$,
$i=1,\dots,q$. Then integral priors for $\\{M_{1},\dots,M_{q},\dots,M_{2q}\\}$
do exist, they are unique and proper priors, and the integral prior for
$\theta_{i}$ satisfies

$\pi_{i}(\theta_{i})=Q_{i}^{2}(\theta_{i}\mid\theta_{i}^{\prime})=\sum_{j=1}^{q}\int\pi_{q+j}^{N}(\theta_{q+j})H_{i}^{j}(\theta_{i},\theta_{q+j})\mathrm{d}\theta_{q+j},$

where the function $H_{i}^{j}(\theta_{i},\theta_{q+j})$ is defined for
$\theta_{i}\in\Theta_{i}$ and $\theta_{q+j}\in\Theta_{j}$, $i,j=$ $1,\dots,q$.

###  3.1 About the compact sets

Without loss of generality consider two models under comparison, $M_{1}$ and
$M_{2}$, with their copies over compact parametric spaces, $M_{3}$ and
$M_{4}$, respectively. It is quite sensible to seek compact sets such that
that the posterior probabilities

$P_{x}^{3}=\int_{\Theta_{3}}\pi_{1}^{N}(\theta_{1}\mid
x)\mathrm{d}\theta_{1}\,\,\,\mathrm{and}\,\,\,P_{x}^{4}=\int_{\Theta_{4}}\pi_{2}^{N}(\theta_{2}\mid
x)\mathrm{d}\theta_{2}\,$

are close to one. While there are a variety of ways to choose these compact
sets, we opt for the simple solution where both $\Theta_{3}$ and $\Theta_{4}$
are a Cartesian products of compact intervals. More specifically, our proposal
is made of credible intervals with posterior probability $1-\alpha$ and
$\alpha\in(0,1)$ small, typically $\alpha\leq 0.05$. The following example
illustrates the behavior of this method.

######  Example 1.

For the data $x=(x_{1},\dots,x_{n})$, consider the models

$M_{1}:f_{1}(x_{j}\mid\theta_{1})=\frac{1}{\theta_{1}}\exp(-x_{j}/\theta_{1}),\,\,\,j=1,\dots,n,\,\,\,\theta_{1}\in(0,1),$

and

$M_{2}:f_{2}(x_{j}\mid\theta_{2})=\frac{1}{\theta_{2}}\exp(-x_{j}/\theta_{2}),\,\,\,j=1,\dots,n,\,\,\,\theta_{2}\in(1,+\infty),$

with $\pi_{i}^{N}(\theta_{i})=1/\theta_{i}$, $i=1,2$. For a sample of size
$n=15$ and mean $0.7$, we implement the above methodology of integral priors,
thus defining the collection of models $\\{M_{1},M_{2},M_{3},M_{4}\\}$, where
$\Theta_{i+2}=[a_{i+2},b_{i+2}]$ is a $100(1-\alpha)$% credible interval based
on the quantiles $\alpha/2$ and $1-\alpha/2$ of $\pi_{i}^{N}(\theta_{i}\mid
x)$, $i=1,2$. We also opt for using training samples $z$ of size one for both
$\theta_{1}$ or for $\theta_{2}$. Simulation from the posterior distributions
is straightforward:

$\theta_{1}=1/(1-z^{-1}\log u_{1}),\,\,\,u_{1}\sim\mathcal{U}(0,1),$

and

$\theta_{2}=-z/\log(u_{2}(1-e^{-z})+e^{-z}),\,\,\,u_{2}\sim\mathcal{U}(0,1).$

For $\alpha\in\\{0.05,0.01,0.005,0.001\\}$, Table (1) provides
$\hat{B}_{12}(x)$, the approximate values of $B_{12}(x)$, as well as the
estimated posterior probability of model $M_{1}$ if the prior probabilities of
$M_{1}$, $M_{2}$, $M_{3}$ and $M_{4}$ are $1/2$, $1/2$, $0$, and $0$,
respectively. Hence, $\hat{P}(M_{1}\mid
x)=\hat{B}_{12}(x)/(1+\hat{B}_{12}(x))$. It is worth noting that under the
unconstrained exponential model

$f(x_{j}\mid\theta)=\theta^{-1}\exp(-x_{j}/\theta),\,\,\,j=1\dots,n,\,\,\,\theta\in(0,+\infty),$

and $\pi^{N}(\theta)\propto 1/\theta$, the posterior probability of
$\Theta_{1}=(0,1)$ is 0.89.

Table 1: Estimation of the Bayes factor $B_{12}(x)$ and posterior probability $P(M_{1}\mid x)$, based on $10^{6}$ simulations from the integral priors. $\alpha$ | $a_{3}$ | $b_{3}$ | $a_{4}$ | $b_{4}$ | $\hat{B}_{12}(x)$ | $\hat{P}(M_{1}\mid x)$
---|---|---|---|---|---|---
0.05 | 0.44 | 0.97 | 1.00 | 1.62 | 8.98 | 0.90
0.01 | 0.39 | 0.99 | 1.00 | 1.90 | 9.10 | 0.90
0.005 | 0.37 | 1.00 | 1.00 | 2.02 | 9.09 | 0.90
0.001 | 0.34 | 1.00 | 1.00 | 2.28 | 9.16 | 0.90

The histograms corresponding to the chains $(\log\theta_{i}^{t})$,
$i=1,2,3,4$, using the associated Markov chains with the compact sets for
$\alpha=0.001$ (last row in Table (1)), are shown in Figure 8.

Figure 8: Histograms based on $10^{6}$ simulations from the integral priors
for $\log\theta_{i}$, $i=1,2$ (first row: $\log\theta_{1}$ left,
$\log\theta_{2}$ right), and $i=3,4$ (second row: $\log\theta_{3}$ left,
$\log\theta_{4}$ right).

##  4 Testing the mean under dependent samples

In the specific context of paired samples in which each subject is measured
twice, resulting in pairs of observations, there are three natural hypotheses
regarding the mean difference, $\mu$, namely, $\mu<0$, $\mu=0$, and $\mu>0$.
This multiple non-nested hypothesis testing problem has been previously
studied in Berger and Mortera (1999), assuming known variance and comparing
the encompassing arithmetic intrinsic Bayes factor (EIBF), the encompassing
expected intrinsic Bayes factor (EEIBF), the median intrinsic Bayes factor,
and the fractional Bayes factor. Unfortunately, as pointed out by the authors,
these approaches do not correspond to Bayes factors with respect to genuine
prior distributions. In this section, we therefore apply the theory of
integral priors for this testing problem.

Under model $M_{i}$ the difference data $x=(x_{1},\dots,x_{n})$ are
independently drawn from the Normal distribution
$\mathcal{N}(\mu_{i},\sigma_{i}^{2})$, $i=1,2,3$, with $\mu_{1}<0$,
$\mu_{2}=0$, and $\mu_{3}>0$, that is,

$f_{i}(x\mid\mu_{i},\sigma_{i})=(2\pi)^{-n/2}\sigma_{i}^{-n}\exp\left(-\frac{1}{2\sigma_{i}^{2}}\sum_{j=1}^{n}(x_{j}-\mu_{i})^{2}\right),\,\,\,i=1,3,$
$f_{2}(x\mid\sigma_{2})=(2\pi)^{-n/2}\sigma_{2}^{-n}\exp\left(-\frac{1}{2\sigma_{2}^{2}}\sum_{j=1}^{n}x_{j}^{2}\right),$

and we consider the default estimation priors
$\pi_{i}^{N}(\theta_{i})=1/\sigma_{i}$, $i=1,2,3$, where
$\theta_{1}=(\mu_{1},\sigma_{1})$, $\theta_{2}=\sigma_{2}$, and
$\theta_{3}=(\mu_{3},\sigma_{3})$.

To investigate this setting, we produced a simulated sample $x$ of size $n=15$
from the Normal distribution with mean $-1$ and variance $4$. The observed
sample mean and standard deviation were $-0.927$ and $2.530$, respectively. We
thus consider copies $M_{4}$, $M_{5}$, and $M_{6}$, with compact parametric
spaces $\Theta_{i}$, $i=4,5,6$, of the form

$\theta_{4}=(\mu_{4},\sigma_{4})\in\Theta_{4}=[a_{4},b_{4}]\times[c_{4},d_{4}],$
$\theta_{5}=\sigma_{5}\in\Theta_{5}=[c_{5},d_{5}],$

and

$\theta_{6}=(\mu_{6},\sigma_{6})\in\Theta_{6}=[a_{6},b_{6}]\times[c_{6},d_{6}]$

respectively.

For models $M_{1}$ and $M_{3}$, imaginary training samples are samples, $z$,
of size 2, and the corresponding posterior distributions are

$\pi_{i}^{N}(\mu_{i},\sigma_{i}\mid
z)\propto\sigma_{i}^{-3}\exp\left(-\frac{s_{z}^{2}}{2\sigma_{i}^{2}}\right)\exp\left(-\frac{(\bar{z}-\mu_{i})^{2}}{\sigma_{i}^{2}}\right),\,\,\,i=1,3,$

where $\bar{z}$ and $s_{z}\neq 0$ are the sample mean and standard deviation
of $z$, respectively. The posterior distribution $\pi_{i}^{N}(\mu_{i}\mid z)$
is proportional to

$\int_{0}^{+\infty}\sigma_{i}^{-3}\exp\left(-\frac{s_{z}^{2}}{2\sigma_{i}^{2}}\right)\exp\left(-\frac{(\bar{z}-\mu_{i})^{2}}{\sigma_{i}^{2}}\right)\mathrm{d}\sigma_{i}=\frac{1}{s_{z}^{2}+2(\bar{z}-\mu_{i})^{2}}\,,$

an half-Cauchy density, and therefore, to simulate from
$\pi_{i}^{N}(\mu_{i}\mid z)$ we can generate $u,v\sim\mathcal{U}(0,1)$ and
take

$\mu_{1}=\bar{z}-\frac{s_{z}}{\sqrt{2}}\cot\left(\frac{u}{2}(\pi-2\arctan(\sqrt{2}\bar{z}/s_{z}))\right)$

and

$\mu_{3}=\bar{z}-\frac{s_{z}}{\sqrt{2}}\cot\left(\frac{v-1}{2}(\pi+2\arctan(\sqrt{2}\bar{z}/s_{z}))\right).$

To simulate from $\pi_{i}^{N}(\sigma_{i}\mid\mu_{i},z)$ we can generate
$\xi_{i}$ from the exponential distribution $\mathcal{E}xp(1)$ and take

$\sigma_{i}=\sqrt{\frac{s_{z}^{2}+2(\bar{z}-\mu_{i})^{2}}{2\xi_{i}}},\,\,\,i=1,3.$

For model $M_{2}$, imaginary training samples are samples of size 1, $z\neq
0$, and the posterior distribution is the inverse Gamma distribution

$\pi_{2}^{N}(\sigma_{2}\mid
z)\propto\sigma_{2}^{-2}\exp\left(-\frac{z^{2}}{2\sigma_{2}^{2}}\right).$

Therefore, to simulate from $\pi_{2}^{N}(\sigma_{2}\mid z)$ we can generate
$\xi_{2}\sim\mathcal{G}a(1/2,1)$ and take $\sigma_{2}=|z|/\sqrt{2\xi_{2}}$.

For $\alpha\in\\{0.05,0.01,0.005\\}$ we again select the compact intervals as
the $100(1-\alpha)$% credible intervals based on the quantiles $\alpha/2$ and
$1-\alpha/2$ of the marginal posterior distribution of
$\pi_{i}^{N}(\theta_{i}\mid y)$, $i=1,2,3$, that is

$1-\alpha=\int_{a_{i+3}}^{b_{i+3}}\pi_{i}^{N}(\mu_{i}\mid
x)\mathrm{d}\mu_{i},\,\,\,i=1,3,$

and

$1-\alpha=\int_{c_{i+3}}^{d_{i+3}}\pi_{i}^{N}(\sigma_{i}\mid
x)\mathrm{d}\sigma_{i},\,\,\,i=1,2,3.$

Based on such compact sets, we apply the methodology of integral priors with
$500,000$ simulations in order to approximate the posterior probabilities of
models $M_{1}$, $M_{2}$ and $M_{3}$ when selecting as set of prior
probabilities $\\{1/3,1/3,1/3,0,0,0\\}$. The etimation of the posterior
probabilities of models $M_{1}$, $M_{2}$, and $M_{3}$ are shown in Table 2,
again with very little variability when $\alpha$ decreases. Figure 9 shows the
histograms based on the simulations from the integral prior for
$\log(-\mu_{1})$ and $\log(\mu_{3})$, corresponding to $\alpha=0.005$.

Table 2: Estimation of the posterior probabilities of models $M_{1}$, $M_{2}$, and $M_{3}$ in the Normal mean difference example. $\alpha$ | $\hat{P}(M_{1}\mid x)$ | $\hat{P}(M_{2}\mid x)$ | $\hat{P}(M_{3}\mid x)$ | $P_{x}^{4}$ | $P_{x}^{5}$ | $P_{x}^{6}$
---|---|---|---|---|---|---
0.05 | 0.51 | 0.42 | 0.06 | 0.90 | 0.95 | 0.91
0.01 | 0.50 | 0.43 | 0.06 | 0.98 | 0.99 | 0.98
0.005 | 0.50 | 0.44 | 0.06 | 0.99 | 0.99 | 0.99
Figure 9: Histograms based on simulations from the integral priors for
$\log(-\mu_{1})$ and $\log(\mu_{3})$ in the Normal mean difference example.

##  5 Integral Priors for ANOVA

In the classical ANOVA setting, we now consider a sample $x_{i}$ from the
Normal population $N(\mu_{i},\sigma^{2})$, $i=1,\dots,k$, in order to handle
on the hypothesis test

$\begin{array}[]{l}H_{0}:\,\mu_{1}=\dots=\mu_{k}\\\
H_{1}:\,\mathrm{all\,the}\,\mu_{i}\,\mathrm{are\,not\,equal}.\end{array}$

Hence we are again dealing with a comparison between two models. Under model
$M_{1}$ the $k$ populations have mean $\mu$ and variance $\sigma_{1}^{2}$, and
under model $M_{2}$ the population $i$ has mean $\mu_{i}$ and variance
$\sigma_{2}^{2}$, $i=1,\dots,k$. The default (estimation) prior distributions
in this setting are $\pi_{1}^{N}(\mu,\sigma_{1})=1/\sigma_{1}$ and
$\pi_{2}^{N}(\mu_{1},\dots,\mu_{k},\sigma_{2})=1/\sigma_{2}$. Unfortunately,
these priors are improper and therefore their use is not justified for model
selection. In order to construct integral priors for ANOVA, we again consider
copies $M_{3}$ and $M_{4}$ with compact parametric spaces.

The imaginary training samples for the transitions
$\theta_{1}\rightarrow\theta_{2}$ and $\theta_{2}\rightarrow\theta_{1}$ are
generated as follows:

ITS for $\theta_{1}\rightarrow\theta_{2}$. Given $\theta_{1}=(\mu,\sigma_{1})$

  1. 1.

$j\sim\mathcal{U}\\{1,\dots,k\\}$

  2. 2.

$z_{i}\sim N(\mu,\sigma_{1}^{2})$, $i\neq j$, and
$z_{j}=(z_{j}^{1},z_{j}^{2})$, where $z_{j}^{1}$, $z_{j}^{2}\sim
N(\mu,\sigma_{1}^{2})$

  3. 3.

$z=(z_{1},\dots,z_{k})$ is the ITS

ITS for $\theta_{2}\rightarrow\theta_{1}$. Given
$\theta_{2}=(\mu_{1},\dots,\mu_{k},\sigma_{2})$

  1. 1.

$j_{1},j_{2}\sim\mathcal{U}\\{1,\dots,k\\}$

  2. 2.

$z_{j_{1}}\sim N(\mu_{j_{1}},\sigma_{2}^{2})$, $z_{j_{2}}\sim
N(\mu_{j_{2}},\sigma_{2}^{2})$

  3. 3.

$z=(z_{j_{1}},z_{j_{2}})$ is the ITS

For model $M_{1}$, when the sample size is $2$, that is $x=(z_{1},z_{2})$, the
simulation from $\pi_{1}^{N}(\mu,\sigma_{1}\mid x)$ can be carried out as
follows:

  1. 1.

$\xi\sim\mathcal{G}a(1/2,1)$, $\sigma_{1}=\frac{\mid
z_{1}-z_{2}\mid}{2\sqrt{\xi}}$

  2. 2.

$\mu\sim N(\overline{x},\sigma_{1}^{2}/2)$.

For model $M_{2}$, when the sample size is $n_{j}=2$ for some
$j\in\\{i,\dots,k\\}$ and $n_{i}=1$ for $i\neq j$, then

$\pi_{2}^{N}(\sigma_{2}\mid
x)\propto\sigma_{2}^{-2}\exp\left(-\frac{s_{j}^{2}}{2\sigma_{2}^{2}}\right),$

and therefore the simulation from
$\pi_{2}^{N}(\mu_{1},\dots,\mu_{k},\sigma_{2}\mid x)$ can be carried out as
follows:

  1. 1.

$\xi\sim\mathcal{G}a(1/2,1)$, $\sigma_{2}=\frac{\mid
z_{1}-z_{2}\mid}{2\sqrt{\xi}}$, where $x_{j}=(z_{1},z_{2})$

  2. 2.

$\mu_{i}\sim N(\overline{x}_{i},\sigma_{2}^{2}/n_{i})$, $i=1,\dots,k$.

To implement our methodology we consider $k=3$ populations. The compact
parametric spaces are chosen as product

$\Theta_{3}=[a_{3},b_{3}]\times[c_{3},d_{3}]$

and

$\Theta_{4}=[a_{41},b_{41}]\times[a_{42},b_{42}]\times[a_{43},b_{43}]\times[c_{4},d_{4}].$

For $\alpha=0.05$ and $\alpha=0.01$, we select the intervals as the
$100(1-\alpha)\%$ credible intervals. For sample sizes $n_{1}=n_{2}=n_{3}=10$,
the realizations of simulated sample means are
$\overline{x}_{1}=\overline{x}_{2}=0$,
$\overline{x}_{3}\in\\{0,0.75,1,1.5\\}$, while the realizations of the sample
standard deviations are $s_{1}=s_{2}=s_{3}=1$. Table 3 produces estimated
posterior probability $\hat{P}(M_{1}\mid x)$ under the prior probabilities
$P(M_{1})=P(M_{2})=1/2$ and $P(M_{3})=P(M_{4})=0$.

Table 3: Estimations of the posterior probabilities of model $M_{1}$ in the ANOVA example. | $\overline{x}_{3}$ | $P_{x}^{3}$ | $P_{x}^{4}$ | $\hat{P}(M_{1}\mid x)$
---|---|---|---|---
$\alpha=0.05$ | 0 | 0.90 | 0.82 | 0.91
0.75 | 0.90 | 0.82 | 0.67
1 | 0.90 | 0.82 | 0.37
1.5 | 0.90 | 0.82 | 0.03
$\alpha=0.01$ | 0 | 0.98 | 0.96 | 0.91
0.75 | 0.98 | 0.96 | 0.68
1 | 0.98 | 0.96 | 0.39
1.5 | 0.98 | 0.96 | 0.03

##  6 Conclusions and oncoming research

We advocated in this paper a principled generalization of the methodology of
integral priors (Cano et al., 2008) when two or more models are to be
compared. Our approach proves to be quite automated and generic, both in the
nested and in the nonnested cases, with the only requirement being the
existence of training samples of finite size for all models under comparison.
Since we did not compare here the outcome of more approximate solutions like
intrinsic and fractional pseudo Bayes factors (Robert, 2001), further
experiments could be conducted in that spirit. The examples we considered so
far allow to approximate the Bayes factor with the Monte Carlo estimator, but
extensions to more general settings are directly available when using evidence
approximations (Friel and Wyse, 2011), like bridge sampling, inverse logistic
regression or even harmonic means since some models (Marin and Robert, 2011)
are associated with compact parameter spaces.

Acknowledgements The research of D Salmerón was supported by the Fundación
SéNeCa-Agencia de Ciencia y Tecnología de la Región de Murcia Program for
Excellence in Scientific Research (project 20862/PI/18). CP Robert is
partially funded by the European Union (ERC-2022-SYG-OCEAN-101071601) and by a
Prairie chair from the Agence Nationale de la Recherche (ANR-19-P3IA-0001).

## References

  * Bayarri and García-Donato (2008) Bayarri MJ, García-Donato G (2008). Generalization of Jeffreys Divergence-Based Priors for Bayesian Hypothesis Testing. Journal of the Royal Statistical Society. Series B (Statistical Methodology), 70, 981-1003.
  * Bayarri et al. (2012) Bayarri MJ, Berger JO, Forte A, García-Donato G (2012). Criteria for Bayesian model choice with application to variable selection. The Annals of Statistics, 40, 1550-1577.
  * Berger and Bernardo (1992) Berger JO, Bernardo JM (1992). On the development of reference priors (with discussion). In Bayesian Statistics 4, 35-60. Oxford University Press, London.
  * Berger and Mortera (1999) Berger JO, Mortera J (1999). Default Bayes Factors for Nonnested Hypothesis Testing. Journal of the American Statistical Association, 94(446), 542–554.
  * Berger and Pericchi (1996) Berger JO, Pericchi LR (1996). The intrinsic Bayes factor for model selection and prediction. Journal of the American Statistical Association, 91, 109-122.
  * Bernardo (1979) Bernardo JM (1979). Reference posterior distributions for Bayesian inference (with discussion). Journal of the Royal Statistical Society, Series B, 41: 113-147.
  * Cano et al. (2018) Cano JA, Iniesta M, Salmerón D (2018). Integral priors for Bayesian model selection: how they operate from simple to complex cases. TEST 27, 968–987.
  * Cano et al. (2007a) Cano JA, Kessler M, Salmerón D (2007a). Integral priors for the one way random effects model. Bayesian Anal 2(1):59–68.
  * Cano et al. (2007b) Cano JA, Kessler M, Salmerón D (2007b). A synopsis of integral priors for the one way random effects model. In: Bernardo JM et al (eds) Bayesian statistics 8. Oxford University Press, Oxford, pp 577–582.
  * Cano et al. (2008) Cano JA, Salmerón D, Robert CP (2008). Integral equation solutions as prior distributions for Bayesian model selection. Test 17(3):493–504.
  * Cano and Salmerón (2013) Cano JA, Salmerón D (2013). Integral priors and constrained imaginary training samples for nested and non-nested Bayesian model comparison. Bayesian Anal 8(2):361–380
  * Cano and Salmerón (2016) Cano JA, Salmerón D (2016). A review of the developments on integral priors for Bayesian model selection. Beio 32(2):96–111.
  * Consonni et al. (2018) Consonni G, Fouskakis D, Liseo B, Ntzoufras I (2018). Prior distributions for objective Bayesian analysis. Bayesian Analysis, 13, 627-679.
  * Fouskakis et al. (2015) Fouskakis D, Ntzoufras I, Draper D (2015). Power-Expected-Posterior Priors for Variable Selection in Gaussian Linear Models. Bayesian Analysis, 10, 75-107.
  * Friel and Wyse (2011) Friel N, Wyse J (2011). Estimating the evidence – a review, arXiv 1111.1957.
  * Jeffreys (1961) Jeffreys H (1961). Theory of Probability. Oxford University Press, London.
  * Marin and Robert (2011) Marin JM, Robert CP (2011). Importance sampling methods for Bayesian discrimination between embedded models. in Frontiers of Statistical Decision Making and Bayesian Analysis (M.-H. Chen, D.K. Dey, P. Müller and D. Sun, ed.), Springer, New York.
  * Moreno (1997) Moreno E (1997). Bayes Factors for Intrinsic and Fractional Priors in Nested Models. Bayesian Robustness. Lecture Notes-Monograph Series, 31, 257–270.
  * Moreno et al. (1998) Moreno E, Bertolino F, Racugno W (1998). An intrinsic limiting procedure for model selection and hypotheses testing. Journal of the American Statistical Association, 93, 1451-1460.
  * Moreno (2005) Moreno E (2005). Objective Bayesian methods for one-sided testing. Test 14, 181–198.
  * Pérez and Berger (2002) Pérez JM, Berger JO (2002). Expected posterior priors for model selection. Biometrika, 89, 491-511.
  * Robert (1995) Robert CP (1995). Simulation of truncated normal variables. Statistics and Computing, 5, 121–125.
  * Robert (2001) Robert CP (2001). The Bayesian Choice (second edition). Springer, New York.
  * Salmerón et al. (2015) Salmerón D, Cano JA, Robert CP (2015). Objective Bayesian hypothesis testing in binomial regression models with integral prior distributions. Stat Sin 25(3):1009–1023.