A general measure of the impact of priors in Bayesian statistics via Stein's Method

TL;DR

This paper introduces a new measure based on Wasserstein distance to quantify how different prior choices affect posterior distributions in Bayesian models, demonstrated on common statistical models.

Contribution

It proposes a general, quantifiable measure of prior impact in Bayesian analysis using Stein's Method and Wasserstein distance, applicable across various models.

Findings

Effective in measuring prior influence on posteriors.

Applicable to normal, Binomial, and Poisson models.

Provides a unified framework for prior impact assessment.

Abstract

We propose a measure of the impact of any two choices of prior distributions by quantifying the Wasserstein distance between the respective resulting posterior distributions at any fixed sample size. We illustrate this measure on the normal, Binomial and Poisson models.

Equations

$$d_{\mathcal{W}}(P_1,P_2)=\underset{h\in \mathcal{H}}{\sup }|\mathrm{E}[h(X_1)]-\mathrm{E}[h(X_2)]|$$

$$p_i(\theta ;x)={\kappa }_i(x)p_i(\theta )\ell (x;\theta ),i=1,2,$$

$${\tau }_i(\theta ;x)=\frac{1}{p_i(\theta ;x)}{\int }_{a_i}^{\theta }({\mu }_i-y)p_i(y;x)dy,i=1,2,$$

$$\rho (\theta )=\frac{p_2(\theta )}{p_1(\theta )}.$$

$$|{\mu }_1-{\mu }_2|=\frac{|\mathrm{E}[{\tau }_1({\Theta }_1;x){\rho }^{\prime }({\Theta }_1)]|}{\mathrm{E}[\rho ({\Theta }_1)]}\le d_{\mathcal{W}}(P_1,P_2)\le \frac{\mathrm{E}[{\tau }_1({\Theta }_1;x)|{\rho }^{\prime }({\Theta }_1)|]}{\mathrm{E}[\rho ({\Theta }_1)]}$$

$$|{\mu }_1-{\mu }_2|\le d_{\mathcal{W}}(P_1,P_2)\le ||{\rho }^{\prime }|{|}_{\infty }\frac{\mathrm{Var}[{\Theta }_1]}{\mathrm{E}[\rho ({\Theta }_1)]}$$

$$\mathrm{E}[{\tau }_1({\Theta }_1;x)|{\rho }^{\prime }({\Theta }_1)|]\frac{{\kappa }_2(x)}{{\kappa }_1(x)}.$$

$$1={\int }_{\Theta }{p}_2(\theta ;x)d\theta ={\int }_{\Theta }\frac{{\kappa }_2(x)}{{\kappa }_1(x)}\rho (\theta )p_1(\theta ;x)d\theta =\frac{{\kappa }_2(x)}{{\kappa }_1(x)}\mathrm{E}[\rho ({\Theta }_1].$$

$$\frac{|\mathrm{E}[{\tau }_1({\Theta }_1;x){\rho }^{\prime }({\Theta }_1)]|}{\mathrm{E}[\rho ({\Theta }_1)]}$$

$$d_{\mathcal{W}}(P_1,P_2)=\frac{\mathrm{E}[{\tau }_1({\Theta }_1;x)|{\rho }^{\prime }({\Theta }_1)|]}{\mathrm{E}[\rho ({\Theta }_1)]}.$$

$$x\mapsto \frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(x-\mu {)}^2}{2{\sigma }^2}\right),x\in \mathbb{R},$$

$$\ell(x;\sigma^{2})=(2\pi\sigma^{2})^{-\frac{n}{2}}\exp\Bigg{\{}-\frac{1}{2}\sum_{i=1}^{n}{\frac{(x_{i}-\mu)^{2}}{\sigma^{2}}}\Bigg{\}}.$$

$$\sigma^{2}\mapsto(\sigma^{2})^{-\frac{n}{2}-1}\exp\Bigg{\{}-\frac{1}{2}\sum_{i=1}^{n}{\frac{(x_{i}-\mu)^{2}}{\sigma^{2}}}\Bigg{\}}.$$

