A First Approach to Closeness Distributions

Abstract

Probabilistic graphical models allow us to encode a large probability distribution as a composition of smaller ones. It is oftentimes the case that we are interested in incorporating in the model the idea that some of these smaller distributions are likely to be similar to one another. In this paper we provide an information geometric approach on how to incorporate this information and see that it allows us to reinterpret some already existing models. Our proposal relies on providing a formal definition of what it means to be close. We provide an example on how this definition can be actioned for multinomial distributions. We use the results on multinomial distributions to reinterpret two already existing hierarchical models in terms of closeness distributions.

1. Introduction

Bayesian modeling [1] builds on our ability to describe a given process in probabilistic terms, known as probabilistic modeling. As stated in [2]: “Statistical methods and models commonly involve multiple parameters that can be regarded as related or connected in such a way that the problem implies dependence of the joint probability model for these parameters”. Hierarchical modeling [3] is widely used for that purpose in areas such as epidemiological modeling [4] or to model oil or gas production [5]. The motivation for this paper comes from realizing that many hierarchical models can be understood, from a high level perspective, as defining a distribution over the multiple parameters that establishes that distributions which are closer to each other, are more likely. Thus, the main motivation is to start providing the mathematical tools that allow a probabilistic modeler to build hierarchical (and non-hierarchical) models starting from that geometrical concepts.

We start by introducing a simple example to illustrate the kind of problems we are interested in solving. Consider the problem of estimating a parameter $\theta$ using data from a small experiment and a prior distribution constructed from similar previous experiments. The specific problem description is borrowed from [2]:

In the evaluation of drugs for possible clinical application, studies are routinely performed on rodents. For a particular study drawn from the statistical literature, suppose the immediate aim is to estimate θ, the probability of a tumor in a population of female laboratory rats of type ‘F344’ that receive a zero dose of the drug (a control group). The data show that 4 out of 14 rats developed endometrial stromal polyps (a kind of tumor). (...) Typically, the mean and standard deviation of underlying tumor risks are not available. Rather, historical data are available on previous experiments on similar groups of rats. In the rat tumor example, the historical data were in fact a set of observations of tumor incidence in 70 groups of rats (

Table 1. Tumor incidence in 70 historical groups of rats and in the current group of rats (from [6]). The table displays the values of: (number of rats with tumors)/(total number of rats).

Previous experiments:
0/20	0/20	0/20	0/20	0/20	0/20	0/20	0/19	0/19	0/19
0/19	0/18	0/18	0/17	1/20	1/20	1/20	1/20	1/19	1/19
1/18	1/18	2/25	2/24	2/23	2/20	2/20	2/20	2/20	2/20
2/20	1/10	5/49	2/19	5/46	3/27	2/17	7/49	7/47	3/20
3/20	2/13	9/48	10/50	4/20	4/20	4/20	4/20	4/20	4/20
4/20	10/48	4/19	4/19	4/19	5/22	11/46	12/49	5/20	5/20
6/23	5/19	6/22	6/20	6/20	6/20	16/52	15/47	15/46	9/24
Current experiment: 4/14

). In the ith historical experiment, let the number of rats with tumors be $y_i$ and the total number of rats be $n_i$. We model the $y_i$’s as independent binomial data, given sample sizes $n_i$ and study-specific means ${\theta }_i$.

Example. Estimating the risk of tumor in a group of rats.

We can depict our graphical model (for more information on the interpretation of the graphical models in this paper the reader can consult [7,8]) as shown in Figure 1, where current and historical experiments are a random sample from a common population, having $\mathbf{h}$ as hyperparameters, which follow f as prior distribution. Equationally, our model can be described as:

$$\begin{array}{cc}\mathbf{h}\sim f\hfill & \end{array}$$

(1)

$$\begin{array}{ccc}{\theta }_i\sim g\left(\mathbf{h}\right)\hfill & & \hfill \forall i\in [1:71]\end{array}$$

