Bayesian Networks - Learning from Data

Learning a Bayesian network can be split into two problems:

1. Parameter learning: Given a set of data samples and a DAG that captures the dependencies between the variables, estimate the (conditional) probability distributions of the individual variables.

2. Structure learning: Given a set of data samples, estimate a DAG that captures the dependencies between the variables.

Conditional Independence¶

In a Bayesian Network, we represent observed variables by shaded nodes (i.e. $\textcolor{dodgerblue}{Ⓑ}$, all other variables are understood as unobserved (hidden, i.e. $Ⓐ$).

Three elementary patterns appear repeatedly in Bayesian networks:

Chain (serial connection)¶

Fork (diverging connection)¶

Collider (converging connection)¶

d-Separation (directional separation)¶

Think of paths as pipes and statistical dependence as water:

• If at least one path is open, dependence can flow.

• If all paths are blocked, then $A$ and $C$ are conditionally independent given $Z$.

Testing Conditional Independence in Data¶

Given a candidate causal DAG $G$ that could have generated a dataset $D$:

1. Use d-Separation to list all implied conditional independencies (e.g. $X \perp Z \mid Y$).

2. For each implied independence, test on data:

   • Use regression or other statistical tools,

   • Check whether the conditional dependence is approximately zero.

For example for the following graph $G$:

• Suppose $G$ implies $W \perp Z_1 \mid X$.

• Regress $W$ on $X$ and $Z_1$ using $w = \beta_X x + \beta_{Z_1} z_1$

• If $\beta_{Z_1} \approx 0$ (statistically insignificant), this is consistent with the model.

• If $\beta_{Z_1}$ is clearly non-zero, the independence fails in data, so $G$ must be rejected.

Markov equivalence and equivalence classes¶

Many different DAGs can induce the same set of conditional independencies. These DAGs are Markov equivalent.

Factorization examples¶

Consider three variables $X,Y,Z$ where we know that

Several factorizations of $P(X,Y,Z)$ are consistent with this:

• One possible factorization:

  Loading...

  If $P(y \mid x,z)$ does not depend on $x$, we can write this as

  Loading...

• Another factorization:

  Loading...

Different DAGs correspond to different factorization orders, but if they encode the same conditional independencies ($X \perp Y \mid Z$), they belong to the same Markov equivalence class.

Parameter sharing¶

In many models, multiple nodes share the same local conditional distribution. This is called parameter sharing.

Consider a chain of four binary variables:

Estimating the full joint distribution

• Full joint distribution $P(X,Y,Z,W)$ has $2^4 = 16$ entries (15 degrees of freedom)

• With limited data, many combinations $(x,y,z,w)$ are rare or absent

• Frequency estimates for these cells are unreliable

Using factorization

• Each factor involves at most 2 variables

• We estimate each conditional with more data per parameter

• This yields more stable estimates from the same dataset

$\implies$ Exploiting the graph structure dramatically reduces the number of parameters and improves estimation in high dimensions.

Hidden Markov Model (HMM)¶

Probabilistic Inference¶

Input:

• A Bayesian network $P(X_1,\dots,X_n)$,

• Evidence $E=e$ for some variables $E \subseteq X$,

• Query variables $Q \subseteq X$,

Output:

The output is a table that specifies a probability for each assignment of values to Q.

Parameter learning with fully observed data¶

Assume:

• The BN structure (graph) is known,

• Training data $D_{\text{train}}$ consists of full assignments to all variables $X_1,\dots,X_n$,

• We want to estimate parameters $\vartheta$ (local conditional probabilities).

The maximum likelihood objective is

Taking logs converts the product into a sum, and the problem decouples across variables and parent configurations.

For each node $X_i$ and parent assignment $x_{\text{Parents}(i)}$, we solve:

subject to $\sum_{x_i} p(x_i \mid x_{\text{Parents}(i)}) = 1$.

Solution:

Maximum likelihood via count-and-normalize¶

Using pgmpy for MLE¶

Suppose we have the following data:

We know that the variables relate as follows:

Parameter learning is the task to estimate the values of the conditional probability distributions (CPDs), for the variables fruit, size, and tasty.

State counts¶

To make sense of the given data, we can start by counting how often each state of the variable occurs. If the variable is dependent on parents, the counts are done conditionally on the parents states, i.e. separately for each parent configuration:

We can see, for example, that as many apples as bananas were observed and that 5 large bananas were tasty, while only 1 was not.

Maximum Likelihood Estimation¶

A natural estimate for the CPDs is to simply use the relative frequencies, with which the variable states have occured. We observed 7 apples among a total of 14 fruits, so we might guess that about 50% of fruits are apples.

This approach is Maximum Likelihood Estimation (MLE). According to MLE, we should fill the CPDs in such a way, that $P(\text{data}|\text{model})$ is maximal. This is achieved when using the relative frequencies. pgmpy supports MLE as follows:

• mle.estimate_cpd(variable) computes the state counts and divides each cell by the (conditional) sample size.

• The mle.get_parameters()-method returns a list of CPDs for all variables of the model.

• The built-in fit()-method of BayesianNetwork provides more convenient access to parameter estimators:

While very straightforward, the ML estimator has the problem of overfitting to the data. In above CPD, the probability of a large banana being tasty is estimated at 0.833, because 5 out of 6 observed large bananas were tasty. Fine. But note that the probability of a small banana being tasty is estimated at 0.0, because we observed only one small banana and it happened to be not tasty. But that should hardly make us certain that small bananas aren’t tasty!

We simply do not have enough observations to rely on the observed frequencies. If the observed data is not representative for the underlying distribution, ML estimations will be extremly far off.

Laplace smoothing and Dirichlet priors¶

Pure maximum likelihood can assign zero probability to unseen events:

• If some $(x_i,x_{\text{Parents}(i)})$ combination never appears in training data, then

  $$\hat{p}(x_i \mid x_{\text{Parents}(i)}) = 0.$$

  (14)

To avoid this, we use Laplace smoothing:

• Add $\lambda > 0$ to each count before normalization.

For discrete $X_i$:

Interpretation

• Equivalent to a Dirichlet prior on the conditional distribution and doing MAP estimation.

• Larger $\lambda$ produces stronger smoothing (probabilities closer to uniform).

• As we collect more data, the posterior concentrates around the empirical frequencies.

Smoothing in Python¶

The Bayesian Parameter Estimator starts with already existing prior CPDs, that express our beliefs about the variables before the data was observed. Those “priors” are then updated, using the state counts from the observed data.

One can think of the priors as consisting in pseudo state counts, that are added to the actual counts before normalization. Unless one wants to encode specific beliefs about the distributions of the variables, one commonly chooses uniform priors, i.e. ones that deem all states equiprobable.

A very simple prior is the so-called K2 prior (Laplace smoothing), which simply adds 1 to the count of every single state. A somewhat more sensible choice of prior is BDeu (Bayesian Dirichlet equivalent uniform prior). For BDeu we need to specify an equivalent sample size N and then the pseudo-counts are the equivalent of having observed N uniform samples of each variable (and each parent configuration):

The estimated values in the CPDs are now more conservative. In particular, the estimate for a small banana being not tasty is now around 0.64 rather than 1.0. Setting equivalent_sample_size to 10 means that for each parent configuration, we add the equivalent of 10 uniform samples (here: +5 small bananas that are tasty and +5 that aren’t).

BayesianEstimator, too, can be used via the fit()-method:

Learning with hidden or missing variables¶

In many applications, some variables are unobserved in the training data.

Example:

• In a movie-rating model with variables $(G,R_1,R_2)$, the genre $G$ might be hidden:

  • Observed: ratings $(R_1,R_2)$,

  • Hidden: $G$.

We then aim to maximize the marginal likelihood over observed variables:

where

sums over hidden variables $H$.

Direct optimization is hard because:

• The sum over $h$ is inside the product over examples,

• There is no simple closed-form solution like in the fully observed case.

This motivates the EM algorithm.

Expectation–Maximization (EM) algorithm¶

EM solves maximum marginal likelihood by iteratively constructing pseudo-complete datasets. It comes with the following guarantees:

• The marginal log-likelihood does not decrease at each iteration

• Converges to a local optimum (not necessarily global)

EM in Python¶

For this example, we simulate some data from the alarm model and use it to learn back the model parameters. In this example, we simply use the structure to the alarm model.

The Expectation Maximization (EM) estimator can work in the case when latent variables are present in the model. To simulate this scenario, we will specify some of the variables in our new_model as latents and drop those variables from samples to simulate missing data.

Markov networks and factor graphs¶

A Markov network (also called a Markov random field) is an undirected graphical model that represents a joint distribution over $X = (X_1,\dots,X_n)$ via factors.

We define:

• A set of factors $f_1,\dots,f_m$, each depending on a subset of variables.

• For any assignment $x$, a weight

  $$\text{Weight}(x) = \prod_{j=1}^m f_j(x).$$

  (20)

The normalized joint distribution is

where $Z$ is the partition function.

A factor graph is a bipartite graph with:

• Variable nodes $X_i$,

• Factor nodes $f_j$,

• Edges connecting factors to the variables they depend on.

Factor graphs generalize Bayesian networks and are convenient for describing models like object tracking, image models, and general constraint systems.

Object tracking: MAP trajectory vs marginal probabilities¶

In the object tracking example:

• $X_1,X_2,X_3$ are discrete positions,

• Observation factors $o_i(x_i)$ reflect sensor likelihoods,

• Transition factors $t_i(x_i,x_{i+1})$ enforce smooth motion.

The factor graph defines the weight

Two different tasks:

1. MAP assignment (most likely trajectory):

   $$(x_1^{*},x_2^{*},x_3^{*}) = \arg\max_{x_1,x_2,x_3} \text{Weight}(x_1,x_2,x_3).$$

   (23)

2. Marginal probabilities:

   $$P(X_2 = v) = \sum_{x_1,x_3} P(X_1=x_1,X_2=v,X_3=x_3) = \frac{1}{Z} \sum_{x_1,x_3} \text{Weight}(x_1,v,x_3).$$

   (24)

Lesson:

• The most probable trajectory may place the object at one position at time 2,

• While the marginal $P(X_2=v)$ may be higher for another position when summing over all trajectories.

Thus, returning only the MAP trajectory does not convey uncertainty about individual positions.

Gibbs sampling for approximate marginals¶

Computing exact marginals

is often infeasible in large models. Gibbs sampling approximates these marginals via a Markov chain.

Algorithm (high level):

1. Initialize a complete assignment $x^{(0)}$.

2. For $t = 0,1,2,\dots$:

   • Choose a variable index $i$,

   • Sample a new value for $X_i$ from the conditional

     $$X_i^{(t+1)} \sim P\big(X_i \mid X_{-i} = x_{-i}^{(t)}\big),$$

     (26)

     where $X_{-i}$ denotes all variables except $X_i$.

3. After a burn-in period, collect samples $x^{(t)}$ and estimate marginals by frequencies:

   $$\hat{P}(X_i = v) = \frac{\#\{t: X_i^{(t)} = v\}}{\text{number of collected samples}}.$$

   (27)

Key properties:

• Gibbs sampling is a form of Markov Chain Monte Carlo (MCMC).

• Under mild conditions (irreducibility, aperiodicity), the chain converges to the true joint distribution.

• It only needs local computations (looking at factors involving $X_i$), making it scalable to large models.

The slides illustrate Gibbs sampling on the object tracking / Markov network example.