$$\sigma^{2}\mapsto\frac{\beta^{\alpha}}{\Gamma(\alpha)}(\sigma^{2})^{-\alpha-1}\exp\Bigg{\{}-\frac{\beta}{\sigma^{2}}\Bigg{\}}$$

$$\rho ({\sigma }^2)=\frac{p_2({\sigma }^2)}{p_1({\sigma }^2)}\propto \frac{({\sigma }^2{)}^{-\alpha -1}\exp \left(-\frac{\beta }{{\sigma }^2}\right)}{({\sigma }^2{)}^{-1}}=({\sigma }^2{)}^{-\alpha }\exp \left(-\frac{\beta }{{\sigma }^2}\right).$$

$$d_{\mathcal{W}}(P_1,P_2)$$

$${\rho }^{\prime }({\sigma }^2)$$

$$\mathrm{E}[\rho ({\Theta }_1)]$$

$$\mathrm{E}[{\tau }_1({\Theta }_1)|{\rho }^{\prime }({\Theta }_1)|]$$

$$\frac{\mathrm{E}[{\tau }_1({\Theta }_1)|{\rho }^{\prime }({\Theta }_1)|]}{\mathrm{E}[\rho ({\Theta }_1)]}$$

$$\frac{\left|\alpha {\sum }_{i=1}^n(x_i-\mu {)}^2-(\frac{n}{2}-1)\beta \right|}{(\frac{n}{2}+\alpha -1)(\frac{n}{2}-1)}\le d_{\mathcal{W}}(P_1,P_2)\le \frac{\left(\frac{1}{2}{\sum }_{i=1}^n(x_i-\mu {)}^2\alpha +\frac{n}{2}\beta +(2\alpha -1)\beta \right)}{(\frac{n}{2}+\alpha -1)(\frac{n}{2}-1)}.$$

$$x\mapsto \left(\genfrac{}{}{0ex}{}{n}{x}\right){\theta }^x(1-\theta {)}^{n-x}$$

$$p_1(\theta )=\frac{1}{\theta (1-\theta )}.$$

$$\theta \mapsto \left(\genfrac{}{}{0ex}{}{n}{x}\right){\theta }^{x-1}(1-\theta {)}^{n-x-1}$$

$$\rho (\theta )=\frac{p_2(\theta )}{p_1(\theta )}\propto {\theta }^{\alpha }(1-\theta {)}^{\beta }.$$

$$d_{\mathcal{W}}(P_1,P_2)\ge \left|\frac{x}{n}-\frac{\alpha +x}{n+\alpha +\beta }\right|=\frac{\left|n\alpha -(\alpha +\beta )x\right|}{n(n+\alpha +\beta )}.$$

$$d_{\mathcal{W}}(P_1,P_2)\le \frac{1}{n}\left(\alpha +(\beta -\alpha )\frac{\alpha +x}{\alpha +\beta +n}\right)$$

$$x\mapsto \frac{\exp (-\theta ){\theta }^x}{x!}$$

$$\ell (x;\theta )=\frac{\exp (-n\theta ){\theta }^{{\sum }_{i=1}^n{x}_i}}{{\prod }_{i=1}^n{x}_i!}.$$

$$\rho (\theta )=\frac{p_2(\theta )}{p_1(\theta )}\propto {\theta }^{{\alpha }_2-{\alpha }_1}\exp \left(-({\beta }_2-{\beta }_1)\theta \right).$$

Full text

††journal: Journal of LaTeX Templates

1 Introduction

A key question in Bayesian analysis is the effect of the prior on the posterior, and how we can measure this effect. Will the posterior distributions derived with distinct priors become very similar if more and more data are gathered? It has been proved formally in Diaconis and Freedman [1] and Diaconis and Freedman [2] that, under certain regularity conditions, the impact of the prior is waning as the sample size increases. From a practical viewpoint it is more important to know what happens at finite sample size $n$. Recently, Ley et al. [3] have provided a partial answer to this question by investigating the Wasserstein distance between the posterior distribution based on a given prior of interest and the no-prior data-only based posterior, which allows detecting at fixed sample size $n$ the effect of the prior of interest. This distance being mostly impossible to calculate exactly, they have provided sharp upper and lower bounds. This work is, to the best of the authors’ knowledge, the first to quantify at any sample size the prior effect. However, it strongly relies on the assumption that one prior is the flat uniform prior (or data-only prior), and hence it does not allow a direct comparison between two general priors.