(2)

$$\begin{array}{cccc}\hfill y_i& \sim Binomial(n_i,{\theta }_i)\hfill & & \hfill \forall i\in [1:71].\end{array}$$

(3)

The model used for this problem in [2] is the Beta-Binomial model, where g is taken to be the Beta distribution, hence $\mathbf{h}=(\alpha ,\beta )$ (see Figure 2). Furthermore, in [2] the prior f over $\alpha ,\beta$ is taken to be proportional to ${(\alpha +\beta )}^{-5/2}$, giving the model

$$\begin{array}{cc}\hfill p(\alpha ,\beta )& \propto {(\alpha +\beta )}^{-5/2}\hfill \end{array}$$

(4)

$$\begin{array}{ccc}\hfill {\theta }_i& \sim Beta(\alpha ,\beta )\hfill & \hfill \forall i\in [1:71]\end{array}$$

(5)

$$\begin{array}{ccc}\hfill y_i& \sim Binomial(n_i,{\theta }_i)\hfill & \hfill \forall i\in [1:71].\end{array}$$

(6)

The presentation of the model in [2] simply introduces the assumption that “the Beta prior distribution with parameters $(\alpha ,\beta )$ is a good description of the population distribution of the ${\theta }_i$’s in the historical experiments” without further justification. In this paper we would like to show that a large part of this model can be obtained from the intuitive idea that the probability distributions for rats with tumors in each group are similar. To do that we develop a framework for encoding as a probability distribution the assumption that two probability distributions are close to each other, and rely on information geometric concepts to model the idea of closeness.

We start by introducing the general concept of closeness distribution in Section 2. Then, we analyze the particular case in which we choose to measure remoteness between distributions by means of the Kullback Leibler divergence in the family of multinomial distributions in Section 3. The results from Section 3 are used in Section 4 to reinterpret the Beta Binomial model proposed in [2] for the rodents example, and in Section 5 to reinterpret the Hierarchical Dirichlet Multinomial model proposed by Azzimonti et al. in [9,10,11]. We are convinced that closeness distributions could play a relevant role in probabilistic modeling, allowing for more explicitly geometrically inspired probabilistic models. This paper is just a first step towards a proper definition and understanding of closeness distributions.

2. Closeness Distributions

We start by introducing the formal framework required to discuss the probability distributions. Then, we formalize what we mean by remoteness through a remoteness function, and we introduce closeness distributions as those that implement a remoteness function.

2.1. Probabilities over Probabilities

Information geometry [12] has shown us that most families of probability distributions can be understood as a Riemannian manifold. Thus, we can work with probabilities over probabilities by defining random variables which take values in a Riemannian manifold. Here, we only introduce some fundamental definitions. For a more detailed overview of measures and probability see [13], of Riemannian manifolds see [14]. Finally, Pennec provides a good overview of the probability on Riemannian manifolds in [15].

We start by noting that each manifold $\mathcal{M}$, has an associated $\sigma$-algebra, ${\mathcal{L}}_{\mathcal{M}}$, the Lebesgue $\sigma$-algebra of $\mathcal{M}$ (see Section 1, chapter XII in [16]). Furthermore, the existence of a metric g induces a measure ${\eta }_g$ (see Section 1.3 in [15]). The volume of $\mathcal{M}$ is defined as $Vol\left(\mathcal{M}\right)={\int }_{\mathcal{M}}1d{\eta }_g.$

Definition 1.

Let $(\mathsf{\Omega },\mathcal{F},P)$ be a probability space and $(\mathcal{M},g)$ be a Riemannian manifold. A random variable (referred to as a random primitive in [15]) $\mathbf{x}$ taking values in $\mathcal{M}$ is a measurable function from Ω to $\mathcal{M}.$ Furthermore, we say that $\mathbf{x}$ has a probability density function (p.d.f.) $p_{\mathbf{x}}$ (real, positive, and integrable function) if:

$$\forall \mathcal{X}\in {\mathcal{L}}_{\mathcal{M}}P(\mathbf{x}\in \mathcal{X})={\int }_{\mathcal{X}}{p}_{\mathbf{x}}d{\eta }_g,\hspace{1em}\mathit{and}\hspace{1em}P\left(\mathcal{M}\right)=1.$$

We would like to highlight that the density function $p_{\mathbf{x}}$ is intrinsic to the manifold. If $x^{\prime }=\pi \left(x\right)$ is a chart of the manifold defined almost everywhere, we obtain a random vector ${\mathbf{x}}^{\prime }=\pi \left(\mathbf{x}\right)$. The expression of $p_{\mathbf{x}}$ in this parametrization is

$$p_{{\mathbf{x}}^{\prime }}\left(x^{\prime }\right)=p_{\mathbf{x}}\left({\pi }^{-1}\left(x^{\prime }\right)\right).$$

Let $f:\mathcal{M}\to \mathbb{R}$ be a real function on $\mathcal{M}$. We define the expectation of f under $\mathbf{x}$ as

$$\mathbb{E}\left[f\left(\mathbf{x}\right)\right]=\int f\left(x\right)p_{\mathbf{x}}\left(x\right){\eta }_g\left(dx\right)$$

We have to be careful when computing $\mathbb{E}\left[f\right(\mathbf{x}\left)\right]$ so that we do it independently of the parametrization. We have to use the fact that

$$\int f\left(x\right)p_{\mathbf{x}}\left(x\right){\eta }_g\left(dx\right)=\int f\left({\pi }^{-1}\left(x^{\prime }\right)\right)p_{{\mathbf{x}}^{\prime }}\left(x^{\prime }\right)\sqrt{\mid G\left(x^{\prime }\right)\mid }d{x}^{\prime },$$

where $G\left(x^{\prime }\right)$ is the Fisher matrix at $x^{\prime }$ in the parametrization $\pi$. Hence,

$$\mathbb{E}\left[f\left(\mathbf{x}\right)\right]=\int f\left({\pi }^{-1}\left(x^{\prime }\right)\right){\rho }_{{\mathbf{x}}^{\prime }}\left(x^{\prime }\right)d{x}^{\prime }.$$

where ${\rho }_{{\mathbf{x}}^{\prime }}\left(x^{\prime }\right)=p_{{\mathbf{x}}^{\prime }}\left(x^{\prime }\right)\sqrt{\mid G\left(x^{\prime }\right)\mid }=p_{\mathbf{x}}\left({\pi }^{-1}\left(x^{\prime }\right)\right)\sqrt{\mid G\left(x^{\prime }\right)\mid }$ is the expression of $p_{\mathbf{x}}$ in the parametrization for integration purposes, that is, its expression with respect to the Lebesgue measure $d{x}^{\prime }$ instead of $d{\eta }_g.$

We note that ${\rho }_{{\mathbf{x}}^{\prime }}$ depends on the chart used whereas $p_{\mathbf{x}}$ is intrinsic to the manifold.

2.2. Formalizing Remoteness and Closeness

Intuitively, the objective of this section is to create a probability distribution over pairs of probability distributions that assigns higher probability to those pairs of probability distributions which are “close”.

We assume that we measure how distant two points are in $\mathcal{M}$ by means of a remoteness function $r:\mathcal{M}\times \mathcal{M}\to \mathbb{R}$, such that $r(x,y)\ge 0$ for each $x,y\in \mathcal{M}$. Note that r does not need to be transitive, symmetric, or reflexive.

As can be seen in Appendix A, r induces a total order ${\le }_r$ in $\mathcal{M}\times \mathcal{M}.$ We say that two remoteness functions $r,s$ are order-equivalent if ${\le }_r={\le }_s$.

Proposition 1.

Let $\gamma ,\beta \in \mathbb{R},\gamma ,\beta >0$. Then, $\gamma ·r+\beta$ is order-equivalent to $r.$

Proof. 

$a{\le }_{r}b$ iff $r\left(a\right)\le r\left(b\right)$ iff $\gamma ·r\left(a\right)\le \gamma ·r\left(b\right)$ iff $\gamma ·r\left(a\right)+\beta \le \gamma ·r\left(b\right)+\beta$ iff $a{\le }_{\gamma ·r+\beta }b$. □

We say that a probability density function $p:\mathcal{M}\times \mathcal{M}\to \mathbb{R}$ implements a remoteness function r if ${\ge }_p={\le }_r.$ This is equivalent to stating that for each $x,y,z,t\in \mathcal{M}$ we have that $p(x,y)\ge p(z,t)$ iff $r(x,y)\le r(z,t)$. That is, a density function implements a remoteness function r if it assigns higher probability density to those pairs of points which are closer according to r.

Once we have clarified what it means for a probability to implement a remoteness function, we introduce a specific way of creating probabilities to that.

Definition 2.

Let $f_r:\mathcal{M}\times \mathcal{M}\to \mathbb{R}$ be $f_r(x,y)=exp(-r(x,y))$. If $Z_r={\int }_{\mathcal{M}}{\int }_{\mathcal{M}}{f}_{r}d{\eta }_{g}d{\eta }_g$ is finite, we define the density function

$$p_r(x,y)=\frac{f_r(x,y)}{Z_r}=\frac{exp(-r(x,y\left)\right)}{Z_r}.$$

(7)

We refer to the corresponding probability distribution as a closeness distribution.

Note that $p_r$ is defined intrinsically. Following the explanation in the previous section, let $\pi$ be a chart of $\mathcal{M}$ defined almost everywhere. The representation of this pdf in the parametrization $x^{\prime },y^{\prime }=(\pi \left(x\right),\pi \left(y\right))$ is simply

$$p_{r^{\prime }}(x^{\prime },y^{\prime })=p_r({\pi }^{-1}\left(x^{\prime }\right),{\pi }^{-1}\left(y^{\prime }\right))$$

(8)

and its representation for integration purposes is

$${\rho }_{r^{\prime }}(x^{\prime },y^{\prime })=p_r({\pi }^{-1}\left(x^{\prime }\right),{\pi }^{-1}\left(y^{\prime }\right))\sqrt{\left|G\right(x^{\prime }\left)\right|}\sqrt{\left|G\right(y^{\prime }\left)\right|}$$

(9)

Proposition 2.

It it exists, $p_r$ implements r.

Proof. 

The exponential is a monotonous function and the minus sign in the exponent is used to revert the order. □

Proposition 3.

If r is measurable and $\mathcal{M}$ has finite volume, then $Z_r$ is finite, and hence $p_r$ implements r.

Proof. 

Note that since $r(x,y)\ge 0$, we have that $f_r(x,y)\le 1$, and hence $f_r$ is bounded. Furthermore, $f_r$ is measurable since it is a composition of measurable functions. Now, since any bounded measurable function in a finite volume space is integrable, $Z_r$ is finite. □

Obviously, once we have established a closeness distribution $p_r$ we can define its marginal and conditional distributions in the usual way. We note $p_r\left(x\right)$ (resp. $p_r\left(y\right)$) as the marginal over x (resp. y). We note $p_r\left(x\right|y)$ (resp. $p_r\left(y\right|x)$) as the conditional density of x given y (resp. y given x).

3. KL-Closeness Distributions for Multinomials

In this section we study closeness distributions on $M_n$ (the family of multinomial distributions of dimension n, or the family of finite discrete distributions over $n+1$ atoms). To do that, first we need to establish the remoteness function. It is well known that there is an isometry between $M_n$ and the positive orthant of the n dimensional sphere ($S_n$) of radius 2 (see Section 7.4.2. in [17]). This isometry allows us to compute the volume of the manifold as the area of the sphere of radius 2 on the positive orthant.

Proposition 4.

The volume of $M_n$ is $Vol\left(M_n\right)=\frac{{\pi }^{\frac{n+1}{2}}}{\Gamma \left(\frac{n+1}{2}\right)}.$

Proof. 

The area of a sphere $S_n$ of radius r is $A_{n,r}=\frac{2{\pi }^{\frac{n+1}{2}}{r}^n}{\Gamma \left(\frac{n+1}{2}\right)}.$ Taking $r=2$, $A_{n,2}=\frac{{\pi }^{\frac{n+1}{2}}{2}^{n+1}}{\Gamma \left(\frac{n+1}{2}\right)}.$ Now, there are $2^{n+1}$ orthants, so the positive orthant amounts for $\frac{1}{2^{n+1}}$ of that area, as stated. □

Figure 3 shows that the volume of the space of multinomial distributions over $n+1$ atoms reaches its maximum at $n=7$. The main takeover of Proposition 4 is that the volume of $M_n$ is finite, because this allows us to prove the following result:

Proposition 5.

For any measurable remoteness function r on $M_n$ there is a closeness distribution $p_r$ implementing it.

Proof. 

Directly from Proposition 3 and the fact that $M_n$ has finite volume. □

A reasonable choice of remoteness function for a statistical manifold is the Kullback-Leibler (KL) divergence. The next section analyzes the closeness distributions that implement KL in ${\mathcal{M}}_n$.

3.1. Closeness Distributions for KL as Remoteness Function

Let $\theta \in {\mathcal{M}}_n.$ Thus, $\theta$ is a discrete distribution over $n+1$ atoms. We write ${\theta }_i$ to represent $p(x=i|\theta ).$ Note that each ${\theta }_i$ is independent of the parametrization an thus it is an intrinsic quantity of the distribution.

Let $\theta ,\mu \in {\mathcal{M}}_n$. The KL divergence between $\theta$ and $\mu$ is

$$D(\mu ,\theta )=\sum _{i=1}^{n+1}{\mu }_{i}log\frac{{\mu }_i}{{\theta }_i}.$$

We want to study the closeness distributions that implement KL in ${\mathcal{M}}_n$. The detailed derivation of these results can be found in Appendix B. The closeness pdf according to Equation (7) is

$$\begin{array}{c}\hfill p_D(\mu ,\theta )=\frac{1}{Z_D}\prod _{i=1}^{n+1}{{\theta }_i}^{{\mu }_i}\prod _{i=1}^{n+1}{{\mu }_i}^{-{\mu }_i}\end{array}$$

The marginal for $\mu$ is

$$\begin{array}{c}\hfill p_D\left(\mu \right)=\frac{1}{Z_D}\prod _{i=1}^{n+1}{{\mu }_i}^{-{\mu }_i}B(\mu +\frac{1}{2})\end{array}$$

where $B\left(\alpha \right)=\frac{{\prod }_{i=1}^k\Gamma \left({\alpha }_i\right)}{\Gamma \left({\sum }_{i=1}^k{\alpha }_i\right)}$ is the multivariate Beta function.

The conditional for $\theta$ given $\mu$:

$$\begin{array}{c}\hfill p_D(\theta \mid \mu )=\frac{{\prod }_{i=1}^{n+1}{{\theta }_i}^{{\mu }_i}}{B(\mu +\frac{1}{2})}\end{array}$$

(10)

Equation (10) is very similar to the expression of a Dirichlet distribution. In fact, the expression of $p_D(\theta \mid \mu )$ for integration purposes in the expectation parameterization, namely ${\rho }_D(\theta \mid \mu )$, is that of a Dirichlet distribution:

$$\begin{array}{c}\hfill {\rho }_D(\theta \mid \mu )=Dirichlet(\theta ;\mu +\frac{1}{2})\end{array}$$

(11)

Equation (11) deserves some attention. We have defined the joint density $p_D(\mu ,\theta )$ so that pairs of distributions $(\mu ,\theta )$ that are close in terms of KL divergence are assigned a higher probability than pairs of distributions $({\mu }^{*},{\theta }^{*})$ which are further away in terms of KL. Hence, the conditional $p_D(\theta \mid \mu )$ assigns a larger probability to those distributions $\theta$ which are close in terms of KL to $\mu .$ This means that whenever we have a probabilistic model which encodes two multinomial distributions $\theta$ and $\mu$, and we are interested in introducing that $\theta$ should be close to $\mu$, we can introduce the assumption that $\theta \sim Dirichlet(\mu +\frac{1}{2}).$

Interesting as it is for modeling purposes, the use of Equation (11) however does not allow the modeler to convey information regarding the strength of the link. That is, $\theta$’s in the KL-surrounding of $\mu$ will be more probable, but there is no way to establish how much more probable. We know by Proposition 1 that for any remoteness function r, we can select $\gamma ,\beta >0$, and $\gamma ·r+\beta$ is order-equivalent to r. We can take advantage of that fact and use $\gamma$ to encode the strength of the probabilistic link between $\theta$ and $\mu$. If instead of using the KL (D) as the remoteness function, we opt for $\gamma ·D$, following a parallel development to the one above we will find that

$${\rho }_{\gamma ·D}(\theta \mid \mu )=Dirichlet(\theta ;\gamma \mu +\frac{1}{2}).$$

(12)

Now, Equation (12) allows the modeler to fix a large value of $\gamma$ to encode that it is extremely unlikely that $\theta$ separates from $\mu$, or a value of $\gamma$ close to 0 to encode that the link between $\theta$ and $\mu$ is highly loose. Furthermore it is important to realize that Equation (12) allows us to interpret any already existing model which incorporates Dirichlet (or Beta) distributions with the only requirement that each of its concentration parameters is larger than $\frac{1}{2}.$ Say we have a model in which $\theta \sim Dirichlet\left(\alpha \right)$. Then, defining $\mu$ by coordinates as ${\mu }_i=\frac{{\alpha }_i-\frac{1}{2}}{-\frac{n+1}{2}+{\sum }_{i=1}^{n+1}{\alpha }_i}$, we can interpret the model as imposing $\theta$ to be close to $\mu$ with intensity $\gamma =\frac{{\alpha }_1-\frac{1}{2}}{{\mu }_1}.$ Note that, extending this interpretation a bit to the extreme, since the strength of the link reduces as $\gamma \to 0$, a “free” Dirichlet will have all of its weights set to $\frac{1}{2}.$ This coincides with the classical prior suggested by Jeffreys [18,19] for this very same problem. This is reasonable since Jeffreys’ prior was constructed to be independent of the parametrization, that is, to be intrinsic to the manifold, similarly to what we are doing.

3.2. Visualizing the Distributions

In the previous section we have seen an expression for $p_{\gamma ·D}(\theta \mid \mu ).$ Since the KL divergence is not symmetric, we have that $p_{\gamma ·D}(\mu \mid \theta )$ is different from $p_{\gamma ·D}(\theta \mid \mu ).$ Unfortunately, we have not been able to provide a closed form expression for $p_{\gamma ·D}(\mu \mid \theta ).$ However, it is possible to compute it numerically in order to compare both conditionals.

Figure 4 shows a comparison of $p_{\gamma ·D}(\mu \mid \theta )$ and $p_{\gamma ·D}(\theta \mid \mu ).$ According to what is suggested in [20], for a proper interpretation of the densities we show its density function, which is intrinsic to the manifold, instead of its expression in the parametrization, as is commonly done. Note that from Equation (12), the value of $p_{\gamma ·D}(\theta \mid \mu )$ is 0 at $\theta =0$ and $\theta =1.$ In Figure 4, we can see that this is not the case for $p_{\gamma ·D}(\mu \mid \theta )$ neither at $\mu =0$ nor at $\mu =1$. In fact we see that $p_{\gamma ·D}(\theta \mid \mu )$ always starts below $p_{\gamma ·D}(\mu \mid \theta )$ at $\theta =0$ (resp. $\mu =0)$. Then, as $\theta$ (resp. $\mu$) grows, it is always the case that $p_{\gamma ·D}(\theta \mid \mu )$ goes over $p_{\gamma ·D}(\mu \mid \theta )$, to end decreasing again below it when $\theta$ (resp $\mu$) approaches to 1.

4. Reinterpreting the Beta-Binomial Model

We are now ready to go back to the rodents example provided in the introduction. The main idea we would like this hierarchical model to capture is that the ${\theta }_i$’s are somewhat similar. We do this by introducing a new random variable $\mu$ to which we would like each ${\theta }_i$ to be close to (see Figure 5). Furthermore, we introduce another variable $\gamma$ that controls how tightly coupled the ${\theta }_i$’s are to $\mu$. Now, $\mu$ represents a proportion, and priors for proportions have been well studied, including the “Bayes-Laplace rule” [21] which recommends $Beta(1,1)$, the Haldane prior [22] which is an improper prior ${lim}_{\alpha \to 0^{+}}Beta(\alpha ,\alpha )$, and the Jeffreys’ prior [18,19] $Beta(\frac{1}{2},\frac{1}{2}).$ Following the arguments in the previous section, here we stick with the Jeffreys’ prior. A more difficult problem is the selection of the prior for $\gamma$, where we still do not have a well founded choice. Note that taking a look at Equation (12), $\gamma$’s role acts similarly (although not exactly equal) to an equivalent sample size. Thus, the prior over $\gamma$ could be thought as a prior over the equivalent sample size with which $\mu$ will be incorporated as prior into the determination of each of the ${\theta }_i$’s. In case the size of each sample ($n_i$) is large, there will not be much difference between a hierarchical model and modeling of each of the 71 experiments as independent experiments. So, it makes sense for the prior over $\gamma$ to concentrate on a relatively small equivalent sample sizes. Following this line of thought we propose $\gamma$ to follow a $Gamma(\alpha =1,\beta =0.1).$

To summarize, the hierarchical model we obtain based on closeness probability distributions is:

$$\begin{array}{cc}\hfill \mu & \sim Beta(\frac{1}{2},\frac{1}{2})\hfill \end{array}$$

(13)

$$\begin{array}{cc}\hfill \gamma & \sim Gamma(1,0.1)\hfill \end{array}$$

(14)

$$\begin{array}{ccc}\hfill {\theta }_i& \sim Beta(\gamma \mu +\frac{1}{2},\gamma (1-\mu )+\frac{1}{2})\hfill & \hfill \forall i\in [1:71]\end{array}$$

(15)

$$\begin{array}{ccc}\hfill y_i& \sim Binomial(n_i,{\theta }_i)\hfill & \hfill \forall i\in [1:71].\end{array}$$

(16)

Figure 6 shows that the posteriors generated by both models are similar, and put the parameter $\mu$ (the pooled average) between 0.08 and 0.15, and the parameter $\gamma$ (the intensity of the link between $\mu$ and each of the ${\theta }_i$’s) between 5 and 25. Furthermore, the model is relatively insensitive to the parameters of the prior for $\gamma$ as long as they do create a sparse prior. Thus, we see that selecting the prior as $\Gamma (1,0.5)$ creates a prior that is too concentrated on the low values of $\gamma$ (that is, it imposes a relatively mild closeness link between $\mu$ and each of the ${\theta }_i$’s). This changes the estimation a lot. However, $\Gamma (1,0.01)$ creates a posterior similar to that of $\Gamma (1,0.1)$, despite being more spread.

5. Hierarchical Dirichlet Multinomial Model

