# The diameter of a stochastic matrix: A new measure for sensitivity analysis
in Bayesian networks

Manuele Leonelli  Jim Q. Smith  Sophia K. Wright

###### Abstract

Bayesian networks are one of the most widely used classes of probabilistic
models for risk management and decision support because of their
interpretability and flexibility in including heterogeneous pieces of
information. In any applied modelling, it is critical to assess how robust the
inferences on certain target variables are to changes in the model. In
Bayesian networks, these analyses fall under the umbrella of sensitivity
analysis, which is most commonly carried out by quantifying dissimilarities
using Kullback-Leibler information measures. In this paper, we argue that
robustness methods based instead on the familiar total variation distance
provide simple and more valuable bounds on robustness to misspecification,
which are both formally justifiable and transparent. We introduce a novel
measure of dependence in conditional probability tables called the diameter to
derive such bounds. This measure quantifies the strength of dependence between
a variable and its parents. We demonstrate how such formal robustness
considerations can be embedded in building a Bayesian network.

###### keywords:

Bayesian networks , Edge strength , Sensitivity analysis , Total variation
distance

[label1]organization=School of Science and Technology, IE University, Madrid,
country = Spain

[label2]organization=Department of Statistics, University of Warwick,
Coventry, country = UK [label3]organization=The Alan Turing Institute, London,
country = UK

##  1 Introduction

Bayesian networks (BNs) (e.g Koller and Friedman, 2009; Smith, 2010) are now a
well-established and widely-used AI modelling approach for a wide range of
risk management applications. They support decision makers by providing an
intuitive graphical framework to reason about the dependence of various risk
factors and an efficient platform to perform inferential queries, scenario and
sensitivity analyses. Their use as a decision support tool in business and OR
has been increasing over the years, including case studies in project
management (van Dorp, 2020), supply chain (Garvey et al., 2015), marketing
(Hosseini, 2021), and logistics (Qazi, 2022), among others.

BNs are defined by two components: a directed acyclic graph (DAG) where each
node is a variable of interest and edges represent the, possibly causal,
relationship between them; a conditional probability table (CPT) for each node
of the DAG reporting the probability distribution of the associated variable
conditional on its parents. BNs are highly interpretable due to their
graphical nature, representing the probabilistic relationships between
variables, making it easy for users to understand and trace the influence of
one variable on another. With explainability now recognized as critical for
the use of AI in applied research (Rudin, 2019), including in OR (De Bock et
al., 2023), BNs stand out by providing transparent and intuitive explanations,
thereby enhancing trust and clarity in decision-making processes.

The underlying DAG and the associated CPTs can be learned from data using
machine learning algorithms or elicited using experts’ opinions and knowledge.
There is now a vast amount of algorithms to learn BN from data (e.g. Scutari
et al., 2019, for a review), as well as various software implementing such
routines (most notably the bnlearn R package, Scutari, 2010). Furthermore,
there are now protocols to guide the construction of BNs from the input of
experts (French, 2011; Renooij, 2001; Werner et al., 2017; Wilkerson and
Smith, 2021).

No matter how the BN was constructed, in any real-world modelling with BNs, it
is critical to assess the importance of various risk factors and evaluate the
overall robustness of the model to misspecifications of its inputs. Such a
step is usually referred to as _sensitivity analysis_ (Borgonovo, 2017; Razavi
et al., 2021) and is a fundamental ingredient of applied mathematical
modelling in any area of science (Borgonovo and Plischke, 2016; Saltelli et
al., 2000, 2021). A variety of sensitivity methods for BNs have been
introduced (see e.g. Ballester-Ripoll and Leonelli, 2023; Chan and Darwiche,
2002; Rohmer, 2020; van der Gaag et al., 2007) and implemented in various
pieces of software, (e.g. the bnmonitor R package, Leonelli et al., 2023).

Most sensitivity methods in BNs either use the KL-divergence or the Chan-
Darwiche distance to quantify the dissimilarity between two BNs with different
input parameters (Chan and Darwiche, 2005; Leonelli and Riccomagno, 2022;
Renooij, 2014). However, these measures depend very heavily on the accurate
specification of small probabilities since they are specified in log-
probabilities and ratio of probabilities, respectively. In many applied
scenarios, mainly when BNs are used as a decision support tool, the
misspecification of improbable events has only a small impact on the required
outputs of a decision analysis. For this reason, in this paper, we propose the
use of a much more conventional distance measure (widely used in probability
theory and stochastic analysis, Sason and Verdú, 2016), namely the total
variation distance. Although it is often difficult to derive explicit formulae
for the impacts of deviation in variation, it is nevertheless straightforward
to tightly bound such deviations, as we demonstrate below.

In particular, we introduce a novel measure, called diameter, which summarizes
the information carried by a CPT of the BN in variation distance. By amount of
information, we broadly mean strength of dependence between a variable and its
parents in the DAG, although the exact intuition varies depending on the type
of analysis conducted, as we illustrate in the following sections. The
diameter will be essential for constructing tight bounds to the effect of the
model’s perturbations derived from a suite of sensitivity studies.

Unlike almost all sensitivity methods available to BNs that require a complete
specification of the whole BN, the methods proposed here can rely on a partial
specification of the CPTs. Despite the critical role of sensitivity analysis
during the knowledge engineering process of constructing and eliciting a BN
(Coupé et al., 2000; Laskey and Mahoney, 2000), methods based on a partial
model specification are scarce (see e.g. Albrecht et al., 2014; Boneh et al.,
2006, for two exceptions).

An implementation of the developed routines is freely available via the
bnmonitor R package (Leonelli et al., 2023). Proofs of the main results are
collated in the appendix.

##  2 Bayesian networks and sensitivity analysis

###  2.1 Bayesian networks

Let $G=([n],E)$ be a DAG with vertex set $[n]=\\{1,\dots,n\\}$ and edge set
$E$. Let $\bm{X}=(X_{i})_{i\in[n]}$ be categorical random variables with joint
probability mass function (pmf) $p$ and sample space
$\mathbb{X}=\times_{i\in[n]}\mathbb{X}_{i}$. For $A\subset[n]$, we let
$\bm{X}_{A}=(X_{i})_{i\in A}$ and $\bm{x}_{A}=(x_{i})_{i\in A}$ where
$\bm{x}_{A}\in\mathbb{X}_{A}=\times_{i\in A}\mathbb{X}_{i}$. We say that $p$
is Markov to $G$ if, for $\bm{x}\in\mathbb{X}$,

$p(\bm{x})=\prod_{i\in[n]}p(x_{i}\mid\bm{x}_{\Pi_{i}}),$ (1)

where $\Pi_{k}$ is the parent set of $k$ in $G$ and
$p(x_{k}|\bm{x}_{\Pi_{k}})$ is a shorthand for
$p(X_{k}=x_{k}|\bm{X}_{\Pi_{k}}=\bm{x}_{\Pi_{k}})$. The factorization in
Equation (1) is equivalent to a set of conditional independences as formalized
by the local Markov property. In a BN with DAG $G$ it holds that, for every
$i\in[n]$,

$X_{i}\mathchoice{\mathrel{\hbox to0.0pt{$\displaystyle\perp$\hss}\mkern
2.0mu{\displaystyle\perp}}}{\mathrel{\hbox
to0.0pt{$\textstyle\perp$\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptstyle\perp$\hss}\mkern
2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptscriptstyle\perp$\hss}\mkern
2.0mu{\scriptscriptstyle\perp}}}X_{ND_{i}}|X_{\Pi_{i}},$

where $ND_{i}$ is the set of non-descendants of $i$ in $G$. For categorical
variables, each term of the factorization in Equation (1) is taken from the
CPT of the relevant vertex. The CPT is a stochastic matrix (i.e. non-negative
rows summing to one) where each row correspond to a combination of the
parents.

The DAG associated to a BN provides an intuitive overview of the relationships
between variables of interest. However, it also provides a framework to assess
if any generic conditional independence holds for a specific subset of the
variables via the so-called _d-separation_ criterion (see e.g. Koller and
Friedman, 2009). Efficient propagation of probabilities and evidence to
compute any (conditional) probability involving a specific subset of the
variables can be carried out exploiting the underlying DAG structure, which we
discuss in Section 2.3 below. Next, we introduce two examples to illustrate
BNs.

###  2.2 Two examples

(a) The asia BN

(b) The triangulated asia BN.

Figure 1: Illustration of the triangulation of a BN.

We first consider the well-known asia BN reported in Figure 1(a), representing
the possible causes of shortness of breath, also called dyspnea (dysp). A
visit to Asia (asia) may cause tuberculosis (tub), while smoking may cause
bronchitis (bronc) and lung cancer (lung). The vertex either is an or gate,
which is equal to yes if either lung or tub is equal to yes. The output of a
chest x-ray (xray) only depends on the presence of lung cancer or
tuberculosis, while dyspnea depends on all three diseases (bronchitis, lung
cancer, tuberculosis). The factorization of the pmf for this BN is

$p(\texttt{xray}|\texttt{either})p(\texttt{dysp}|\texttt{bronc},\texttt{either})p(\texttt{either}|\texttt{lung},\texttt{tub})\\\
p(\texttt{bronc}|\texttt{smoke})p(\texttt{lung}|\texttt{smoke})p(\texttt{tub}|\texttt{asia})p(\texttt{smoke})p(\texttt{asia}).$

bronc |  either | $P(\texttt{dysp}|\texttt{either},\texttt{bronc})$
---|---|---
|  |  yes |  no
yes |  yes |  0.9 |  0.1
yes |  no |  0.8 |  0.2
no |  yes |  0.7 |  0.3
no |  no |  0.1 |  0.9
lung |  tub | $P(\texttt{either}|\texttt{tub},\texttt{lung})$
---|---|---
|  |  yes |  no
yes |  yes |  1.0 |  0.0
yes |  no |  1.0 |  0.0
no |  yes |  1.0 |  0.0
no |  no |  0.0 |  1.0

Table 1: CPTs from the asia BN associated to the vertices dysp (left) and
either (right).

Table 1 reports two CPTs from the asia BN. Each row represents a combination
of the parent variables and the associated numeric entries are non-negative
and sum to one. Therefore, they are stochastic matrices. The variable either
is deterministically defined by its parents.

The second example is a BN introduced in Varando et al. (2024), learned over
data from the 2012 Italian enterprise innovation survey collected by ISTAT
(ISTAT, 2015), the Italian national statistical institute. The survey reports
information about medium-sized Italian companies and their involvement with
innovation in 2010-2012. The analysis aims to assess which factors related to
innovation are connected with changes in the company revenue. The considered
variables are reported in Table 2 and details about data pre-processing can be
found in Varando et al. (2024). Figure 2 reports the learned BN using the tabu
algorithm of the bnlearn R package (Scutari, 2010), comprising of 15 vertices
and 38 edges. The output variable (GROWTH) is directly affected by the number
of employees of the company (EMP12) and whether or not the company carried out
other innovation activities in 2010-2012 (INPD). Its CPT is reported in Table
3.