Our aim in the present paper is to fill this gap and to extend the methodology of Ley et al. [3] to the setting where we can quantify the distributional distance between posteriors resulting from any two priors. As in Ley et al. [3], we opt here for the Wasserstein-1 distance defined as

[TABLE]

for $X_{1}$ and $X_{2}$ random variables with respective distribution functions $P_{1}$ and $P_{2}$, and where $\mathcal{H}$ stands for the class of Lipschitz-1 functions. We shall provide sharp lower and upper bounds on this distance for any two posteriors. For the upper bound case, our approach relies on a variant of Stein’s density approach Stein et al. [4], see also Ley et al. [5] for a recent account.

The present paper is organized as follows. In Section 2 we provide the notations and terminology used throughout the paper, state and prove our main theorem. Then in Section 3 we illustrate the strength of our new general measure of the difference between two priors by considering several examples of well-known distributions. In each case, we compare the effects of distinct priors on the posterior distribution of the parameter we are interested in. The choices of our priors are motivated by research papers that have discussed various choices of viable priors for a certain parameter.

2 A general method to compare the effects of two priors

We start by fixing our notations. We consider independent and identically distributed observations $X_{1},\ldots,X_{n}$ from a parametric model with parameter of interest $\theta\in\Theta\subseteq\mathbb{R}$. Our framework allows the common distribution of the $X_{i}$’s to be continuous as well as discrete. For the sake of simplicity, we shall denote the sampling density or likelihood of $X_{1},\ldots,X_{n}$ by $\ell(x;\theta)$ where $x=(x_{1},\ldots,x_{n})$ are the observed values. Now consider two distinct (possibly improper) prior densities $p_{1}(\theta)$ and $p_{2}(\theta)$ for $\theta$. The two resulting posterior densities for $\theta$ then take on the guise

[TABLE]

where $\kappa_{1}(x),\kappa_{2}(x)$ are normalizing constants. We write $(\Theta_{1},P_{1})$ and $(\Theta_{2},P_{2})$ the couples of random variables and cumulative distribution functions corresponding respectively to the densities $p_{1}(\theta;x)$ and $p_{2}(\theta;x)$.

Suppose that both posterior distributions have finite means $\mu_{1}$ and $\mu_{2}$. For our general result we borrow a quantity from the Stein Method literature, namely the Stein kernel $\tau_{i}$ of $P_{i}$ defined as

[TABLE]

where $a_{i}$ is the lower bound of the support $I_{i}=(a_{i},b_{i})$ of $p_{i}$. One can readily see that this function is always positive and vanishes at the boundaries of the support. For more information about Stein kernels, we refer the reader to Ley et al. [3]; for the present paper no further knowledge than the definition (1) is necessary.

The key element in our developments is the fact that the densities $p_{1}(\theta;x)$ and $p_{2}(\theta;x)$ are nested, meaning that one support is included in the other. We here suppose that ${I}_{2}\subseteq{I}_{1}$ which allows us to write $p_{2}(\theta;x)=\frac{\kappa_{2}(x)}{\kappa_{1}(x)}\rho(\theta)p_{1}(\theta;x)$ with

[TABLE]

With these notations in hand, we obtain the following theorem which is a general way of comparing the two priors $p_{1}(\theta)$ and $p_{2}(\theta)$ and follows from Stein’s Method for nested densities Ley et al. [3].

Theorem 1**.**

Consider $\mathcal{H}$ the set of Lipschitz-1 functions on $\mathbb{R}$. Assume that $\theta\mapsto\rho(\theta)$ is differentiable on ${I}_{2}$ and satisfies (i) ${\rm E}[|\Theta_{1}-\mu_{1}|\rho(\Theta_{1})]<\infty$, (ii) $\left(\rho(\theta)\int_{a_{1}}^{\theta}(h(y)-{\rm E}[h(\Theta_{1})])p_{1}(y;x)dy\right)^{\prime}$ is integrable for all $h\in\mathcal{H}$ and (iii) $\lim_{\theta\rightarrow a_{2},b_{2}}\rho(\theta)\int_{a_{1}}^{\theta}(h(y)-{\rm E}[h(\Theta_{1})])p_{1}(y;x)dy=0$ for all $h\in\mathcal{H}$. Then

