# Imprecise Bayesian Neural Networks

Michele Caprio, Souradeep Dutta, Kuk Jang, Vivian Lin, Radoslav Ivanov, Oleg
Sokolsky, Insup Lee  PRECISE Center, Dept. of Computer and Information
Science, University of Pennsylvania, 3330 Walnut Street, Philadelphia, PA
19104<EMAIL_ADDRESS>Dept. of
Computer Science, Rensselaer Polytechnic Institute, 110 Eighth Street, Troy,
NY 12180<EMAIL_ADDRESS>
###### Abstract.

Uncertainty quantification and robustness to distribution shifts are important
goals in machine learning and artificial intelligence. Although Bayesian
Neural Networks (BNNs) allow for uncertainty in the predictions to be
assessed, different sources of uncertainty are indistinguishable. We present
Imprecise Bayesian Neural Networks (IBNNs); they generalize and overcome some
of the drawbacks of standard BNNs. These latter are trained using a single
prior and likelihood distributions, whereas IBNNs are trained using credal
prior and likelihood sets. They allow to distinguish between aleatoric and
epistemic uncertainties, and to quantify them. In addition, IBNNs are more
robust than BNNs to prior and likelihood misspecification, and to distribution
shift. They can also be used to compute sets of outcomes that enjoy
probabilistic guarantees. We apply IBNNs to two case studies. One, for motion
prediction in autonomous driving scenarios, and two, to model blood glucose
and insulin dynamics for artificial pancreas control. We show that IBNNs
performs better when compared to an ensemble of BNNs benchmark.

###### Key words and phrases:

Bayesian deep learning; imprecise probabilities; credal sets; epistemic and
aleatory uncertainties; uncertainty quantification; machine learning
robustness

###### 2010 Mathematics Subject Classification:

Primary: 68T37; Secondary: 68T05, 68W25

##  1\. Introduction

One of the greatest virtues an individual can have is arguably being aware of
their own ignorance, and acting cautiously as a consequence. Similarly, an
autonomous system using neural networks (NNs) would greatly benefit from
understanding the probabilistic properties of the NN’s output (e.g., variance,
robustness to distribution shift), in order to incorporate them into any
further decision-making. In this paper, we present a generalization of
Bayesian neural networks that allows us to give a machine such a desirable
quality.

In the last few years, there has been a proliferation of work on calibrating
(classification) NNs, in order to estimate the confidence in their outputs
[31] or to produce conformal sets that are guaranteed to contain the true
label, in a probably approximately correct (PAC) sense [58]. While such
methods are a promising first step, they require a calibration set (in
addition to the original training set) and cannot be directly used on out-of-
distribution data without further examples.

Bayesian neural networks (BNNs) offer one approach to overcome the above
limitations. The Bayesian paradigm provides a rigorous framework to analyze
and train uncertainty-aware neural networks, and more generally to support the
development of learning algorithms [39]. In addition, it overcomes some of the
drawbacks of deep learning models, namely that they are prone to overfitting,
which adversely affects their generalization capabilities, and that they tend
to be overconfident about their predictions when they provide a confidence
interval. BNNs, though, are trained using a single prior, which may still
suffer from miscalibration and robustness issues [50].

In this work we introduce imprecise Bayesian neural networks (IBNNs). Unlike
other techniques in the fields of artificial intelligence (AI) and machine
learning (ML) involving imprecise probabilities – that typically only focus on
classification problems – IBNNs can be used for classification, prediction,
and regression. They capture the ambiguity the designer faces when selecting
which prior to choose for the parameters of a neural network and which
likelihood distribution to choose for the training data at hand. An IBNN can
be defined as a NN trained using credal prior and likelihood sets. Credal sets
are convex sets of probability measures (see Remark 1); we use them to train
IBNNs in order to overcome some of the drawbacks of BNNs. In particular, they
allow to counter the criticism to the practice in (standard) Bayesian
statistics of (i) using a single, arbitrary prior to represent the initial
state of ignorance of the agent, (ii) using non-informative priors to model
ignorance, and (iii) using a single, arbitrary likelihood to represent the
agent’s knowledge about the sampling model.111Criticisms (i) and (iii) are
also pointed out in [52, Section 2.2]. As a consequence, credal sets make the
analysis more robust to prior and likelihood misspecification. In addition,
they make it possible to quantify and distinguish between epistemic (i.e.
reducible) and aleatoric (i.e. irreducible) uncertainties. This is desirable
in light of several areas of recent ML research, such as Bayesian deep
learning [19, 41], adversarial example detection [64], and data augmentation
in Bayesian classification [40].

The motivation for training an artificial NN using credal sets is twofold: (i)
a single probability distribution does not suffice to represent ignorance in
the sense of lack of knowledge [35], and (ii) working with credal sets allows
to be robust against prior and likelihood misspecification. A more in-depth
discussion can be found in Appendix A. In addition, we point out that despite
a hierarchical Bayesian model (HBM) approach seems to be a viable alternative
to one based on credal sets, these latter are better justified
philosophically, and do not suffer from the same theoretical shortcomings of
HBM approaches [6, 35, 38, 71]. A more detailed explanation can be found in
Appendix B.

We summarize our contributions next: (1) We present IBNNs, and develop the
theoretical tools and the procedure required to use them in practice. (2) We
present theoretical results to show that IBNNs are more robust than BNNs to
prior and likelihood misspecification, and to distribution shifts. We also
prove how IBNNs can be used to specify sets of outcomes that enjoy
probabilistic guarantees. (3) We apply IBNNs to model two safety critical
systems. One, motion prediction for autonomous driving, and two, the human
insulin and blood glucose dynamics for artificial pancreas control. We
demonstrate improvements in both these settings with respect to ensemble of
BNNs methods.

Before moving on, let us point out that while IBNNs pay a computational price
coming from the use of credal sets, they are able to quantify epistemic
uncertainty, unlike BNNs, and to do so in a principled manner, unlike ensemble
of BNNs. In addition, they require less stringent assumptions on the nature of
the prior and likelihood ambiguity faced by the agent than other imprecise-
probabilities-based techniques, as explained in Appendices B and M. We also
stress that if the scholar prioritizes computational efficiency, they should
use backpropagation-based methods, as they are much faster than Bayesian
techniques. In safety-critical situations instead – where using uncertainty-
informed methods is crucial for the downstream task performance – IBNNs are
the natural choice.

Structure of the paper. Section 2 presents the needed preliminary concepts,
followed by section 3 that introduces some important theoretical results. We
discuss the applied aspects of IBNNs in in section 4. We present our
experimental results in section 5, and we examine the related work in section
6. Section 7 concludes our work. In the appendices, we give further
theoretical and philosophical arguments and we prove our claims.

##  2\. Background and Preliminaries

In this section, we present the background notions that are needed to
understand our main results. In section 2.1 we introduce Bayesian neural
networks, and section 2.2 discusses (finitely generated) credal sets, upper
probabilities and lower probabilities. The reader that is familiar with these
concepts can skip to section 3.

###  2.1. Bayesian neural networks

In line with the recent survey on BNNs by [39], Bayes’ theorem can be stated
as $P(H\mid D)=[P(D\mid H)P(H)]/P(D)=P(D,H)/\int
P(D,H^{\prime})\text{d}H^{\prime}$, where $H$ is a hypothesis about which the
agent holds some prior beliefs, and $D$ is the data the agent uses to update
their initial opinion. Probability distribution $P(D\mid H)$ represents how
likely it is to observe data $D$ if hypothesis $H$ were to be true, and is
called likelihood, while probability distribution $P(H)$ represents the
agent’s initial opinion around the plausibility of hypothesis $H$, and is
called prior. The evidence available is encoded in $P(D)=\int
P(D,H^{\prime})\text{d}H^{\prime}$, while posterior probability $P(H\mid D)$
represents the agent’s updated opinion. Using Bayes’ theorem to train a
predictor can be understood as learning from data $D$: the Bayesian paradigm
offers an established way of quantifying uncertainty in deep learning models.

BNNs are stochastic artificial neural networks (ANNs) trained using a Bayesian
approach [29, 39, 47, 66, 72]. The goal of ANNs is to represent an arbitrary
function $y=\Phi(x)$. Let $\theta$ represent the parameters of the network,
and call $\Theta$ the space $\theta$ belongs to. Stochastic neural networks
are a type of ANN built by introducing stochastic components to the network.
This is achieved by giving the network either a stochastic activation or
stochastic weights to simulate multiple possible models with their associated
probability distribution. This can be summarized as $\theta\sim p(\theta)$,
$y=\Phi_{\theta}(x)+\varepsilon$, where $\Phi$ depends on $\theta$ to
highlight the stochastic nature of the neural network, $p$ is the density of a
probability measure $P$ on $\Theta$,222We can write $p$ as the Radon-Nikodym
derivative of $P$ with respect to some $\sigma$-finite dominating measure
$\mu$, that is, $p=\text{d}P/\text{d}\mu$. and $\varepsilon$ represents random
noise to account for the fact that function $\Phi_{\theta}$ is just an
approximation.

To design a BNN, the first step is to choose a deep neural network
architecture, that is, functional model $\Phi_{\theta}$. Then, the agent
specifies the stochastic model, that is, a prior distribution over the
possible model parametrization $p(\theta)$, and a prior confidence in the
predictive power of the model $p(y\mid x,\theta)$. Given the usual assumption
that multiple data points from the training set are independent, the product
$\prod_{(x,y)\in D}p(y\mid x,\theta)$ represents the likelihood of outputs
$y\in D_{\mathbf{y}}$ given inputs $x\in D_{\mathbf{x}}$ and parameter
$\theta$, where (a) $D=D_{\mathbf{x}}\times D_{\mathbf{y}}$ is the training
set; (b) $D_{\mathbf{x}}=\\{x_{i}\\}_{i=1}^{n}$ is the collection of training
inputs, which is a subset of the space $\mathcal{X}$ of inputs; (c)
$D_{\mathbf{y}}=\\{y_{i}\\}_{i=1}^{n}$ is the collection of training outputs,
which is a subset of the space $\mathcal{Y}$ of outputs. In the paper, we call
“likelihood” both $p(y\mid x,\theta)$ and $\prod_{(x,y)\in D}p(y\mid
x,\theta)$, as no confusion arises.

The model parametrization can be considered to be hypothesis $H$. Following
[39], we assume independence between model parameters $\theta$ and training
inputs $D_{\mathbf{x}}$, in formulas
$D_{\mathbf{x}}\perp\\!\\!\\!\perp\theta$. Hence, Bayes’ theorem can be
rewritten as

$p(\theta\mid D)=\frac{p(D_{\mathbf{y}}\mid
D_{\mathbf{x}},\theta)p(\theta)}{\int_{\Theta}p(D_{\mathbf{y}}\mid
D_{\mathbf{x}},\theta^{\prime})p(\theta^{\prime})\text{d}\theta^{\prime}}\propto
p(D_{\mathbf{y}}\mid D_{\mathbf{x}},\theta)p(\theta).$

Notice that the equality comes from having assumed
$D_{\mathbf{x}}\perp\\!\\!\\!\perp\theta$. Posterior density $p(\theta\mid D)$
is high dimensional and highly nonconvex [37, 39], so computing it and
sampling from it is a difficult task. The first issue is tackled using
Variational Inference (VI) procedures, while Markov Chain Monte Carlo (MCMC)
methods address the second challenge. Both are reviewed – in the context of
machine learning – in [39, Section V], where the authors also inspect their
limitations. BNNs can be used for prediction, regression, and classification
[39, Section II]; besides having a solid theoretical justification, there are
practical benefits from using BNNs, as presented in [39, Section III].

###  2.2. Imprecise probabilities

As IBNNs are rooted in the theory of imprecise probabilities (IPs), in this
section we give a gentle introduction to the IP concepts we will use
throughout the paper.

IBNNs are based on the Bayesian sensitivity analysis (BSA) approach to IPs,
that in turn is grounded in the dogma of ideal precision (DIP) [5], [71,
Section 5.9]. The DIP posits that in any problem there is an ideal probability
model which is precise, but which may not be precisely known. We call this
condition ambiguity [23, 28].

Facing ambiguity can be represented mathematically by a set $\mathcal{P}$ of
priors and a set $\mathcal{L}$ of likelihoods that seem “plausible” or “fit”
to express the agent’s beliefs on the parameters of interest and their
knowledge of the data generating process (DGP). Generally speaking, the
farther apart the “boundary elements” of the sets (i.e. their infimum and
supremum), the higher the agent’s ambiguity. Of course, if $\mathcal{P}$ and
$\mathcal{L}$ are singletons we go back to the usual Bayesian paradigm.

A procedure based on sets $\mathcal{P}$ and $\mathcal{L}$ yields results that
are more robust to prior and likelihood misspecification than a regular
Bayesian method. In the presence of prior ignorance and indecisiveness about
the sampling model, it is better to give answers in the form of intervals or
sets, rather than arbitrarily select a prior and a likelihood, and then
update. Sets $\mathcal{P}$ and $\mathcal{L}$ allow to represent indecision,
thus leading to less informative but more robust conclusions.

######  Remark 1.

Throughout the paper, we denote by $\Pi=\\{P_{1},\ldots,P_{k}\\}$,
$k\in\mathbb{N}$, a finite set of probabilities on a generic space $\Omega$,
such that for all $j\in\\{1,\ldots,k\\}$, $P_{j}$ cannot be written as a
convex combination of the other $k-1$ components of $\Pi$. We denote by
$\Pi^{\prime}$ its convex hull $\Pi^{\prime}\equiv\text{Conv}\Pi$, i.e., the
set of probabilities $Q$ on $\Omega$ that can be written as
$Q(A)=\sum_{j=1}^{k}\alpha_{j}P_{j}(A)$, for all $A\subset\Omega$, where the
$\alpha_{j}$’s are elements of $[0,1]$ that sum up to $1$. In the literature,
it is referred to as a finitely generated credal set [13, 51]. Notice then
that the extreme elements of $\Pi^{\prime}$ correspond to the elements of
$\Pi$, that is, $\text{ex}\Pi^{\prime}=\Pi$. Simple graphical representations
of finitely generated credal sets are given in Figures 1 and 2.

Figure 1. Suppose we are in a $3$-class classification setting, so
$\Omega=\\{\omega_{1},\omega_{2},\omega_{3}\\}$. Then, any probability measure
$P$ on $\Omega$ can be seen as a probability vector. For example, suppose
$P(\\{\omega_{1}\\})=0.6$, $P(\\{\omega_{2}\\})=0.3$, and
$P(\\{\omega_{3}\\})=0.1$. We have that $P\equiv(0.6,0.3,0.1)^{\top}$. Since
its elements are positive and sum up to $1$, probability vector $P$ belongs to
the unit simplex, the purple triangle in the figure. Then, we can specify
$\Pi=\\{P_{1},\ldots,P_{5}\\}$, and obtain as a consequence that
$\Pi^{\prime}=\text{Conv}\Pi$ is the orange pentagon. It is a convex polygon
with finitely many extreme elements, and it is the geometric representation of
a finitely generated credal set.

Let us now introduce the concepts of lower and upper probabilities. The lower
probability $\underline{P}$ associated with $\Pi$ is given by
$\underline{P}(A)=\inf_{P\in\Pi}P(A)$, for all $A\subset\Omega$. The upper
probability $\overline{P}$ associated with $\Pi$ is defined as the conjugate
to $\underline{P}$, that is,
$\overline{P}(A):=1-\underline{P}(A^{c})=\sup_{P^{\prime}\in\Pi}P^{\prime}(A)$,
for all $A\subset\Omega$. These definitions hold even if $\Pi$ is not finite.
Then, we have the following important result.

######  Proposition 2.

$\overline{P}$ is the upper probability for $\Pi$ if and only if it is also
the upper probability for $\Pi^{\prime}$. That is,
$\overline{P}(A)=\sup_{P\in\Pi}P(A)=\sup_{P^{\prime}\in\Pi^{\prime}}P^{\prime}(A)$,
for all $A\subset\Omega$.

A version of Proposition 2 for finitely additive probability measures can be
found in [71, Section 3.6]. A simple graphical representation of upper and
lower probabilities for a set $A$ is given in Figure 2.

Figure 2. In this figure, a replica of [24, Figure 1],
$\Pi=\\{P_{1},P_{2}\\}$, where $P_{1}$ and $P_{2}$ are two Normal
distributions whose pdf’s $p_{1}$ and $p_{2}$ are given by the dashed blue and
brown curves, respectively. Their convex hull is
$\Pi^{\prime}=\text{Conv}\Pi=\\{Q:Q=\beta P_{1}+(1-\beta)P_{2}\text{, for all
}\beta\in[0,1]\\}$. The pdf $q$ of an element $Q$ of $\Pi^{\prime}$ is
depicted by a solid black curve. In addition, let $A=[-0.8,-0.4]$. Then,
$\underline{P}(A)=\int_{-0.8}^{-0.4}p_{2}(\omega)\text{d}\omega\approx 0$,
while $\overline{P}(A)$ is given by the red shaded area under $p_{1}$, that
is, $\overline{P}(A)=\int_{-0.8}^{-0.4}p_{1}(\omega)\text{d}\omega$.

##  3\. Theoretical results

In this section, we provide the procedure to follow in order to compute the
posterior credal set in the context of IBNNs, and we show how IBNNs are able
to capture both aleatoric and epistemic uncertainties associated with the
analysis. The training method is presented in Section 3.1. We show that IBNNs
are more robust to distribution shifts than regular BNNs as a result of
Proposition 4. We first give the formal definition of an IBNN.

######  Definition 3 (IBNN).

An IBNN is a stochastic artificial neural network trained using finitely
generated credal prior and likelihood sets.

###  3.1. IBNN procedure

Recall that $D=D_{\mathbf{x}}\times D_{\mathbf{y}}$ denotes the training set,
where $D_{\mathbf{x}}=\\{x_{i}\\}_{i=1}^{n}\subset\mathcal{X}$ is the
collection of training inputs,
$D_{\mathbf{y}}=\\{y_{i}\\}_{i=1}^{n}\subset\mathcal{Y}$ is the collection of
training outputs, and $\mathcal{X}$ and $\mathcal{Y}$ denote the input and
output spaces, respectively. We then denote by $P$ a generic prior on the
parameters $\theta\in\Theta$ of a BNN having pdf/pmf $p$, and by
$P_{x,\theta}$ a generic likelihood on the space $\mathcal{Y}$ of outputs
having pdf/pmf $p_{x,\theta}$. The act of computing the posterior from prior
$P$ and likelihood $P_{x,\theta}$ using a BNN is designated by
$\mathfrak{post}(P,P_{x,\theta})$. The IBNN procedure follows.

Procedure Imprecise Bayesian Neural Network. S1 Specify a finite set
$\mathcal{P}$ of plausible prior probabilities on the parameters $\theta$ of
the neural network, and a finite set $\mathcal{L}_{x,\theta}$ of plausible
likelihoods. S2 Compute posterior $P(\cdot\mid D)\equiv
P_{D}=\mathfrak{post}(P,P_{x,\theta})$ on the parameters $\theta$ of the
neural network, for all $P\in\mathcal{P}$ and all
$P_{x,\theta}\in\mathcal{L}_{x,\theta}$.

Step S2 performs an element-wise application of Bayes’ rule for all the
elements of $\mathcal{P}$ and $\mathcal{L}_{x,\theta}$. We obtain a finite set
$\mathcal{P}_{D}$ of posteriors whose cardinality $\\#\mathcal{P}_{D}$ is
given by $\\#\mathcal{P}\times\\#\mathcal{L}_{x,\theta}$. Its convex hull
$\text{Conv}\mathcal{P}_{D}$ is the credal posterior set. By Lemma 2 we have
that the upper and lower probabilities of $\mathcal{P}_{D}$ and
$\text{Conv}\mathcal{P}_{D}$ coincide. In applications, every element of
$\mathcal{P}$ should not be obtainable as a convex combination of other
elements. This to avoid carrying out redundant operations in Step S2. The same
holds for the elements of $\mathcal{L}_{x,\theta}$. Furthermore, choosing $2$
to $5$ priors and likelihoods is usually enough to safely hedge against prior
and likelihood misspecification.

Our procedure prescribes the scholar to specify a finite set $\mathcal{P}$ of
priors on the parameters of the neural network, and to specify a finite set
$\mathcal{L}_{x,\theta}$ of likelihoods (that can be used to capture different
architectures of the neural network) that approximate the true data generating
process. For every pair of prior $P$ and likelihood $P_{x,\theta}$, the
associated posterior $P(\cdot\mid D)\equiv P_{D}$ on the parameters $\theta$
of the network captures the distribution of these latter, revised according to
the available information encapsulated in the likelihood. Eventually, as
explained in section 4, we use the posteriors to derive the predictive
distributions $\hat{P}$ on the output space $\mathcal{Y}$, which are the
distributions of interests for the downstream tasks.

Notice that in the case that $\mathcal{P}$ and $\mathcal{L}_{x,\theta}$ are
both not singletons, for all $A\subset\Theta$ the interval
$[\underline{P}_{D}(A),\overline{P}_{D}(A)]$ is wider than the case when one
or the other is a singleton. In the limiting case where both are singletons,
we retrieve the usual Bayesian updating, so the interval shrinks down to a
point.

Robustness Properties of IBNNs. We follow the Procedure above to compute the
posterior credal set for the parameters of our IBNN. Being trained using
credal sets makes IBNNs more robust to distribution shifts than BNNs. To see
this, we present the following general result, and then we apply it to our
case.

Let $\mathscr{P}$ be a generic set of probabilities, and consider a
probability measure $P^{\prime}$ such that $P^{\prime}\not\in\mathscr{P}$.

######  Proposition 4.

Call $d$ any metric and $div$ any divergence on the space of probability
measures of interest. Let
$d(\mathscr{P},P^{\prime}):=\inf_{P\in\mathscr{P}}d(P,P^{\prime})$ and
$div(\mathscr{P}\|P^{\prime}):=\inf_{P\in\mathscr{P}}div(P\|P^{\prime})$.
Then, for all $P\in\mathscr{P}$, $d(\mathscr{P},P^{\prime})\leq
d(P,P^{\prime})$ and $div(\mathscr{P}\|P^{\prime})\leq div(P\|P^{\prime})$.

Proposition 4 holds if $\mathscr{P}$ is any set of probabilities, not just a
credal sets. In Appendix J, we show that the above result still holds if the
elements of $\mathscr{P}$ and $P^{\prime}$ are defined on Euclidean spaces
having different dimensions [7].

Let us now apply Proposition 4 to IBNNs. Suppose that, when designing a
regular BNN, an agent chooses likelihood $\mathbf{P}_{x,\theta}$, while when
designing a more general IBNN, they start by specifying a finite set of
likelihoods
$\mathcal{L}_{x,\theta}=\\{P^{1}_{x,\theta},\ldots,P^{k}_{x,\theta}:k\in\mathbb{N}\text{,
}\theta\in\Theta\text{, }x\in\mathcal{X}\\}$, and then let the induced credal
set $\text{Conv}\mathcal{L}_{x,\theta}$ represent their uncertainty around the
sampling model. Assume that $\mathbf{P}_{x,\theta}\in\mathcal{L}_{x,\theta}$
(this means that when designing the regular BNN, the agent chooses arbitrarily
which of the elements of $\mathcal{L}_{x,\theta}$ to use) and that the
“oracle” data generating process $P^{o}_{x,\theta}$ is different from
$\mathbf{P}_{x,\theta}$, $P^{o}_{x,\theta}\neq\mathbf{P}_{x,\theta}$, so that
we are actually in the presence of distribution shift. Then, we have two
cases. (1) If the true sampling model $P^{o}_{x,\theta}$ belongs to
$\text{Conv}\mathcal{L}_{x,\theta}$, then the distance – measured via a metric
or a divergence – between $\text{Conv}\mathcal{L}_{x,\theta}$ and
$P^{o}_{x,\theta}$ is $0$ while that between $\mathbf{P}_{x,\theta}$ and
$P^{o}_{x,\theta}$ is positive. (2) If
$P^{o}_{x,\theta}\not\in\text{Conv}\mathcal{L}_{x,\theta}$, then the distance
between $\text{Conv}\mathcal{L}_{x,\theta}$ and $P^{o}_{x,\theta}$ is no
larger than the distance between $\mathbf{P}_{x,\theta}$ and
$P^{o}_{x,\theta}$, no matter (i) which metric or distance we use (Proposition
4), (ii) whether or not $P^{o}_{x,\theta}$ and the elements of
$\text{Conv}\mathcal{L}_{x,\theta}$ are defined on the same Euclidean space
(Appendix J, Lemma 16). A visual representation is given in Figure 3.

Figure 3. IBNNs are more robust to distribution shifts than regular BNNs. Here
$\text{Conv}\mathcal{L}_{x,\theta}$ is the convex hull of five plausible
likelihoods, and $d$ denotes a generic metric on the space
$\Delta(\mathcal{Y},\mathcal{A}_{\mathbf{y}})$ of probabilities on
$\mathcal{Y}$. We see how
$d(\text{Conv}\mathcal{L}_{x,\theta},P^{o}_{x,\theta})<d(\mathbf{P}_{x,\theta},P^{o}_{x,\theta})$;
if we replace metric $d$ by a generic divergence $div$, the inequality would
still hold.

Computational aspects of IBNNs. A computational bottleneck of the IBNN
procedure appears to be step S2, which is a combinatorial task. We have to
calculate $\\#\mathcal{P}\times\\#\mathcal{L}_{x,\theta}$ many posteriors, but
this procedure allows us to forego any additional assumptions on the nature of
the lower and upper probabilities that are often times required by other
imprecise-probabilities-based techniques.333If we are willing to make such
assumptions, Theorem 9 in Appendix D shows how to compute the upper posterior
using only upper prior and upper likelihood. Clearly, the procedure is
simplified if either $\mathcal{P}$ or $\mathcal{L}_{x,\theta}$ are singletons.

The posteriors in $\mathcal{P}_{D}$ are typically very high-dimensional and
highly non-convex; we use variational inference (VI), outlined in [39, Section
V], to approximate them in an efficient fashion. That is, we project every
$P_{D}\in\mathcal{P}_{D}$ onto a set $\mathbb{S}$ of “well-behaved”
distributions (e.g. Normals) using the KL divergence. In formulas,
$\breve{P}_{D}=\operatorname*{arg\,min}_{Q\in\mathbb{S}}KL(Q\|P_{D})$. By
“well-behaved”, we mean that they have to satisfy the conditions in [75,
Sections 2 and 3].444We assume that the conditions on the priors and the
likelihoods given in [75] are satisfied. This ensures that as the sample size
goes to infinity, the approximated posteriors converge to the true data
generating process. We also point out how despite using a VI approximation for
the exact posteriors, by Proposition 4 the credal set of approximated
posteriors is closer to the “oracle” posterior $P^{o}_{D}$ than any of its
elements. As a consequence, working with credal sets leads to VI posterior
approximations that are better than the ones resulting from a single BNN, or
an ensemble of BNNs, where several BNNs are combined into one.

######  Remark 5.

Although highly unlikely in practice, it is theoretically possible that the VI
approximation of the (finite) set $\mathcal{P}_{D}$ of posteriors is a
singleton, see Figure 4. While the conditions in [75] guarantee that
asymptotically the approximated posteriors coincide with the true data
generating process, this typically does not happen with finite-dimensional
datasets. As a consequence, obtaining a singleton when projecting
$\mathcal{P}_{D}$ onto $\mathbb{S}$ may result in an underestimation of the
uncertainties faced by the scholar. In that case, we either consider a
different set – whose elements still satisfy the conditions in [75, Sections 2
and 3] – on which to project $\mathcal{P}_{D}$ according to the KL divergence,
or we use a different “projection operator”, that is, a divergence different
from the KL. For example, [75] suggest Rényi and $\chi^{2}$ divergences, or
Hellinger and total variation metrics. Alternatively, we can consider a
different approximation strategy altogether, for instance the Laplace
approximation [60].

######  Remark 6.

Assume that the “oracle” prior $P^{o}$ is in the prior credal set
$\text{Conv}\mathcal{P}$ and that the “oracle” likelihood $P^{o}_{x,\theta}$
is in the likelihood credal set $\text{Conv}\mathcal{L}_{x,\theta}$. Then, it
is immediate to see that the “oracle” posterior $P^{o}_{D}$ belongs to the
posterior credal set $\text{Conv}\mathcal{P}_{D}$. Naturally, this does not
imply that the posterior credal set collapses to $P^{o}_{D}$. In general, it
is unlikely that a finite amount of data is able to completely annihilate all
the epistemic uncertainty faced by the agent. What may happen is that if the
training set is large enough, $\text{Conv}\mathcal{P}_{D}$ may be inscribed in
a ball of small radius around $P^{o}_{D}$. This does not mean that we suffer
from under-confidence due to larger-than-necessary epistemic uncertainty.
Rather, the relative epistemic uncertainty (measured by the difference between
prior and posterior uncertainty, divided by the prior uncertainty) drops
significantly. In addition, working with sets of prior and likelihoods allows
us to hedge against prior and likelihood misspecification: the statements
$P^{o}\in\text{Conv}\mathcal{P}$ and
$P^{o}_{x,\theta}\in\text{Conv}\mathcal{L}_{x,\theta}$ are either strong
assumptions, or very hard to verify.

Figure 4. Let $\Delta(\Omega,\mathcal{B})$ denote the space of probability
measures on $\Theta$. Let
$\mathcal{P}_{D}=\\{P_{D}^{1},P_{D}^{2},P_{D}^{3}\\}$, so that
$\text{Conv}\mathcal{P}_{D}$ is the black segment. Then, if we project the
elements of $\mathcal{P}_{D}$ onto $\mathbb{S}_{1}$ via the KL divergence, we
obtain the same distribution $\breve{\mathbf{P}}_{D}$. This is detrimental to
the analysis because such an approximation underestimates epistemic (and
possibly also aleatoric) uncertainty faced by the agent. Then, the scholar
should specify a different set $\mathbb{S}_{2}$ of “well-behaved”
distributions onto which project the elements of $\mathcal{P}_{D}$. In the
figure, we see that $P_{D}^{1}$, $P_{D}^{2}$, and $P_{D}^{3}$ are projected
onto $\mathbb{S}_{2}$ via the KL divergence to obtain $\breve{P}_{D}^{1}$,
$\breve{P}_{D}^{2}$, and $\breve{P}_{D}^{3}$. The convex hull of these latter,
captured by the red shaded triangle, represents the variational approximation
of $\text{Conv}\mathcal{P}_{D}$.

###  3.2. Aleatoric and epistemic uncertainties

In [35], the authors study uncertainty in the context of supervised learning.
In particular, they extensively review the existing approaches to quantify
aleatoric and epistemic uncertainty (AU and EU, respectively). The former
refers to irreducible uncertainty, the variability in the outcome of an
experiment which is due to inherently random effects. An example of AU is coin
flipping: the data generating process in this experiment has a stochastic
component that cannot be reduced by any additional source of information. EU,
instead, corresponds to reducible uncertainty, caused by a lack of knowledge
about the best model. Continuing the coin example, this corresponds to holding
a belief that the coin is biased. After a few tosses, we realize whether or
not the coin is in fact biased, or if instead it is fair. Notice how, in the
precise case – that is, when the agent specifies a single distribution – EU is
absent. A single probability measure is only able to capture the idea of
aleatoric uncertainty, since it represents a case in which the agent knows
exactly the true data generating process. This is a well-studied property of
credal sets [35, Page 458]. EU should not be confused with the concept of
epistemic probability [15, 16, 71]. In the subjective probability literature,
epistemic probability can be captured by a single distribution. Its best
definition can be found in [71, Sections 1.3.2 and 2.11.2]. There, the author
specifies how epistemic probabilities model logical or psychological degrees
of partial belief of the agent. We remark, though, how de Finetti and Walley
work with finitely additive probabilities, while in this paper we use
countably additive probabilities.

Recall that, given a probability measure $P$ on a generic space $\Omega$, the
(Shannon) entropy of $P$ is defined as $H(P):=\mathbb{E}[-\log
p]=-\int_{\Omega}\log[p(\omega)]P(\text{d}\omega)$ if $\Omega$ is uncountable,
where $p$ denotes the pdf of $P$. If $\Omega$ is at most countable, we have
that $H(P)=-\sum_{\omega\in\Omega}P(\\{\omega\\})\log[P(\\{\omega\\})]$. The
entropy primarily captures the shape of the distribution, namely its
“peakedness” or non-uniformity [20, 35], and hence informs about the
predictability of the outcome of a random experiment: the higher its value,
the lower the predictability. Now, consider a generic set of probabilities
$\mathscr{P}$ on $\Omega$. Then, we can define the imprecise versions of the
Shannon entropy as proposed by [1, 35],
$\overline{H}(P):=\sup_{P\in\mathscr{P}}H(P)$ and
$\underline{H}(P):=\inf_{P\in\mathscr{P}}H(P)$, called the upper and lower
Shannon entropy, respectively.555In Appendix E, we provide bounds to the
values of upper and lower entropy. The upper entropy is a measure of total
uncertainty since it represents the minimum level of predictability associated
with the elements of $\mathscr{P}$. In [1, 35], the authors posit that it can
be decomposed as a sum of aleatoric and epistemic uncertainties, and that this
latter can be specified as the difference between upper and lower entropy,
thus obtaining

$\underbrace{\overline{H}(P)}_{\text{total
uncertainty}}=\underbrace{\underline{H}(P)}_{\text{aleatoric
uncertainty}}+\underbrace{\left[\overline{H}(P)-\underline{H}(P)\right]}_{\text{epistemic
uncertainty}}.$

We have the following proposition.

######  Proposition 7.

Let $\Pi,\Pi^{\prime}$ be sets of probability measures as the ones considered
in Remark 1. Then,
$\sup_{P\in\Pi}H(P)=\overline{H}(P)\leq\overline{H}(P^{\prime})=\sup_{P\in\Pi^{\prime}}H(P)$.

Proposition 7 tells us that the upper entropy of the extreme elements in
$\Pi=\text{ex}\Pi^{\prime}$ is a lower bound for the upper entropy of the
whole credal set $\Pi^{\prime}$. An immediate consequence of Proposition 7 is
that the lower entropy of the extreme elements in $\Pi$ is an upper bound for
the lower entropy of the whole credal set $\Pi^{\prime}$.

In the context of IBNNs, we are especially interested in EU and AU associated
with the predictive credal set, introduced in the next section. It is the set
of distributions on the output set $\mathcal{Y}$ induced by the posterior
distributions in $\mathcal{P}_{D}$. We also point out a salient feature of
this approach: the EU and AU computed from the predictive credal set embed
uncertainty comparable to an uncountably infinite ensemble of BNNs, i.e. an
ensemble of BNNs of cardinality $\aleph_{1}$, despite the simple and intuitive
mathematics over a finite set. They are not merely pessimistic bounds to the
uncertainties associated with a finite ensemble of BNNs.

##  4\. Practical aspects

In this section, we first illustrate how to elicit credal prior and likelihood
sets that are needed for step S1 of the procedure in section 3.1. Then, we
describe how IBNN are instrumental in specifying a set of outputs that enjoys
probabilistic guarantees.

Outlined in [39, Sections IV-B and IV-C1], for classification, the standard
process for BNNs involves

  * •

A Normal prior with zero mean and diagonal covariance $\sigma^{2}I$ on the
coefficients of the network, that is, $p(\theta)=\mathcal{N}(0,\sigma^{2}I)$.
In the context of IBNNs, we could specify e.g.
$\mathcal{P}=\\{P:p(\theta)=\mathcal{N}(\mu,\sigma^{2}I)\text{,
}\mu\in\\{\mu_{-},\mathbf{0},\mu_{+}\\}\text{, }\sigma^{2}\in\\{3,7\\}\\}$.
That is, the extreme elements of the prior credal set are five independent
Normals having different levels of “fatness” of the tails, and centered at a
vector $\mu_{+}$ having positive entries, a vector $\mu_{-}$ having negative
entries, and a vector $\mathbf{0}$ having entries equal to $0$. They capture
the ideas of positive bias, negative bias, and no bias of the coefficients,
respectively. This is done to hedge against possible prior misspecification.

  * •

A Categorical likelihood, $p(y\mid x,\theta)=\text{Cat}(\Phi_{\theta}(x))$,
whose parameter is given by the output of a functional model $\Phi_{\theta}$.
In the context of IBNNs, we could specify the set of extreme elements of the
likelihood credal set as
$\mathcal{L}_{x,\theta}=\\{P_{x,\theta}:p_{x,\theta}(y)=\text{Cat}(\Phi_{s,\theta}(x)),s\in\mathcal{S}\\}$,
where $\mathcal{S}=\\{1,\ldots,S\\}\subset\mathbb{N}$ is a generic index set.
Specifying set $\mathcal{L}_{x,\theta}$, then, corresponds to eliciting a
finite number $S$ of possible (parametrized) architectures for the neural
network, $\Phi_{s,\theta}$, $s\in\mathcal{S}=\\{1,\ldots,S\\}$, and obtain, as
a consequence, $S$ categorical distributions
$\\{\text{Cat}(\Phi_{s,\theta}(x))\\}_{s\in\mathcal{S}}$. This captures the
ambiguity around the true data generating process faced by the agent, and
allows them to hedge against likelihood misspecification.

More in general, we can use the priors and likelihoods that better fit the
type of analysis we are performing. For example, for the choice of the priors
we refer to [26], where the authors study the problem of selecting the right
type of prior for BNNs. We compute the
$N:=\\#\mathcal{P}\times\\#\mathcal{L}_{x,\theta}$ posteriors
$\mathfrak{post}(P,P_{x,\theta})$, for all $P\in\mathcal{P}$ and all
$P_{x,\theta}\in\mathcal{L}_{x,\theta}$. They are the elements of set
$\mathcal{P}_{D}$, whose convex hull $\text{Conv}\mathcal{P}_{D}$ represents
the posterior credal set. The posteriors $P_{D}\in\mathcal{P}_{D}$ are
distributions on the parameter space $\Theta$.666Assume for simplicity that
they all have density with respect to some ($\sigma$-finite) dominating
measure. Every such posterior induces a distribution on the output space
$\mathcal{Y}$ via the predictive distribution

$p(\tilde{y}\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})=\int_{\Theta}p(\tilde{y}\mid\theta,\tilde{x})\cdot
p(\theta\mid x_{1},y_{1},\ldots,x_{n},y_{n})\text{ d}\theta,$

where $\tilde{x}$ is a new input and $\tilde{y}$ is the associated output (see
Appendix F for more details). In symbols, we can write
$P_{D}\rightsquigarrow\hat{P}$, for all $P_{D}\in\mathcal{P}_{D}$, and let
$\hat{\mathcal{P}}:=\\{\hat{P}_{1},\ldots,\hat{P}_{N}\\}$. Its convex hull
$\text{Conv}\hat{\mathcal{P}}$ is the predictive credal set. To save notation,
we write $\hat{p}_{k}(\cdot)\equiv
p_{k}(\cdot\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})$, for all
$k\in\\{1,\ldots,N\\}$. Now, for every such distribution we compute the
$\alpha$-level Highest Density Region (HDR) $R(\hat{p}_{k}^{\alpha})$,
$\alpha\in(0,1)$. As defined in [11, Section 1], it is a subset of the output
space $\mathcal{Y}$ such that

$\int_{R(\hat{p}_{k}^{\alpha})}\hat{p}_{k}(y)\text{d}y\geq 1-\alpha\quad\text{
and }\quad\int_{R(\hat{p}_{k}^{\alpha})}\text{d}y\text{ is a minimum.}$

We need $\int_{R(\hat{p}_{k}^{\alpha})}\text{d}y$ to be a minimum because we
want $R(\hat{p}_{k}^{\alpha})$ to be the smallest possible region that gives
us the desired probabilistic coverage. Equivalently, we can write that
$R(\hat{p}_{k}^{\alpha})=\\{y\in\mathcal{Y}:\hat{p}_{k}(y)\geq\hat{p}_{k}^{\alpha}\\}$,
where $\hat{p}_{k}^{\alpha}$ is a constant value. In particular, it is the
largest constant such that $\hat{P}_{k}[y\in R(\hat{p}_{k}^{\alpha})]\geq
1-\alpha$ [34]. In dimension $1$, $R(\hat{p}_{k}^{\alpha})$ can be interpreted
as the narrowest interval – or union of intervals – in which the value of the
(true) output falls with probability of at least $1-\alpha$ according to
distribution $\hat{P}_{k}$. As we can see, HDR’s are a Bayesian counterpart of
confidence intervals. We give a simple visual example in Figure 5.

Figure 5. The $0.25$-HDR from a Normal Mixture density. This picture is a
replica of [34, Figure 1]. The geometric representation of “$75\%$ probability
according to $\hat{P}$” is the area between the pdf curve $\hat{p}(y)$ and the
horizontal bar corresponding to $\hat{p}^{0.25}$. A higher probability
coverage (according to $\hat{P}$) would correspond to a lower constant, so
$\hat{p}^{\alpha}<\hat{p}^{0.25}$, for all $\alpha<0.25$. In the limit, we
recover $100\%$ coverage at $\hat{p}^{0}=0$.

Finally, compute the $\alpha$-level Imprecise Highest Density Region (IHDR)
$IR_{\alpha}:=\cup_{k=1}^{N}R(\hat{p}_{k}^{\alpha})$. By taking the union of
the HDR’s, we ensure that all the probability measures in the predictive
credal set $\text{Conv}\hat{\mathcal{P}}$ assign probability of at least
$1-\alpha$ to the event $\\{y\in IR_{\alpha}\\}$; this is a consequence of
Lemma 2. In turn, this implies that $\underline{\hat{P}}(y\in
IR_{\alpha})=\inf_{k\in\\{1,\ldots,N\\}}\hat{P}_{k}(y\in IR_{\alpha})\geq
1-\alpha$. This can be interpreted as the event $A=\\{\text{The (true) output
belongs to }IR_{\alpha}\\}$ having lower probability of at least $1-\alpha$, a
probabilistic guarantee for the set of outputs $IR_{\alpha}$ generated by our
procedure.

This method allows to quantify and control the difference between the upper
and lower probabilities of $IR_{\alpha}$. To see this, notice that
$\overline{\hat{P}}(y\in IR_{\alpha})\leq 1$, so $\overline{\hat{P}}(y\in
IR_{\alpha})-\underline{\hat{P}}(y\in IR_{\alpha})\leq\alpha$. Notice also
that the size of $IR_{\alpha}$ is an increasing function of both AU and EU. As
a consequence, it is related, but it is not equal, to the AU the agent faces.
If we want to avoid to perform the procedure only to discover that
$IR_{\alpha}$ is “too big”, then we can add an “AU check” at the beginning.
This, together with computing $IR_{\alpha}$ in a classification setting, is
explored in Appendices G and H.

We conclude with a remark. We first point out that in applications the
predictive density $\hat{p}$ is computed using the VI approximation
$\breve{p}$ of the posterior density, that is, by solving
$\int_{\Theta}p(\tilde{y}\mid\theta,\tilde{x})\cdot\breve{p}(\theta\mid
x_{1},y_{1},\ldots,x_{n},y_{n})\text{ d}\theta$. This allows to obtain
$\hat{p}$ even when the exact posterior is intractable. Second, we observe
that the scholar is mostly interested in the AU and EU associated with the
credal set $\text{Conv}\hat{\mathcal{P}}$ of predictive distributions. This
because they are ultimately interested in reporting the uncertainty around the
predicted outputs given a new input in the problem at hand, more than the
uncertainty on the parameters of the BNNs and on the likelihoods. Finally, a
tacit assumption is that the EU and AU associated with the predictive credal
set $\text{Conv}\hat{\mathcal{P}}$ are similar whether we compute its elements
using the exact posteriors or their VI approximation. This is needed because
in real-world applications the agent would not be able to make the ambiguity
around the correct prior and likelihood percolate towards downstream tasks
without using variational inference (or other types of approximations).

##  5\. Experiments

In this section, we focus on the downstream task performance of IBNNs. We show
that they perform better than an ensemble of BNNs (EBNN) – introduced formally
in the next paragraph. To demonstrate the utility of the proposed method, we
study the behavior of certain safety-critical settings under distribution
shifts and its ramifications. One, for motion prediction in autonomous driving
scenarios, and two, to model blood glucose and insulin dynamics for artificial
pancreas control.

In [22], the authors pursue the same heuristic endeavor as we do, but take an
ensemble route. They too consider different BNNs, but instead of keeping them
separate and use them to build a predictive credal set, they average them out.
Similar to theirs, we elicit the following procedure, that we call ensemble of
BNNs (EBNN). Consider $k$ different BNNs, and compute the posterior
distribution on the parameters. They induce $k$ predictive distributions on
the output space $\mathcal{Y}$, each having mean $\mu_{j}$ and variance
$\sigma^{2}_{j}$. We call EBNN distribution $P_{\text{ens}}$ a Normal having
mean $\mu_{\text{ens}}=1/k\sum_{j=1}^{k}\mu_{j}$ and covariance matrix
$\sigma^{2}_{\text{ens}}I$, where
$\sigma^{2}_{\text{ens}}=1/k\sum_{j=1}^{k}\sigma^{2}_{j}+1/(k-1)\sum_{j=1}^{k}(\mu_{j}-\mu_{\text{ens}})^{2}$.
We use the $\alpha$-level HDR associated with $P_{\text{ens}}$ as a benchmark
for our $IR_{\alpha}$. Being rooted in IP theory, IBNNs are a better justified
procedure than EBNN from a theoretical standpoint. In this section, we show
with two applications that their implementation performs better than EBNN. We
do not consider Bayesian Model Averaging (BMA) as a baseline for IBNNs because
BMA suffers from pitfalls when used in the context of Bayesian neural networks
[36]. In addition, we do not compare IBNNs against belief tracking techniques
because these latter require extra assumptions that IBNNs do not, see Appendix
L.

Uncertainty Quantification Benchmark Experiments. In Appendix I, we include
details of experiments where we train IBNNs for image classification tasks for
standard datasets like CIFAR10 [44], SVHN [56], Fashion-MNIST [74], and MNIST
[49]. There, we discuss how IBNNs are better than EBNN at disentangling AU and
EU, and at quantifying them. While for IBNNs an increase in the corruption of
the dataset corresponds to an increase in both EU and AU, for EBNN the AU
either stays constant, or it decreases. This counterintuitive – and erroneous
– behavior makes us conclude that we are able to capture the uncertainties
better than the baseline.

###  5.1. Motion Prediction for Autonomous Racing

In this case study, we demonstrate the utility of IBNNs for motion prediction
in autonomous driving scenarios. An important challenge in autonomous driving
is understanding the intent of other agents and predicting their future
trajectories to allow for safety-aware planning. In autonomous racing, where
control is pushed to the dynamical limits, accurate and robust predictions are
even more essential for outperforming opponent agents while assuring safety.
Again, IBNNs provide a straightforward method for quantifying uncertainty and
deriving robust prediction regions for anticipating an agent’s behavior.

We use the problem settings in [69] to define the problem of obtaining
prediction sets for future positions of an autonomous racing agent. Our
results show that the prediction regions have improved coverage when compared
to EBNN. These results hold in both in-distribution and out-of-distribution
settings, which are described below.

#### Problem.

Let $O^{i}(t,l)\equiv O^{i}=\\{\pi^{i}_{t-l},\ldots,\pi^{i}_{t}\\}$ denote the
$i$-th trajectory instance of an agent at time $t$, consisting of the observed
positions from time $t-l$ up to time $t$. Let then $C^{i}$ be a time-invariant
context variable. Let also $F^{i}(t,h)\equiv
F^{i}=\\{\pi^{i}_{t+1},\dots,\pi^{i}_{t+h}\\}$ be the collection of the next
$h$ future positions. We wish to obtain a model $M$ that predicts region
$\mathcal{R}_{\alpha}$ with probabilistic guarantees. In particular, for EBNN
$\mathcal{R}_{\alpha}$ is the $\alpha$-level HDR of $P_{\text{ens}}$, so that
$P_{\text{ens}}(F^{i}\in\mathcal{R}_{\alpha})\geq 1-\alpha$, while for IBNN
$\mathcal{R}_{\alpha}=IR_{\alpha}$, so that
$\underline{\hat{P}}(F^{i}\in\mathcal{R}_{\alpha})\geq 1-\alpha$.

The dataset consists of instances of $(O^{i},F^{i})$ divided into a training
set $\mathcal{D}_{\text{train}}$ and a testing set
$\mathcal{D}_{\text{test}}$. We train an uncertainty-aware model on
$\mathcal{D}_{\text{train}}$ that computes the triplet
$(F^{i}_{l},F^{i}_{m},F^{i}_{u})=M(O^{i},C^{i})$ where $F^{i}_{l}$,
$F^{i}_{u}$, $F^{i}_{m}$ are the lower, upper, and mean predictions of the
future positions, respectively.

The dataset $\mathcal{D}_{\text{all}}$ is created by collecting simulated
trajectories of autonomous race cars in the F1Tenth-Gym [57] (details in
[69]). As shown in Figure 6, different racing lines were utilized including
the center, right, left, and optimal racing line for the Spielberg track.

Figure 6. Motion Prediction for F1Tenth-Gym Environment [57]. Data is
collected by simulating various racing lines on the Spielberg Track.

We denote these by $\mathcal{D}_{\text{center}}$,
$\mathcal{D}_{\text{right}}$, $\mathcal{D}_{\text{left}}$, and
$\mathcal{D}_{\text{race}}$, respectively. Position $\pi$ is a vector
$\pi=(a,b,\vartheta,v)^{\top}$, where $a$ and $b$ are coordinates in a
$2$-dimensional Euclidean space, and $\vartheta$ and $v$ are the heading and
speed, respectively. In total, the $\mathcal{D}_{\text{all}}$ consists of
$34686$ train instances, $4336$ validation instances, and $4336$ test
instances.

#### In-distribution vs. Out-of-distribution.

We consider the prediction task to be in-distribution when
$\mathcal{D}_{\text{train}},\mathcal{D}_{\text{test}}\subset\mathcal{D}_{\text{all}}$.
It is out-of-distribution (OOD) when
$\mathcal{D}_{\text{train}}\subset\mathcal{D}_{\text{center}}\cup\mathcal{D}_{\text{right}}\cup\mathcal{D}_{\text{left}}$
and $\mathcal{D}_{\text{test}}\subset\mathcal{D}_{\text{race}}$.

#### Metrics.

We train the ensemble of BNNs and the IBNN models, $M_{\text{ens}}$ and
$M_{\text{IBNN}}$ respectively, using the same architecture and different
seeds. We compare the performance with respect to the test set by computing
the single-step coverage, where each prediction time-step is treated
independently, and the multi-step coverage, which considers the entire
$h$-step prediction.

Figure 7 (A) depicts a sample of the in-distribution evaluation for each of
the models. For a given trajectory, the red boxes indicate when the prediction
region did not cover the actual trajectory at that time-step. Qualitatively,
$M_{\text{IBNN}}$ has less missed timesteps when compared to $M_{\text{ens}}$.
Table 1 shows that IBNNs perform better in terms of both one-step and multi-
step coverage. Similar results can be observed for the OOD scenario. There,
all models were trained on racing lines which are predominantly parallel to
the track curvature. As a consequence, when the test set consists of instances
with higher curvatures, the overall coverage of all models degrades. This can
be seen in Figure 7 (B), where the prediction of the models (orange) tends to
be straight while the actual trajectory is more curved (green). Despite this,
the figure and the coverage metrics in Table 1 show how IBNN exhibits a more
robust behavior.

In-Distribution Results
---
| Ensemble | IBNN
$1-\alpha$ | $0.9$ | $0.95$ | $0.99$ | $0.9$ | $0.95$ | $0.99$
One-step | 0.962 | 0.980 | 0.992 | 0.992 | 0.995 | 0.997
Multi-step | 0.638 | 0.826 | 0.937 | 0.914 | 0.948 | 0.979
Out-of-Distribution Results
| Ensemble | IBNN
$1-\alpha$ | $0.9$ | $0.95$ | $0.99$ | $0.9$ | $0.95$ | $0.99$
One-step | 0.919 | 0.950 | 0.980 | 0.979 | 0.988 | 0.995
Multi-step | 0.532 | 0.703 | 0.860 | 0.825 | 0.884 | 0.943
Table 1. F1Tenth coverage results. We report one-step coverage and multi-step
coverage across 3 different values of $\alpha$. IBNNs exceed coverage of EBNNs
in all settings.

(a)

(b)

Figure 7. Left: F1Tenth in-distribution Results. Given an input of past
observations, IBNNs exhibit better coverage of the future target trajectory.
Predictions which do not cover the target within the desired $1-\alpha$ level
are indicated in red. Right: F1Tenth Out-of-distribution (OOD) Results. Robust
performance is exhibited by IBNNs when compared to EBNN in OOD settings.

###  5.2. Artificial Pancreas Control

Overall Setup. In this next case study we consider the problem of data-driven
control of human blood glucose-insulin dynamics, using an artificial pancreas
system, Figure 8.

Figure 8. The Bayesian neural networks predict a future blood glucose value.
These individual predictions are combined to get a robust estimate of the true
value as an interval. This is used by the MPC control algorithm to recommend
insulin dosage for the patient. The patient block in our experiment is
simulated using the virtual patient models from the UVa-Padova simulator.

External insulin delivery is accomplished by using an insulin pump controlled
by the artificial pancreas software, which attempts to regulate the blood-
glucose (BG) level of the patient within the euglycemic range of
$[70,180]mg/dl$ [46]. Levels below $70mg/dl$ lead to hypoglycemia, which can
lead to loss of consciousness, coma or even death. On the other hand, levels
above $300mg/dl$ lead to a condition called the ketoacidosis, where the body
can break down fat due to lack of insulin, and lead to build up of ketones. In
order to treat this situation, patients receive external insulin delivery
through insulin pumps. Artificial Pancreas (AP) systems can remedy this
situation by measuring the blood glucose level, and automatically injecting
insulin into the blood stream. Thus, we define the _unsafe regions_ of the
space as $G(t)\in(-\infty,70)\cup(300,\infty)$, where $G(t)$ is the BG value
at time $t$. This is the shaded region in Figure 9.

Figure 9. Starting from an initial glucose value, the task of the artificial
pancreas controller is to maintain the blood glucose value within safe
operating limits using insulin as a mode of control.

Neural Network Models and Controller. Deep neural networks are effective in
capturing the BG-insulin dynamics for personalized medical devices [46]. This
allows for improved device performance. Even though standard feedforward
neural networks can be used, Bayesian neural networks (BNN), and especially a
collection of multiple BNNs, offer a better alternative towards uncertainty
aware predictions. Here, we test the ramifications of these prediction sets,
when used inside an online receding horizon control scheme for insulin
delivery. We use the standard MPC control scheme for this purpose, well-known
in the literature [21]. More formally, let $G(t)$ and $I(t)\equiv I_{t}$ be
the blood-glucose and insulin values at time $t$, respectively. We denote the
finite length trajectory of length $H$ as
$\overleftarrow{G}_{H}(t):=[G(t-H+1),\ldots,G(t)]$, and
$\overleftarrow{I}_{H}(t):=[I(t-H+1),\ldots,I(t)]$. An uncertainty aware model
$M$ computes the triplet
$(G_{l}(t+l),G_{m}(t+l),G_{u}(t+l))=M(\overleftarrow{G}_{H}(t),\overleftarrow{I}_{H}(t))$,
where $G_{m}$ is the mean prediction output, and $G_{l},G_{u}$ are the lower
and upper predictions of the glucose value, respectively. By design, it is
true that $G_{l}\leq G_{m}\leq G_{u}$. A model predictive control algorithm –
whose cost function we denote by $J$ – solves
$\operatorname*{arg\,min}_{I_{0},I_{1},\dots,I_{k-1}}\sum_{i=0}^{k-1}J(M(\overleftarrow{G}_{H}(t+i),\overleftarrow{I}_{H}(t+i)))$.

After every time step, the control algorithm picks the first insulin input
$I_{0}$ as the insulin bolus for the patient, and discards the rest. Cost
function $J$ takes into account three factors, (i) Distance of the mean
prediction level $G_{m}$ at each time step from a target value of $120mg/dl$,
(ii) Distance of upper and lower predictions ($G_{u}$ and $G_{l}$) from the
unsafe regions of the state space $G(t)>300$ and $G(t)<70$, and (iii) Total
insulin injected $\sum_{t=0}^{k-1}I_{t}$.

Starting with some initial glucose value $G(0)$, we measure the performance of
the artificial pancreas controller as the fraction of time it spends in the
unsafe regions,

$t_{\text{unsafe}}=\frac{1}{L}\sum_{t=1}^{L}\mathbbm{1}\left\\{{G(t)\in(-\infty,70)\cup(300,\infty)}\right\\},$

where $\mathbbm{1}\\{\cdot\\}$ denotes the indicator function. A lower value
is more desirable. We compare EBNN and IBNNs as different realizations of the
model $M$.

Distribution Shift using Meals. A well known problem with learnt models is
distribution shift. Bayesian neural networks can address this issue by
apprising the end user of the increased uncertainty. For regression models of
the type described above, this appears as larger prediction intervals
$[G_{l},G_{u}]$. The artificial pancreas controller can run into this
situation in the following way: the insulin-glucose time series data collected
for training the data-driven model $M$ can be without meals, while at test
time the patient can have meals. This creates a distribution shift between the
training and test time data. Fortunately, the UVa-Padova simulator [14] allows
us to create datasets with and without meal inputs. In this case study, the
training data was obtained by randomly initializing the BG value in the range
$[120,190]$, and simulating the patient for $720$ minutes. The controller was
executed at $5$ minutes intervals. At test time the patient was supplied meals
at specific time intervals (for details, see Appendix K). This creates a
significant distribution shift since meals are effectively an unknown variable
which can affect the system state. However, from the controller’s perspective
this is practical, since patients can have unannounced meals.

Results and Discussion. To capture the difference in performance between EBNN
and IBNN, we compute
$\text{Perf}_{\text{diff}}:=(t_{\text{unsafe}}^{\text{EBNN}}-t_{\text{unsafe}}^{\text{IBNN}})/t_{\text{unsafe}}^{\text{EBNN}}$.
Both $t_{\text{unsafe}}^{\text{EBNN}}$ and $t_{\text{unsafe}}^{\text{IBNN}}$
depend on interval $[G_{l},G_{u}]$; for EBNN, this latter corresponds to the
$\alpha$-level HDR associated with EBNN distribution $P_{\text{ens}}$, while
for IBNN it corresponds to the IHDR $IR_{\alpha}$. We consider one case in
which an IBNN is trained using a credal prior set and only one likelihood (we
choose different seeds which initialize the prior distributions but we keep
the same architecture for the BNNs), and another case in which we do the
opposite (we use the same seed and different architectures).

We report $\text{Perf}_{\text{diff}}$, across different choices in Table 2. We
observe that a more conservative estimate, as informed by the IBNN model as
compared to the EBNN framework, results in controllers which respect the
safety limits better. To see that IBNNs are more conservative than EBNNs,
notice that when combining predictive distributions from multiple BNNs, an
IBNN combines the predictions via a finitely generated credal set (FGCS) whose
extrema are the individual distributions. On the contrary, EBNNs take an
average of the individual distributions to compute the ensemble distribution
$P_{\text{ens}}$. The union of the HDRs of the predictive distributions is
more conservative (i.e., broader) than the HDR of the single ensemble
distribution. While the approach by EBNNs seems like a reasonable choice on
the surface, it falls short in capturing the uncertainty necessary for the
downstream task. For more details on this case study, see Appendix K.

$1-\alpha$ |  $0.9$ |  $0.95$ |  $0.99$
---|---|---|---
Varying Seeds |  2.3% |  3.5% |  5.2%
Varying |  |  |
Architectures |  0.5% |  -3.8% |  4.4%
Table 2. We report the performance improvements when using IBNNs as compared
to EBNNs across $3$ different values of $\alpha$. Row $1$ corresponds to the
case where the individual BNNs are trained with different seeds for the prior
distribution; and Row $2$ is the case when the BNNs have different
architectures.

##  6\. Related Work

In [12], the authors introduce credal classifiers (CCs) as a generalization of
classifiers based on Bayesian networks. Unlike CCs, IBNNs do not require
independence assumptions between non-descendant, non-parent variables. In
addition, IBNNs avoid NP-hard complexity issues of searching for optimal
structure in the space of Bayesian networks [10]. In [52], an epistemic
convolutional neural network (ECNN) is developed that explicitly models the
epistemic uncertainty induced by training data of limited size and quality. A
clear distinction is that ECNNs measure uncertainty in target-level
representations whereas IBNNs identify the uncertainty measure on the output
space. Despite the merit of their work, we believe IBNNs achieve greater
generality, since they are able to quantify aleatoric and epistemic
uncertainty and are applicable to problems beyond classification. For a review
of the state of the art concerning the distinction between EU and AU we refer
to [35] and to [52]. We also point out how IBNNs have been recently used to
solve prior-likelihood conflicts in Bayesian statistics [54]. Further
references can be found in Appendix M.

##  7\. Conclusion

We presented IBNNs, a generalization of BNNs that allows to distinguish
between AU and EU, and to quantify them. We showed how they can be used to
specify a set of outputs that enjoys probabilistic guarantees, and we applied
them to safety-critical settings.

We point out how a region that improves $IR_{\alpha}$, meaning that it would
be tighter, is the following:
$IR^{\prime}_{\alpha}=\\{y\in\mathcal{Y}:\underline{\hat{p}}(y)\geq\underline{\hat{p}}^{\alpha}\\}$,
where $\underline{\hat{p}}$ is the lower pdf/pmf evaluated at $y$, and
$\underline{\hat{p}}^{\alpha}$ is the largest constant such that
$\underline{\hat{P}}(y\in IR^{\prime}_{\alpha})\geq 1-\alpha$. The problem
with $IR^{\prime}_{\alpha}$ is that, while the highest density regions
associated with the predictive distributions can be computed using off-the-
shelf tools, calculating $\underline{\hat{p}}$ and
$\underline{\hat{p}}^{\alpha}$ would have been much more computationally
expensive. In addition, it would have required to come up with a new technique
to find $\underline{\hat{p}}$ and $\underline{\hat{p}}^{\alpha}$. We defer
studying this to future work.

We also plan to apply IBNNs to continual learning (CL) to overcome the curse
of dimensionality and to capture an agent’s preference over the tasks to
perform.

Furthermore, we intend to relate IBNNs to Bayesian model selection (BMS) [27].
This latter suffers from the same problem as “regular” Bayesian inference.
That is, while it tries to come up with a sophisticate prior that induces
shrinkage, it still relies on the “correctness” of that prior, i.e. on
correctly specifying the prior’s parameters. In the future, an interesting way
of combining IBNNs with BMS will be to use a finite number of regularized
horseshoe priors, as suggested by [27, Section 3.2], as extreme elements of
the prior credal set.

We also call attention to the fact that an IBNN is a model-based approach. The
relationship with model-free approaches such as conformal prediction will be
the object of future studies. In particular, we are interested in finding in
which cases IHDRs are narrower than conformal regions, and vice versa.

Finally, we point out how one possible way of easing the burden of the
combinatorial task in step S2 of the IBNN procedure is to specify a prior
credal set whose size strikes the perfect balance between being “vague enough”
so that we do not underestimate the EU, and being “small enough” so that the
IBNN is actually implementable. We suspect conjugacy of the priors may play a
key role in this endeavor. Because of its centrality, we defer the study of
“optimal prior credal sets” to future work.