Name | Explanation | Levels
---|---|---
GP (P) | Belongs to an industrial group | Yes/No
LARMAR (L) | Main market | Regional/National/International
INPDGD (D) | Product innovation 2010-2012 | Yes/No
INPDSV (V) | Service innovation 2010-2012 | Yes/No
INPD (N) | Other innovations 2010-2012 | Yes/No
INABA (A) | Abandoned innovation 2008-2010 | Yes/No
INONG (I) | Ongoing innovation from 2008-2010 | Yes/No
CO (C) | Cooperation agreements for innovation | Yes/No
ORG (O) | New organization practices | Yes/No
MKT (M) | New marketing practices | Yes/No
PUB (B) | Contracts with public institutions | Yes/No
EMP12 (2) | Number of employees in 2012 | 10-49/50-249/$>$250
EMPUD (E) | Employees with degree | 0%/1-10%/$>$10%
RR (R) | Research & development | Yes/No
GROWTH (G) | Increased revenue 2012/2010 | Yes/No

Table 2: Variables from the 2012 ISTAT enterprise innovation survey. Figure 2:
The istat BN learned over the 2012 ISTAT enterprise innovation data.

EMP12 |  INPD | $P(\texttt{GROWTH}|\texttt{INPD},\texttt{EMP12})$
---|---|---
|  |  yes |  no
10-49 |  yes |  0.563 |  0.437
10-49 |  no |  0.469 |  0.531
50-249 |  yes |  0.608 |  0.392
50-249 |  no |  0.557 |  0.443
$>250$ |  yes |  0.636 |  0.364
$>250$ |  no |  0.590 |  0.410

Table 3: CPT associated to the vertex GROWTH in the istat BN in Figure 2.

###  2.3 Computing probabilities in BNs

Probabilistic inference in BNs is known to be NP-hard (Cooper, 1990). However,
algorithms that use the underlying DAG to speed up computations have been
defined. One of the most famous algorithms transforms the original DAG into a
junction tree (see e.g. Koller and Friedman, 2009). We will show in Section 5
that the junction tree can also be used for sensitivity investigations (also
used in Kjærulff and van der Gaag, 2000).

The junction tree algorithm first transforms the original DAG into an
undirected graph, by first applying moralization (the addition of an edge
between any two parents of the same vertex not joined by an edge), then
dropping the directionality of the edges and lastly triangulating the graph
(adding undirected edges until the resulting graph is such that every cycle of
length strictly greater than 3 possesses a chord, that is, an edge joining two
nonconsecutive vertices of the cycle). The result of this process over the
asia BN is shown in Figure 1(b).

The cliques of the triangulated graph, its maximal fully connected subgraphs,
$C_{1},\dots,C_{m}$, can be totally ordered starting from any clique including
a root of the original graph. Let $S_{i}=C_{i}\cap\cup_{j=1}^{i-1}C_{j}$ be
the separator of $C_{i}$ from the preceding cliques. The cliques can always be
ordered to respect the running intersection property, meaning that there is at
least one $j<i$ such that $S_{i}\subset C_{j}$ for any
$i\in[m]\setminus\\{1\\}$. This implies that the result of intersecting a
clique with all previous cliques is contained within one or more earlier
cliques. Any clique ordering and choice of separator containment can be
depicted by a so-called junction tree: an undirected tree graph with vertices
$C_{1},\dots,C_{m}$ and an undirected edge between $C_{i}$ and $C_{j}$ if
$C_{i}$ is the chosen clique respecting $S_{j}\subset C_{i}$. The
factorization in Equation (1) can be equivalently written in terms of the
cliques as

$p(\bm{x})=\frac{\prod_{i\in[m]}p(\bm{x}_{C_{i}})}{\prod_{i\in[m]\setminus\\{1\\}}p(\bm{x}_{S_{i}})}=\prod_{i\in[m]\setminus\\{1\\}}p(\bm{x}_{C_{i}}|\bm{x}_{S_{i}})p(\bm{x}_{C_{1}}).$

For the asia BN, a possible ordering of its cliques is
$C_{1}=\\{\texttt{asia},\texttt{tub}\\}$,
$C_{2}=\\{\texttt{tub},\texttt{lung},\texttt{either}\\}$,
$C_{3}=\\{\texttt{either},\texttt{xray}\\}$,
$C_{4}=\\{\texttt{bronc},\texttt{lung},\texttt{either}\\}$,
$C_{5}=\\{\texttt{bronc},\texttt{either},\texttt{dysp}\\}$, and
$C_{6}=\\{\texttt{smoke},\texttt{bronc},\texttt{lung}\\}$, with separators
$S_{2}=\\{\texttt{tub}\\}$, $S_{3}=\\{\texttt{either}\\}$,
$S_{4}=\\{\texttt{lung},\texttt{either}\\}$,
$S_{5}=\\{\texttt{bronc},\texttt{either}\\}$, and
$S_{6}=\\{\texttt{bronc},\texttt{lung}\\}$. The resulting junction tree is
reported in Figure 3, where, as customary, we labeled the edges with the
separators.

Figure 3: The junction tree representing the asia BN. Initial of variable
names are used. Cliques are circled, while separators are squared.

A straightforward consequence of the junction tree structure is that any two
cliques, $C_{i}$ and $C_{j}$ say, are connected by a simple path. Assume that
$C_{j}$ includes descendants of $C_{i}$ in the original DAG. This implies a
sequence of separators in the unique path between $C_{i}$ and $C_{j}$, say
$S_{i+1}$,…$S_{j}$. Letting $S_{k}^{*}=S_{k}\setminus\cup_{l=k+1}^{j}S_{l}$,
for $k=i+1,\dots,j$, we have that

$\displaystyle p(\bm{x}_{C_{i}},\bm{x}_{C_{j}})$
$\displaystyle=\sum_{\bm{x}_{T_{ij}}\in\mathbb{X}_{T_{ij}}}p(\bm{x}_{C_{i}},\bm{x}_{C_{j}},\bm{x}_{T_{ij}})$
$\displaystyle=\sum_{\bm{x}_{T_{ij}}\in\mathbb{X}_{T_{ij}}}p(\bm{x}_{C_{j}}|\bm{x}_{S^{*}_{j}})p(\bm{x}_{S^{*}_{j}}|\bm{x}_{S^{*}_{j-1}})\cdots
p(\bm{x}_{S^{*}_{i+2}}|\bm{x}_{S^{*}_{i+1}})p(\bm{x}_{C_{i}}),$ (2)

where $T_{ij}=\cup_{k=i+3}^{j-1}S^{*}_{k}\setminus\left\\{C_{j}\cup
C_{i}\right\\}$. Equation (2) can be thought of as the pmf of a “donating”
clique $C_{i}$ and a “target” clique $C_{j}$, expressed in terms of a sequence
of transitions in a non-homogeneous Markov chain. It means that standard
results from non-homogeneous Markov chain theory can be used to measure the
extent of the diminishing effect of information as it passes along the simple
path from $C_{i}$ to $C_{j}$. In particular, it is well-known that variation
distance in an ergodic, acyclic Markov chain contracts as information is
propagated through the system (e.g. Roberts and Rosenthal, 2004). This
observation will be critical for the developments of sensitivity bounds we
introduce later on in this paper.

###  2.4 Sensitivity analysis in BNs

Sensitivity methods for BNs have been widely studied (see Rohmer, 2020, for
the most comprehensive yet only partial review). Here, we provide an overview
of BN models’ most traditional sensitivity investigations. Notice that most
methods require a complete specification of the whole model: its DAG and all
entries of its CPTs.

The most widely used sensitivity analysis studies the effect of perturbations
of the CPT entries on outputs of interest. Sensitivity functions
mathematically model the relationship between inputs and outputs (Castillo et
al., 1997; Leonelli et al., 2017; van der Gaag et al., 2007). Because of the
computational complexity of deriving sensitivity functions for multiple
parameter variations (Chan and Darwiche, 2004; Kwisthout and van der Gaag,
2008), in practice, most often only perturbations of one CPT entry at a time
are performed, although recent methodological advances matching BNs to more
flexible structures have made more complex investigations possible (Ballester-
Ripoll and Leonelli, 2022b; Salmani and Katoen, 2023).

Another type of sensitivity investigation studies the overall effect of a node
on an output. This is often quantified by the mutual information between the
associated variables (Kjaerulff and Madsen, 2008). Albrecht et al. (2014)
approached this problem by only considering the DAG of the BN, thus without
requiring the CPTs to be defined, and introduced the distance weighted
influence between two variables $X_{j}$ and $X_{i}$ of a DAG. Let $S_{ji}$ be
the set of active, simple trails from $j$ to $i$ (see e.g. Koller and
Friedman, 2009) and $w\in(0,1]$. The distance weighted influence of $X_{j}$ on
$X_{i}$ is

$DWI(X_{j},X_{i},w)=\sum_{s\in S_{ji}}w^{|s|},$ (3)

where $|s|$ is the length of the trail $s$. So $DWI$ measures how connected
two vertices of a DAG are, where longer trails have a smaller contribution to
the influence. If $w=1$, DWI simply counts the number of active, simple
trails. Another way to measure the influence of a node on an output is to
quantify how valuable it would be to observe the associated variable, the so-
called sensitivity to evidence approach (Ballester-Ripoll and Leonelli, 2022a;
Gómez-Villegas et al., 2014).

Another class of sensitivity methods assesses the relevance of the edges of a
BN. One standard way to do this for data-learned BNs is to use a non-
parametric bootstrap approach and learn a BN for each replication: edges that
do not appear frequently are deemed to have less strength (Scutari and
Nagarajan, 2013). Another possibility is to use sensitivity functions to
quantify the effect of an edge removal (Renooij, 2010). Given the complexity
of probabilistic inference, deleting edges having a small impact on the
inferences made by the model (Choi et al., 2005; Choi and Darwiche, 2006) is
often desirable.

In Section 4 we introduce new methods for the node relevance and edge strength
problems mentioned above, but we also address investigations that have
received less attention in the literature.

##  3 Total variation and the diameter

We next introduce our new metric, the diameter, based on the total variation
distance between pmfs.

######  Definition 1.

Let $p$ and $p^{\prime}$ be two pmfs over the same sample space $\mathbb{X}$.
The _total variation distance_ between $p$ and $p^{\prime}$ is