[TABLE]

and, if the variance of $\Theta_{1}$ exists,

[TABLE]

where $||\cdot||_{\infty}$ stands for the infinity norm.

Proof. Our proof relies on Theorem 3.1 of Ley et al. [3] about Stein’s Method for nested densities. Let us start by establishing the upper bounds in (2) and (3). The function $\pi_{0}$ from Ley et al. [3] coincides in our case with $\frac{\kappa_{2}(x)}{\kappa_{1}(x)}\rho(\theta)$, yielding our conditions on $\rho(\theta)$ and the upper bound

[TABLE]

The final form of the upper bound in (2) is obtained by noticing that

[TABLE]

The upper bound in (3) follows from ${\rm E}[\tau_{1}(\Theta_{1};x)|\rho^{\prime}(\Theta_{1})|]\leq||\rho^{\prime}||_{\infty}{\rm E}[\tau_{1}(\Theta_{1};x)]$ and from the property of a Stein kernel to satisfy ${\rm E}[\tau_{1}(\Theta_{1};x)\varphi^{\prime}(\Theta_{1})]={\rm E}[(\Theta_{1}-\mu_{1})\varphi(\Theta_{1})]$ for all differentiable function $\varphi$ for which ${\rm E}[(\Theta_{1}-\mu_{1})\varphi(\Theta_{1})]$ exists, implying in particular that ${\rm E}[\tau_{1}(\Theta_{1};x)]={\rm Var}[\Theta_{1}]$.

The lower bound $\frac{|{\rm E}[\tau_{1}(\Theta_{1};x)\rho^{\prime}(\Theta_{1})]|}{{\rm E}[\rho(\Theta_{1})]}$ is obtained along the same lines as the upper bound. It remains to show that it equals the simplified expression $|\mu_{1}-\mu_{2}|$. Using again the above-mentioned Stein kernel property, we get

[TABLE]

$\square$

We attract the reader’s attention to the fact that the lower bound $|\mu_{1}-\mu_{2}|$ can be readily obtained by noticing that the identity function $h(x)=x$ belongs to the class of Lipschitz-1 functions. For the comparison of lower and upper bounds it is however relevant to provide as well the equivalent expression of the lower bound in terms of the Stein kernel and $\rho$. Indeed, a close examination of both upper and lower bounds reveals that, if $\rho$ is a monotone increasing or decreasing function, the bounds do coincide, leading to

Corollary 1**.**

If in addition to the conditions of Theorem 1 we assume that the ratio $\rho$ is monotone increasing or decreasing, then

[TABLE]

We now briefly discuss the relevance of the bounds from Theorem 1. The impact of a prior can best be measured by considering its effect on the resulting posterior distribution, and our result precisely allows us to quantify the difference in terms of Wasserstein distance between two posteriors resulting from distinct priors. As can be appreciated from (2), our bounds are sharp since they contain the same quantities in both the upper and lower bounds; this fact is further underlined by the equality in Corollary 1. The special case of one prior being the improper prior $p_{1}(\theta)=1$ allows us to retrieve the measure of the influence of the prior given in Ley et al. [3].

3 Comparison of various priors for various distributions

We now illustrate the strength of Theorem 1 by comparing popular choices of priors for parameters of three famous distributions, namely the normal, the Binomial and the Poisson models. We will work out in detail the case of the normal distribution. As we shall see, the bounds that we get allow us to conclude that, in all these cases, the difference between the resulting priors vanishes at the speed of $n^{-1}$ independently of the observations, but that at finite sample size $n$ the observations do play a role.

3.1 Priors for the scale parameter of the normal distribution

Our first example concerns the normal $N(\mu,\sigma^{2})$ distribution with probability density function

[TABLE]