## References

  * [1] Joaquín Abellán, George Jiří Klir, and Serafín Moral.  Disaggregated total uncertainty measure for credal sets.  International Journal of General Systems, 1(35):29–44, 2006.
  * [2] Joaquín Abellán and Serafín Moral.  A non-specificity measure for convex sets of probability distributions.  International Journal of Uncertainty, Fuzziness and Knowledge-Based Systems, 8(3):357–367, 2000.
  * [3] Alexander Amini, Wilko Schwarting, Ava Soleimany, and Daniela Rus.  Deep evidential regression.  In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 14927–14937. Curran Associates, Inc., 2020.
  * [4] David Avis and Komei Fukuda.  Reverse search for enumeration.  Discrete Applied Mathematics, 65:21–46, 1996.
  * [5] James O. Berger.  The robust Bayesian viewpoint.  In Joseph B. Kadane, editor, Robustness of Bayesian Analyses. Amsterdam : North-Holland, 1984.
  * [6] Jose M. Bernardo.  Reference posterior distributions for Bayesian inference.  Journal of the Royal Statistical Society: Series B, 41(2):113–128, 1979.
  * [7] Yuhang Cai and Lek-Heng Lim.  Distances between probability distributions of different dimensions.  IEEE Transactions on Information Theory, 68(6):4020–4031, 2022\.
  * [8] Michele Caprio and Ruobin Gong.  Dynamic precise and imprecise probability kinematics.  In Enrique Miranda, Ignacio Montes, Erik Quaeghebeur, and Barbara Vantaggi, editors, Proceedings of the Thirteenth International Symposium on Imprecise Probability: Theories and Applications, volume 215 of Proceedings of Machine Learning Research, pages 72–83. PMLR, 11–14 Jul 2023\.
  * [9] Simone Cerreia-Vioglio, Fabio Maccheroni, and Massimo Marinacci.  Ergodic theorems for lower probabilities.  Proceedings of the American Mathematical Society, 144:3381–3396, 2015.
  * [10] David M. Chickering, David Heckerman, and Christopher Meek.  Large-sample learning of Bayesian networks is NP-hard.  The Journal of Machine Learning Research, 5:1287–1330, 2004.
  * [11] Frank P. A. Coolen.  Imprecise highest density regions related to intervals of measures.  Memorandum COSOR, 9254, 1992.
  * [12] Giorgio Corani, Alessandro Antonucci, and Marco Zaffalon.  Bayesian Networks with Imprecise Probabilities: Theory and Application to Classification, chapter 4 of Data Mining: Foundations and Intelligent Paradigms: Volume 1: Clustering, Association and Classification, pages 49–93.  Berlin, Germany : Springer, 2012.
  * [13] Fabio Gagliardi Cozman.  Credal networks.  Artificial Intelligence, 120:199–233, 2000.
  * [14] Chiara Dalla Man, Francesco Micheletto, Dayu Lv, Marc D. Breton, Boris Kovatchev, and Claudio Cobelli.  The UVA/PADOVA type 1 diabetes simulator: New features.  Journal of Diabetes Science and Technology, 8(1):26–34, 2013.
  * [15] Bruno de Finetti.  Theory of Probability, volume 1.  New York : Wiley, 1974.
  * [16] Bruno de Finetti.  Theory of Probability, volume 2.  New York : Wiley, 1975.
  * [17] Thierry Denœux.  An evidential neural network model for regression based on random fuzzy numbers.  In Sylvie Le Hégarat-Mascle, Isabelle Bloch, and Emanuel Aldea, editors, Belief Functions: Theory and Applications, pages 57–66, Cham, 2022\. Springer International Publishing.
  * [18] Thierry Denœux.  Quantifying prediction uncertainty in regression using random fuzzy sets: the ENNreg model.  IEEE Transactions on Fuzzy Systems, pages 1–10, 2023.
  * [19] Stefan Depeweg, Jose-Miguel Hernandez-Lobato, Finale Doshi-Velez, and Steffen Udluft.  Decomposition of uncertainty in Bayesian deep learning for efficient and risk-sensitive learning.  In International Conference on Machine Learning, pages 1184–1193. PMLR, 2018.
  * [20] Didier Dubois and Eyke Hüllermeier.  Comparing probability measures using possibility theory: A notion of relative peakedness.  International Journal of Approximate Reasoning, 45(2):364–385, 2007\.  Eighth European Conference on Symbolic and Quantitative Approaches to Reasoning with Uncertainty (ECSQARU 2005).
  * [21] Souradeep Dutta, Taisa Kushner, and Sriram Sankaranarayanan.  Robust data-driven control of artificial pancreas systems using neural networks.  In Milan Češka and David Šafránek, editors, Computational Methods in Systems Biology, pages 183–202, Cham, 2018. Springer International Publishing.
  * [22] Romain Egele, Romit Maulik, Krishnan Raghavan, Prasanna Balaprakash, and Bethany Lusch.  Autodeuq: Automated deep ensemble with uncertainty quantification.  CoRR, abs/2110.13511, 2021.
  * [23] Daniel Ellsberg.  Risk, ambiguity, and the Savage axioms.  The Quarterly Journal of Economics, 75(4):643–669, 1961.
  * [24] Emlyn Flint, Florence Chikurunhe, and Anthony Seymour.  Regime-based tactical allocation for equity factors and balanced portfolios.  SSRN Electronic Journal, 01 2017.
  * [25] Benoît Fortin, Samir Hachour, and François Delmotte.  Multi-target PHD tracking and classification using imprecise likelihoods.  International Journal of Approximate Reasoning, 90:17–36, 2017\.
  * [26] Vincent Fortuin, Adrià Garriga-Alonso, Mark van der Wilk, and Laurence Aitchison.  BNNpriors: A library for bayesian neural network inference with different prior distributions.  Software Impacts, 9:100079, 2021.
  * [27] Soumya Ghosh, Jiayu Yao, and Finale Doshi-Velez.  Model selection in Bayesian neural networks via horseshoe priors.  Journal of Machine Learning Research, 20(182):1–46, 2019.
  * [28] Itzhak Gilboa and Massimo Marinacci.  Ambiguity and the Bayesian paradigm.  In Daron Acemoglu, Manuel Arellano, and Eddie Dekel, editors, Advances in Economics and Econometrics, Tenth World Congress, volume 1. Cambridge : Cambridge University Press, 2013.
  * [29] Ethan Goan and Clinton Fookes.  Bayesian Neural Networks: An Introduction and Survey, pages 45–87.  Cham, Switzerland : Springer International Publishing, 2020.
  * [30] Ruobin Gong and Xiao-Li Meng.  Judicious judgment meets unsettling updating: dilation, sure loss, and Simpson’s paradox.  Statistical Science, 36(2):169–190, 2021.
  * [31] Chuan Guo, Geoff Pleiss, Yu Sun, and Kilian Q Weinberger.  On calibration of modern neural networks.  In International conference on machine learning, pages 1321–1330. PMLR, 2017.
  * [32] Dan Hendrycks and Thomas G. Dietterich.  Benchmarking neural network robustness to common corruptions and perturbations.  CoRR, abs/1903.12261, 2019.
  * [33] Peter J. Huber and Elvezio M. Ronchetti.  Robust statistics.  Wiley Series in Probability and Statistics. Hoboken, New Jersey : Wiley, 2nd edition, 2009.
  * [34] Rob J. Hyndman.  Computing and graphing highest density regions.  The American Statistician, 50(2):120–126, 1996.
  * [35] Eyke Hüllermeier and Willem Waegeman.  Aleatoric and epistemic uncertainty in machine learning: an introduction to concepts and methods.  Machine Learning, 3(110):457–506, 2021.
  * [36] Pavel Izmailov, Patrick Nicholson, Sanae Lotfi, and Andrew G Wilson.  Dangers of bayesian model averaging under covariate shift.  In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 3309–3322. Curran Associates, Inc., 2021.
  * [37] Pavel Izmailov, Sharad Vikram, Matthew D Hoffman, and Andrew Gordon Gordon Wilson.  What are Bayesian neural network posteriors really like?  In Marina Meila and Tong Zhang, editors, Proceedings of the 38th International Conference on Machine Learning, volume 139, pages 4629–4640. PMLR, 2021.
  * [38] Harold Jeffreys.  An invariant form for the prior probability in estimation problems.  Proceedings of the Royal Society A, 186(1007):453–461, 1946.
  * [39] Laurent Valentin Jospin, Hamid Laga, Farid Boussaid, Wray Buntine, and Mohammed Bennamoun.  Hands-on Bayesian neural networks – A tutorial for deep learning users.  IEEE Computational Intelligence Magazine, 17(2):29–48, 2022.
  * [40] Sanyam Kapoor, Wesley J Maddox, Pavel Izmailov, and Andrew Gordon Wilson.  On uncertainty, tempering, and data augmentation in bayesian classification.  arXiv preprint arXiv:2203.16481, 2022.
  * [41] Alex Kendall and Yarin Gal.  What uncertainties do we need in bayesian deep learning for computer vision?  Advances in neural information processing systems, 30, 2017.
  * [42] John Klein, Christèle Lecomte, and Pierre Miché.  Hierarchical and conditional combination of belief functions induced by visual tracking.  International Journal of Approximate Reasoning, 51(4):410–428, 2010\.
  * [43] Ranganath Krishnan, Mahesh Subedar, and Omesh Tickoo.  MOPED: efficient priors for scalable variational inference in bayesian deep neural networks.  CoRR, abs/1906.05323, 2019.
  * [44] Alex Krizhevsky, Vinod Nair, and Geoffrey Hinton.  CIFAR-10 (Canadian Institute for Advanced Research).  Available at http://www.cs.toronto.edu/ kriz/cifar.html, 2009.
  * [45] Meelis Kull and Peter A. Flach.  Reliability Maps: A Tool to Enhance Probability Estimates and Improve Classification Accuracy.  In Toon Calders, Floriana Esposito, Eyke Hüllermeier, and Rosa Meo, editors, Machine Learning and Knowledge Discovery in Databases, Lecture Notes in Computer Science, pages 18–33, Berlin, Heidelberg, 2014. Springer.
  * [46] Taisa Kushner, David Bortz, David M. Maahs, and Sriram Sankaranarayanan.  A data-driven approach to artificial pancreas verification and synthesis.  In Proceedings of the 9th ACM/IEEE International Conference on Cyber-Physical Systems, ICCPS ’18, pages 242–252. IEEE Press, 2018.
  * [47] Jouko Lampinen and Aki Vehtari.  Bayesian approach for neural networks - review and case studies.  Neural Networks, 4:257–274, 2001.
  * [48] Daniel Lassiter.  Representing credal imprecision: from sets of measures to hierarchical Bayesian models.  Philosophical Studies, 177(6):1463–1485, 2020.
  * [49] Yann Lecun, Léon Bottou, Yoshua Bengio, and Patrick Haffner.  Gradient-based learning applied to document recognition.  Proceedings of the IEEE, 86(11):2278–2324, 1998.
  * [50] Peter Lenk and Bryan Orme.  The value of informative priors in Bayesian inference with sparse data.  Journal of Marketing Research, 46(6):832–845, 2009.
  * [51] Isaac Levi.  The Enterprise of Knowledge.  London, UK : MIT Press, 1980.
  * [52] Shireen Kudukkil Manchingal and Fabio Cuzzolin.  Epistemic deep learning.  Available at arxiv:2206.07609, 2022.
  * [53] Massimo Marinacci and Luigi Montrucchio.  Introduction to the mathematics of ambiguity.  In Itzhak Gilboa, editor, Uncertainty in economic theory: a collection of essays in honor of David Schmeidler’s 65th birthday. London : Routledge, 2004.
  * [54] Alexander Marquardt, Julian Rodemann, and Thomas Augustin.  An empirical study of prior-data conflicts in Bayesian neural networks.  Poster presented at 13th International Symposium on Imprecise Probabilities: Theories and Applications (ISIPTA) 2023, 2023.
  * [55] Norman Mu and Justin Gilmer.  MNIST-C: A robustness benchmark for computer vision.  CoRR, abs/1906.02337, 2019.
  * [56] Yuval Netzer, Tao Wang, Adam Coates, Alessandro Bissacco, Bo Wu, and Andrew Y. Ng.  Reading digits in natural images with unsupervised feature learning.  In NIPS Workshop on Deep Learning and Unsupervised Feature Learning 2011, 2011.
  * [57] Matthew O’Kelly, Hongrui Zheng, Dhruv Karthik, and Rahul Mangharam.  F1TENTH: An Open-source Evaluation Environment for Continuous Control and Reinforcement Learning.  In Proceedings of the NeurIPS 2019 Competition and Demonstration Track, pages 77–89. PMLR, 2020.
  * [58] Sangdon Park, Osbert Bastani, Nikolai Matni, and Insup Lee.  PAC confidence sets for deep neural networks via calibrated prediction.  In 8th International Conference on Learning Representations, 2020\.
  * [59] Luis Raul Pericchi.  Sets of prior probabilites and bayesian robustness.  Documentation Section on the website of the Society for Imprecise Probability Theory and Applications (SIPTA), 1998.
  * [60] Hippolyt Ritter, Aleksandar Botev, and David Barber.  A scalable laplace approximation for neural networks.  In International Conference on Learning Representations, 2018.
  * [61] Robin Senge, Stefan Bösner, Krzysztof Dembczyński, Jörg Haasenritter, Oliver Hirsch, Norbert Donner-Banzhoff, and Eyke Hüllermeier.  Reliable classification: Learning classifiers that distinguish aleatoric and epistemic uncertainty.  Information Sciences, 255:16–29, January 2014.
  * [62] Murat Sensoy, Lance Kaplan, and Melih Kandemir.  Evidential deep learning to quantify classification uncertainty.  In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
  * [63] Marek Smieja and Jacek Tabor.  Entropy of the mixture of sources and entropy dimension.  IEEE Transactions on Information Theory, 58(5):2719–2728, 2012\.
  * [64] Lewis Smith and Yarin Gal.  Understanding measures of uncertainty for adversarial example detection.  arXiv preprint arXiv:1803.08533, 2018.
  * [65] Bjørnar Tessem.  Interval probability propagation.  International Journal of Approximate Reasoning, 7(3):95–120, 1992\.
  * [66] D. Michael Titterington.  Bayesian methods for neural networks and related models.  Statistical Science, 19(1):128–139, 2004.
  * [67] Krasymyr Tretiak, Georg Schollmeyer, and Scott Ferson.  Neural network model for imprecise regression with interval dependent variables.  Available at arXiv:2206.02467, 2022.
  * [68] Matthias C.M. Troffaes and Gert de Cooman.  Lower Previsions.  Chichester, United Kingdom : John Wiley and Sons, 2014.
  * [69] Renukanandan Tumu, Lars Lindemann, Truong Nghiem, and Rahul Mangharam.  Physics Constrained Motion Prediction with Uncertainty Quantification, 2023.  arXiv:2302.01060.
  * [70] Peter Walley.  Coherent lower (and upper) probabilities.  Technical report, University of Warwick, Coventry, 1981.
  * [71] Peter Walley.  Statistical Reasoning with Imprecise Probabilities, volume 42 of Monographs on Statistics and Applied Probability.  London : Chapman and Hall, 1991.
  * [72] Hao Wang and Dit-Yan Yeung.  A survey on Bayesian deep learning.  ACM Computing Surveys, 53(5):1–37, 2021.
  * [73] Larry A. Wasserman and Joseph B. Kadane.  Bayes’ theorem for Choquet capacities.  The Annals of Statistics, 18(3):1328–1339, 1990.
  * [74] Han Xiao, Kashif Rasul, and Roland Vollgraf.  Fashion-MNIST: a novel image dataset for benchmarking machine learning algorithms.  CoRR, abs/1708.07747, 2017.
  * [75] Fengshuo Zhang and Chao Gao.  Convergence rates of variational posterior distributions.  The Annals of Statistics, 48(4):2180 – 2207, 2020.

##  Appendix A Why do we need IPs?

The main motivations for training an artificial stochastic neural network
using credal sets are two. Let $(\Omega,\mathcal{F})$ be the measurable space
of interest.

  * (i)

A single probability distribution does not suffice to represent ignorance in
the sense of lack of knowledge; this is well documented in the literature, see
e.g. [35] and references therein. Consider the example of complete ignorance
(CI) in the case of a finite state space $\Omega$ [35, Section 3.3]. In
standard Bayesian analysis, CI is modeled in terms of the uniform distribution
$\text{Unif}(\Omega)$; this is justified by Laplace’s “principle of
indifference”. Then, however, it is not possible to distinguish between
precise probabilistic knowledge about a random event – called prior
indifference; think of the tossing of a fair coin – and a complete lack of
knowledge due to an incomplete description of the experiment – called prior
ignorance. Another problem is given by the additive nature of probability
distributions. Consider again the example of a uniform distribution. First,
let us observe that it is not invariant under reparametrization. In addition,
if we model the ignorance about the length $x$ of the side of a cube in
$\mathbb{R}^{3}$ via a uniform measure on the interval
$[l,u]\subset\mathbb{R}$, then this does not yield a uniform distribution of
$x^{3}$ on $[l^{3},u^{3}]$, which suggests some degree of informedness about
the cube’s volume. Finally, as pointed out in [71], if we ask a subject – even
an expert – about their opinion regarding some events, it is much more likely
that they will report interval of probabilities rather than single values.

  * (ii)

Working with credal sets allows to achieve prior and likelihood robustness:
realistically large sets $\mathcal{P}$ of priors and $\mathcal{L}$ of
likelihoods are elicited. Using credal sets, the agent recognizes that prior
beliefs and knowledge about the sampling model are limited and imprecise.
Combining each pair of functions in $\mathcal{P}$ and $\mathcal{L}$ using
Bayes’ rule, a class of posterior distributions – reflecting the updated state
of uncertainty – is formed. If the available information is not sufficient to
identify a unique posterior distribution, or a set of posteriors whose
diameter is small, credal sets allow to represent indecision, thus leading to
a less informative but more robust conclusions.777Here “diameter” has to be
understood as the distance between upper and lower probability of event $A$,
for all $A\in\mathcal{F}$.

##  Appendix B On the use of credal sets

Let us address a critique raised against the use of credal sets. [48] argues
against the use of sets of probabilities to model an agent’s prior beliefs and
their knowledge of the sampling model, while debating in favor of using
hierarchical Bayesian models. As reported in [35, Secton 4.6.2], the argument
against credal sets that is more cogent for the machine learning literature is
that modeling a lack of knowledge in a set-based manner may hamper the
possibility of inductive inference, up to a point where learning from
empirical data is not possible any more. With this, we mean the following. As
[59] points out, the natural candidate for a class of priors to represent
complete ignorance is the class $\mathcal{P}_{\text{all}}$ of all
distributions. When this class leads to non-vacuous and useful conclusions,
these are quite compelling and uncontroversial. It turns out that the
posterior probabilities obtained from this class are vacuous, that is, their
lower and upper bounds are $0$ and $1$: no finite sample is enough to
annihilate a sufficiently extreme prior belief. There is then a compromise to
be made, and this is the compromise of near-ignorance. The near-ignorance
class should be vacuous a priori in some respects, typically the ones that are
the most important for the analysis at hand. This way of proceeding is labeled
as arbitrary by [48], who instead advocates for the use of hierarchical
Bayesian procedures. We find this critique not compelling, as during the
analysis the job of the agent is to model reality: as pointed out in [35,
Secton 5], statistical inference is not possible without underlying
assumptions, and conclusions drawn from data are always conditional on those
assumptions. If we were to work every time with the maximum level of
generality, we would hardly be able to reach any conclusions. For example, in
a statistical analysis we never consider the state $\Omega$ of apparently
possible states [71, section 2.1.2], that is, the one that contains all the
states $\omega$ that are logically consistent with the available information.
If we consider a coin toss, we let the state space be $\Omega=\\{\text{heads,
tails}\\}$, certainly not $\Omega=\\{$heads, tails, coin landing on its edge,
coin braking into pieces on landing, coin disappearing down a crack in the
floor$\\}$. The same holds for sets of probabilities: it makes much more sense
to work with near-ignorance credal sets than to work with
$\mathcal{P}_{\text{all}}$. A final reason to rebut the point in [48] is that
the problems indicated in (i) in section A that make the use of the uniform
prior distribution – often interpreted as representing epistemic uncertainty
in standard Bayesian inference – at least debatable are inherited by
hierarchical Bayesian modeling, as specified in [6, 35, 38]. Furthermore a
single distribution that is the result of a hierarchical Bayesian procedure is
unable to gauge epistemic uncertainty [35].

###  B.1. A note on how to specify credal sets

