An Empirical Evaluation of the Approximation of Subjective Logic Operators Using Monte Carlo Simulations

Let $\Omega$ be a discrete collection of $M$ mutually exclusive and exhaustive events. A subjective opinion $\omega$ is a triple:

such that

where $\mathbf{b}\in\mathbb{R}^{M}$, with $b_{i}\in\mathbb{R}_{\geq 0}$, is the belief vector expressing the probability mass that the modeler places on each event $x_{i}$ in $\Omega$, $u\in\mathbb{R}_{\geq 0}$ is the uncertainty scalar quantifying the uncertainty of the modeler in its definition of $\mathbf{b}$, and $\mathbf{a}\in\mathbb{R}^{M}$ is the prior vector encoding a prior probability distribution over the events in $\Omega$. This subjective opinion is called a multinomial opinion.
Notice that the constraint in Equation 2 limits the degrees of freedom of $\mathbf{b}$ and $u$ to $M$ and, consequently, defines a $M$-dimensional simplex on which subjective opinions may be represented.

The limit-case multinomial opinion is the binomial opinion for $M=2$. In this case $\Omega=\{x,\overline{x}\}$ and the subjective opinion in Equation 1 may be re-written for simplicity as:

such that

where $b\in\mathbb{R}_{\geq 0}$ is the belief scalar expressing the probability of $x$, $d\in\mathbb{R}_{\geq 0}$ is the disbelief scalar expressing the probability of $\overline{x}$, $u\in\mathbb{R}_{\geq 0}$ is the uncertainty scalar and $a\in\mathbb{R}_{\geq 0}$ is a scalar expressing the prior probability of $x$.
Having only two degrees of freedom, binomial opinions in the form $\left(b,d,u,a\right)$ belong to a two-dimensional simplex and may be visualized together with $a$ in a barycentric coordinate system¹¹1See http://folk.uio.no/josang/sl/BV.html for an illustration..

In order to ground SL, a mapping has been defined between multinomial opinions and Dirichlet pdfs and between binomial opinions and Beta pdfs.

Given a mapping constant $W\in\mathbb{R}_{\geq 0}$, it is possible to define a unique mapping from opinions to pdfs. Let $\omega=\left(\mathbf{b},u,\mathbf{a}\right)$ be a multinomial opinion with $u\neq 0$; $\omega$ can be mapped to a Dirichlet pdf $p$ with distribution $\mathtt{Dir}\left(\boldsymbol{\alpha}\right)$, where the vector of parameters $\boldsymbol{\alpha}$ is defined as:

For binomial opinions, specifically, given a mapping constant $W\in\mathbb{R}_{\geq 0}$, it is possible to define a unique mapping from opinions to pdfs. Let $\omega=\left(b,d,u,a\right)$ be a binomial opinion with $u\neq 0$; $\omega$ can be mapped to a Beta pdf $p$ with distribution $\mathtt{Beta}\left(\alpha,\beta\right)$, where $\alpha$ and $\beta$ are parameters defined as:

Notice that, for reasons of consistency, $W$ is usually fixed to $2$ [4]. We then have a mapping $s$ from opinion $\omega$ to pdf $p$:

Vice versa, given a mapping constant $W\in\mathbb{R}_{\geq 0}$ and a fixed prior distribution $\mathbf{a}$, it is possible to define a unique mapping from pdfs to opinions. Let $p$ be a Dirichlet pdf with distribution $\mathtt{Dir}\left(\boldsymbol{\alpha}\right)$ with $\alpha_{i}>1$; $p$ can be mapped to a multinomial opinion $\omega=\left(\mathbf{b},u,\mathbf{a}\right)$, where the parameters are computed as:

Again, for a binomial opinion, given a mapping constant $W\in\mathbb{R}_{\geq 0}$ and a fixed prior distribution $a$, it is possible to define a unique mapping from pdfs to opinions. Let $p$ be a Beta pdf with distribution $\mathtt{Beta}\left(\alpha,\beta\right)$ with $\alpha,\beta>1$; $p$ can be mapped to a binomial opinion $\omega=\left(b,d,u,a\right)$, where the parameters are computed as:

For reasons of consistency, $W$ is usually fixed to $2$ [4]. Given a prior distribution $a$, this generates the mapping $t$ from pdf $p$ to opinion $\omega$:

SL defines several operators over subjective opinions, such as addition, product or fusion [4]. In general, these operators are computed over the parameters of subjective opinions. Let $\omega_{X}=\left(\mathbf{b}_{X},u_{X},\mathbf{a}_{X}\right)$ and $\omega_{Y}=\left(\mathbf{b}_{Y},u_{Y},\mathbf{a}_{Y}\right)$ be two subjective opinions and let $\circ_{SL}:\mathcal{S}\times\mathcal{S}\rightarrow\mathcal{S}$ be a generic operator over the space of subjective opinions $\mathcal{S}$. Then, $\omega_{Z}=\left(\mathbf{b}_{Z},u_{Z},\mathbf{a}_{Z}\right)$ resulting from the application of the operator to $\omega_{X}$ and $\omega_{Y}$ is given as:

where $f_{b},f_{u},f_{a}:\mathcal{S}\times\mathcal{S}\rightarrow\mathbb{R}_{\geq 0}$ are operator-specific functions returning the values of belief, uncertainty and prior for the opinion $\omega_{Z}$.

When properly defined, SL operators can be used to approximate operations over probability distribution functions. Suppose we are given two pdfs, $p_{X}$ and $p_{Y}$, and we want to compute a generic operation over them, $\circ_{P}:\mathcal{P}\times\mathcal{P}\rightarrow\mathcal{P}$ over the space of probability distributions $\mathcal{P}$. Computing this operation over probability distributions may be very complex. However, if we have an SL operator $\circ_{SL}:\mathcal{S}\times\mathcal{S}\rightarrow\mathcal{S}$ that approximates $\circ_{P}$, we may find a workaround computing $p_{X}\circ_{P}p_{Y}$ by projecting the two distribution onto the opinions $\omega_{X}$ and $\omega_{Y}$, computing the resulting opinion $\omega_{Z}=\omega_{X}\circ_{SL}\omega_{Y}$, and then mapping the result back onto a probability distribution function $p_{Z}^{SL}$. In this way, the resulting pdf $p_{Z}^{SL}$ provides an easy-to-compute approximation of the real pdf $p_{Z}$ (see Figure 1).

MC simulations are stochastic numerical algorithms designed to find approximate solutions through repeated random sampling. This paradigm has been applied in many areas of research to solve problems whose exact analytical solution is impossible or too difficult to derive. In statistics, MC simulations are widely used to evaluate probability distributions whose analytical form can not be explicitly expressed. Let $X$ be a random variable with a probability distribution $p_{X}$ on the support $\Omega$; let us also assume that the analytical form of $p_{X}$ is unknown but that we can sample realizations $x_{i}$ of the random variable $X$; then, MC simulations allow us to draw a large number of independent samples $x_{i}$ and use them to (i) compute useful empirical statistical descriptors $\hat{S}_{X}$ of the probability distribution $p_{X}$, or, eventually, (ii) reconstruct the approximate shape of the probability distribution $p_{X}$.

In order to compute useful empirical statistical descriptors $S_{X}$ of the probability distribution $p_{X}$, MC simulations rely on integration and on the law of large numbers. Let $S_{X}$ be a statistics of the probability distribution $p_{X}$ that can be computed from a function $f\left(\cdot\right)$ applied to the samples $x_{i}$. The statistics $S_{X}$ is then defined as:

By the law of large numbers, an estimator of $S_{X}$ can be computed using $N$ samples of $x_{i}$ as:

It is immediate to see that using Equation 11 and choosing an appropriate function $f(\cdot)$ we can directly estimate useful statistics of the distribution $p_{X}$, such as moments and quantiles. Thus, through a MC simulation we can sample points from $p_{X}$ and compute informative estimator statistics $\hat{S}_{X}$.

