Continuous-Time Bayesian Networks

Creative Commons License

aGrUM

interactive online version

In [1]:
import pyAgrum as gum
import pyAgrum.lib.notebook as gnb
import pyAgrum.ctbn as ct
import pyAgrum.ctbn.notebook as ctbnnb

A CTBN(Continuous-Time Bayesian Networks) is a graphical model that allows a Bayesian Network to evolve over continuous time.

(For more details see https://pyagrum.readthedocs.io/en/latest/ctbn.html)

Creating a CTBN

In [2]:
# 1. Creating the ctbn
ctbn = ct.CTBN()

# 2. Adding the random variables
ctbn.add(gum.LabelizedVariable("cloudy?", "cloudy?"))
ctbn.add(gum.LabelizedVariable("rain?", "rain?"))

# 3. Adding the arcs
ctbn.addArc("cloudy?", "rain?")
ctbn.addArc("rain?", "cloudy?")

# 4. Filling the CIMS (Conditional Intensity Matrixes)
ctbn.CIM("cloudy?")[{"rain?": "0"}] = [[-.1, .1],
                                        [1, -1]]
ctbn.CIM("cloudy?")[{"rain?": "1"}] = [[-100, 100],
                                        [.5, -.5]]

ctbn.CIM("rain?")[{"cloudy?": "0"}] = [[-0.5, 0.5],
                                        [1000, -1000]]
ctbn.CIM("rain?")[{"cloudy?": "1"}] = [[-2, 2],
                                        [2, -2]]

# 5. Plotting the CTBN (whith/whithout the CIMs)
ctbnnb.showCtbnGraph(ctbn)
ctbnnb.showCIMs(ctbn)

../_images/notebooks_71-PyModels_CTBN_5_0.svg
cloudy?#j
rain?
cloudy?#i
0
1
0
0
-0.10000.1000
1
1.0000-1.0000
1
0
-100.0000100.0000
1
0.5000-0.5000

CIM(cloudy?)
rain?#j
cloudy?
rain?#i
0
1
0
0
-0.50000.5000
1
1000.0000-1000.0000
1
0
-2.00002.0000
1
2.0000-2.0000

CIM(rain?)

cloudy?#i gives the state cloudy? is currently in and cloudy?#j gives its next state.

The CIM object, built around a pyAgrum.Potential is used to store the parameters of the variables, describing the waiting time before switching to another state.

A useful tool for CIMs, mainly used for simple inference, is amalgamation. It allows to merge many CIMs into one.

In [3]:
cimA = ctbn.CIM("cloudy?")
cimB = ctbn.CIM("rain?")
amal = cimA.amalgamate(cimB) # or cimA * cimB

gnb.showPotential(amal.getPotential()) # as a Potential
cloudy?#i
rain?#i
rain?#j
cloudy?#j
0
1
0
0
0
-0.60001.0000
1
0.1000-3.0000
1
0
0.50000.0000
1
0.00002.0000
1
0
0
1000.00000.0000
1
0.00002.0000
1
0
-1100.00000.5000
1
100.0000-2.5000
In [4]:
print("\nJoint intensity matrix as numpy array\n")
print(amal.toMatrix()) # as a numpy array

Joint intensity matrix as numpy array

[[-6.0e-01  1.0e-01  5.0e-01  0.0e+00]
 [ 1.0e+00 -3.0e+00  0.0e+00  2.0e+00]
 [ 1.0e+03  0.0e+00 -1.1e+03  1.0e+02]
 [ 0.0e+00  2.0e+00  5.0e-01 -2.5e+00]]

A random variable of a ctbn can be identified by either its name or id.

In any case, most class functions parameters allow for a variable to be identified by either its name or id (type NameOrId).

In [5]:
labels_with_id = ctbn.labels(0)
labels_with_name = ctbn.labels("cloudy?")
print(f"labels using id : {labels_with_id}, using name : {labels_with_name}")

variable_with_id = ctbn.variable(0)
variable_with_name = ctbn.variable("cloudy?")
print(f"The variables 0 and 'cloudy?' are the same ? {variable_with_id==variable_with_name}")
labels using id : ('0', '1'), using name : ('0', '1')
The variables 0 and 'cloudy?' are the same ? True

All functions with NameOrId parameters work the same way

Random CTBN generator

The CtbnGenerator module can be used to quickly generate CTBNs. The parameters value are drawn at random from a value range. The diagonal coefficients are just the negative sum of the ones on the same line (in a CIM). It is also possible to choose the number of variables, the maximum number of parents a variable can have, and the domain size of the variables.

In [6]:
randomCtbn = ct.randomCTBN(valueRange=(1,10), n=8, parMax=2, modal=2)

ctbnnb.showCtbnGraph(randomCtbn)
ctbnnb.showCIMs(randomCtbn)
../_images/notebooks_71-PyModels_CTBN_15_0.svg
V1#j
V2
V5
V1#i
v1_1
v1_2
v2_1
v5_1
v1_1
-4.31724.3172
v1_2
4.8171-4.8171
v5_2
v1_1
-8.47468.4746
v1_2
2.6312-2.6312
v2_2
v5_1
v1_1
-5.84185.8418
v1_2
4.7838-4.7838
v5_2
v1_1
-6.24306.2430
v1_2
2.0341-2.0341

CIM(V1)
V2#j
V7
V8
V2#i
v2_1
v2_2
v7_1
v8_1
v2_1
-6.38786.3878
v2_2
8.9027-8.9027
v8_2
v2_1
-4.07414.0741
v2_2
6.3407-6.3407
v7_2
v8_1
v2_1
-9.50649.5064
v2_2
5.8454-5.8454
v8_2
v2_1
-8.68408.6840
v2_2
6.5747-6.5747

CIM(V2)
V3#j
V4
V2
V3#i
v3_1
v3_2
v4_1
v2_1
v3_1
-7.40487.4048
v3_2
3.9423-3.9423
v2_2
v3_1
-8.09098.0909
v3_2
6.0274-6.0274
v4_2
v2_1
v3_1
-6.42396.4239
v3_2
3.9850-3.9850
v2_2
v3_1
-8.06078.0607
v3_2
1.6419-1.6419

CIM(V3)
V4#j
V2
V1
V4#i
v4_1
v4_2
v2_1
v1_1
v4_1
-1.13421.1342
v4_2
3.7340-3.7340
v1_2
v4_1
-8.80928.8092
v4_2
2.7529-2.7529
v2_2
v1_1
v4_1
-9.46009.4600
v4_2
5.6839-5.6839
v1_2
v4_1
-6.79176.7917
v4_2
9.8443-9.8443

CIM(V4)
V5#j
V1
V4
V5#i
v5_1
v5_2
v1_1
v4_1
v5_1
-3.30463.3046
v5_2
8.3787-8.3787
v4_2
v5_1
-4.09244.0924
v5_2
7.1730-7.1730
v1_2
v4_1
v5_1
-5.90625.9062
v5_2
5.7211-5.7211
v4_2
v5_1
-1.34301.3430
v5_2
1.8059-1.8059

CIM(V5)
V6#j
V8
V7
V6#i
v6_1
v6_2
v8_1
v7_1
v6_1
-3.76143.7614
v6_2
2.4733-2.4733
v7_2
v6_1
-1.87821.8782
v6_2
4.0112-4.0112
v8_2
v7_1
v6_1
-7.59137.5913
v6_2
4.7263-4.7263
v7_2
v6_1
-8.15588.1558
v6_2
7.8676-7.8676

CIM(V6)
V7#j
V7#i
v7_1
v7_2
v7_1
-8.48168.4816
v7_2
7.6979-7.6979

CIM(V7)
V8#j
V5
V4
V8#i
v8_1
v8_2
v5_1
v4_1
v8_1
-5.75925.7592
v8_2
9.6744-9.6744
v4_2
v8_1
-8.75358.7535
v8_2
8.5555-8.5555
v5_2
v4_1
v8_1
-2.26682.2668
v8_2
2.5436-2.5436
v4_2
v8_1
-8.93728.9372
v8_2
5.2977-5.2977

CIM(V8)

Inference with CTBNs

Exact inference & Forward Sampling inference

  • Amalgamation is used to compute exact inference. It consists in merging all the CIMs of a CTBN into one CIM. Then we use the properties of markov process : \(P(X_t) = P(X_0)*exp(Q_Xt)\) for enough time to notice convergence.

  • Forward sampling works this way:

    • At time \(t\) :

    • Draw a “time until next change” \(dt\) value for each variable’s exponential distribution of their waiting time.

    • Select the variable with smallest transition time : \(next time = t + dt\).

    • Draw the variable next state uniformly.

    • loop until \(t = timeHorizon\).

When doing forward sampling, we first iterate this process over a given amount of iterations (called burnIn). We consider the last values as the initial configuration for sampling. The reason is that the initial state of a CTBN is unknown without real data.

In [7]:
# 1. Initialiaze Inference object
ie = ct.SimpleInference(ctbn)
ie_forward = ct.ForwardSamplingInference(ctbn)
gnb.sideBySide(ctbnnb.getCtbnGraph(ctbn), ctbnnb.getCtbnInference(ctbn, ie), ctbnnb.getCtbnInference(ctbn, ie_forward))
ctbn cloudy? cloudy? rain? rain? cloudy?->rain? rain?->cloudy?
structs Inference in   2.82ms cloudy? 2024-07-27T17:53:25.628657 image/svg+xml Matplotlib v3.9.1, https://matplotlib.org/ rain? 2024-07-27T17:53:25.641730 image/svg+xml Matplotlib v3.9.1, https://matplotlib.org/ cloudy?->rain? rain?->cloudy?
structs Inference in 127.02ms cloudy? 2024-07-27T17:53:25.943402 image/svg+xml Matplotlib v3.9.1, https://matplotlib.org/ rain? 2024-07-27T17:53:25.956043 image/svg+xml Matplotlib v3.9.1, https://matplotlib.org/ cloudy?->rain? rain?->cloudy?

In any case, inference can be made using makeInference() function from every Inference object.

In [8]:
import time
# Simple inference
time1 = time.time()
ie.makeInference(t=5000)
pot1 = ie.posterior("cloudy?")
elapsed1 = time.time() - time1

# Sampling inference
time2 = time.time()
ie_forward.makeInference(timeHorizon=5000, burnIn=100)
pot2 = ie_forward.posterior("cloudy?")
elapsed2 = time.time() - time2

gnb.sideBySide(gnb.getPotential(pot1), gnb.getPotential(pot2), captions=[f"Simple inference {round(elapsed1*1000,2)} ms", f"Sampling inference {round(elapsed2*1000, 2)} ms"])
cloudy?
0
1
0.83340.1666

Simple inference 0.59 ms
cloudy?
0
1
0.82350.1765

Sampling inference 176.49 ms

Plotting a sample

Here is an example of a trajectory plot

In [9]:
# 1_1. Find trajectory stored in the inference class
ie_forward.makeInference()
traj = ie_forward.trajectory # a single trajectory List[Tuple[float, str, str]]

# 1_2. Or find it in a csv file
ie_forward.writeTrajectoryCSV('out/trajectory.csv', n=1, timeHorizon=1000, burnIn=100) # to save a trajectory
traj = ct.readTrajectoryCSV('out/trajectory.csv') # a dict of trajectories

# 2. Plot the trajectory
ct.plotTrajectory(ctbn.variable("cloudy?"), traj[0], timeHorizon=40, plotname="plot of cloudy?")
../_images/notebooks_71-PyModels_CTBN_22_0.svg

It is also possible to plot the proportions of the states a variable is in over time.

Let’s take the exemple of a variable with domain size=5.

In [10]:
# 1. Create CTBN
new_ctbn = ct.CTBN()
X = gum.LabelizedVariable("X", "X", ["a1", "a2", "a3", "a4", "a5"])
new_ctbn.add(X)
new_ctbn.CIM("X")[:] = [[-10, 2, 4, 0, 4],
                        [0.5 , -1, 0.5, 0, 0],
                        [5, 0.5, -10, 3.5, 1],
                        [50, 0.5, 50, -200.5, 100],
                        [1, 0, 2, 7, -10]]
# 2. Sampling
ieX = ct.ForwardSamplingInference(new_ctbn)
ieX.writeTrajectoryCSV("out/traj_quinary.csv", n=3, timeHorizon=100)
trajX = ct.readTrajectoryCSV("out/traj_quinary.csv")

# 3. Plotting
ct.plotFollowVar(X, trajX)
../_images/notebooks_71-PyModels_CTBN_24_0.svg

Learning a CTBN

The ctbn library has tools to learn the graph of a CTBN as well as the CIMs of its variables

Learning the variables

From a sample, we can find the variables name and their labels. However, if the sample is too short, it is possible that some information may be missed.

In [11]:
newCtbn = ct.CTBNFromData(traj)

print("variables found and their labels")
for name in newCtbn.names():
    print(f"variable {name}, labels : {newCtbn.labels(name)}")
variables found and their labels
variable cloudy?, labels : ('0', '1')
variable rain?, labels : ('0', '1')

Learning the graph

In [12]:
# 1. Create a Learner.
learner = ct.Learner('out/trajectory.csv')

# 2. Run learning algorithm
learned_ctbn = learner.learnCTBN()

# 3. Plot the learned CTBN
gnb.show(learned_ctbn)
../_images/notebooks_71-PyModels_CTBN_29_0.svg

Let’s note that having a test object is necessary since we want to be able to test independency between variables. For now only 2 tests are available : (Fisher and Chi2) test and Oracle test. Oracle test is only used for benchmarking since it uses the initial ctbn as a reference.

Fitting parameters

It is also possible to approximate the transition time distributions of each variable, i.e the CIMs.

Note that conditioning variables are taken from the ctbn.

In [13]:
cims_before = ctbnnb.getCIMs(ctbn)
learner.fitParameters(ctbn)
cims_after = ctbnnb.getCIMs(ctbn)

gnb.sideBySide(cims_before, cims_after, captions=["original CIMs", "learned CIMs"])
cloudy?#j
rain?
cloudy?#i
0
1
0
0
-0.10000.1000
1
1.0000-1.0000
1
0
-100.0000100.0000
1
0.5000-0.5000

CIM(cloudy?)
rain?#j
cloudy?
rain?#i
0
1
0
0
-0.50000.5000
1
1000.0000-1000.0000
1
0
-2.00002.0000
1
2.0000-2.0000

CIM(rain?)

original CIMs
cloudy?#j
rain?
cloudy?#i
0
1
0
0
-0.09400.0940
1
1.0970-1.0970
1
0
-83.146083.1460
1
0.4500-0.4500

CIM(cloudy?)
rain?#j
cloudy?
rain?#i
0
1
0
0
-0.50300.5030
1
1083.4120-1083.4120
1
0
-1.98501.9850
1
1.9350-1.9350

CIM(rain?)

learned CIMs

Example from [Nodelman et al, 2003]

This example comes from :

U. Nodelman, C.R. Shelton, and D. Koller (2003). "Learning Continuous Time Bayesian Networks." Proc. Nineteenth Conference on Uncertainty in Artificial Intelligence (UAI) (pp. 451-458).
In [14]:
# 1. Creating the ctbn
drugEffect = ct.CTBN()

drugEffect.add(gum.LabelizedVariable("Eating", ""))
drugEffect.add(gum.LabelizedVariable("Hungry", ""))
drugEffect.add(gum.LabelizedVariable("Full stomach", ""))
drugEffect.add(gum.LabelizedVariable("Uptake", ""))
drugEffect.add(gum.LabelizedVariable("Concentration", "?"))
drugEffect.add(gum.LabelizedVariable("Barometer", ""))
drugEffect.add(gum.LabelizedVariable("Joint pain", ""))
drugEffect.add(gum.LabelizedVariable("Drowsy", ""))

drugEffect.addArc("Hungry","Eating")
drugEffect.addArc("Eating","Full stomach")
drugEffect.addArc("Full stomach","Hungry")
drugEffect.addArc("Full stomach","Concentration")
drugEffect.addArc("Uptake","Concentration")
drugEffect.addArc("Concentration","Drowsy")
drugEffect.addArc("Barometer","Joint pain")
drugEffect.addArc("Concentration","Joint pain")

gnb.show(drugEffect)
../_images/notebooks_71-PyModels_CTBN_34_0.svg
In [15]:
from pyAgrum.ctbn.CTBNGenerator import randomCIMs

randomCIMs(drugEffect,(0.1,10))
drugEffect.CIM("Eating")[{"Hungry":0}]=[[-0.1,0.1],
                                        [10,-10]]
drugEffect.CIM("Eating")[{"Hungry":1}]=[[-2,2],
                                        [0.1,-0.1]]
drugEffect.CIM("Full stomach")[{"Eating":0}]=[[-0.1,0.1],
                                        [10,-10]]
drugEffect.CIM("Full stomach")[{"Eating":1}]=[[-2,2],
                                        [0.1,-0.1]]
drugEffect.CIM("Hungry")[{"Full stomach":0}]=[[-0.1,0.1],
                                        [10,-10]]
drugEffect.CIM("Hungry")[{"Full stomach":1}]=[[-2,2],
                                        [0.1,-0.1]]

ieX = ct.ForwardSamplingInference(drugEffect)
ieX.writeTrajectoryCSV("out/drugEffect.csv", n=3, timeHorizon=100)
trajX = ct.readTrajectoryCSV("out/drugEffect.csv")

ct.plotFollowVar(drugEffect.variable("Hungry"), trajX)
ct.plotFollowVar(drugEffect.variable("Full stomach"), trajX)
ct.plotFollowVar(drugEffect.variable("Eating"), trajX)
ct.plotFollowVar(drugEffect.variable("Concentration"), trajX)
../_images/notebooks_71-PyModels_CTBN_35_0.svg
../_images/notebooks_71-PyModels_CTBN_35_1.svg
../_images/notebooks_71-PyModels_CTBN_35_2.svg
../_images/notebooks_71-PyModels_CTBN_35_3.svg

Input/output with pickle

In [16]:
import pickle
with open("test.pkl","bw") as f:
  pickle.dump(drugEffect,f)
drugEffect
Out[16]:
ctbn Eating Eating Full stomach Full stomach Eating->Full stomach Hungry Hungry Hungry->Eating Full stomach->Hungry Concentration Concentration Full stomach->Concentration Uptake Uptake Uptake->Concentration Joint pain Joint pain Concentration->Joint pain Drowsy Drowsy Concentration->Drowsy Barometer Barometer Barometer->Joint pain
In [17]:
with open("test.pkl","br") as f:
  copyDrug=pickle.load(f)
copyDrug
Out[17]:
ctbn Eating Eating Full stomach Full stomach Eating->Full stomach Hungry Hungry Hungry->Eating Full stomach->Hungry Concentration Concentration Full stomach->Concentration Uptake Uptake Uptake->Concentration Joint pain Joint pain Concentration->Joint pain Drowsy Drowsy Concentration->Drowsy Barometer Barometer Barometer->Joint pain
In [ ]: