10. Dynamic Bayesian Network, Markov Chain

Let’s see how we can represent a Markov Chain (MC) as a Dynamic Bayesian Network (DBN). We will verify our results with the stationary distribution or steady state.

10.1. Model

The model, \(P\), has two states, sunny and rainy. Our initial state, \(s\), will be sunny.

import numpy as np

P = np.array([[0.9, 0.1], [0.5, 0.5]])
s = np.array([1.0, 0.0])

10.2. Power Method

We can use the Power Method to find the steady state.

\(s^n = s^0P^n\)

We will continue to compute \(s^n\) until \(||P^n - P^{n-1}|| \lt \delta\).

s.dot(P)

array([0.9, 0.1])

s.dot(P).dot(P)

array([0.86, 0.14])

import pandas as pd

P_curr = s.dot(P)
trace_df = [list(P_curr) + [np.nan]]

for it in range(100):
    P_next = P_curr.dot(P)
    d = np.linalg.norm(P_next - P_curr, 1)

    if d < 0.0001:
        break

    trace_df.append(list(P_next) + [d])
    P_curr = P_next

trace_df = pd.DataFrame(trace_df, columns=['sunny', 'rainy', 'delta'])
P_next

array([0.83335081, 0.16664919])

Here is the trace of \(s^n\) until convergence.

trace_df

	sunny	rainy	delta
0	0.900000	0.100000	NaN
1	0.860000	0.140000	0.080000
2	0.844000	0.156000	0.032000
3	0.837600	0.162400	0.012800
4	0.835040	0.164960	0.005120
5	0.834016	0.165984	0.002048
6	0.833606	0.166394	0.000819
7	0.833443	0.166557	0.000328
8	0.833377	0.166623	0.000131

10.3. Eigen decomposition

We can also find the steady state directly through eigen decomposition. The steady state is the eigenvector associated with the eigenvalue closest to 1.

def get_stationary_p(X):
    S, U = np.linalg.eig(X.T)
    r = (U[:,np.isclose(S, 1)][:,0] / U[:,np.isclose(S, 1)][:,0].sum()).real
    return r

get_stationary_p(P)

array([0.83333333, 0.16666667])

10.4. Dynamic Bayesian Network

Another way to represent a MC is with a DBN. In general, the DBN structure consists of an endless number of nodes corresponding to time slices. However, in a DBN, we can unroll the DBN structure to just 2 consecutive time slices. It is assumed that the intra- (within a time slice) and inter-relationships (between two consecutive time slices) do not change over time.

from pybbn.graph.dag import Bbn
from pybbn.graph.edge import Edge, EdgeType
from pybbn.graph.jointree import EvidenceBuilder
from pybbn.graph.node import BbnNode
from pybbn.graph.variable import Variable
from pybbn.pptc.inferencecontroller import InferenceController

# create the nodes
w_curr = BbnNode(Variable(0, 'weather', ['rainy', 'sunny']), [1.0, 0.0])
w_next = BbnNode(Variable(2, 'weather_next', ['rainy', 'sunny']), [0.9, 0.1, 0.5, 0.5])

# create the network structure
bbn = Bbn() \
    .add_node(w_curr) \
    .add_node(w_next) \
    .add_edge(Edge(w_curr, w_next, EdgeType.DIRECTED))

# convert the BBN to a join tree
join_tree = InferenceController.apply(bbn)

for node in [w_curr, w_next]:
    potential = join_tree.get_bbn_potential(node)
    print(node)
    print(potential)
    print('-' * 15)

0|weather|rainy,sunny
0=rainy|1.00000
0=sunny|0.00000
---------------
2|weather_next|rainy,sunny
2=rainy|0.90000
2=sunny|0.10000
---------------

We can perform inference on the DBN and use the probability vector of the next state as input to the following time slice.

import pandas as pd

def get_probs(jt, node):
    entries = jt.get_bbn_potential(node).entries
    return np.array([e.value for e in entries])

P_curr = get_probs(join_tree, w_curr)
jt = InferenceController.reapply(join_tree, {w_curr.id: list(P_curr)})

trace_df = [list(P_curr) + [np.nan]]

for _ in range(20):
    P_next = get_probs(jt, w_next)
    d = np.linalg.norm(P_next - P_curr, 1)

    if d < 0.0001:
        break

    trace_df.append(list(P_next) + [d])
    P_curr = P_next

    jt = InferenceController.reapply(jt, {w_curr.id: list(P_curr)})

trace_df = pd.DataFrame(trace_df, columns=['sunny', 'rainy', 'delta'])
P_next

array([0.83335081, 0.16664919])

trace_df

	sunny	rainy	delta
0	1.000000	0.000000	NaN
1	0.900000	0.100000	0.200000
2	0.860000	0.140000	0.080000
3	0.844000	0.156000	0.032000
4	0.837600	0.162400	0.012800
5	0.835040	0.164960	0.005120
6	0.834016	0.165984	0.002048
7	0.833606	0.166394	0.000819
8	0.833443	0.166557	0.000328
9	0.833377	0.166623	0.000131