$d_{V}(p,p^{\prime})=\frac{1}{2}\sum_{x\in\mathbb{X}}|p(x)-p^{\prime}(x)|.$

As quantified through the absolute differences, deviation in variation
corresponds much more closely to the types of error we would envisage
experiencing within either an elicitation exercise or through misestimation.

We next define a measure of the overall total variation distance between the
rows of a CPT. As already noticed, CPTs are stochastic matrices where each row
is a pmf.

######  Definition 2.

The (upper) _diameter_ of a $n\times m$ stochastic matrix $P$ with rows
$p_{1},\dots,p_{n}$, denoted as $d^{+}(P)$, is

$d^{+}(P)=\max_{i,j\in[n]}d_{V}(p_{i},p_{j}).$

The lower diameter is

$d^{-}(P)=\min_{i,j\in[n]}d_{V}(p_{i},p_{j}).$

Our main focus is on the upper diameter which we henceforth simply refer to as
the diameter, while the lower diameter will become relevant in Section 4.4
only. The diameter is the largest variation distance between any two rows of a
stochastic matrix. The use of the maximum distance is motivated by
establishing bounds on the effect of probability perturbations and their
propagation throughout the DAG.

We next characterize properties of the diameter and its relationship with
conditional independence and marginalization, which will further shed light on
its interpretation. We therefore now focus specifically on stochastic matrices
representing CPTs. For a CPT we use the notation $P_{\cdot|\cdot}$, where the
elements on the rhs of the subscript are the conditioning variables, while
those on the lhs of the subscript are those of which the pmf is evaluated. For
ease of interpretation, we use $X,Y,Z$ to denote generic categorical random
variables (although all the results below can be written for categorical
random vectors).

######  Proposition 1.

For two categorical random variables $X$ and $Y$ it holds

$d^{+}(P_{Y|X})=0\Leftrightarrow Y\mathchoice{\mathrel{\hbox
to0.0pt{$\displaystyle\perp$\hss}\mkern
2.0mu{\displaystyle\perp}}}{\mathrel{\hbox
to0.0pt{$\textstyle\perp$\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptstyle\perp$\hss}\mkern
2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptscriptstyle\perp$\hss}\mkern
2.0mu{\scriptscriptstyle\perp}}}X.$

The closer the diameter is to zero, the less dependent two variables are. To
see this, it is easy to check that whenever some non-trivial function of $Y$
can be written as a deterministic function of $X$ then $d^{+}(P_{Y|X})=1$, its
maximum value. So when changing the levels of $X$ has a minimum impact on the
pmf of $Y$ then $d^{+}(P_{Y|X})\approx 0$. Note that unless $P_{Y|X}$ is
symmetric, $d^{+}(P_{Y|X})\neq d^{+}(P_{X|Y})$, in fact the difference between
these can be arbitrarily close to 1 (Wright, 2018).

tub | lung | bronc | either | xray | dysp
---|---|---|---|---|---
0.04 | 0.09 | 0.30 | 1.00 | 0.93 | 0.80

Table 4: Diameter of the non-root nodes of the asia BN.

GP | LARMAR | INPDGD | INPDSV | INPD | INABA | INONG | CO | ORG | MKT | PUB | RR | EMPUD | GROWTH
---|---|---|---|---|---|---|---|---|---|---|---|---|---
0.666 | 0.611 | 0.857 | 0.859 | 0.648 | 0.500 | 0.761 | 0.484 | 0.546 | 0.697 | 0.344 | 0.941 | 0.521 | 0.167

Table 5: Diameter of the non-root nodes of the istat BN. Figure 4: Heatmap of
the diameters of the istat BN. Darker colors represent higher diameter values.

Tables 4 and 5 report the diameters of the non-root nodes of the asia and
istat BNs. For the asia BN it can be noticed highly different values of the
diameter, with the variables tub and lung almost independent of their parents,
while the diameter of either is one since it is a deterministic function of
its parents (or gate). For the istat BN, the lowest diameter is the one of the
target variable GROWTH, demonstrating that a company’s revenue weakly depends
on its parents. The CPT of GROWTH in Table 3 includes very similar rows.
Figure 4 gives a visualization of the istat BN diameters, where nodes with a
darker color have CPTs whose diameter is larger.

Next we formalize how the diameter behaves under marginalization.

######  Proposition 2.

For three categorical random variables $X,Y,Z$ it holds

$d^{+}(P_{Y|X})\leq d^{+}(P_{Y|XZ}),$

and equality holds if and only if $Y\mathchoice{\mathrel{\hbox
to0.0pt{$\displaystyle\perp$\hss}\mkern
2.0mu{\displaystyle\perp}}}{\mathrel{\hbox
to0.0pt{$\textstyle\perp$\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptstyle\perp$\hss}\mkern
2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptscriptstyle\perp$\hss}\mkern
2.0mu{\scriptscriptstyle\perp}}}Z|X$.

The intuition behind this result is that each row of $P_{Y|X}$ is a weighted
average of some rows of $P_{Y|XZ}$ and therefore the rows of $P_{Y|X}$ are
necessarily “closer” to each other. This implies that the diameter of a CPT
based on a subset of the variables in the parent set can always be bounded by
the diameter of the CPT based on the complete set of parents.

The last result writes the diameter of a random vector as the sum of the
diameter of individual CPTs.

######  Proposition 3.

For three categorical random variables $X,Y,Z$ it holds

$d^{+}(P_{YZ|X})\leq\min\\{d^{+}(P_{Y|XZ})+d^{+}(P_{Z|X}),1\\}.$

These results are useful to construct bounds over the CPTs of vertices of a
junction tree, by using the diameters of the original CPTs of the BN, thus not
having to compute any new information. We showcase their usefulness in Section
5.

In the next section we discuss how the diameter can be used for various
sensitivity investigations in BNs which are fully defined (all CPTs have been
fully learned or elicited). However, we envisage that the diameter can also be
directly elicited from experts given a DAG. The sensitivity methods we develop
could then drive the complete elicitation of the BN by focusing on the CPTs
which have the biggest impact on the output of interest. A more comprehensive
discussion of the use of the diameter within the elicitation of a BN is given
in Section 7.

##  4 Sensitivity analysis using the diameter

###  4.1 Edge strength

Figure 5: The asia BN with edge strengths as edges’ labels and widths.

The first problem we consider is the quantification of the strength of an edge
in a BN. Let $P_{i}$ be the CPT of $Y_{i}$ and with $P_{i|\bm{x}}$ we denote
the sub-CPT of $P_{i}$ including only the rows specified by $\bm{x}$.

######  Definition 3.

The strength of edge $(j,i)$ in a BN is defined as

$\delta_{ji}=\max_{\bm{x}\in\mathbb{X}_{\Pi_{i}\setminus
j}}d^{+}(P_{i|\bm{x}}).$

So $\delta_{ji}$ is the largest diameter out of all CPTs for every combination
of all parents of $i$ excluding $j$. To illustrate this, consider the CPT of
GROWTH in Table 3. The strength of the edge (INPD,GROWTH) is the largest of
the diameters of the the three CPTS where EMP12 is fixed to its three levels
10-49/50-249/$>250$.

Figures 5 and 6 report the edge strengths as edge labels and widths in the
asia and istat BN, respectively. They vary from almost zero to exactly one in
the case of a functional relationship. Interestingly, the two edges with the
lowest strength in the istat BN are those into the output variable GROWTH,
again indicating how the revenue of an enterprise is almost independent of all
other factors.

Figure 6: The istat BN with edge strengths as edges’ labels and widths.

Just as the diameter represents marginal independence, edge strength denotes
conditional independence.

######  Proposition 4.

$X_{i}\mathchoice{\mathrel{\hbox to0.0pt{$\displaystyle\perp$\hss}\mkern
2.0mu{\displaystyle\perp}}}{\mathrel{\hbox
to0.0pt{$\textstyle\perp$\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptstyle\perp$\hss}\mkern
2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptscriptstyle\perp$\hss}\mkern
2.0mu{\scriptscriptstyle\perp}}}X_{j}|\bm{X}_{\Pi_{i}\setminus
j}\Leftrightarrow\delta_{ji}=0.$

Thus, in a formal sense, $\delta_{ji}$ is a measure of the extent by which
this conditional independence is violated and the merit of knowing the value
of $X_{j}$ to predict $X_{i}$ once we know the value of
$\bm{X}_{\Pi_{i}\setminus j}$.

The following proposition links edge strength to the diameter of a CPT.

######  Proposition 5.

It holds that $\delta_{ji}\leq d^{+}(P_{i})$ and
$d^{+}(P_{i})\leq\sum_{j\in\Pi_{i}}\delta_{ji}$. Also if $|\Pi_{i}|=1$, then
$\delta_{ji}=d^{+}(P_{i})$.

Notice that to derive the edge strength $\delta_{ji}$ from elicitation only,
$|\mathbb{X}_{\Pi_{i}\setminus j}|$ diameters must be defined. Edges that
appear to be weak could then be discarded before attempting a full
quantitative elicitation of the complete CPT, since the size of CPTs increase
quadratically with the number of parents.

###  4.2 Edge weigthed influence

Given a DAG and edges labeled with their strength, we may be interested in
quantifying the effect of a node on another. For this task, we define a novel
measure we henceforth call edge weighted influence. Recall that $S_{ji}$ is
the set of active, simple trails from $j$ to $i$.

######  Definition 4.

The edge weigthed influence of $X_{j}$ on $X_{i}$, $EWI(X_{j},X_{i})$, is
defined as:

$EWI(X_{j},X_{i})=\sum_{s\in S_{ji}}\left(\prod_{(k,l)\in
s}\delta_{kl}\right)^{|s|}.$ (4)

The definition of the edge weighted influence is inspired by the one of the
distance weighted influence where, instead of giving a weight $w$ to every
edge of the BN, we consider the edge strengths $\delta_{ji}$. The edge
weighted influence sits inbetween mutual information, which requires a full BN
definition, and the distance weighted influence, only requiring the DAG, since
it is based on the DAG together with some measure of edge strength.

In our examples, the edge weighted distance is derived using the edge
strengths computed from the full CPTs, but any measure of edge strength could
be considered. For instance, it could be the proportion of times an edge has
appeared in non-parametric bootstrap structural learning. An alternative is if
these strengths were somehow elicited directly from experts.