A probability distribution $p_{X}$ is completely characterized by the collection of all its moments; if we know the parametric form of the function $p_{X}$ from which we are sampling from, but we ignore the exact value of its parameters, we can compute an estimate $\hat{p}_{X}$ by setting the moments to the estimated values $\hat{M}_{i}\left[X\right]$. In several scenarios of interest it may actually be possible to compute analytically the value of few lower moments $M_{i}\left[X\right]$ of interest (such as, mean and variance); this approach is well-known and it has been used in the study of SL operators as well (see, for instance, [7]). In general, though, MC simulation and integration provide an empirical and robust way to compute estimators of the $i$-th moments $\hat{M}_{i}\left[X\right]$ of a probability distribution $p_{X}$, even when no exact analytical formula for computing the moments $M_{i}\left[X\right]$ of interest is available.

Beyond computing statistics, it is possible to use samples $x_{i}$ generated in a MC simulation to reconstruct the actual probability distribution $p_{X}$. A standard approach to reconstruct a continuous function $p_{X}$ from a set of finite points $x_{i}$ is kernel density estimation (KDE). Any function may be expressed as a convolution with a kernel function $\kappa(\cdot)$:

Practically, it is possible to get an empirical approximation using only a finite set of points $x_{i}$:

where the kernel $\kappa(\cdot)$ is a symmetric function, like a triangular function or a Gaussian, and $w$ denotes the width of the kernel; empirical rules are available to select an optimal value for this parameter in relation to the number of samples available [11]. Thus, using the same MC simulation procedure to sample points from $p_{X}$ it is possible also to estimate an approximate probability distribution $\hat{p}_{X}$.

Suppose we are given two probability distributions, $p_{X}$ and $p_{Y}$, and suppose we want to compute the distribution $p_{Z}$ determined by the application of operation $\circ_{P}:\mathcal{P}\times\mathcal{P}\rightarrow\mathcal{P}$, that is, $p_{Z}=p_{X}\circ_{P}p_{Y}$. If the pdf $p_{Z}$ can not be computed analytically, MC simulations may be used to sample from $p_{Z}$ and to estimate a pdf that approximates $p_{Z}$. As a first solution, we could rely on the samples $\{z_{1},z_{2}\dots z_{N}\}$ obtained by sampling from $p_{X}$ and $p_{Y}$ to estimate the moments $\hat{M}_{i}\left[Z\right]$ and then instantiate a moment-matching approximation $\hat{p}_{Z}^{\mathrm{MM}}$ (see Figure 2). Alternatively, we could use the same samples $\{z_{1},z_{2}\dots z_{N}\}$ from $p_{Z}$ to perform a kernel-density estimation and compute the KDE approximation $\hat{p}_{Z}^{\mathrm{KDE}}$ (see Figure 3). Notice that, differently from the SL approximation $p_{Z}^{SL}$, we decorate the approximations computed via MC simulations $\hat{p}_{Z}^{\mathrm{MM}}$ and $\hat{p}_{Z}^{\mathrm{KDE}}$ with a hat to underline that they are empirical statistics.

SL operators are defined to be extremely efficient. Indeed, given two opinions $\omega_{X}$ and $\omega_{Y}$ and the generic operator $\circ_{SL}:\mathcal{S}\times\mathcal{S}\rightarrow\mathcal{S}$, the computation of $\omega_{Z}=\omega_{X}\circ_{SL}\omega_{Y}$ usually requires only a limited number of function evaluations, as shown in Equation 9. The number of evaluations is $2\cdot M+1$, where $M$ is the number of events over which the opinions are defined. Thus, the overall complexity is $\mathcal{O}\left(M\right)$: it depends only on the number of events considered, and it is independent of the actual form of the mapped distributions. This makes SL operators an attractive choice especially when working in lower dimensions.

The MC approach is, by definition, computationally intensive. The computational complexity of a MC simulation scales as a function of the number $N$ of samples that must be produced. Each iteration requires random sampling and the execution of all the operations necessary to sample from $p_{Z}$. Overall, the computational complexity of the MC simulation is $\mathcal{O}\left(N\right)$. If we are using MC simulations to estimate a moment-matching approximation $\hat{p}_{Z}^{\mathrm{MM}}$, MC integration allows us to compute statistics from the samples generated during the MC simulation with no additional overall computational complexity. However, if we want to estimate the actual pdf via KDE we have to take into account an increase in the overall computational complexity from the linear order to the quadratic order $\mathcal{O}\left(N^{2}\right)$. Computing the pdf $\hat{p}_{Z}^{\mathrm{KDE}}$ is then a significantly computationally expensive procedure.