with location parameter $\mu\in\mathbb{R}$ and dispersion parameter $\sigma^{2}>0$. Our parameter of interest here is $\sigma^{2}$, and we consider $\mu$ to be fixed. Our motivation for studying various priors for $\sigma^{2}$ comes from Kavetski et al. [6] who consider a storm depth multiplier model to represent rainfall uncertainty, where the errors appear under multiplicative form and are assumed to be normal. They fix the mean $\mu$ but state that “less is understood about the degree of rainfall uncertainty, i.e., the multiplier variance” and therefore study various priors for $\sigma^{2}$. The two priors that they investigate are Jeffreys’ prior and an Inverse Gamma prior, and we shall here show how our methods allow to directly measure the difference in impacts between these two priors.

Consider the setting where $x=(x_{1},\ldots,x_{n})$ is a random sample from the $N(\mu,\sigma^{2})$ population with fixed location $\mu$. The likelihood function $\ell(x;\sigma^{2})$ of the normal model can be factorized into

[TABLE]

The Jeffreys’ prior is a popular prior as it is invariant under reparametrization and is proportional to the square-root of the Fisher information associated with the parameter of interest. In the present setting it is proportional to $\frac{1}{\sigma^{2}}$ which leads to the posterior density being proportional to

[TABLE]

The latter density can be recognized to be an Inverse Gamma distribution with parameters $(\frac{n}{2},\frac{1}{2}\sum_{i=1}^{n}(x_{i}-\mu)^{2})$. Now, it is easy to see that the Inverse Gamma prior

[TABLE]

with positive real parameters $(\alpha,\beta)$ yields, in combination with the likelihood, as posterior also an Inverse Gamma distribution with parameters $(\frac{n}{2}+\alpha,\frac{1}{2}\sum_{i=1}^{n}(x_{i}-\mu)^{2}+\beta)$. We consider Jeffreys’ prior as first prior $P_{1}$ and the Inverse Gamma prior as second prior $P_{2}$, leading to the ratio

[TABLE]

With this in hand, Conditions (ii) and (iii) of Theorem 1 are readily fulfilled. Condition (i) is equivalent to requiring that $n/2+\alpha-1>0$ which is trivially the case as soon as we have at least two observations. Therefore, denoting by $\mu_{1}$ and $\mu_{2}$ the means of $P_{1}$ and $P_{2}$, respectively, we can calculate the lower bound as follows:

[TABLE]

In order to obtain the upper bound, we need to calculate

[TABLE]

and, writing $\Theta_{1}$ the random variable associated with $P_{1}$,

[TABLE]

The Stein kernel for the Inverse Gamma distribution with parameters $(\frac{n}{2},\frac{1}{2}\sum_{i=1}^{n}(x_{i}-\mu)^{2})$ corresponds to $\tau_{1}(\sigma^{2})=\frac{(\sigma^{2})^{2}}{\frac{n}{2}-1}$. Therefore we have

[TABLE]

Putting the ends together, we thus get as upper bound

[TABLE]

We attract here the reader’s attention to the fact that the same proportionality constant appears on the numerator and denominator, which is why we we left it out from our calculations for the sake of readability.

Finally we find that the Wasserstein distance between the posteriors based on the Jeffreys’ and the Inverse Gamma priors for $\sigma^{2}$ is bounded as follows:

[TABLE]

We see that both the lower and upper bounds are of the order $O(n^{-1})$, which is a quantification of the well-known fact that the impact of the choice of the prior wanes asymptotically. Clearly, when $\alpha=\beta=0$, the distance equals zero which is natural as then the two priors coincide. It is interesting to remark that the observations are centered around $\mu$ in the bounds, as this variation is precisely what the parameter $\sigma^{2}$ measures. Finally we remark that in the calculation of ${\rm E}[\tau_{1}(\Theta_{1})|\rho^{\prime}(\Theta_{1})|]$ we have used the triangular inequality for $\left|-\alpha\Theta_{1}+\beta\right|$, while alternative upper boundings or direct calculations with absolute values would lead to a sharper upper bound. However, those will not be as readable as our bounds, which anyway are of the desired order of $n^{-1}$.

3.2 Priors for the success parameter in a Binomial model

Many random phenomena worth studying have binary outcomes and therefore can be modelled using the famous Binomial distribution $Bin(n,\theta)$ with probability mass function

[TABLE]