|  | DWI |
---|---|---|---
Name | Mutual Information | $w=0.1$ | $w=0.2$ | $w=0.5$ | $w=1$ | EWI
GP | 0.00132 (6) | 0.01121 (11) | 0.05152 (13) | 0.53125 (14) | 5 (12) | 0.00645 (4)
LARMAR | 0.00094 (10) | 0.01363 (9) | 0.07683 (9) | 1.70898 (7) | 65 (5) | 0.00145 (8)
INPDGD | 0.00166 (4) | 0.01397 (7) | 0.08270 (7) | 1.91797 (5) | 55 (7) | 0.00075 (10)
INPDSV | 0.00104 (9) | 0.01448 (6) | 0.09372 (5) | 3.05078 (2) | 136 (3) | 0.00351 (6)
INPD | 0.00410 (2) | 0.11111 (1) | 0.24992 (1) | 0.96875 (10) | 5 (12) | 0.09434 (2)
INABA | 0.00029 (13) | 0.00337 (14) | 0.04374 (14) | 2.56641 (4) | 131 (4) | 0.00005 (13)
INONG | 0.00128 (7) | 0.00508 (12) | 0.05262 (12) | 1.88672 (6) | 57 (6) | 0.00009 (12)
CO | 0.00076 (11) | 0.01390 (8) | 0.08017 (8) | 1.53906 (8) | 31 (8) | 0.00089 (9)
ORG | 0.00147 (5) | 0.02221 (3) | 0.09952 (4) | 0.90625 (12) | 7 (10) | 0.00352 (5)
MKT | 0.00074 (12) | 0.01618 (4) | 0.10517 (3) | 3.42188 (1) | 203 (2) | 0.00032 (11)
PUB | 0.00001 (14) | 0.00503 (13) | 0.05401 (11) | 2.83301 (3) | 208 (1) | 0.00000 (14)
EMP12 | 0.00503 (1) | 0.11111 (1) | 0.24992 (1) | 0.96875 (10) | 5 (12) | 0.12116 (1)
EMPUD | 0.00105 (8) | 0.01321 (10) | 0.06752 (10) | 0.78125 (13) | 7 (10) | 0.00232 (7)
RR | 0.00277 (3) | 0.01465 (5) | 0.08339 (6) | 1.32812 (9) | 19 (9) | 0.00774 (3)
|  | 0.65566 | 0.52805 | -0.67472 | -0.80443 | 0.81538

Table 6: Node influence measures in the istat BN considering GROWTH as output
node. In parenthesis the node ranking for each column. The last row reports
the Spearman correlation between mutual information and the other influence
measures.

To illustrate the use of the edge weigthed distance we consider the istat BN,
since asia has a too simple topology. Table 6 reports the mutual information,
distance weighted influence for various choices of $w$, and the edge weighted
distance between GROWTH and every other node. As already noticed by Albrecht
et al. (2014) the $DWI$ measures varies greatly with $w$, with little
intuition behind what an optimal $w$ value might actually be. On the other
hand, $EWI$ does not require the choice of $w$, of course at the cost of
having to specify edge strengths, which on the other hand have a much more
straightforward meaning. The highest Spearman correlation with mutual
information, reported in the last row of Table 6, is attained by $EWI$
showcasing how the underlying DAG together with edge strengths gives a very
precise approximation to the full BN model.

Figure 7: Heatmap of $EWI$ in the istat BN with GROWTH as output (reported in
light blue).

The two variables with the greatest effect on GROWTH are, as expected, its
parents EMP12 and INPD. However, because of the weak dependence between GROWTH
and its parents and the actual DAG topology, the EWI of the other nodes
quickly vanishes. Figure 7 reports a heatmap of the EWI which clearly provide
a visualization of this. Such a heatmap can be highly valuable during
elicitation of the entries of the CPTs to convince practitioners on focusing
on the specification of the probabilities of the most important factors first.

###  4.3 Level amalgamation

One practical issue found by discrete BN modellers is the number of levels
each random variable within the system should be assigned. Obviously there is
a trade-off here. The finer the division of levels, the more nuanced the BN
can be. On the other hand, the fewer the number of levels, the easier it will
be to faithfully elicit or efficiently estimate the probabilities within a BN.
The technology of the diameter can be adapted to guide possible amalgamations
of the levels of categorical variables, just as when considering whether or
not to keep a weak edge in the system.

A first consideration to take into account is that the interpretation of the
states can still be understood and quantified by experts. For ordinal
variables is therefore recommended to only consider merging consecutive levels
(for instance for the variable EMP12 with levels 10-49/50-249/$>50$, one
should consider merging only the first two or the last two levels).

The second step in level amalgamation is deciding how to combine the
probabilities associated with those levels that are to be amalgamated. Suppose
the levels of variable $X_{i}$ are being merged. In the CPT $P_{i}$ then the
probabilities of the two associated levels are simply summed up. However, also
the CPTs of a child of $Y_{i}$ need to be adapted because of amalgamation. We
recommend taking a simple row average, because the convexity of variation
distance (Lemma 3) guarantees that the diameter of the original CPT does not
increase, and, more importantly, this method does not require additional
information. In practice we can calculate the diameter of the CPTs with levels
amalgamated and combine the closest states first, then find the next closest
states and add to the amalgamation iteratively until the combination appears
to induce a significant change from the original diameter.

|  |  GP |  LARMAR |  INPDGD |  INPD |  CO |  ORG |  EMPUD |  GROWTH
---|---|---|---|---|---|---|---|---|---
EMP12 |  10-49/50-249 |  0.461 |  0.541 |  0.857 |  0.607 |  0.484 |  0.495 |  0.404 |  0.123
|  50-249/$>$250 |  0.538 |  0.597 |  0.788 |  0.577 |  0.376 |  0.456 |  0.500 |  0.153
|  Original |  0.666 |  0.611 |  0.857 |  0.648 |  0.484 |  0.546 |  0.521 |  0.167

|  |  CO |  ORG |  RR
---|---|---|---|---
EMPUD |  0%/$<$10% |  0.484 |  0.483 |  0.907
|  $<$10%/$>$10% |  0.403 |  0.458 |  0.941
|  Original |  0.484 |  0.546 |  0.941

|  |  PUB
---|---|---
LARMAR |  Regional/National |  0.111
|  National/International |  0.186
|  Original |  0.344

Table 7: Diameters of the CPTs in the istat BN resulting from the amalgamation
of levels for the variables EMP12, EMPUD, and LARMAR.

Table 7 reports the diameters resulting from the amalgamation of levels for
variables with more than two levels in the istat BN. It can be seen that
merging the levels of LARMAR leads to a strong decrease in the diameter. On
the other hand, overall the diameter of the children CPTs of the variables
EMP12 and EMPUD does not show a strong decrease, highlighting that the
amalgamation of a pair of levels would have a small effect on the BN.

###  4.4 Asymmetry strength

There is now an increasingly amount of evidence that standard conditional
independence may be too restrictive to faithfully and fully represent
dependence patterns in data (Eggeling et al., 2019; Leonelli and Varando,
2024; Pensar et al., 2015). For BNs this means that there are equalities
between probability distributions within the CPTs, which therefore have no
graphical counterpart in the underlying DAG.

The simplest class of constraints that could appear in a CPT of a BN is the
so-called context-specific conditional independence (Boutilier et al., 1996).
We say that $X_{i}$ is conditionally independent of $X_{j}$ given context
$X_{k}=x_{k}$ if $p(x_{i}|x_{j},x_{k})=p(x_{i}|x_{k})$ for all
$x_{i}\in\mathbb{X}_{i},x_{j}\in\mathbb{X}_{j}$ and a specific
$x_{k}\in\mathbb{X}_{k}$. A context-specific independence reduces to a
standard, symmetric independence if it holds for all $x_{k}\in\mathbb{X}_{k}$.
Consider the CPT in Table 8 of a random variable $X_{i}$ conditional on
$X_{j}$ and $X_{k}$, all taking levels high, medium, and low. It can be seen
that $X_{i}$ is conditionally independent of $X_{j}$ when $X_{k}=\text{high}$
since all rows where $X_{k}=\text{high}$ have the same pmf.

$X_{j}$ |  $X_{k}$ | $P(X_{i}|X_{j},X_{k})$
---|---|---
|  |  high |  medium |  low
high |  high |  0.500 |  0.300 |  0.200
high |  medium |  0.400 |  0.200 |  0.500
high |  low |  0.300 |  0.100 |  0.600
medium |  high |  0.500 |  0.300 |  0.200
medium |  medium |  0.300 |  0.500 |  0.200
medium |  low |  0.200 |  0.400 |  0.400
low |  high |  0.500 |  0.300 |  0.200
low |  medium |  0.200 |  0.200 |  0.600
low |  low |  0.200 |  0.200 |  0.600

Table 8: Example of a CPT embedding context-specific and partial conditional
indepenences.

Pensar et al. (2016) introduced a more general extension to conditional
independence called partial conditional independence. We say that $X_{i}$ is
partially conditionally independent of $X_{j}$ in the domain
$\mathcal{D}_{j}\subseteq\mathbb{X}_{j}$ given context $X_{k}=x_{k}$ if
$p(x_{i}|x_{j},x_{k})=p(x_{i}|\tilde{x}_{j},x_{k})$ holds for all
$(x_{i},x_{j}),(x_{i},\tilde{x}_{j})\in\mathbb{X}_{i}\times\mathcal{D}_{j}$.
Partial and context-specific independence coincides if
$\mathcal{D}_{j}=\mathbb{X}_{j}$. Furthermore, the sample space
$\mathbb{X}_{j}$ must contain more than two elements for a non-trivial partial
conditional independence to hold. The CPT in Table 8 embeds the partial
conditional independence between $X_{i}$ and $X_{k}$ in the domain
$\\{\text{medium},\text{low}\\}$ in the context $X_{j}=\text{low}$.

Currently, there are no methods to assess whether additional equalities are
present in the CPTs of a BN, beyond visual investigation. Here we demonstrate
that the diameter can be also be used for this task, just as when considering
whether or not to keep a weak edge in the system.

######  Definition 5.

The index of context-specific independence between $X_{i}$ and $X_{j}$, for
$j\in\Pi_{i}$, in the context $\bm{x}\in\mathbb{X}_{\Pi_{i}\setminus j}$ is

$\delta_{\bm{x}i}^{+}=\ d^{+}(P_{i|\bm{x}}),$

and the index of partial independence is

$\delta_{\bm{x}i}^{-}=\ d^{-}(P_{i|\bm{x}}).$

The index of context-specific independence computes the upper diameter of a
CPT where all parents but one are fixed, while the index of partial
independence computes the lower diameter for the same CPT. If the lower
diameter is close to zero it means that there are at least two rows of the CPT
which are very similar to each other, thus implying a partial conditional
independence.

GP
---
EMP12 | Context
10-49 | 0.255
50-249 | 0.224
$>250$ | 0.284

EMP12
---
GP | Context | Partial
No | 0.398 | 0.007
Yes | 0.265 | 0.131

Table 9: Context-specific and partial independence indices for the CPT of the
variable EMPUD in the istat BN.

Propositions 4 and 5 could be straightforwardly extended to link the context-
specific index to context-specific independence and bound the index with the
diameter of the CPT.

As an illustration Table 9 reports the context-specific and partial indices
for the CPT of the variable EMPUD in the istat BN. The context-specific
independence indices are quite close to each other, with a very strong
dependence between EMP12 and EMPUD in the context GP = No. The partial
independence index is only reported for the EMP12 variables since it is
ternary. It can be seen that in the context GP = No it is quite plausible the
presence of a partial independence: the probability distribution of EMPUD is
very similar for at least two levels of EMP12 in the context GP = No. C shows
that a model selection algorithm for a class of models that formally embeds
also non-symmetric types of independence does indeed learn the mentioned
partial independence.

##  5 Error propagation

It is well known that when using standard propagation algorithms on updating
one of the clique margins, say $C_{i}$, the knock-on effect on the other
clique margins becomes weaker and weaker as the updated cliques become
progressively more remote from $C_{i}$. This property is exploited by Albrecht
et al. (2014) when defining the DWI influence, for instance, and thus also in
the EWI distance. It is also known that if the underlying DAG is a polytree,
the mutual influence between two nodes decreases with the distance (number of
edges) between them.

The extent of the deviation can be bounded using variation distance, providing
an upper limit to the potential error in the distributions of focus variables
induced by the misspecification of various CPTs in the BN. This is
particularly useful when we elicit a large BN and want to know how far away
from target nodes we need to elicit the corresponding CPTs accurately. The
following result provides the basis for the derivation of such bounds.

######  Theorem 1.

Consider two categorical random variables $X$ and $Y$. Let $p$ and
$p^{\prime}$ two pmfs over $(X,Y)$ such that $p(y|x)=p^{\prime}(y|x)$ for all
$x\in\mathbb{X}$ and $y\in\mathbb{Y}$. Then

$d_{V}(p(y),p^{\prime}(y))\leq d^{+}(P_{Y|X})d_{V}(p(x),p^{\prime}(x)).$

The interpretation of this result is as follows. Suppose the CPT of $Y$ given
$X$ has been specified accurately, but that the margin probability of $Y$ is
uncertain. Then our marginal beliefs about $Y$ are no more uncertain than
those about $X$, because by definition $d^{+}(P_{Y|X})\leq 1$. More
importantly, we have a bound on how much the uncertainty, quantified in terms
of total variation, reduces in terms of $d^{+}(P_{Y|X})$ \- a measure of how
far away $Y$ is from independence of $X$.

We can now use Theorem 1 to provide a bound of the effect of perturbation of
an output variable on the cliques of the junction tree of a BN.

######  Theorem 2.

Consider an output variable $X_{j}$ which has no children in the BN and assume
it is in clique $C_{j}$ with no loss of generality. Let $C_{i}$ be another
clique and $S_{i+1},\dots,S_{j}$ be the separators along the unique path
between $C_{i}$ and $C_{j}$. Let
$S_{k}^{*}=S_{k}\setminus\cup_{l=k+1}^{j}S_{l}$, for $k=i+1,\dots,j$,
$P_{k}^{*}$ be the CPT representing
$p(\bm{x}_{S^{*}_{k}}|\bm{x}_{S^{*}_{k-1}})$ for $k=i+2,\dots,j$ and $P_{j}$
be the CPT of $p(\bm{x}_{j}|\bm{x}_{S^{*}_{j}})$. Then for any two pmfs $p$
and $p^{\prime}$ we have that

$d_{V}(p(x_{j}),p^{\prime}(x_{j}))\leq
d^{+}(P_{j})\prod_{k=i+2}^{j}d^{+}\left(P^{*}_{k}\right)d_{V}(p(\bm{x}_{C_{i}}),p^{\prime}(\bm{x}_{C_{i}})).$

The proof of this result follows easily from Equation (2) and successive
application of Theorem 1. Notice that the quantities $d^{+}(P^{*}_{k})$ may or
may not be readily available from a fully defined BN. In the best case
scenario, their computation may simply require the use of Propositions 2 and 3
which take advantage of simple properties of the diameter. In other cases,
when a separator $S_{k+1}$ includes parents of variables in $S_{k}$, the
required CPT must be constructed using inferential queries over the BN.

Next we define the impact of one clique upon a variable of interest in order
to ascertain the diminishing effect of errors downstream in the chain of a
junction tree.

######  Definition 6.

The impact of $C_{i}$ on an output variable $Y_{j}$ in clique $C_{j}$ is
$d^{+}(P_{j})\prod_{k=i+2}^{j}d^{+}\left(P^{*}_{k}\right)$.

The impact of one clique on an output is a simple measure of the maximum
possible influence the misspecification of one set of probabilities could have
on another as measured by a bound on the variation distance. As an
illustration, Figure 8 shows the junction tree of the asia BN together with
the clique impact when the variable xray is the output of interest. The impact
of the clique either, lung, tub is 0.93, which is simply the diameter of the
CPT of xray. The impact of the cliques asia, tub and lung, either, bronc is
also 0.93 since this is equal to the diameter of the CPT of xray times the
diameter of the CPT of either, which is equal to 1. The remaining cliques have
impact $0.93\cdot 1\cdot 0.09$: that is, the impact of bronc, lung, either
times the diameter of the CPT of lung.

Figure 8: The junction tree representing the asia BN where cliques are colored
by their impact on xray.

Notice that because the diameter is bounded by one, we have the following nice
property, confirming that cliques further away from the output have a smaller
effect.

######  Proposition 6.

The impact of a clique $C_{i}$ on an output variable $Y_{j}$ in clique $C_{j}$
is smaller than the impact of $C_{k}$ on $C_{j}$ for any $C_{k}$ along the
unique path between $C_{i}$ and $C_{j}$.

These influences provide a very useful tool for prioritisation of the
elicitation in a BN. If we can obtain estimates of influence across a junction
tree (either from direct elicitation or alternatively after having performed a
preliminary coarse elicitation of the corresponding CPTs) then we can use
these influences to identify which of those CPTs to refine. A further
discussion of how the results of this paper can be used for elicitation is
next.

##  6 The diameter as an elicitation tool

We have hinted throughout this paper at possible ways in which the diameter
can be used to evaluate the state of a BN elicitation from experts. These
evaluations can be embedded into a formal protocol. However, there are many
considerations that a user needs to make before undertaking model
construction: transparency of the model, computational issues, elicitation
constraints and so forth, which vary in importance depending on the context of
the model. So, setting a bound on any effects or perturbations against
differing approaches is often best undertaken informally. We acknowledge that
the framework we have presented here is sufficiently formal to admit
generalisation, and this is work that we plan to undertake in the future.

To implement our techniques as efficiently as possible, we recommend two
differing approaches tailored to the specific circumstances of the modeller.
Firstly, there are occasions in which we have obtained provisional information
from one expert who can recommend nodes, levels, interactions and provisional
CPTs before undertaking a more formal elicitation conference with multiple
experts (see e.g. Barons et al., 2018). In this particular scenario, we can
begin to design the analysis by using the bounds discussed earlier on the
preliminary values stated by the expert. We recommend starting by eliciting
attributes and nodes of interest before working systematically backwards along
the chain of inference to discover parent nodes and conditional independences,
performing variation measures on preliminary CPT values to determine the
efficacy of including variables in the model. Of course, after the complete
elicitation has taken place, the robustness analyses suggested above can be
repeated for a final sensitivity analysis.

In situations where we are starting the model with no such preliminary
information, it may be wisest to attempt to elicit the value of the diameter
of the CPT directly before eliciting the entire matrix so that complete
elicitation is not undertaken before we can derive concrete bounds on the
usefulness of this data harvesting exercise. To elicit the diameter directly,
we need to ascertain the largest differences between rows of a CPT, which
corresponds to requesting the “best case scenario” and “worst case scenario”
probabilities before calculating the variation distance between the two.

By following this procedure, we, therefore, continuously appraise and compare
each possible simplification against the potential accuracy of an analysis,
weighted against the issues provided by a simpler model representation. We
have demonstrated above that in many cases, the effects of various
simplifications are often minimal, and approximations based on these
simplifications are justified from a pragmatic point of view. We also note
that some of the best approximations frequently differ from those currently
undertaken in practice. For example, using an approximation that deletes an
edge can cause significant changes whilst allowing dependence only on subsets
of levels, which performs much better.

##  7 Discussion

We have demonstrated here how the properties of variation distance can be
harnessed to study the robustness of a discrete BN, if certain target
variables are known a priori to be those of primary interest. We took
advantage of the fact that the variation distance naturally embeds conditional
independence relationships between variables. This observation enabled us to
devise a seamless way of looking at perturbed versions of a BN for various
sensitivity problems, including edge strength, node influence, level
amalgamation and asymmetric dependence. An implementation of the developed
methods is freely available to practitioners through the bnmonitor R package.

Much work is still to be undertaken, starting with refining the bounds we have
developed here. Similarly, within this paper, we have had no space to consider
the robustness of the choice of probability distribution on the entries of the
CPTs of a BN. In Smith and Daneshkhah (2010), BN robustness associated with
the inputs of the distribution in terms of the local DeRobertis distance
(DeRoberts and Hartigan, 1981) is studied, while Smith and Rigat (2012)
provided bounds on posterior variation distances. Therefore, a reasonably
straightforward extension to the variation bounds we have presented here can
be developed by carefully combining our results with the local DeRobertis
distance to provide a comprehensive robustness analysis.

Our ideas also apply directly to dynamic BNs where the system’s sensitivity
can be even more critical because the dynamic nature of the problem makes the
model much more complex. Throughout this paper, for simplicity, we have
considered only robustness as it applies to finite categorical BNs. However,
using the approach developed here, the variation distance on these highly
structured and complex Markov processes can help us determine the robustness
of dynamic BNs to dynamic effects. Furthermore, the diameter could be extended
to perform sensitivity studies in mixed BNs, where both categorical and
continuous variables are considered.