It is evident that, taking into account computational complexity only, SL operators dominate MC simulations, with or without KDE, especially considering that the number $N$ of samples in a MC simulation is required to grow large in order to return reliable results even in low dimensions.

The bias of SL operators is dependent on the definition of the specific operator, and a generic theoretical treatment is not possible. Moreover, an analytic study of the bias is not always available for all possible SL operators. In Section 7 we will consider the case study of the binomial operator for subjective logic and we will analyze more in detail its specific bias.

MC simulations are known to provide asymptotically unbiased estimators. If we estimate a statistics $S_{X}$ of the pdf $p_{X}$ using a MC integration as in Equation 11, then $\hat{S}_{X}$ is an asymptotically unbiased estimator, that is, in the limit of infinite samples, it converges to the true quantity it approximates:

where we made explicit the dependence of $\hat{S}_{X}$ on the number of samples $N$.

If we use MC integration to estimate the moments $\hat{M}\left[Z\right]$ for a moment-matching approximation $\hat{p}_{Z}^{\mathrm{MM}}$, the MC simulation provides us with unbiased estimators of the moments; this means that, by increasing the number of samples generated in a MC simulation, we can get arbitrarily close to the true value of the estimated quantity. However, notice that while the estimated moments $\hat{M}\left[Z\right]$ are asymptotically unbiased, the $\hat{p}_{Z}^{\mathrm{MM}}$ is biased; this bias is due to the limited set of moments $\hat{M}\left[Z\right]$ used to approximate $p_{Z}$.

If we use a MC simulation to estimate the true pdf directly via KDE, the empirical pdf $\hat{p}_{Z}^{\mathrm{KDE}}$ is biased. In this case, it is known that the width parameter $w$ of KDE regulates the trade-off between bias and variance. In general, the bias can be shown to be proportional to the width $w$ of the kernel $\kappa(\cdot)$:

under the constraint that $w$ can not be reduced to zero, for statistical and computational reasons [11]. When using a Gaussian kernel, the widely-adopted Silverman rule suggests the adoption of a kernel width of the following size:

where $\hat{\sigma}$ is the empirical standard deviation computed from the samples:

where $\hat{\mu}$ is the empirical mean. It follows, then, that the bias of the KDE approximation is proportional to:

As said, this bias can never be reduced to zero. However, in specific computational setting, this bias may be bounded by finding an optimal trade-off between the number of samples $N$ and the empirical standard deviation $\hat{\sigma}$. In particular, if the domain of $p_{Z}$ is a discrete domain, as in the case of multinomial opinions and Dirichlet pdfs which underlie subjective opinions, then the empirical standard deviation $\hat{\sigma}$ may be bounded and it may be possible to estimate the magnitude of the bias as a function of the number of samples $N$.

In summary, from the point of view of approximation, MC simulations represents a safer choice than SL operators, as they are grounded in solid theory and they allow us to quantify and to control the bias. The lack of any bound for SL operators may be seen as an obstacle in adopting them when working in critical domains where precise approximations are required. In the next section, we will introduce our protocol to solve this problem and estimate the degree of approximation of SL operators.

In the previous sections we illustrated two methodologies for finding an approximation of the pdf $p_{Z}$, one based on SL operators ($p_{Z}^{SL}$) and one relying on MC simulations ($\hat{p}_{Z}^{\mathrm{KDE}},\hat{p}_{Z}^{\mathrm{MM}}$). Figure 4 merges the graphs in Figure 1, 2 and 3 to illustrate the alternative computational paths that are offered to compute an approximation of $p_{Z}$; starting from the distributions $\left(p_{X},p_{Y}\right)$, the upper path represents the SL approach to finding an approximation of $p_{Z}=p_{X}\circ_{P}p_{Y}$, while the lower paths represent MC approaches to finding an approximation of the same quantity $p_{Z}$.