Recently [9,10,11], Azzimonti et al. have proposed a hierarchical Dirichlet multinomial model to estimate conditional probability tables (CPTs) in Bayesian networks. Given two discrete finite random variables X (over domain $\mathcal{X}$) and Y over domain ($\mathcal{Y}$) which are part of a Bayesian network, and such that Y is the only parent of X in the network, the CPT for X is responsible of storing $p\left(X\right|Y).$ The usual CPT model (the so called Multinomial-Dirichlet) adheres to parameter independence and stores $\left|\mathcal{Y}\right|$ different independent Dirichlet distributions over each of the ${\theta }_{X|y}$. Instead, Azzimonti et al. propose the hierarchical Multinomial-Dirichlet model, where “the parameters of different conditional distributions belonging to the same CPT are drawn from a common higher-level distribution”. Their model can be summarized equationally as

$$\begin{array}{cc}\hfill \alpha & \sim s·Dirichlet\left({\alpha }_0\right)\hfill \\ \hfill {\theta }_{X|y}& \sim Dirichlet\left(\alpha \right)\hfill & \hfill \forall y\in \mathcal{Y}\\ \hfill X_y& \sim Categorical\left({\theta }_{X|y}\right)\hfill & \hfill \forall y\in \mathcal{Y}\end{array}$$

and graphically as shown in Figure 7.

The fact that the Dirichlet distribution is the conditional of a closeness distribution allows us to think about this model as a generalization of the model presented for the rat example. Thus, the Hierarchical Dirichlet Multinomial model can be understood as introducing the assumption that there is a probability distribution with parameter $\mu$, that is close in terms of its KL divergence to each of the $y\in \mathcal{Y}$ different distributions, each of them parameterized by ${\theta }_{X|y}.$ Thus, in equational terms, we have that the model can be rewritten as

$$\begin{array}{cc}\hfill \mu & \sim Dirichlet(\frac{1}{2},\dots ,\frac{1}{2})\hfill \end{array}$$

(17)

$$\begin{array}{cc}\hfill \gamma & \sim Gamma(1,0.1)\hfill \end{array}$$

(18)

$$\begin{array}{ccc}\hfill {\theta }_{X|y}& \sim Dirichlet(\gamma \mu +\frac{1}{2})\hfill & \hfill \forall y\in \mathcal{Y}\end{array}$$

(19)

$$\begin{array}{ccc}\hfill X_y& \sim Categorical\left({\theta }_{X|y}\right)\hfill & \hfill \forall y\in \mathcal{Y}\end{array}$$

(20)

and depicted as shown in Figure 8. Note that $\gamma$ in our reinterpreted model plays a role quite similar to the one that s played on Azzimonti’s model. To maintain the parallel with the model developed for the rodents example, here we have also assumed a $Gamma(1,0.1)$ as prior over $\gamma$, instead of the punctual distribution assumed in [10], but we could easily mimic their approach and specify a single value for $\gamma .$

Note that we are not claiming that we are improving the Hierarchical Dirichlet Multinomial model, we are just reinterpreting it in a way that is easier to understand conceptually.

6. Conclusions and Future Work

We have introduced the idea and the formalization remoteness functions and closeness distributions in Section 2. We have proven that any remoteness function induces a closeness distribution, provided that the volume of the space of distributions is finite. We have particularized closeness distributions for multinomials in Section 3, taking the remoteness function as the KL-divergence. By analyzing two examples, we have shown that closeness distributions can be a useful tool to the probabilistic model builder. We have seen that they can provide additional rationale and geometric intuitions for some commonly used hierarchical models.

In the Crowd4SDG and Humane-AI-net projects, these mathematical tools could prove useful in the understanding and development of consensus models for citizen science, improving the ones presented in [23]. Our plan is to study this in future work.

In this paper we have concentrated on discrete closeness distributions. The study of continuous closeness distributions remains for future work.