## References

  * Albrecht et al. (2014) Albrecht, D., Nicholson, A.E., Whittle, C., 2014.  Structural sensitivity for the knowledge engineering of Bayesian networks, in: European Workshop on Probabilistic Graphical Models, Springer. pp. 1–16.
  * Ballester-Ripoll and Leonelli (2022a) Ballester-Ripoll, R., Leonelli, M., 2022a.  Computing Sobol indices in probabilistic graphical models.  Reliability Engineering & System Safety , 108573.
  * Ballester-Ripoll and Leonelli (2022b) Ballester-Ripoll, R., Leonelli, M., 2022b.  You only derive once (YODO): Automatic differentiation for efficient sensitivity analysis in Bayesian networks, in: International Conference on Probabilistic Graphical Models, PMLR. pp. 169–180.
  * Ballester-Ripoll and Leonelli (2023) Ballester-Ripoll, R., Leonelli, M., 2023.  The YODO algorithm: An efficient computational framework for sensitivity analysis in Bayesian networks.  International Journal of Approximate Reasoning 159, 108929.
  * Barons et al. (2018) Barons, M.J., Hanea, A.M., Wright, S.K., Baldock, K.C., Wilfert, L., Chandler, D., Datta, S., Fannon, J., Hartfield, C., Lucas, A., Ollerton, J., Potts, S.G., Carreck, N.L., 2018.  Assessment of the response of pollinator abundance to environmental pressures using structured expert elicitation.  Journal of Apicultural Research 57, 593–604.
  * Boneh et al. (2006) Boneh, T., Nicholson, A.E., Sonenberg, E.A., 2006.  Matilda: A visual tool for modeling with Bayesian networks.  International Journal of Intelligent Systems 21, 1127–1150.
  * Borgonovo (2017) Borgonovo, E., 2017.  Sensitivity analysis: An introduction for the management scientist.  Springer.
  * Borgonovo and Plischke (2016) Borgonovo, E., Plischke, E., 2016.  Sensitivity analysis: A review of recent advances.  European Journal of Operational Research 248, 869–887.
  * Boutilier et al. (1996) Boutilier, C., Friedman, N., Goldszmidt, M., Koller, D., 1996.  Context-specific independence in Bayesian networks, in: Proceedings of the 12th Conference on Uncertainty in Artificial Intelligence, pp. 115–123.
  * Carli et al. (2022) Carli, F., Leonelli, M., Riccomagno, E., Varando, G., 2022.  The R package stagedtrees for structural learning of stratified staged trees.  Journal of Statistical Software 102, 1–30.
  * Castillo et al. (1997) Castillo, E., Gutiérrez, J.M., Hadi, A.S., 1997.  Sensitivity analysis in discrete Bayesian networks.  IEEE Transactions on Systems, Man, and Cybernetics-Part A: Systems and Humans 27, 412–423.
  * Chan and Darwiche (2002) Chan, H., Darwiche, A., 2002.  When do numbers really matter?  Journal of Artificial Intelligence Research 17, 265–287.
  * Chan and Darwiche (2004) Chan, H., Darwiche, A., 2004.  Sensitivity analysis in Bayesian networks: From single to multiple parameters, in: Proceedings of the 20th Conference on Uncertainty in Artificial Intelligence, pp. 67–75.
  * Chan and Darwiche (2005) Chan, H., Darwiche, A., 2005.  A distance measure for bounding probabilistic belief change.  International Journal of Approximate Reasoning 38, 149–174.
  * Choi et al. (2005) Choi, A., Chan, H., Darwiche, A., 2005.  On Bayesian network approximation by edge deletion, in: Proceedings of the 21st Conference on Uncertainty in Artificial Intelligence, pp. 128–135.
  * Choi and Darwiche (2006) Choi, A., Darwiche, A., 2006.  An edge deletion semantics for belief propagation and its practical impact on approximation quality, in: Proceedings of the National Conference on Artificial Intelligence, p. 1107.
  * Collazo et al. (2018) Collazo, R.A., Görgen, C., Smith, J.Q., 2018.  Chain event graphs.  CRC Press.
  * Cooper (1990) Cooper, G.F., 1990.  The computational complexity of probabilistic inference using Bayesian belief networks.  Artificial Intelligence 42, 393–405.
  * Coupé et al. (2000) Coupé, V.M., Van Der Gaag, L.C., Habbema, J.D.F., 2000.  Sensitivity analysis: An aid for belief-network quantification.  The Knowledge Engineering Review 15, 215–232.
  * De Bock et al. (2023) De Bock, K.W., Coussement, K., De Caigny, A., Slowiński, R., Baesens, B., Boute, R.N., Choi, T.M., Delen, D., Kraus, M., Lessmann, S., et al., 2023.  Explainable AI for operational research: A defining framework, methods, applications, and a research agenda.  European Journal of Operational Research .
  * DeRoberts and Hartigan (1981) DeRoberts, L., Hartigan, J.A., 1981.  Bayesian inference using intervals of measures.  The Annals of Statistics , 235–244.
  * van Dorp (2020) van Dorp, J.R., 2020.  A dependent project evaluation and review technique: A Bayesian network approach.  European Journal of Operational Research 280, 689–706.
  * Eggeling et al. (2019) Eggeling, R., Grosse, I., Koivisto, M., 2019.  Algorithms for learning parsimonious context trees.  Machine Learning 108, 879–911.
  * French (2011) French, S., 2011.  Aggregating expert judgement.  Revista de la Real Academia de Ciencias Exactas, Fisicas y Naturales. Serie A. Matematicas 105, 181–206.
  * van der Gaag et al. (2007) van der Gaag, L.C., Renooij, S., Coupé, V.M., 2007.  Sensitivity analysis of probabilistic networks, in: Advances in Probabilistic Graphical Models. Springer, pp. 103–124.
  * Garvey et al. (2015) Garvey, M.D., Carnovale, S., Yeniyurt, S., 2015.  An analytical framework for supply network risk propagation: A Bayesian network approach.  European Journal of Operational Research 243, 618–627.
  * Gómez-Villegas et al. (2014) Gómez-Villegas, M.A., Main, P., Viviani, P., 2014.  Sensitivity to evidence in Gaussian Bayesian networks using mutual information.  Information Sciences 275, 115–126.
  * Hosseini (2021) Hosseini, S., 2021.  A decision support system based on machined learned Bayesian network for predicting successful direct sales marketing.  Journal of Management Analytics 8, 295–315.
  * ISTAT (2015) ISTAT, 2015.  Italian innovation survey 2010-2012.  , http://www.istat.it/en/archive/87787.
  * Kjærulff and van der Gaag (2000) Kjærulff, U., van der Gaag, L.C., 2000.  Making sensitivity analysis computationally efficient, in: Proceedings of the 16th conference on Uncertainty in Artificial Intelligence, pp. 317–325.
  * Kjaerulff and Madsen (2008) Kjaerulff, U.B., Madsen, A.L., 2008.  Bayesian networks and influence diagrams.  Springer Science+ Business Media 200, 114.
  * Koller and Friedman (2009) Koller, D., Friedman, N., 2009.  Probabilistic graphical models: principles and techniques.  MIT press.
  * Kwisthout and van der Gaag (2008) Kwisthout, J., van der Gaag, L.C., 2008.  The computational complexity of sensitivity analysis and parameter tuning, in: Proceedings of the 24th Conference on Uncertainty in Artificial Intelligence, pp. 349–356.
  * Laskey and Mahoney (2000) Laskey, K.B., Mahoney, S.M., 2000.  Network engineering for agile belief network models.  IEEE Transactions on Knowledge and Data Engineering 12, 487–498.
  * Leonelli et al. (2017) Leonelli, M., Görgen, C., Smith, J.Q., 2017.  Sensitivity analysis in multilinear probabilistic models.  Information Sciences 411, 84–97.
  * Leonelli et al. (2023) Leonelli, M., Ramanathan, R., Wilkerson, R.L., 2023.  Sensitivity and robustness analysis in Bayesian networks with the bnmonitor R package.  Knowledge-Based Systems 278, 110882.
  * Leonelli and Riccomagno (2022) Leonelli, M., Riccomagno, E., 2022.  A geometric characterization of sensitivity analysis in monomial models.  International Journal of Approximate Reasoning 151, 64–84.
  * Leonelli and Varando (2024) Leonelli, M., Varando, G., 2024.  Structural learning of simple staged trees.  Data Mining and Knowledge Discovery 38, 1520–1544.
  * Pensar et al. (2015) Pensar, J., Nyman, H., Koski, T., Corander, J., 2015.  Labeled directed acyclic graphs: A generalization of context-specific independence in directed graphical models.  Data Mining and Knowledge Discovery 29, 503–533.
  * Pensar et al. (2016) Pensar, J., Nyman, H., Lintusaari, J., Corander, J., 2016.  The role of local partial independence in learning of Bayesian networks.  International Journal of Approximate Reasoning 69, 91–105.
  * Qazi (2022) Qazi, A., 2022.  Adoption of a probabilistic network model investigating country risk drivers that influence logistics performance indicators.  Environmental Impact Assessment Review 94, 106760.
  * Razavi et al. (2021) Razavi, S., Jakeman, A., Saltelli, A., Prieur, C., Iooss, B., Borgonovo, E., Plischke, E., Piano, S.L., Iwanaga, T., Becker, W., et al., 2021.  The future of sensitivity analysis: An essential discipline for systems modeling and policy support.  Environmental Modelling & Software 137, 104954.
  * Renooij (2001) Renooij, S., 2001.  Probability elicitation for belief networks: Issues to consider.  The Knowledge Engineering Review 16, 255–269.
  * Renooij (2010) Renooij, S., 2010.  Bayesian network sensitivity to arc-removal, in: Proceedings of the 5th European Workshop on Probabilistic Graphical Models, pp. 233–240.
  * Renooij (2014) Renooij, S., 2014.  Co-variation for sensitivity analysis in Bayesian networks: Properties, consequences and alternatives.  International Journal of Approximate Reasoning 55, 1022–1042.
  * Roberts and Rosenthal (2004) Roberts, G.O., Rosenthal, J.S., 2004.  General state space Markov chains and MCMC algorithms.  Probability Surveys 1, 20–71.
  * Rohmer (2020) Rohmer, J., 2020.  Uncertainties in conditional probability tables of discrete Bayesian belief networks: A comprehensive review.  Engineering Applications of Artificial Intelligence 88, 103384.
  * Rudin (2019) Rudin, C., 2019.  Stop explaining black box machine learning models for high stakes decisions and use interpretable models instead.  Nature Machine Intelligence 1, 206–215.
  * Salmani and Katoen (2023) Salmani, B., Katoen, J.P., 2023.  Automatically finding the right probabilities in Bayesian networks.  Journal of Artificial Intelligence Research 77, 1637–1696.
  * Saltelli et al. (2021) Saltelli, A., Jakeman, A., Razavi, S., Wu, Q., 2021.  Sensitivity analysis: A discipline coming of age.  Environmental Modelling & Software 146, 105226.
  * Saltelli et al. (2000) Saltelli, A., Tarantola, S., Campolongo, F., 2000.  Sensitivity analysis as an ingredient of modeling.  Statistical Science 15, 377–395.
  * Sason and Verdú (2016) Sason, I., Verdú, S., 2016.  $f$-divergence inequalities.  IEEE Transactions on Information Theory 62, 5973–6006.
  * Scutari (2010) Scutari, M., 2010.  Learning Bayesian networks with the bnlearn R package.  Journal of Statistical Software 35, 1–22.
  * Scutari et al. (2019) Scutari, M., Graafland, C.E., Gutiérrez, J.M., 2019.  Who learns better Bayesian network structures: Accuracy and speed of structure learning algorithms.  International Journal of Approximate Reasoning 115, 235–253.
  * Scutari and Nagarajan (2013) Scutari, M., Nagarajan, R., 2013.  Identifying significant edges in graphical models of molecular networks.  Artificial intelligence in medicine 57, 207–217.
  * Smith (2010) Smith, J.Q., 2010.  Bayesian decision analysis: Principles and practice.  Cambridge University Press.
  * Smith and Anderson (2008) Smith, J.Q., Anderson, P.E., 2008.  Conditional independence and chain event graphs.  Artificial Intelligence 172, 42–68.
  * Smith and Daneshkhah (2010) Smith, J.Q., Daneshkhah, A., 2010.  On the robustness of Bayesian networks to learning from non-conjugate sampling.  International Journal of Approximate Reasoning 51, 558–572.
  * Smith and Rigat (2012) Smith, J.Q., Rigat, F., 2012.  Isoseparation and robustness in parametric Bayesian inference.  Annals of the Institute of Statistical Mathematics 64, 495–519.
  * Varando et al. (2024) Varando, G., Carli, F., Leonelli, M., 2024.  Staged trees and asymmetry-labeled DAGs.  Metrika , 1–28.
  * Werner et al. (2017) Werner, C., Bedford, T., Cooke, R.M., Hanea, A.M., Morales-Napoles, O., 2017.  Expert judgement for dependence in probabilistic modelling: A systematic literature review and future research directions.  European Journal of Operational Research 258, 801–819.
  * Wilkerson and Smith (2021) Wilkerson, R.L., Smith, J.Q., 2021.  Customized structural elicitation, in: Expert Judgement in Risk and Decision Analysis. Springer, pp. 83–113.
  * Wright (2018) Wright, S., 2018.  Robustness in Bayesian Networks.  Ph.D. thesis. Department of Statistics, University of Warwick.