Now, approximate methods trade off precision in the results for simplicity in computation. In order to make a grounded decision on which approximation path in Figure 4 to use, it is necessary to quantify the trade-off between computational complexity and bias. As discussed in Section 4 and 5, in the case of KDE approximation via MC simulations, both complexity and bias are known. However, in the case of SL operators, we may easily derive their computational complexity, but we have no simple way of evaluating their bias. Exploiting the properties of MC integration and the idea of distance between pdfs, it is possible to assess the degree of approximation of SL operators in a computational fashion by relating them to MC simulations.

A simple way to evaluate how well a pdf $p$ approximates another pdf $q$ is to estimate the distance between them, $D\left[p,q\right]$, where $D\left[\cdot,\cdot\right]$ is a measure of distance or divergence between pdfs [8]. The degree of approximation of $p_{Z}^{SL}$ could then be obtained by measuring the distance from the true pdf $p_{Z}$:

However, since the true pdf $p_{Z}$ is taken to be unknown or hard to compute, it is challenging to get a direct estimate of these quantities. Since we can not rely directly on $p_{Z}$, we can instead exploit MC simulations and its properties.

From Equation 15 in Section 5, we know that the KDE estimation $\hat{p}_{Z}^{\mathrm{KDE}}$ is biased and we know how to evaluate it. Moreover, from Equation 18 in Section 5, we see that this bias depends on the number of samples $N$ and the standard deviation $\hat{\sigma}$. Now, if the domain of $p_{Z}$ is a discrete domain, as in the case of multinomial opinions and Dirichlet pdfs, then the empirical standard deviation $\hat{\sigma}$ may be bounded and it may be possible to estimate the magnitude of the bias as a function of the number of samples $N$. It may be possible to select a number of samples $N$ that shrinks the bias to a negligible quantity; in such case, we can then accept the KDE estimation $\hat{p}_{Z}^{\mathrm{KDE}}$ as a close approximation of the true pdf $p_{Z}$:

We underline that this approximation holds only under the assumption that, for an increasing number of samples $N$, the bias of the estimate $\hat{p}_{Z}^{\mathrm{KDE}}(N)$ tends, if not to zero, to a quantity whose order of magnitude is negligible with respect to further analysis; in other words, the validity of the approximation in Equation 20 is conditional on the pdf $p_{Z}$ we are considering, the analysis we will be carrying out, and the number of samples we can produce (for an example of an evaluation of these conditions, see the application to the case study of the product of Beta pdfs in Section 7 and Section 7.1).

The approximation in Equation 20 is extremely useful because it means that while we can not evaluate absolute distances with respect to the true distribution $p_{Z}$, we can still evaluate the relative distance between the KDE approximation and the SL approximation, and use it as a proxy for the distance between the SL approximation $p_{Z}^{SL}$ and the true distribution $p_{Z}$:

Thus, given only a finite set of samples $N$ we can obtain an empirical statistic of the distance as:

If the condition in Equation 20 holds, we expect the distance $D\left[\hat{p}_{Z}^{\mathrm{KDE}}(N),p_{Z}^{SL}\right]$ to be orders of magnitudes greater than $D\left[p_{Z},\hat{p}_{Z}^{\mathrm{KDE}}(N)\right]$; this would indeed confirm that the bias of $\hat{p}_{Z}^{\mathrm{KDE}}(N)$ is negligible and that the computation of $\hat{D}\left[p_{Z},p_{Z}^{SL}\right]$ (using a finite number of samples) provides a good estimate of the degree of approximation of the SL approximation.

So far, we have discussed distance measures in abstract terms. The quantity $\hat{D}\left[p_{Z},p_{Z}^{SL}\right]$ may indeed be computed using different pdf distance, such as $\phi$-divergences or integral probability metrics [13].

In this paper, we will rely on computing a simple integral distance, defined as:

This distance $D_{I}\left[p,q\right]$ is the same as the total variation distance except for the scaling constant:

The constant $\frac{1}{2}$ rescales the distance on the interval $[0,1]$. However, in order to get an absolute evaluation of how the mass of the two distributions $p$ and $q$ overlaps, we drop the scaling constant.

The choice of an integral distance $D_{I}\left[\cdot,\cdot\right]$ is justified for three reasons. First, from a conceptual point of view, an integral distance allows us to get a complete picture of the difference between two pdfs. While measures based on the evaluation of a limited set of synthetic statistics such as moments would provide us with a rough evaluation of the difference between two distributions, an integral distance provides a more precise way to assess the distribution of the mass of probability, taking into account, for instance, the potential presence of multiple modes or how mass subtly distributes on the tails.

Second, from a computational point of view, the integral distance $D_{I}\left[\cdot,\cdot\right]$ allows us, once again, to exploit MC integration. Recall that we want to get an estimation of $\hat{D}\left[p_{Z},p_{Z}^{SL}\right]$ through the approximation $D\left[\hat{p}_{Z}^{\mathrm{KDE}}(N),p_{Z}^{SL}\right]$. Now, if we reconstruct $\hat{p}_{Z}^{\mathrm{KDE}}$ via KDE, we can estimate the integral distance via MC integration over the domain $\Omega$ of the events as:

Third, from a theoretical point of view, the integral in Equation 25 is related to the bias:

Thus, using the integral distance $D_{I}\left[\hat{p}_{Z}^{\mathrm{KDE}},p_{Z}^{SL}\right]$ we can obtain an estimation of the distance $\hat{D}\left[p_{Z},p_{Z}^{SL}\right]$ as well as an upper bound on the bias of $p_{Z}^{SL}$. Notice that the absolute value in the integral distance provides a more honest evaluation of the absolute difference between pdfs, avoiding an averaging effect in absence of the absolute value operator.

Figure 5 summarizes our overall framework to evaluate the degree of approximation of the SL approximation as the integral distance $\int\left|p_{Z}^{SL}-\hat{p}_{Z}^{\mathrm{KDE}}\right|$, under the assumption that $D\left[p_{Z},\hat{p}_{Z}^{\mathrm{KDE}}(N)\right]\approx 0$. This approach is generic and it is not tied to the SL approximation. If the condition in Equation 20 can be guaranteed, the same approach may be used to get an estimation of the distance between the true pdf $p_{Z}$ and other potential approximation. For instance, Figure 5 shows our methodology applied also to the problem of estimating the distance from the true pdf of the moment-matching approximation $D\left[p_{Z},\hat{p}_{Z}^{\mathrm{MM}}(N)\right]$ by computing the distance $\int\left|\hat{p}_{Z}^{\mathrm{MM}}-\hat{p}_{Z}^{\mathrm{KDE}}\right|$.

Let $X$ be a random variable on the support $[0,1]$; we say that $X$ follows a Beta distribution $X\sim\mathtt{Beta}\left(\alpha,\beta\right)$ with parameters $\alpha\in\mathbb{R}_{\geq 0}$ and $\beta\in\mathbb{R}_{\geq 0}$ when its probability density function $p_{X}$ has the following form:

where $B\left(\alpha,\beta\right)$ is the Beta function.

Let $X\sim\mathtt{Beta}\left(\alpha_{X},\beta_{X}\right)$ and $Y\sim\mathtt{Beta}\left(\alpha_{Y},\beta_{Y}\right)$ be two independent Beta random variables with associated pdfs $p_{X}$ and $p_{Y}$. Let us define a third random variable $Z$ as the product of the two Beta random variables $Z=X\cdot Y$. The probability density function $p_{Z}$ of $Z$ does not follow a Beta distribution anymore, and its precise analytical form can not be easily expressed using elementary functions [9]. An analytical solution to the evaluation of the pdf of the product of two Beta distributions has been offered in [10]²²2This paper actually presents the more generic solution to the problem of multiplying two general Beta distribution, which subsume the multiplication of two simple Beta distributions as defined above.:

where $F_{D}^{(3)}$ is the Lauricella D hyper-geometric series. While this formula provides an elegant solution to the problem of finding the pdf of the product of two Beta pdfs, its straightforward evaluation is challenging as the Lauricella function requires the computation of factorial products and series.

Other analytical approaches to evaluate the product of two or more Beta distributions have been proposed, including methods relying on high-order functions, such as the Meijer G-function or Fox’s H function, or modeling the pdf of the product using an infinite mixture of simpler distributions [14, 1]. These approaches also present computational challenges, despite more efficient solutions have been investigated [14, 1].

Finally, common approaches rely on MC simulations to sample points from the probability distribution of $Z$ and to compute statistics of the pdf by matching the moments or the quantiles of $Z$ [9], as we reviewed in Section 3. Notice that when considering the product of two independent random variables $Z=X\cdot Y$, it is straightforward to compute the mean $E\left[Z\right]$ and the variance $Var\left[Z\right]$ of $p_{Z}$ analytically as:

Thus, if we were to perform a moment-matching approximation considering only the first two moments, we could instantiate such an approximation without running any MC simulation with a constant computational complexity of $\mathcal{O}\left(1\right)$.

An alternative solution to compute the product of Beta distributions is based on the use of the SL operator for binomial multiplication. Given two binomial opinions $\omega_{X}$ and $\omega_{Y}$ defined on different domains, the binomial opinion $\omega_{Z}$ resulting from the multiplication $\omega_{X}\cdot\omega_{Y}$ is computed as [4]:

Practically, a binomial product operator allows us to evaluate the combination of two opinions over two different facts. In the domain of probability distributions, the multiplication of opinions $\omega_{Z}=\omega_{X}\cdot\omega_{Y}$ translates into the multiplication of the mapped pdfs $p_{Z}=p_{X}\cdot p_{Y}$.

Now, assume we are interested in computing the product $Z=X\cdot Y$, where $X$ and $Y$ are two independent Beta random variables. Since an analytic solution is hard to compute, we may decide to rely either on the SL approximation or on a MC approximation.

Concerning moment-matching approximations we may consider a Gaussian pdf and a Beta pdf. Using a Gaussian pdf is a choice motivated by the simplicity and the ubiquity of this distribution; however, this is clearly a naive choice, as a Gaussian pdf has an unbounded support, is symmetrical and it assumes that all the statistical moments greater than the second are zero. Using a Beta distribution is a more prudent choice: even if it is known that the product of two Betas is not, in general, a Beta distribution, a Beta pdf still fits the right support and it may have other moments different from zero. In order to evaluate the parameters of our Gaussian and Beta approximation, we may rely on MC simulations or on an analytic evaluation. If we opt for MC simulations, we can use MC integration to estimate the mean $\hat{\mu}$ and the variance $\hat{\sigma}^{2}$ of $p_{Z}$; the Gaussian approximation $\hat{p}_{Z}^{\mathrm{GAUSS}}$ is then instantiated as $\mathtt{N}\left(\hat{\mu},\hat{\sigma}^{2}\right)$, while the Beta approximation $\hat{p}_{Z}^{\mathrm{BETA}}$ is defined as $\mathtt{Beta}\left(\frac{-\hat{\mu}\left(\hat{\sigma}^{2}+\hat{\mu}^{2}-\hat{\mu}\right)}{\hat{\sigma}^{2}},\frac{\left(\hat{\mu}-1\right)\left(\hat{\sigma}^{2}+\hat{\mu}^{2}-\hat{\mu}\right)}{\hat{\sigma}^{2}}\right)$, thus guaranteeing that $p_{Z}$, $\hat{p}_{Z}^{\mathrm{GAUSS}}$ and $\hat{p}_{Z}^{\mathrm{BETA}}$ have the same mean and variance. If we rely on an analytic approach, we can easily compute the mean $\mu$ and the variance $\sigma^{2}$ of $p_{Z}$ using Equation 31; as before, the Gaussian approximation ${p}_{Z}^{\mathrm{GAUSS}}$ is then instantiated as $\mathtt{N}\left({\mu},{\sigma}^{2}\right)$, while the Beta approximation ${p}_{Z}^{\mathrm{BETA}}$ is defined as $\mathtt{Beta}\left(\frac{-{\mu}\left({\sigma}^{2}+{\mu}^{2}-{\mu}\right)}{{\sigma}^{2}},\frac{\left({\mu}-1\right)\left({\sigma}^{2}+{\mu}^{2}-{\mu}\right)}{{\sigma}^{2}}\right).$