As pointed out in [12, Section 3.3], there is a way of obtaining credal sets
starting from sets of probability intervals; in addition, standard algorithms
can compute the extreme elements of a credal set for which a probability
interval has been provided [4]. However, the resulting number of extrema is
exponential in the size of the possibility space [65];888Recall that the
possibility space of a random variable is the space of the values it can take
on. for this reason we prefer to specify prior and likelihood finitely
generated credal sets instead.

##  Appendix C A further IP concept: the core

Because of the conjugacy property of upper and lower probabilities, let us
focus on upper probabilities only. We say that upper probability
$\overline{P}$ is concave if $\overline{P}(A\cup
B)\leq\overline{P}(A)+\overline{P}(B)-\overline{P}(A\cap B)$, for all
$A,B\in\mathcal{F}$. Recall that $\Delta(\Omega,\mathcal{F})$ denotes the set
of all probability measures on $(\Omega,\mathcal{F})$. Upper probability
$\overline{P}$ is compatible with the convex set [30]

$\displaystyle\text{core}(\overline{P}):$
$\displaystyle=\\{P\in\Delta(\Omega,\mathcal{F}):P(A)\leq\overline{P}(A)\text{,
}\forall A\in\mathcal{F}\\}$
$\displaystyle=\\{P\in\Delta(\Omega,\mathcal{F}):\underline{P}(A)\leq
P(A)\leq\overline{P}(A)\text{, }\forall A\in\mathcal{F}\\}$

where the second equality is a characterization [9, Page 3389]. Notice that
the core is convex [53, Section 2.2]. We assume it is nonempty and
weak⋆-closed.999Recall that in the weak⋆ topology, a net
$(P_{\alpha})_{\alpha\in I}$ converges to $P$ if and only if
$P_{\alpha}(A)\rightarrow P(A)$, for all $A\in\mathcal{F}$. Then, it is
weak⋆-compact as a result of [53, Proposition 3].

Since the core is convex, the set $\text{ex}[\text{core}(\overline{P})]$ of
extreme points of the core is well defined. It contains all the elements of
the core that cannot be written as a convex combination of one another. The
following important result is a consequence of [71, Theorem 3.6.2].

######  Theorem 8.

Suppose $\text{core}(\overline{P})$ is nonempty and weak⋆-closed. Then, the
following holds.

  * (a)

$\text{ex}[\text{core}(\overline{P})]\neq\emptyset$.

  * (b)

$\text{core}(\overline{P})$ is the closure in the weak⋆ topology of the convex
hull of $\text{ex}[\text{core}(\overline{P})]$.

  * (c)

If $\overline{P}(A)=\sup_{P\in\text{core}(\overline{P})}P(A)$, for all
$A\in\mathcal{F}$, then
$\overline{P}(A)=\sup_{P\in\text{ex}[\text{core}(\overline{P})]}P(A)$, for all
$A\in\mathcal{F}$.

So in order to define an upper probability $\overline{P}$ that setwise
dominates the elements of $\text{core}(\overline{P})$ it is enough to specify
the extreme points of the core.

##  Appendix D A new Bayes’ theorem for IPs

We present Theorem 9, a result that – although appealing – does not lend
itself well to be used to train an IBNN.

Call $\Theta$ the parameter space of interest and assume it is Polish, that
is, the topology for $\Theta$ is complete, separable, and metrizable. This
ensures that the set $\Delta(\Theta,\mathcal{B})$ of probability measures on
$\Theta$ is Polish as well, where $\mathcal{B}$ denotes the Borel
$\sigma$-algebra for $\Theta$. Let $\mathscr{X}$ be the set of all bounded,
non-negative, $\mathcal{B}$-measurable functionals on $\Theta$. Call
$\mathscr{D}=\mathcal{X}\times\mathcal{Y}$ the sample space endowed with the
product $\sigma$-algebra
$\mathcal{A}=\mathcal{A}_{\mathbf{x}}\times\mathcal{A}_{\mathbf{y}}$, where
$\mathcal{A}_{\mathbf{x}}$ is the $\sigma$-algebra endowed to $\mathcal{X}$
and $\mathcal{A}_{\mathbf{y}}$ is the $\sigma$-algebra endowed to
$\mathcal{Y}$. Let the agent elicit
$\mathcal{L}_{\theta}:=\\{P_{\theta}\in\Delta(\mathscr{D},\mathcal{A}):\theta\in\Theta\\}$.
Assume that each $P_{\theta}\in\mathcal{L}_{\theta}$ has density
$L(\theta)=p(D\mid\theta)$ with respect to some $\sigma$-finite dominating
measure $\nu$ on $(\mathscr{D},\mathcal{A})$; this represents the likelihood
function for $\theta$ having observed data $D\subset\mathscr{D}$. We assume
for now that $L\in\mathscr{X}$, for all $D\subset\mathscr{D}$.

Let the agent specify a set $\mathcal{P}$ of probabilities on
$(\Theta,\mathcal{B})$. Then, compute $\overline{P}$, and consider
$\mathcal{P}^{\text{co}}:=\text{core}(\overline{P})$; it represents the
agent’s initial beliefs.101010Superscript “co” stands for convex and core. We
assume that every $P\in\mathcal{P}^{\text{co}}$ has density $p$ with respect
to some $\sigma$-finite dominating measure $\mu$ on $(\Theta,\mathcal{B})$,
that is, $p=\frac{\text{d}P}{\text{d}\mu}$. We require the agent’s beliefs to
be represented by the core for two main reasons. The first, mathematical, one
is to ensure that the belief set is compatible with the upper probability. The
second, philosophical, one is the following [8]. A criticism brought forward
by [71, Section 2.10.4.(c)] is that, given an upper probability
$\overline{P}$, there is no cogent reason for which the agent should choose a
specific $P_{T}$ that is dominated by $\overline{P}$, or – for that matter – a
collection of “plausible” probabilities. Because the core considers all
(countably additive) probability measures that are dominated by
$\overline{P}$, it is the perfect instrument to reconcile Walley’s behavioral
and sensitivity analysis interpretations.

Let the agent compute $\overline{P}_{\theta}$, and consider
$\mathcal{L}^{\text{co}}_{\theta}:=\text{core}(\overline{P}_{\theta})$; it
represents the set of plausible likelihoods. Let

$\mathscr{L}:=\left\\{{L=\frac{\text{d}P_{\theta}}{\text{d}\nu}\text{,
}P_{\theta}\in\mathcal{L}_{\theta}^{\text{co}}}\right\\},$ (1)

and denote by $\overline{L}(\theta):=\sup_{L\in\mathscr{L}}L(\theta)$ and by
$\underline{L}(\theta):=\inf_{L\in\mathscr{L}}L(\theta)$, for all
$\theta\in\Theta$. Call

$\displaystyle\mathcal{P}^{\text{co}}_{D}:=\Bigg{\\{}P_{D}\in\Delta(\Theta,\mathcal{B}):\frac{\text{d}P_{D}}{\text{d}\mu}=\text{
}$ $\displaystyle p(\theta\mid
D)=\frac{L(\theta)p(\theta)}{\int_{\Theta}L(\theta)p(\theta)\text{d}\theta}\text{,
}$ $\displaystyle p=\frac{\text{d}P}{\text{d}\mu}\text{,
}P\in\mathcal{P}^{\text{co}}\text{,
}L=\frac{\text{d}P_{\theta}}{\text{d}\nu}\text{,
}P_{\theta}\in\mathcal{L}_{\theta}^{\text{co}}\Bigg{\\}}$

the class of posterior probabilities when the prior is in
$\mathcal{P}^{\text{co}}$ and the likelihood is in
$\mathcal{L}_{\theta}^{\text{co}}$, and let
$\overline{P}_{D}(A)=\sup_{P_{D}\in\mathcal{P}^{\text{co}}_{D}}P_{D}(A)$, for
all $A\in\mathcal{B}$. Then, the following is a generalization of Bayes’
theorem in [73].

######  Theorem 9.

Suppose $\mathcal{P}^{\text{co}},\mathcal{L}_{\theta}^{\text{co}}$ are
nonempty and weak⋆-closed. Then for all $A\in\mathcal{B}$,

$\overline{P}_{D}(A)\leq\frac{\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}{\mathbf{c}},$
(2)

where
$\mathbf{c}:=\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)+\inf_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\underline{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)$,
provided that the ratio is well defined, where $\mathbbm{1}_{A}$ denotes the
indicator function for $A\in\mathcal{B}$. In addition, if $\overline{P}$ is
concave, then the inequality in (2) is an equality for all $A\in\mathcal{B}$.

This result is particularly appealing because, given some assumptions, it
allows to perform a (generalized) Bayesian update of a prior upper probability
(PUP) by carrying out only one operation, even when the likelihood is ill
specified so that a set of likelihoods is needed. We also have the following.

######  Lemma 10.

Suppose $\mathcal{P}^{\text{co}},\mathcal{L}_{\theta}^{\text{co}}$ are
nonempty and weak⋆-closed. Then, if $\overline{P}$ is concave, we have that
$\overline{P}_{D}$ is concave as well.

This lemma is important because it tells us that the generalized Bayesian
update of Theorem 9 preserves concavity, and so it can be applied to
successive iterations. If at time $t$ the PUP is concave, then the PUP at time
$t+1$ – that is, the posterior upper probability at time $t$ – will be concave
too. Necessary and sufficient conditions for a generic upper probability to be
concave are given in [53, Section 5].

In the future, these results can be generalized to the case in which the
elements of $\mathscr{X}$ are unbounded using techniques in [68], and to the
case in which the elements of $\mathscr{X}$ are $\mathbb{R}^{d}$-valued, for
some $d\in\mathbb{N}$, since we never used specific properties of $\mathbb{R}$
in our proofs.

Despite being attractive, the generalized Bayesian update of Theorem 9 hinges
upon three assumptions, namely that $\mathcal{P}^{\text{co}}$ and
$\mathcal{L}^{\text{co}}_{\theta}$ are both cores of an upper probability,
that they are nonempty and weak⋆-closed, and that the prior upper probability
$\overline{P}$ is concave. As the proverb goes, there is no free lunch. Having
to check these assumptions, together with computing a supremum, an infimum,
and the integrals in (2), makes Theorem 9 inadequate to be applied in the
context of IBNNs.

##  Appendix E Bounds on upper and lower entropy

In this section we find an upper bound for $\overline{H}(P)$ and a lower bound
for $\underline{H}(P)$ that are extremely interesting. We first need to
introduce three new concepts.

######  Definition 11.

Consider a set $\mathscr{P}$ of probabilities on a generic measurable space
$(\Omega,\mathcal{F})$. We say that lower probability $\underline{P}$ is
convex if $\underline{P}(A\cup
B)\geq\underline{P}(A)+\underline{P}(B)-\underline{P}(A\cap B)$, for all
$A,B\in\mathcal{F}$.

Then, let $\mathfrak{P}$ be either an upper or a lower probability, and
consider a generic bounded function $f$ on $(\Omega,\mathcal{F})$, that is,
$f\in B(\Omega)$. We define the Choquet integral of $f$ with respect to
$\mathfrak{P}$ as follows

$\displaystyle\int_{\Omega}f(\omega)\mathfrak{P}(\text{d}\omega):=\int_{0}^{\infty}\mathfrak{P}\left(\\{\omega\in\Omega:f(\omega)\geq
t\\}\right)\text{d}t+\int_{-\infty}^{0}\left[\mathfrak{P}\left(\\{\omega\in\Omega:f(\omega)\geq
t\\}\right)-\mathfrak{P}(\Omega)\right]\text{d}t,$

where the right hand side integrals are (improper) Riemann integrals. If
$\mathfrak{P}$ is additive, then the Choquet integral reduces to the standard
additive integral.

Finally, if $\Omega$ is uncountable, define for all $\omega\in\Omega$

$\underline{\pi}(\omega):=\inf_{P\in\mathscr{P}}\frac{\text{d}P}{\text{d}\mu}(\omega)\quad\text{and}\quad\overline{\pi}(\omega):=\sup_{P\in\mathscr{P}}\frac{\text{d}P}{\text{d}\mu}(\omega).$

We call them lower and upper densities, respectively.

The following theorem gives the desired bounds.

######  Theorem 12.

Consider a set $\mathscr{P}$ of probabilities on a generic measurable space
$(\Omega,\mathcal{F})$. If $\Omega$ is uncountable, assume that every
$P\in\mathscr{P}$ is dominated by a $\sigma$-finite measure $\mu$ and that the
Radon-Nikodym derivatives $\frac{\text{d}P}{\text{d}\mu}$ are continuous and
bounded, for all $P\in\mathscr{P}$. Define

$H(\underline{P}):=\begin{cases}-\int_{\Omega}\log\left[\overline{\pi}(\omega)\right]\underline{P}(\text{d}\omega)&\text{
if }\Omega\text{ is uncountable}\\\
-\sum_{\omega\in\Omega}\underline{P}(\\{\omega\\})\log[\overline{P}(\\{\omega\\})]&\text{
if }\Omega\text{ is at most countable}\end{cases}$

and similarly

$H(\overline{P}):=\begin{cases}-\int_{\Omega}\log\left[\underline{\pi}(\omega)\right]\overline{P}(\text{d}\omega)&\text{
if }\Omega\text{ is uncountable}\\\
-\sum_{\omega\in\Omega}\overline{P}(\\{\omega\\})\log[\underline{P}(\\{\omega\\})]&\text{
if }\Omega\text{ is at most countable}\end{cases},$

Then, $\overline{H}(P)\leq H(\overline{P})$ and $\underline{H}(P)\geq
H(\underline{P})$. In addition, if $\overline{P}$ is concave, the first bound
is tighter, and if $\underline{P}$ is convex, the second bound is tighter.

######  Remark 13.

In [35, Section 4.6.1], the authors point out that [2] presents a
generalization of the Hartley measure, called generalized Hartley measure
$GH(P)$, that can be used to disaggregate the total uncertainty captured by
$\overline{H}(P)$ into aleatoric and epistemic uncertainties. We prefer not to
introduce it in the present work because $GH(P)$ is defined based on the mass
function of a belief function [30, Definition 2.4].111111A belief function is
a mathematical concept that should not be confused with the term “belief” we
used throughout the paper to address the agent’s knowledge. This entails that
the authors assume that lower probability $\underline{P}$ associated with the
set of probabilities of interest is a belief function, and so that for every
collection $\\{A,A_{1},\ldots,A_{k}\\}$ such that $A\subset A_{i}$, the
following holds

$\underline{P}(A)\geq\sum_{\emptyset\neq
I\subset\\{1,\ldots,k\\}}(-1)^{\\#I-1}\underline{P}(\cap_{i\in I}A_{i}),$

for all $k\in\mathbb{N}$. As it is immediate to see, this is assumption is not
needed in the context of IBNNs.

##  Appendix F How to derive a predictive distribution

Suppose we performed a Bayesian updating procedure so to obtain posterior pdf
$p(\theta\mid x_{1},y_{1},\ldots,x_{n},y_{n})$. Recall that
$\\{(x_{i},y_{i})\\}_{i=1}^{n}\in(\mathcal{X}\times\mathcal{Y})^{n}$ denotes
the training set. We obtain the predictive distribution
$p(\tilde{y}\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})$ on $\mathcal{Y}$ as
follows

$\displaystyle p(\tilde{y}\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})$
$\displaystyle=\int_{\Theta}p(\tilde{y},\theta\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})\text{
d}\theta$
$\displaystyle=\int_{\Theta}p(\tilde{y}\mid\theta,\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})\cdot
p(\theta\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})\text{ d}\theta$
$\displaystyle=\int_{\Theta}p(\tilde{y}\mid\theta,\tilde{x})\cdot p(\theta\mid
x_{1},y_{1},\ldots,x_{n},y_{n})\text{ d}\theta,$

where $p(\tilde{y}\mid\theta,\tilde{x})$ is the likelihood used to derive the
posterior. Notice that the last equality comes from output $\tilde{y}$ only
depending on input $\tilde{x}$ and parameter $\theta$, and from having assumed
$D_{\mathbf{x}}\perp\\!\\!\\!\perp\theta$ (see section 2.1). From an applied
point of view, a sample from
$p(\tilde{y}\mid\tilde{x},x_{1},y_{1},\ldots,x_{n},y_{n})$ is obtained as
follows:

  1. (1)

specify input $\tilde{x}$;

  2. (2)

sample a parameter $\tilde{\theta}$ from the posterior, $\tilde{\theta}\sim
p(\theta\mid x_{1},y_{1},\ldots,x_{n},y_{n})$;

  3. (3)

plug $\tilde{\theta}$ in the likelihood and sample $\tilde{y}\sim
p({y}\mid\tilde{\theta},\tilde{x})$.

##  Appendix G Aleatoric uncertainty check for $\alpha$-level IHDR