##  Appendix A Auxiliary results

In this appendix we present results that will be used in the proofs of the
main propositions and theorems of the paper.

######  Lemma 1.

Let $p$ and $p^{\prime}$ be two pmfs over the same sample space $\mathbb{X}$.
Then

$\sum_{x\in\mathbb{X}}\min\\{p(x),p^{\prime}(x)\\}=1-d_{V}(p,p^{\prime})$ (5)

###### Proof.

By definition of total variation distance we have

$\sum_{x\in\mathbb{X}}\min\\{p(x),p^{\prime}(x)\\}+d_{V}(p,p^{\prime})=\sum_{x\in\mathbb{X}}\min\\{p(x),p^{\prime}(x)\\}\\\
+\frac{1}{2}\sum_{x\in\mathbb{X}}|p(x)-p^{\prime}(x)|.$

Let $\mathbb{X}_{1}$ be the subset of $\mathbb{X}$ including those $x$ such
that $p(x)\leq p^{\prime}(x)$ and
$\mathbb{X}_{2}=\mathbb{X}\setminus\mathbb{X}_{1}$. Therefore

$\displaystyle\sum_{x\in\mathbb{X}}\min\\{p(x),p^{\prime}(x)\\}+d_{V}(p,p^{\prime})=$
$\displaystyle\sum_{x\in\mathbb{X}_{1}}p(x)+\sum_{x\in\mathbb{X}_{2}}p^{\prime}(x)+$
$\displaystyle\frac{1}{2}\sum_{x\in\mathbb{X}_{1}}(p^{\prime}(x)-p(x))+\frac{1}{2}\sum_{x\in\mathbb{X}_{2}}(p(x)-p^{\prime}(x))$
$\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{x\in\mathbb{X}_{1}}p(x)+\frac{1}{2}\sum_{x\in\mathbb{X}_{1}}p^{\prime}(x)+$
$\displaystyle\frac{1}{2}\sum_{x\in\mathbb{X}_{2}}p^{\prime}(x)+\frac{1}{2}\sum_{x\in\mathbb{X}_{2}}p(x)$
$\displaystyle=$ $\displaystyle\frac{1}{2}+\frac{1}{2}=1$

∎

######  Lemma 2.

Let $p$ and $p^{\prime}$ be two pmfs over the same sample space $\mathbb{X}$.
Call

$\displaystyle\beta$ $\displaystyle=1-d_{V}(p,p^{\prime})$ $\displaystyle
p^{*}$ $\displaystyle=\frac{\min\\{p,p^{\prime}\\}}{\beta}$
$\displaystyle\bar{p}$
$\displaystyle=\frac{p-\min\\{p,p^{\prime}\\}}{1-\beta}$
$\displaystyle\bar{p}^{\prime}$
$\displaystyle=\frac{p^{\prime}-\min\\{p,p^{\prime}\\}}{1-\beta}$

Then

$p=\beta p^{*}+(1-\beta)\bar{p}$ (6)

and similarly for $p^{\prime}$

###### Proof.

First notice that from Lemma 1, both $p^{*}$ and $\bar{p}$ are pmfs. For
$p^{*}$ this follows directly from (5). For $\bar{p}$ notice that

$\displaystyle\sum_{x\in\mathbb{X}}\frac{p(x)-\min\\{p(x),p^{\prime}(x)\\}}{1-\beta}$
$\displaystyle=\frac{1}{d_{V}(p,p^{\prime})}\sum_{x\in\mathbb{X}}\left(p(x)-\min\\{p(x),p^{\prime}(x)\\}\right)$
$\displaystyle=\frac{1}{d_{V}(p,p^{\prime})}(1-(1-d_{V}(p,p^{\prime})))=1.$

By substituting the definitions of $p^{*}$ and $\bar{p}$ into (6) the result
follows. ∎

######  Lemma 3.

Let $p,q_{1},\dots,q_{|\mathbb{X}|}$ be pmfs over the same sample space
$\mathbb{X}$. Let
$\pi=(\pi_{1},\dots,\pi_{|\mathbb{X}|})\in\Delta_{|\mathbb{X}|-1}$ and
$q_{\pi}=\sum_{i=1}^{|\mathbb{X}|}\pi_{i}q_{i}$. Then

$d_{V}(p,q_{\pi})\leq\sum_{i=1}^{|\mathbb{X}|}\pi_{i}d_{V}(p,q_{i})$

###### Proof.

$\displaystyle d_{V}(p,q_{\pi})$
$\displaystyle=\frac{1}{2}\sum_{x\in\mathbb{X}}|p(x)-q_{\pi}(x)|=\frac{1}{2}\sum_{x\in\mathbb{X}}\left|p(x)-\sum_{i=1}^{|\mathbb{X}|}\pi_{i}q_{i}(x)\right|$
$\displaystyle=\frac{1}{2}\sum_{x\in\mathbb{X}}\left|\sum_{i=1}^{|\mathbb{X}|}\pi_{i}p(x)-\sum_{i=1}^{|\mathbb{X}|}\pi_{i}q_{i}(x)\right|=\frac{1}{2}\sum_{x\in\mathbb{X}}\left|\sum_{i=1}^{|\mathbb{X}|}\pi_{i}(p(x)-q_{i}(x))\right|$
$\displaystyle\leq\frac{1}{2}\sum_{x\in\mathbb{X}}\sum_{i=1}^{|\mathbb{X}|}\pi_{i}|p(x)-q_{i}(x)|=\sum_{i=1}^{|\mathbb{X}|}\pi_{i}\frac{1}{2}\sum_{x\in\mathbb{X}}|p(x)-q_{i}(x)|$
$\displaystyle=\sum_{i=1}^{|\mathbb{X}|}\pi_{i}d_{V}(p,q_{i})$

∎

For the next lemma we consider two categorical random variables $X$ and $Y$
with sample spaces $\mathbb{X}$, $\mathbb{Y}$, respectively. With $p(X)$,
$p(X,Y)$ and $p(Y|x)$ we denote the marginal, joint and conditional
probability distributions over the corresponding variables. The same notation
is used for another pmf $p^{\prime}$.

######  Lemma 4.

It holds that

$d_{V}(p(X,Y),p^{\prime}(X,Y))\leq\min\\{d_{V}(p(X),p^{\prime}(X))+\max_{x\in\mathbb{X}}d_{V}(p(Y|x),p^{\prime}(Y|x)),1\\}$

###### Proof.

First notice that

$\displaystyle d_{V}(p(X,Y),p^{\prime}(X,Y))$
$\displaystyle=\frac{1}{2}\sum_{x\in\mathbb{X},y\in\mathbb{Y}}|p(x,y)-p^{\prime}(x,y)|$
$\displaystyle=\frac{1}{2}\sum_{x\in\mathbb{X},y\in\mathbb{Y}}|p(y|x)p(x)-p^{\prime}(y|x)p(x)|.$

Calling $r(y|x)=p(y|x)-p^{\prime}(y|x)$, we have

$\displaystyle d_{V}(p(X,Y),p^{\prime}(X,Y))=$ $\displaystyle\hskip
8.5359pt\frac{1}{2}\sum_{x\in\mathbb{X},y\in\mathbb{Y}}|p(y|x)(p(x)-p^{\prime}(x))+r(y|x)p^{\prime}(x)|$
$\displaystyle\leq$ $\displaystyle\hskip
8.5359pt\frac{1}{2}\sum_{\begin{subarray}{c}x\in\mathbb{X}\\\
y\in\mathbb{Y}\end{subarray}}|p(y|x)(p(x)-p^{\prime}(x))|+\frac{1}{2}\sum_{\begin{subarray}{c}x\in\mathbb{X}\\\
y\in\mathbb{Y}\end{subarray}}|r(y|x)p^{\prime}(x)|$ $\displaystyle=$
$\displaystyle\hskip
8.5359pt\frac{1}{2}\sum_{x\in\mathbb{X}}|p(x)-p^{\prime}(x)|\sum_{y\in\mathbb{Y}}p(y|x)+$
$\displaystyle\hskip
8.5359pt\frac{1}{2}\sum_{x\in\mathbb{X}}p^{\prime}(x)\sum_{y\in\mathbb{Y}}|r(y|x)|$
$\displaystyle\leq$ $\displaystyle\hskip
8.5359pt\frac{1}{2}\sum_{x\in\mathbb{X}}|p(x)-p^{\prime}(x)|+\frac{1}{2}\max_{x\in\mathbb{X}}\sum_{y\in\mathbb{Y}}|r(y|x)|$
$\displaystyle=$ $\displaystyle\hskip
8.5359ptd_{V}(p(X),p^{\prime}(X))+\max_{x\in\mathbb{X}}d_{V}(p(Y|x),p^{\prime}(Y|x)).$