Figure 6 provides a concrete instantiation of the diagram in Figure 4, in which the generic operators $\circ_{P}$ and $\circ_{SL}$ have been substituted with multiplication and the generic moment-matching approximation $\hat{p}_{Z}^{\mathrm{MM}}$ has been replaced by the empirical Gaussian approximation $\hat{p}_{Z}^{\mathrm{GAUSS}}$, the empirical Beta approximation $\hat{p}_{Z}^{\mathrm{BETA}}$, the analytic Gaussian approximation ${p}_{Z}^{\mathrm{GAUSS}}$, and the analytic Beta approximation ${p}_{Z}^{\mathrm{BETA}}$.

As discussed earlier, choosing which path to take, whether to follow the SL approximation path in upper part of the graph or opt for one of the MC approximations in the lower part, requires evaluating the trade-off between computational complexity and degree of approximation of the different approaches. As these parameters are known in the case of MC simulations, we will review here the computational complexity and the approximation of the binomial multiplication.

Binomial multiplication is extremely efficient. Given two binomial opinions $\omega_{X}$ and $\omega_{Y}$ it is possible to compute their product $\omega_{Z}$ through a fixed and finite number of arithmetic operations. Independently from the actual form of the mapped distributions, the product is always computed in the same amount of time. As such, the computational complexity of these SL operators is constant $\mathcal{O}\left(1\right)$.

The original paper that introduced the SL operator for binomial multiplication [5] proposed a first qualitative analysis of the degree of approximation of this operator. In particular, it considered the specific instance of the multiplication of two Beta pdfs of the form $X,Y\sim\mathtt{Beta}\left(1,1\right)$; the pdf of $X$ and $Y$ reduces to a uniform distribution over $[0,1]$, which is taken to be a worst-case scenario with maximal variance and entropy. The analytical solution $Z=X\cdot Y$ to this particular case was then computed and graphically compared to the pdf associated with product $\omega_{Z}=\omega_{X}\cdot\omega_{Y}$. This study provided a clear visual appraisal of the difference between the exact pdf and the SL-approximated pdf, but no quantitative estimation were provided for more general cases.

In order to compute a numerical estimation of the degree of approximation of the SL operator for binomial multiplication we want to rely on the framework described in Section 6.

The basic condition expressed in Equation 20 requires the bias of $\hat{p}_{Z}^{\mathrm{KDE}}$ to be bounded and negligible. Recall that this bias, using a Gaussian kernel with width computed using the Silverman rule, is $E_{KDE}\left[p_{Z}-\hat{p}_{Z}^{\mathrm{KDE}}\right]\propto 1.06^{2}\hat{\sigma}^{2}N^{\frac{2}{5}}$, where

Now, notice that on our bounded support $\left[0,1\right]$ we can expect the difference $\left(x_{i}-\hat{\mu}\right)$ to, be at most, in the order of $10^{-1}$. This implies that, in the worst case, the order of magnitude of $\hat{\sigma}$ may be estimated as:

Consequently, relying on Silverman rule in Equation 16, the order of magnitude of largest kernel width $w$ may be bounded as:

As such, from Equation 18, the bias will be proportional to this upper bound:

Thus, for instance, if we were to run our MC simulation sampling $N=10^{5}$ samples, then we can expect the bias of $\hat{p}_{Z}^{\mathrm{KDE}}$ to be in the order of $10^{-4}$. This analysis on the bias allows us to consider the bias negligible if we are comparing it with quantities, such as $D\left[\hat{p}_{Z}^{\mathrm{KDE}}(N),p_{Z}^{SL}\right]$, order of magnitude greater than $10^{-4}$. If this condition is met, then we can estimate the degree of approximation of binomial multiplication adopting the framework illustrated in Figure 5 and instantiated for this specific SL operator as in Figure 7.