The diameter of $IR_{\alpha}$ is a function of the aleatoric uncertainty (AU)
faced by the agent.121212Since the diameter is a metric concept, we assume
that we can find a well-defined metric $d_{\mathbf{y}}$ on $\mathcal{Y}$. If
that is not the case, we substitute the diameter with the notion of
cardinality. If we want to avoid to perform the procedure in section 4 only to
discover that $IR_{\alpha}$ is “too big”, then we can add an “AU check”.

At the beginning of the analysis, compute the lower entropy
$\underline{H}(\hat{P})=\inf_{k\in\\{1,\ldots,N\\}}H(\hat{P}_{k})$ associated
with the set $\hat{\mathcal{P}}=\\{\hat{P}_{1},\ldots,\hat{P}_{N}\\}$ of
extreme elements of the predictive credal set, which gives us the aleatoric
uncertainty within $\hat{\mathcal{P}}$.131313Also, it gives us a bound for the
AU in the convex hull $\text{Conv}\hat{\mathcal{P}}$. We then verify whether
the lower entropy $\underline{H}(\hat{P})$ is “too high”. That is, if
$\underline{H}(\hat{P})>\varphi$, for some $\varphi>0$, we want our procedure
to abstain. This means that if the aleatoric uncertainty in set
$\hat{\mathcal{P}}$ is too high, then our procedure does not return any output
set for input $\tilde{x}$. The value of $\varphi$ can be set equal to the
entropy of the probability measures that are typically used in the context the
agent works in. For example, in medical applications the agent may consider
the entropy of a Normal distribution, while in financial applications the
entropy of a distribution with fatter tails, such as a $t$-distribution or a
Cauchy. We call these reference $\varphi$ values. If we are interested in
verifying whether (a lower bound for) the AU in $\text{Conv}\hat{\mathcal{P}}$
is “too high”, we check whether $\underline{H}(\hat{P})>\varphi+\log(N)$,
where $N=\\#\hat{\mathcal{P}}$. The reason for the $\log(N)$ term is presented
in Appendix I.

If we add this “AU check”, the $2$-tuple that the agent needs to specify to
output an IHDR is $\langle\varphi,\alpha\rangle$.

##  Appendix H $\alpha$-level IHDR in a Classification setting

In classification problems, BNNs compute the probability vector

$\hat{p}:=\frac{1}{\\#\Theta}\sum_{\theta\in\Theta}\Phi_{\theta\mid D}(x),$

where we write $\Phi_{\theta\mid D}$ to highlight the fact that $\theta$ is
sampled from posterior $p(\theta\mid D)$, and then select the most likely
class $\hat{y}:=\operatorname*{arg\,max}_{j}\hat{p}_{j}$, where the
$\hat{p}_{j}$’s are the elements of $\hat{p}$.

When applied to a classification setting, the general procedure introduced in
section 4 becomes the following. Recall that we denote by
$N:=\\#\mathcal{P}\times\\#\mathcal{L}_{x,\theta}$ the number of predictive
distributions on $\mathcal{Y}$ induced by the extreme elements of our
posterior credal set. Assume that $\mathcal{Y}=\\{y_{1},\ldots,y_{J}\\}$, that
is, there are $J\in\mathbb{N}$ possible labels. Then, the predictive
distributions on $\mathcal{Y}$ can be seen as $N$ $J$-dimensional probability
vectors $\hat{p}^{1},\ldots,\hat{p}^{N}$, where
$\hat{p}^{k}=(\hat{p}^{k}_{1}=\hat{p}^{k}(y_{1}),\ldots,\hat{p}^{k}_{J}=\hat{p}^{k}(y_{J}))^{\top}$
for all $k$, $\hat{p}^{k}_{j}\in[0,1]$ for all $k$ and all $j$, and
$\sum_{j=1}^{J}\hat{p}^{k}_{j}=1$ for all $k$.

Now fix any $k\in\\{1,\ldots,N\\}$, and define the partial order $\preceq_{k}$
on $\mathcal{Y}$ as
$y_{l}\preceq_{k}y_{i}\iff\hat{p}^{k}(y_{l})\geq\hat{p}^{k}(y_{i})$ and
$y_{l}\prec_{k}y_{i}\iff\hat{p}^{k}(y_{l})>\hat{p}^{k}(y_{i})$, where
$i,l\in\\{1,\ldots,J\\}$, $i\neq l$. This means that we can order the labels
according to the probability that $\hat{p}^{k}$ assigns to them: the first
label will be the one having highest probability according to $\hat{p}^{k}$,
the second label will have the second-highest probability according to
$\hat{p}^{k}$, and so on.

Now order the label space $\mathcal{Y}$ according to $\preceq_{k}$ so to
obtain

$\mathcal{Y}^{\preceq_{k}}:=\\{y_{1}^{\preceq_{k}},\ldots,y_{J}^{\preceq_{k}}\\}.$

This means that $y_{1}^{\preceq_{k}}\preceq_{k}y_{j}^{\preceq_{k}}$, for all
$j\in\\{2,\ldots,J\\}$, $y_{2}^{\preceq_{k}}\preceq_{k}y_{j}^{\preceq_{k}}$,
for all $j\in\\{3,\ldots,J\\}$ (but
$y_{1}^{\preceq_{k}}\preceq_{k}y_{2}^{\preceq_{k}}$), and so on. That is, we
order the labels from the most to the least likely according to $\hat{p}^{k}$.

Then, we call $\alpha$-level credible set according to $\hat{p}^{k}$,
$\alpha\in(0,1)$, the set

$\displaystyle\begin{split}CS^{k}_{\alpha}:=\bigg{\\{}y_{1}^{\preceq_{k}},\ldots,y_{j}^{\preceq_{k}}:&\sum_{i=1}^{j}\hat{p}^{k}(y_{i}^{\preceq_{k}})\in[1-\alpha,1-\alpha+\varepsilon]\text{,
}j\leq J,\\\ \text{and }&\nexists
j^{\prime}<j:\sum_{i=1}^{j^{\prime}}\hat{p}^{k}(y_{i}^{\preceq_{k}})\in[1-\alpha,1-\alpha+\varepsilon]\bigg{\\}},\end{split}$
(3)

for some $\varepsilon>0$. It corresponds to the $\alpha$-level HDR. Notice
that we require
$\sum_{i=1}^{j}\hat{p}^{k}(y_{i}^{\preceq_{k}})\in[1-\alpha,1-\alpha+\varepsilon]$
because we may need to go slightly above level $1-\alpha$. Just as a toy
example, we may have $7$ labels, $3$ of which would give a $0.945$ coverage,
while $4$ would give a coverage of $0.953$. If we are interested in the
$\alpha=0.05$-level credible set, we ought to include the fourth label, thus
yielding a coverage slightly higher than $1-\alpha=0.95$. The interpretation
to $CS^{k}_{\alpha}$ is the following: it consists of the smallest collection
of labels to which $\hat{p}^{k}$ assigns probability of at least $1-\alpha$
(that is, those having the highest probability of being the correct one given
the input).

Finally, we call $\alpha$-level imprecise credible set, $\alpha\in(0,1)$, the
set

$ICS_{\alpha}:=\bigcup_{k=1}^{N}CS^{k}_{\alpha}.$

In turn, we have that $\underline{\hat{P}}(y\in ICS_{\alpha})\geq 1-\alpha$.

######  Remark 14.

Notice that if a credible set of level $\approx\alpha$ is enough, then we can
replace the left endpoint of the interval in (3) with
$1-(\alpha+\varepsilon_{k})$, for some $\varepsilon_{k}>0$. Strictly speaking,
in this case we obtain an $(\alpha+\varepsilon_{k})$-level credible set, which
we denote by $CS^{k}_{\alpha_{k}}$, where
$\alpha_{k}:=\alpha+\varepsilon_{k}$, $k\in\\{1,\ldots,N\\}$. Going back to
our toy example, we will have a credible set with $3$ labels that yields a
coverage of $1-(\alpha+\varepsilon_{k})=0.945\approx 0.95=1-\alpha$, so
$\varepsilon_{k}=0.005$. In turn, this implies that the imprecise credible set
will have a coverage of
$1-(\alpha+\max_{k\in\\{1,\ldots,N\\}}\varepsilon_{k})$, that is, it will have
level $\alpha+\max_{k\in\\{1,\ldots,N\\}}\varepsilon_{k}$. We denote it by
$ICS_{\tilde{\alpha}}$, where
$\tilde{\alpha}:=\alpha+\max_{k\in\\{1,\ldots,N\\}}\varepsilon_{k}$.

##  Appendix I Experiments: Uncertainty Quantification

Figure 10. The different trends with increasing corruption severity for EBNN
(baseline) versus IBNNs. IBNNs better inform the degree of shift as compared
to EBNN.

Distribution shifts can introduce uncertainties in the system which can render
the predictions meaningless. This can be due to naturally occurring
corruptions as introduced in [32] for image datasets. The authors introduced
$18$ different noise types, which can be varied across $5$ different severity
levels. The intuition is that in the current context, increasing the noise
severity should generally result in higher uncertainty. We evaluate an IBNN on
4 standard image datasets CIFAR10 [44], SVHN [56], Fashion-MNIST [74], and
MNIST [49]. We use a slightly different set of perturbations introduced in
[55] for gray-scale images like MNIST and Fashion-MNIST. Additionally, we
perform cross domain testing for each dataset, where we expect the
uncertainties to be higher. We implement and train a Resnet-20 Bayesian DNN
model inside the library Bayesian-torch [43]. For each dataset, we train $4$
different networks initialized with different seeds on the prior and with the
same architecture. We use a learning-rate of $0.001$, batch-size of $128$, and
train the networks using mean-field variational inference for $200$ epochs.
The inference is carried out by performing multiple forward passes through
parameters drawn from the posterior distribution. We used $20$ Monte-Carlo
samples in the experiments.

Results. Following [22], for EBNN we posit that
$1/k\sum_{j=1}^{k}\sigma^{2}_{j}$ captures the aleatoric uncertainty, and
$1/(k-1)\sum_{j=1}^{k}(\mu_{j}-\mu_{\text{ens}})^{2}$ captures the epistemic
uncertainty; we use these values as baselines.141414Notice that in this case
we retain the assumption that the EBNN distribution on $\mathcal{Y}$ has mean
$\mu_{\text{ens}}$ and covariance matrix $\sigma^{2}_{\text{ens}}I$, but we do
not require that it is a Normal. That is because in this case $\mathcal{Y}$ is
finite. Notice also that for IBNN the AU is given by the lower entropy
$\underline{H}(\hat{P}):=\inf_{\hat{P}\in\hat{\mathcal{P}}}H(\hat{P})$, where
$\hat{\mathcal{P}}$ was defined in section 4, while the EU is given by the
difference $\overline{H}(\hat{P})-\underline{H}(\hat{P})$, and
$\overline{H}(\hat{P}):=\sup_{\hat{P}\in\hat{\mathcal{P}}}H(\hat{P})$. For
IBNN, we look at lower bounds for both EU and AU. We do so because (i) they
are easy to compute and (ii) we expect both EU and AU to increase as the
corruption grows, so increasing lower bounds are a good indicator of the
desired trend. Denote by $\text{ex}\hat{\mathcal{P}}$ the collection of
extreme elements of $\hat{\mathcal{P}}$. Let then
$\sup_{\hat{P}^{\text{ex}}\in\text{ex}\hat{\mathcal{P}}}H(\hat{P}^{\text{ex}})=:\overline{H}(\hat{P}^{\text{ex}})$
and $\sup_{\hat{P}\in\hat{\mathcal{P}}}H(\hat{P})=:\overline{H}(\hat{P})$.
Since by Proposition 7
$\overline{H}(\hat{P}^{\text{ex}})\leq\overline{H}(\hat{P})$ and
$\underline{H}(\hat{P}^{\text{ex}})\geq\underline{H}(\hat{P})$, we have that
$\overline{H}(\hat{P}^{\text{ex}})-\underline{H}(\hat{P}^{\text{ex}})$ is a
lower bound for the predictive EU. Also, by Proposition 7 and [63, Theorem
III.1], we have that the AU is lower bonded by a quantity that can be
approximated by $\underline{H}(\hat{P}^{\text{ex}})-\log(N)$, where
$N=\\#\hat{\mathcal{P}}$. Then, since $\log(N)$ is a fixed quantity, the trend
of $\underline{H}(\hat{P}^{\text{ex}})-\log(N)$ as we increase the corruption
severity is equivalent to the trend of $\underline{H}(\hat{P}^{\text{ex}})$.
We report the latter for ease of computation. We discuss the results for
CIFAR10 presented in Table 3. We report the average uncertainties across all
test samples as measured by the respective procedures. Note that we abstracted
the severity levels to low (severity = $1$), medium (severity = $2,3$), and
high (severity = $4,5$) in this paper. For CIFAR10, we observe that IBNNs
report lower levels of aleatoric uncertainty when subjected to low levels of
corruption, which gradually increases as we move to higher levels. The
epistemic uncertainty grows as well, but the change in the aleatoric component
is more pronounced. However, for EBNN this is not true. Even though the
epistemic uncertainty increases, the aleatoric shows the reverse trend,
contrary to our intuition of aleatoric uncertainty. Using IBNNs, we achieve
the highest uncertainties, both aleatoric and epistemic, when applied to
entirely different test sets, namely MNIST, SVHN, and Fahion-MNIST. We report
these numbers as well.

Summary. We summarize our results in Figure 10. With increasing corruption
severity, the aleatoric uncertainty for IBNNs becomes more pronounced. This is
expected since for higher corruptions, it is hard to reconcile the performance
gap by simple data-augmentation, without any additional knowledge. For EBNN,
even though the epistemic uncertainty grows, the aleatoric part does not. Note
that the absolute uncertainty differs due to the quantities being
fundamentally different: variance versus entropy. However, the relative trends
are consistent for each dataset, demonstrating the utility of IBNNs.

Other datasets. The results for MNIST, Fashion-MNIST, and SVHN are presented
in Tables 4, 5, and 6, respectively. We observe trends similar to CIFAR10’s
ones: the aleatoric uncertainty increases with corruption severity, and is the
highest for completely unknown datasets.

##  Appendix J Distance between a set of distributions $\mathcal{P}$ and a
single distribution $P^{\prime}$ having different dimensions

Many concepts in this section are derived from [7]. Let $m,n\in\mathbb{N}$
such that $m\leq n$, and let $p\in[1,\infty]$. Call $M^{p}(\mathbb{R}^{j})$
the set of probability measures on $\mathbb{R}^{j}$ having finite $p$-th
moment, and $M_{d}(\mathbb{R}^{j})$ the set of probability measures on
$\mathbb{R}^{j}$ having density with respect to some $\sigma$-finite
dominating measure $\mu$, $j\in\\{m,n\\}$. Let
$O(m,n):=\\{V\in\mathbb{R}^{m\times n}:VV^{\top}=I_{m}\\}$, where $I_{m}$
denotes the $m$-dimensional identity matrix, and for any $V\in O(m,n)$ and any
$b\in\mathbb{R}^{m}$, define the following function

$\varphi_{V,b}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{m},\quad
x\mapsto\varphi_{V,b}(x):=Vx+b.$

Let $\mathcal{B}(\mathbb{R}^{n})$ be the Borel $\sigma$-algebra on
$\mathbb{R}^{n}$, and for any
$Q\in\Delta(\mathbb{R}^{n},\mathcal{B}(\mathbb{R}^{n}))$, define
$\varphi_{V,b}(Q):=Q\circ\varphi_{V,b}^{-1}$, the pushforward measure.
Consider then two generic probability measures $Q,S$ such that $Q\in
M^{p}(\mathbb{R}^{m})$ and $S\in M^{p}(\mathbb{R}^{n})$, and call

$\displaystyle\Phi_{p}^{+}(Q,n)$ $\displaystyle:=\\{\alpha\in
M^{p}(\mathbb{R}^{n}):\varphi_{V,b}(\alpha)=Q\text{, for some }V\in
O(m,n),b\in\mathbb{R}^{m}\\},$ $\displaystyle\Phi_{d}^{+}(Q,n)$
$\displaystyle:=\\{\alpha\in
M_{d}(\mathbb{R}^{n}):\varphi_{V,b}(\alpha)=Q\text{, for some }V\in
O(m,n),b\in\mathbb{R}^{m}\\},$ $\displaystyle\Phi^{-}(S,m)$
$\displaystyle:=\\{\beta\in M(\mathbb{R}^{m}):\varphi_{V,b}(S)=\beta\text{,
for some }V\in O(m,n),b\in\mathbb{R}^{m}\\}.$

Recall now the definition of $p$-Wasserstein metric between two generic
distributions defined on the same Euclidean space. Let $P_{1},P_{2}\in
M^{p}(\mathbb{R}^{n})$, for some $n\in\mathbb{N}$ and some $p\in[1,\infty]$.
Then, the $p$-Wasserstein distance between them is defined as

$W_{p}(P_{1},P_{2}):=\left[\inf_{\gamma\in\Gamma(P_{1},P_{2})}\int_{\mathbb{R}^{2n}}\|x-y\|_{2}^{p}\gamma(\text{d}(x,y))\right]^{1/p},$

where $\|\cdot\|_{2}$ denotes the Euclidean distance, $p=\infty$ is
interpreted as the essential supremum, and
$\Gamma(P_{1},P_{2}):=\\{\gamma\in\Delta(\mathbb{R}^{2n},\mathcal{B}(\mathbb{R}^{2n})):\text{proj}_{1}^{n}(\gamma)=P_{2}\text{,
proj}_{2}^{n}(\gamma)=P_{1}\\}$ is the set of couplings between $P_{1}$ and
$P_{2}$, where $\text{proj}_{1}^{n}$ is the projection onto the first $n$
coordinates, and $\text{proj}_{2}^{n}$ is the projection to the last $n$
coordinates.

Recall then the definition of $f$-divergence between two generic distributions
defined on the same Euclidean space. Let
$P_{1},P_{2}\in\Delta(\mathbb{R}^{n},\mathcal{B}(\mathbb{R}^{n}))$, for some
$n\in\mathbb{N}$, and assume $P_{1}\ll P_{2}$. Then, for any convex functional
$f$ on $\mathbb{R}$ such that $f(1)=0$, the $f$-divergence between $P_{1}$ and
$P_{2}$ is defined as

$div_{f}(P_{1}\|P_{2}):=\int_{\mathbb{R}^{n}}f\left(\frac{\text{d}P_{1}}{\text{d}P_{2}}(x)\right)P_{2}(\text{d}x).$

Aside from the Rényi divergence, the $f$-divergence includes just about every
known divergences as special case [7]. The following are the main results of
this section.

######  Lemma 15.

Let $m,n\in\mathbb{N}$ such that $m\leq n$, and let $p\in[1,\infty]$ and $f$
be any convex functional on $\mathbb{R}$ such that $f(0)=1$. Consider a
generic $\mathscr{P}\subset M^{p}(\mathbb{R}^{m})$ and $P^{\prime}\in
M^{p}(\mathbb{R}^{n})$. Let
$\Phi_{p}^{+}(\mathscr{P},n):=\cup_{P\in\mathscr{P}}\Phi_{p}^{+}(P,n)$ and
$\Phi_{d}^{+}(\mathscr{P},n):=\cup_{P\in\mathscr{P}}\Phi_{d}^{+}(P,n)$. Define

  * •

$W^{+}_{p}({P},P^{\prime}):=\inf_{\alpha\in\Phi_{p}^{+}({P},n)}W_{p}(\alpha,P^{\prime})$,
for all $P\in\mathscr{P}$;

  * •

$div^{+}_{f}({P}\|P^{\prime}):=\inf_{\alpha\in\Phi_{d}^{+}({P},n)}div_{f}({P}\|P^{\prime})$,
for all $P\in\mathscr{P}$;

  * •

$W^{+}_{p}(\mathscr{P},P^{\prime}):=\inf_{\alpha\in\Phi_{p}^{+}(\mathscr{P},n)}W_{p}(\alpha,P^{\prime})$;

  * •

$div^{+}_{f}(\mathscr{P}\|P^{\prime}):=\inf_{\alpha\in\Phi_{d}^{+}(\mathscr{P},n)}div_{f}({P}\|P^{\prime})$;

Then, for all $P\in\mathscr{P}$ the following holds

$W^{+}_{p}(\mathscr{P},P^{\prime})\leq
W^{+}_{p}({P},P^{\prime})\quad\text{and}\quad
div^{+}_{f}(\mathscr{P}\|P^{\prime})\leq div^{+}_{f}({P}\|P^{\prime}).$

######  Lemma 16.

Let $m,n\in\mathbb{N}$ such that $m\leq n$, and let $p\in[1,\infty]$ and $f$
be any convex functional on $\mathbb{R}$ such that $f(0)=1$. Consider a
generic $\mathscr{P}\subset M^{p}(\mathbb{R}^{n})$ and $P^{\prime}\in
M^{p}(\mathbb{R}^{m})$. Let
$\Phi^{-}(\mathscr{P},m):=\cup_{P\in\mathscr{P}}\Phi^{-}(P,m)$. Define

  * •

$W^{-}_{p}({P},P^{\prime}):=\inf_{\alpha\in\Phi^{-}({P},m)}W_{p}(\alpha,P^{\prime})$,
for all $P\in\mathscr{P}$;

  * •

$div^{-}_{f}({P}\|P^{\prime}):=\inf_{\alpha\in\Phi^{-}({P},m)}div_{f}({P}\|P^{\prime})$,
for all $P\in\mathscr{P}$;

  * •

$W^{-}_{p}(\mathscr{P},P^{\prime}):=\inf_{\alpha\in\Phi^{-}(\mathscr{P},m)}W_{p}(\alpha,P^{\prime})$;

  * •

$div^{-}_{f}(\mathscr{P}\|P^{\prime}):=\inf_{\alpha\in\Phi^{-}(\mathscr{P},m)}div_{f}({P}\|P^{\prime})$;

Then, for all $P\in\mathscr{P}$ the following holds

$W^{-}_{p}(\mathscr{P},P^{\prime})\leq
W^{-}_{p}({P},P^{\prime})\quad\text{and}\quad
div^{-}_{f}(\mathscr{P}\|P^{\prime})\leq div^{-}_{f}({P}\|P^{\prime}).$

A visual representation of the application of Lemma 16 to IBNNs is given in
Figure 11.151515Notice that in Figure 11 we have that
$\mathscr{L}_{x,\theta}\equiv\mathcal{L}_{x,\theta}$.

Figure 11. We assume that $n>m$ and that the oracle distribution
$P^{o}_{x,\theta}$ belongs to $M^{p}(\mathbb{R}^{m})$, while (finitely
generated) credal set $\text{Conv}\mathcal{L}_{x,\theta}$ is a subset of
$M^{p}(\mathbb{R}^{m})$, for some finite $p\geq 1$. We also assume that
$\mathbf{P}_{x,\theta}$ is one of the extreme elements of
$\text{Conv}\mathcal{L}_{x,\theta}$. We see how
$W^{-}_{p}(\text{Conv}\mathcal{L}_{x,\theta},P^{o}_{x,\theta})<W^{-}_{p}(\mathbf{P}_{x,\theta},P^{o}_{x,\theta})$;
if we replace metric $W^{-}_{p}$ by a generic $f$-divergence $div_{f}$, the
inequality would still hold thanks to Lemma 16.

##  Appendix K Details on Artificial Pancreas Example

Artificial Pancreas Model. An important factor when designing the controller
for an artificial pancreas is to adapt the insulin delivery algorithm to the
particular details of the patient. This is because patients display a wide
range of variability in their response to insulin, depending on age, BMI, and
other physiological parameters. The Bayesian neural network models have $2$
hidden layers, with $10$ neurons each, for the case study with $4$ different
seeds. This choice was informed by the experiments in [21]. For the case study
with different architectures, we trained BNNs with $4$ different widths:
$10,20,30,40$. The horizon length is $H=10$ time steps, and the prediction
horizon is $5$ steps into the future. The neural networks were trained for
$200$ time steps, with a learning rate of $0.001$, and batch size of $128$
using mean-field variational inference. The training dataset consisted of
$28400$ training samples, recorded without meals.

Controller. We implemented a simple model predictive controller, using an off
the shelf implementation of the covariance matrix adaptation evolution
strategy (CMA-ES). The model predictive control planning horizon was $k=5$,
with a fixed seed for the randomized solver.

##  Appendix L Other possible baselines for IBNN

We considered comparing IBNNs to Bayesian model averaging (BMA) in the initial
state of our project. In [36], though, the authors show the pitfalls of BMA in
the context of Bayesian neural networks. As a consequence, we decided not to
use a credible interval resulting from BMA as a baseline for IBNNs.
Furthermore, our experiments may seem like a type of belief tracking [25, 42].
Two comments are in order. First, this line of literature does not use deep
learning techniques. Second, in these works the authors rely on Dempster-
Shafer theory, a field in imprecise probability theory where lower
probabilities are assumed to be belief functions, see Remark 13. We do not
rely on this assumption in our work.

##  Appendix M Further related work

Modeling uncertainty has been a longstanding goal of ML/AI research and a
variety of approaches have been developed for doing so [31, 58, 39]. Recently,
emphasis has been placed on discerning between aleatoric and epistemic
uncertainties [61, 45, 41]. In [67], the authors present an IP-based neural
network which uses a regression technique based on probability intervals.
Contrary to IBNNs, their NN is rooted in the frequentist approach to imprecise
probabilities [33].

In [52, Sections 2.1, 2.3], the authors focus on belief-functions-based
classification methods. IBNNs cannot be directly compared with this
methodologies because (i) they do not require that the scholar expresses their
knowledge via a belief function, but rather through a credal set; (ii) they
can be used for regression, prediction, and classification; (iii) they are
rooted in Bayesian theory, as opposed to Dempster-Shafer theory. Other works
in ML using a belief function approach are those from the field of evidential
machine learning (EML), see e.g. [17, 18] and references therein. Existing
models mainly address clustering, classification, and regression problems. The
reasons why IBNNs cannot be directly compared with methods from the EML
literature are (i) and (iii) above.

We also point out how there exists a field within ML called evidential deep
learning (EDL, not to be confused with EML). It aims at quantifying and
disentangling EU and AU, and its methodologies can be applied to
classification and regression, see e.g. [3, 62] and references therein. IBNNs
are not immediately comparable to models from EDL because these latter are
non-Bayesian. In addition, their definitions of AU and EU are slightly
different from the canonical ones in ensemble deep learning.

##  Appendix N Proofs

###### Proof of Proposition 2.

If $\overline{P}(A)=\sup_{P^{\prime}\in\Pi^{\prime}}P^{\prime}(A)$, for all
$A\in\mathcal{F}$, then it is immediate to see that
$\overline{P}(A)=\sup_{P\in\Pi}P(A)=\sup_{P^{\prime}\in\text{ex}\Pi^{\prime}}P^{\prime}(A)=\sup_{P^{\prime}\in\Pi^{\prime}}P^{\prime}(A)$,
for all $A\in\mathcal{F}$, since $\Pi\subset\Pi^{\prime}$.

Suppose now that $\overline{P}(A)=\sup_{P\in\Pi}P(A)$, for all
$A\in\mathcal{F}$. Then, we have that for all $P\in\Pi$ and all
$A\in\mathcal{F}$, $P(A)\leq\overline{P}(A)$. Pick now any
$P^{\prime}\in\Pi^{\prime}$. We can write it as
$P^{\prime}=\sum_{j=1}^{k}\alpha_{j}P_{j}$, where $\alpha_{j}\in[0,1]$, for
all $j\in\\{1,\ldots,k\\}$, $\sum_{j=1}^{k}\alpha_{j}=1$, and
$\\{P_{j}\\}_{j=1}^{k}=\Pi$. Pick then any $A\in\mathcal{F}$; we have

$P^{\prime}(A)=\sum_{j=1}^{k}\alpha_{j}P_{j}(A)\leq\sum_{j=1}^{k}\alpha_{j}\overline{P}(A)=\overline{P}(A).$

So $\overline{P}(A)\geq P^{\prime}(A)$. Because this holds for all
$P^{\prime}\in\Pi^{\prime}$ and all $A\in\mathcal{F}$, the claim is proven. ∎

###### Proof of Proposition 4.

Pick any metric $d$ on the space $\Delta(\Omega,\mathcal{F})$ of probabilities
on $(\Omega,\mathcal{F})$, and any $\mathbf{P}\in\mathscr{P}$. Because
$\mathbf{P}$ belongs to $\mathscr{P}$, $\inf_{P\in\mathscr{P}}d(P,P^{\prime})$
can only be either equal to or smaller than $d(\mathbf{P},P^{\prime})$. By the
definition of $d(\mathscr{P},P^{\prime})$, if
$\inf_{P\in\mathscr{P}}d(P,P^{\prime})=d(\mathbf{P},P^{\prime})$, then
$d(\mathscr{P},P^{\prime})=d(\mathbf{P},P^{\prime})$. If instead
$\inf_{P\in\mathscr{P}}d(P,P^{\prime})<d(\mathbf{P},P^{\prime})$, then
$d(\mathscr{P},P^{\prime})<d(\mathbf{P},P^{\prime})$. The proof is similar for
a generic divergence $div$ on $\Delta(\Omega,\mathcal{F})$. ∎

###### Proof of Proposition 7.

The result comes immediately from $\Pi^{\prime}$ being a superset of $\Pi$. ∎

###### Proof of Theorem 9.

Assume $\mathcal{P}^{\text{co}},\mathcal{L}_{\theta}^{\text{co}}$ are nonempty
and weak⋆-closed. Then, as a consequence of [53, Proposition 3] they are also
weak⋆-compact. Pick any $A\in\mathcal{B}$. Recall that we can rewrite the
usual Bayes’ updating rule as

$\displaystyle{P}_{D}(A)$
$\displaystyle=\frac{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)+\int_{\Theta}{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)}$
$\displaystyle=\frac{1}{1+\frac{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)}{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}},$

which is maximized when

$\frac{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)}{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}$

is minimized. But

$\frac{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)}{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}\geq\frac{\inf_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\underline{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)}{\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)},$

which proves the inequality in (2). Assume now that $\overline{P}$ is concave.
By [73, Lemma 1], we have that there exists
$\mathbf{P}\in\mathcal{P}^{\text{co}}$ such that

$\displaystyle\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}L(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)=\int_{\Theta}L(\theta)\mathbbm{1}_{A}(\theta)\mathbf{P}(\text{d}\theta),$
(4)

for all $L\in\mathscr{L}$. In addition, by [73, Lemma 4], we have that for all
$X\in\mathscr{X}$ and all $\epsilon>0$, there exists a non-negative, upper
semi-continuous function $h\leq X$ such that

$\displaystyle\begin{split}\left[\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}X(\theta)P(\text{d}\theta)\right]-\epsilon&<\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}h(\theta)P(\text{d}\theta)\\\
&\leq\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}X(\theta)P(\text{d}\theta).\end{split}$
(5)

Let now $X=\overline{L}\mathbbm{1}_{A}$. Notice that since
$\mathcal{L}^{\text{co}}_{\theta}$ is weak⋆-compact, by (1) so is
$\mathscr{L}$. This implies that $\underline{L},\overline{L}\in\mathscr{L}$,
since a compact set always contains its boundary, so
$\overline{L}\in\mathscr{X}$ as well, and in turn
$\overline{L}\mathbbm{1}_{A}\in\mathscr{X}$. Fix then any $L\in\mathscr{L}$
and put $h=L\mathbbm{1}_{A}$. It is immediate to see that $h$ is non-negative
and upper semi-continuous. Then, by (5), we have that for all $\epsilon>0$

$\displaystyle\begin{split}&\left[\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)\right]-\epsilon<\\\
&\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)\leq\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta).\end{split}$
(6)

Combining (4) and(6), we obtain

$\displaystyle\begin{split}&\left[\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)\right]-\epsilon\\\
&<\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)\mathbf{P}(\text{d}\theta)\leq\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta),\end{split}$
(7)

for all $L\in\mathscr{L}$.

Pick now any $\epsilon>0$ and put

$\displaystyle k$
$\displaystyle:=\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)$
$\displaystyle+\inf_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\underline{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)>0.$

Choose any $L\in\mathscr{L}$ and $\delta\in(0,\epsilon k)$. By (7) we have
that
$[\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)]-\delta<\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)\mathbf{P}(\text{d}\theta)$
and that
$[\inf_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\underline{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)]+\delta>\int_{\Theta}{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)\mathbf{P}(\text{d}\theta)$.
Recall that
$\mathbf{c}:=\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)+\inf_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\underline{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)P(\text{d}\theta)$,
and define
$\mathbf{d}:=\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)\mathbf{P}(\text{d}\theta)+\int_{\Theta}{L}(\theta)\mathbbm{1}_{A^{c}}(\theta)\mathbf{P}(\text{d}\theta)$.
Then,

$\displaystyle\mathbf{P}_{D}(A)$
$\displaystyle=\frac{\int_{\Theta}{L}(\theta)\mathbbm{1}_{A}(\theta)\mathbf{P}(\text{d}\theta)}{\mathbf{d}}$
$\displaystyle\geq\frac{\left[\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)\right]-\delta}{\mathbf{c}+\delta-\delta}$
$\displaystyle=\frac{\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}{\mathbf{c}}-\frac{\delta}{k}$
$\displaystyle>\frac{\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}{\mathbf{c}}-\epsilon.$

Since this holds for all $\epsilon>0$, we have that

$\sup_{P_{D}\in\mathcal{P}_{D}^{\text{co}}}P_{D}(A)=\frac{\sup_{P\in\mathcal{P}^{\text{co}}}\int_{\Theta}\overline{L}(\theta)\mathbbm{1}_{A}(\theta)P(\text{d}\theta)}{\mathbf{c}},$

concluding the proof. ∎

###### Proof of Lemma 10.

[70, 73] show that concave upper probabilities are closed with respect to the
generalized Bayes’ rule. In particular, this means that, if we let