Since the maximum value for the total variation distance is one, the result
follows. ∎

##  Appendix B Proofs

###  B.1 Proof of Proposition 1

For the $(\Rightarrow)$ statement, we notice that if $d^{+}(P_{Y|X})=0$ then
any two rows of $P_{Y|X}$ must be the same. Thus $p(Y|x)=p(Y)$ for all
$x\in\mathbb{X}$. Thus $Y\mathchoice{\mathrel{\hbox
to0.0pt{$\displaystyle\perp$\hss}\mkern
2.0mu{\displaystyle\perp}}}{\mathrel{\hbox
to0.0pt{$\textstyle\perp$\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptstyle\perp$\hss}\mkern
2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptscriptstyle\perp$\hss}\mkern
2.0mu{\scriptscriptstyle\perp}}}X$. For the $(\Leftarrow)$ statement, if
$Y\mathchoice{\mathrel{\hbox to0.0pt{$\displaystyle\perp$\hss}\mkern
2.0mu{\displaystyle\perp}}}{\mathrel{\hbox
to0.0pt{$\textstyle\perp$\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptstyle\perp$\hss}\mkern
2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox
to0.0pt{$\scriptscriptstyle\perp$\hss}\mkern
2.0mu{\scriptscriptstyle\perp}}}X$ then $p(Y|x)=p(Y)$ for all
$x\in\mathbb{X}$. Thus every row of $P_{Y|X}$ is equal to $p(Y)$ and its
diameter is zero.

###  B.2 Proof of Proposition 2

Consider any two $x,x^{\prime}\in\mathbb{X}$. Using twice Lemma 3 we have that

$\displaystyle d_{V}(p(Y|x),p(Y|x^{\prime}))$
$\displaystyle=d_{V}\left(\sum_{z\in\mathbb{Z}}p(Y|x,z)p(z|x),\sum_{z^{\prime}\in\mathbb{Z}}p(Y|x^{\prime},z^{\prime})p(z^{\prime}|x^{\prime})\right)$
$\displaystyle\leq\sum_{z^{\prime}\in\mathbb{Z}}p(z^{\prime}|x^{\prime})d_{V}\left(\sum_{z\in\mathbb{Z}}p(Y|x,z)p(z|x),p(Y|x^{\prime},z^{\prime})\right)$
$\displaystyle\leq\sum_{z^{\prime}\in\mathbb{Z}}p(z^{\prime}|x^{\prime})\sum_{z\in\mathbb{Z}}p(z|x)d_{V}\left(p(Y|x,z),p(Y|x^{\prime},z^{\prime})\right)$
$\displaystyle\leq
d^{+}(P_{Y|(x,x^{\prime})Z})\sum_{z^{\prime}\in\mathbb{Z}}p(z^{\prime}|x^{\prime})\sum_{z\in\mathbb{Z}}p(z|x)$
$\displaystyle=d^{+}(P_{Y|(x,x^{\prime})Z}),$

where $P_{Y|(x,x^{\prime})Z}$ is the CPT of $Y$ conditional on all
$z\in\mathbb{Z}$ but only on $x,x^{\prime}\in\mathbb{X}$. The result follows
since

$d^{+}(P_{Y|X})=\max_{x,x^{\prime}\in\mathbb{X}}d_{V}(p(Y|x),p(Y|x^{\prime}))\leq\max_{x,x^{\prime}\in\mathbb{X}}d^{+}(P_{Y|(x,x^{\prime})Z})\leq
d^{+}(P_{Y|XZ}).$

The proof about the equality between diameters follows a similar reasoning to
that of Proposition 1.

###  B.3 Proof of Proposition 3

For any two $x,x^{\prime}\in\mathbb{X}$, Lemma 4 can be rewritten as

$d_{V}(p(Y,Z|x),p^{\prime}(Y,Z|x^{\prime}))\leq
d_{V}(p(Z|x),p^{\prime}(Z|x^{\prime}))+\\\
\max_{z\in\mathbb{Z}}d_{V}(p(Y|z,x),p^{\prime}(Y|z,x^{\prime})).$

Notice that
$\max_{z\in\mathbb{Z}}d_{V}(p(Y|z,x),p^{\prime}(Y|z,x^{\prime}))\leq
d^{+}(P_{Y|ZX})$ and therefore

$d^{+}(P_{YZ|X})=\max_{x,x^{\prime}\in\mathbb{X}}d_{V}(p(Y,Z|x),p^{\prime}(Y,Z|x^{\prime}))\\\
\leq
d^{+}(P_{Y|ZX})+\max_{x,x^{\prime}\in\mathbb{X}}d_{V}(p(Z|x),p^{\prime}(Z|x^{\prime})).$

Noticing that
$\max_{x,x^{\prime}\in\mathbb{X}}d_{V}(p(Z|x),p^{\prime}(Z|x^{\prime}))=d^{+}(P_{Z|X})$,
the result follows.

###  B.4 Proof of Theorem 1

Consider the total variation distance between $p(Y)$ and $p^{\prime}(Y)$.
Using the definition of conditional probability this is equal to

$d_{V}(p(Y),p^{\prime}(Y))=d_{V}(p(X)P_{Y|X},p^{\prime}(X)P_{Y|X})$

Let $P_{Y|X}(x)$ be the $x$-th row of $P_{Y|X}$ and notice that
$p(X)P_{Y|X}=\sum_{x\in\mathbb{X}}p(x)P_{Y|X}(x)$. Using Lemma 2 we now derive
that

$\displaystyle d_{V}(p(Y),p^{\prime}(Y))$ $\displaystyle=d_{V}((\beta
p^{*}+(1-\beta)\bar{p})P_{Y|X},(\beta
p^{*}+(1-\beta)\bar{p}^{\prime})P_{Y|X})$
$\displaystyle=\frac{1}{2}\sum_{x\in\mathbb{X}}|(1-\beta)(\bar{p}(x)-\bar{p}^{\prime}(x))P_{Y|X}(x)|$
$\displaystyle=\frac{1}{2}\sum_{x\in\mathbb{X}}|1-\beta||(\bar{p}(x)-\bar{p}^{\prime}(x))P_{Y|X}(x)|$
$\displaystyle=|1-\beta|d_{V}(\bar{p}P_{Y|X},\bar{p}^{\prime}P_{Y|X})$
$\displaystyle=(1-\beta)d_{V}(\bar{p}P_{Y|X},\bar{p}^{\prime}P_{Y|X}),$

where the last equality follows by the fact that the total variation distance
is at most 1.

Using again the fact that $p(X)P_{Y|X}=\sum_{x\in\mathbb{X}}p(x)P_{Y|X}(x)$
together with Lemma 3, we get that

$\displaystyle d_{V}(p(Y),p^{\prime}(Y))$
$\displaystyle=(1-\beta)d_{V}(\bar{p}P_{Y|X},\bar{p}^{\prime}P_{Y|X})$
$\displaystyle=(1-\beta)d_{V}\left(\sum_{x\in\mathbb{X}}\bar{p}(x)P_{Y|X}(x),\sum_{x^{\prime}\in\mathbb{X}}\bar{p}^{\prime}(x^{\prime})P_{Y|X}(x^{\prime})\right)$
$\displaystyle\leq(1-\beta)\sum_{x^{\prime}\in\mathbb{X}}\bar{p}^{\prime}(x^{\prime})d_{V}\left(\sum_{x\in\mathbb{X}}\bar{p}(x)P_{Y|X}(x),P_{Y|X}(x^{\prime})\right)$
$\displaystyle\leq(1-\beta)\sum_{x^{\prime}\in\mathbb{X}}\sum_{x\in\mathbb{X}}\bar{p}^{\prime}(x^{\prime})\bar{p}(x)d_{V}\left(P_{Y|X}(x^{\prime}),P_{Y|X}(x)\right)$
$\displaystyle\leq(1-\beta)\sum_{x^{\prime}\in\mathbb{X}}\sum_{x\in\mathbb{X}}\bar{p}^{\prime}(x^{\prime})\bar{p}(x)d^{+}_{Y|X},$

where in the last inequality we used the definition of diameter of a CPT.
Using the fact that
$\sum_{x^{\prime}\in\mathbb{X}}\sum_{x\in\mathbb{X}}\bar{p}^{\prime}(x^{\prime})\bar{p}(x)=1$
and that $1-\beta=d_{V}(p(X),p^{\prime}(X))$, the result follows.

##  Appendix C Refining a Bayesian networks with asymmetric independence

In this appendix we showcase the use of a context-specific class of models,
namely staged trees (Collazo et al., 2018; Smith and Anderson, 2008), to learn
conditional independence relationships from data that are not simply
symmetric. Here we simply give an illustration of how a staged tree would
represent the structure of the CPT of the variable EMPUD in the istat BN. With
this aim we construct an event tree with the parents of EMPUD, namely GP and
EMP12, and EMPUD itself, reported in Figure 9.

Each vertex at depth two in this tree represents the conditional distribution
of EMPUD for a specific combination of its parents: it thus represents a row
of the CPT of EMPUD. A staged tree learning algorithm then aims at finding the
best partition of these vertices, where two vertices are in the same set if
their pmfs are assumed to be equal. When this happens two vertices are said to
be in the same stage and are equally colored in the tree graph. We used the
stagedtrees R package (Carli et al., 2022) to learn an optimal partition of
the vertices associated to the pmf of EMPUD, reported in the coloring of the
tree in Figure 9. It can be noticed that the only two vertices that are given
the same color are for the context GP = No and two levels of EMP12. This
correspond to the setup where the index of partial independence was really low
in Table 9. The fact that no other pair of vertices have the same color
suggests that there are no other equalities between the rows of the CPTs of
EMPUD, confirming what suggested by the indices in Table 9.

GP = NoGP = YesEMP12 = 10-49EMP12 = 50-249EMP12 $>$250EMP12 = 10-49EMP12 =
50-249EMP12 $>$250EMPUD = 0%EMPUD $<$10%EMPUD = YesEMPUD = 0%EMPUD $<$10%EMPUD
$>$10%EMPUD = 0%EMPUD $<$10%EMPUD $>$10%EMPUD = 0%EMPUD $<$10%EMPUD
$>$10%EMPUD = 0%EMPUD $<$10%EMPUD $>$10%EMPUD = 0%EMPUD $<$10%EMPUD $>$10%

Figure 9: Staged tree for the variable EMPUD conditional on its parents in the
istat BN.