where $x\in\{1,\ldots,n\}$ is the number of observed successes, the natural number $n$ denotes the number of binary experiments and $\theta\in(0,1)$ stands for the success parameter. In this part, we suppose $n$ to be fixed and we focus on the success parameter. Ghaderinezhad [7] studies the influence of various priors, including the Beta and Haldane priors, on the posterior by bounding the Wasserstein distance between each resulting posterior and the data-only posterior obtained via the uniform prior. We will show how to compare the Haldane and Beta priors directly by means of our Theorem 1.

The Haldane prior is an improper prior representing complete uncertainty and its expression is given by

[TABLE]

The resulting posterior density is of the form

[TABLE]

which represents a Beta distribution with parameters $(x,n-x)$. In order to be a proper distribution it requires at least one success and one failure in the sequence of experiments. Choosing as second prior the conjugate Beta prior with parameters $\alpha,\beta>0$ yields the $Beta(\alpha+x,n-x+\beta)$ posterior. The ratio of prior densities corresponds here to

[TABLE]

One readily sees that the conditions of Theorem 1 are satisfied. The lower bound thus is given by

[TABLE]

The Stein kernel for the Beta distribution with parameters $(x,n-x)$ is $\tau_{1}(\theta)=\frac{\theta(1-\theta)}{n}$. Proceeding as in Section 3.1, we get the upper bound

[TABLE]

It can be seen both the lower and upper bound are of the order of $O(n^{-1})$. This rate of convergence remains even in the extreme cases $x=0$ and $x=n$.

3.3 Priors for the Poisson model

As last example we consider the most popular count distribution, namely the Poisson distribution with probability mass function

[TABLE]

where $x\in N$ is the number of occurrences and the parameter $\lambda>0$ indicates the average number of events in a given time interval. Classical choices of priors are an exponential prior, Jeffreys’ prior or the conjugate Gamma prior. The latter was used in e.g. the study of asthma mortality rates by Gelman et al. [8]. We opt here to compare in all generality two Gamma priors with non-negative real parameters $(\alpha_{1},\beta_{1})$ and $(\alpha_{2},\beta_{2})$, respectively, as this contains as special cases the exponential prior ($\alpha=1$), the uniform prior ($\alpha=1,\beta=0$) and even the Jeffreys’ prior which is proportional to $\theta^{-1/2}$ ($\alpha=1/2,\beta=0$).

The likelihood function $\ell(x;\theta)$ for a set of sampling data $x=(x_{1},\ldots,x_{n})$ coming from a Poisson model reads

[TABLE]

In combination with a conjugate Gamma prior $\frac{\beta_{j}^{\alpha_{j}}}{\Gamma(\alpha_{j})}\theta^{\alpha_{j}-1}\exp(-\beta_{j}\theta)$ this yields a Gamma posterior distribution with updated parameters $(\sum_{i=1}^{n}x_{i}+\alpha_{j},n+\beta_{j})$, $j=1,2$. Now, a quantity of special interest is the ratio

[TABLE]

While it nicely satisfies the conditions of Theorem 1, certain combinations of the prior parameters lead to a monotone increasing or decreasing behavior, and hence the applicability of Corollary 1. Simple manipulations reveal that these parameter combinations are given by $\alpha_{1}<\alpha_{2}\cap\beta_{1}>\beta_{2}$ (increasing) and $\alpha_{1}>\alpha_{2}\cap\beta_{1}<\beta_{2}$ (decreasing), where in each case one of the strict inequalities can turn into a potential equality. The Stein kernel of the Gamma distribution with parameters $(\sum_{i=1}^{n}x_{i}+\alpha_{1},n+\beta_{1})$ being $\tau_{1}(\theta)=\frac{\theta}{n+\beta_{1}}$, straightforward calculations lead to

[TABLE]

and

[TABLE]

and hence

[TABLE]

This exact distance is of the desired order of $O(n^{-1})$. In case of parameter combinations that lead to no monotone behavior of the ratio $\rho$ we can build upper and lower bounds as done in the previous sections; we leave this as exercise for the reader.

ACKNOWLEDGMENTS:

This research is supported by a BOF Starting Grant of Ghent University.