[PYTHON] pgmpy: Discrete Bayesian Network Trial-Inference

Trigger

Interested in meta-learning, graph neural networks, use of knowledge structures including use as prior knowledge, and knowledge storage. Relatedly, I was looking for a library to easily implement a Bayesian network. Since pgmpy looked good, I will record the general flow.

pgmpy：pgmpy is a python library for working with Probabilistic Graphical Models. https://pgmpy.org/

reference

Implement Bayesian network with Titanic data https://qiita.com/YuyaOmori/items/e051f0360d1f9562620b

Bayesian Network: From Introduction to Application to Human Modeling https://staff.aist.go.jp/y.motomura/paper/BSJ0403.pdf

environment

Windows10 Python3.7 Anaconda pgmpy==0.1.9

Installation

pip install pgmpy==0.1.9

If it is nonGPU without pytorch

conda install pytorch torchvision cpuonly -c pytorch

data

Use the following data

python


import pandas as pd
df = pd.DataFrame()
df['t'] = [1, 1, 1, 1, 0, 0, 1, 1, 1, 2, 0, 0, 1, 1, 1, 2, 2, 2, 2, 2]
df['a'] = [2, 2, 2, 2, 1, 1, 1, 1, 2, 1, 1, 2, 0, 0, 0, 1, 1, 2, 2, 2]
df['h'] = [0, 0, 0, 1, 0, 0, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1]

Run

Model structure definition

python


from pgmpy.models import BayesianModel
model = BayesianModel([('t','h'),('a','h')])

A directed acyclic graph of t → h and a → h was used. image.png

Create and check CPD in the model

python


model.fit(df) #Conditions are omitted. Pay particular attention to overfitting by default

print(model.get_cpds('t'))
print(model.get_cpds('a'))
print(model.get_cpds('h'))

output


+------+-----+
| t(0) | 0.2 |
+------+-----+
| t(1) | 0.5 |
+------+-----+
| t(2) | 0.3 |
+------+-----+
+------+------+
| a(0) | 0.15 |
+------+------+
| a(1) | 0.4  |
+------+------+
| a(2) | 0.45 |
+------+------+
+------+------+--------------------+------+------+------+------+------+------+------+
| a    | a(0) | a(0)               | a(0) | a(1) | a(1) | a(1) | a(2) | a(2) | a(2) |
+------+------+--------------------+------+------+------+------+------+------+------+
| t    | t(0) | t(1)               | t(2) | t(0) | t(1) | t(2) | t(0) | t(1) | t(2) |
+------+------+--------------------+------+------+------+------+------+------+------+
| h(0) | 0.5  | 0.3333333333333333 | 0.5  | 1.0  | 0.0  | 0.0  | 1.0  | 0.6  | 0.0  |
+------+------+--------------------+------+------+------+------+------+------+------+
| h(1) | 0.5  | 0.6666666666666666 | 0.5  | 0.0  | 1.0  | 1.0  | 0.0  | 0.4  | 1.0  |
+------+------+--------------------+------+------+------+------+------+------+------+

Reasoning 1

python


from pgmpy.inference import VariableElimination
ve = VariableElimination(model)

#t=1,h=What is a when 1 is set?
print(ve.map_query(variables=['a'], evidence={'t':1, 'h':1}))

output


{'a': 1}

Reasoning 2

python


#t=0,1,When it is 2, a,h Each estimated value? What?
for i in [0,1,2]:
    print(ve.query(variables=['a', 'h'], evidence={'t':i}))

output



+------+------+------------+
| a    | h    |   phi(a,h) |
+======+======+============+
| a(0) | h(0) |     0.0750 |
+------+------+------------+
| a(0) | h(1) |     0.0750 |
+------+------+------------+
| a(1) | h(0) |     0.4000 |
+------+------+------------+
| a(1) | h(1) |     0.0000 |
+------+------+------------+
| a(2) | h(0) |     0.4500 |
+------+------+------------+
| a(2) | h(1) |     0.0000 |
+------+------+------------+

+------+------+------------+
| h    | a    |   phi(h,a) |
+======+======+============+
| h(0) | a(0) |     0.0500 |
+------+------+------------+
| h(0) | a(1) |     0.0000 |
+------+------+------------+
| h(0) | a(2) |     0.2700 |
+------+------+------------+
| h(1) | a(0) |     0.1000 |
+------+------+------------+
| h(1) | a(1) |     0.4000 |
+------+------+------------+
| h(1) | a(2) |     0.1800 |
+------+------+------------+

+------+------+------------+
| a    | h    |   phi(a,h) |
+======+======+============+
| a(0) | h(0) |     0.0750 |
+------+------+------------+
| a(0) | h(1) |     0.0750 |
+------+------+------------+
| a(1) | h(0) |     0.0000 |
+------+------+------------+
| a(1) | h(1) |     0.4000 |
+------+------+------------+
| a(2) | h(0) |     0.0000 |
+------+------+------------+
| a(2) | h(1) |     0.4500 |
+------+------+------------+

Supplement

--model.fit (df) can be split, for example: Note that it may be easier to handle if it is divided.

python


#CPD creation part
from pgmpy.estimators import BayesianEstimator
estimator = BayesianEstimator(model, df)
cpd_ta = estimator.estimate_cpd('t', prior_type='dirichlet', pseudo_counts=[[0],[0],[0]])
cpd_aa = estimator.estimate_cpd('a', prior_type='dirichlet', pseudo_counts=[[0],[0],[0]])
cpd_h = estimator.estimate_cpd('h', prior_type='dirichlet', pseudo_counts=[[0,0,0,0,0,0,0,0,0],[0,0,0,0,0,0,0,0,0]])

#CPD input part
model.add_cpds(cpd_ta, cpd_aa, cpd_h)

--If you want to create CPD arbitrarily, for example:

python


from pgmpy.factors.discrete import TabularCPD
cpd_h = TabularCPD(variable='h', variable_card=2,
                        values=[[1, 0.3, 0.5, 1, 0, 0, 1, 0.6, 0],
                                [0, 0.7, 0.5, 0, 1, 1, 0, 0.4, 1]],
                        evidence=['t', 'a'],
                        evidence_card=[3, 3])

--Structural learning is omitted.

--It seems that there is also a function like Junction tree that supports multiply connected. Checking how much we can do.

――It seems that there is also a method called dynamic Bayesian network that considers the time system.

-[Relay series] My recommendation node-The famous detective "Bayes node" who finds hidden relationships unravels the causal structure between variables https://www.ibm.com/blogs/solutions/jp-ja/spssmodeler-push-node-18/ Tried with the iris dataset for reference. I was in trouble because one of the features did not appear in structural learning, but there was no problem with the others.

--I want to emphasize the extraction of individual pieces. Next, proceed as follows. Hierarchical Bayes (from "The World of Bayes Modeling") https://recruit.cct-inc.co.jp/tecblog/machine-learning/hierarchical-bayesian/ PyMC3. Requires Theano, which only works with python 3.6 or less.

--NumPyro has hierarchical Bayes. PyTorch-based pyro is almost the same. https://pyro.ai/numpyro/bayesian_hierarchical_linear_regression.html#2.-Modelling:-Bayesian-Hierarchical-Linear-Regression-with-Partial-Pooling https://qiita.com/takeajioka/items/ab299d75efa184eb1432

Recommended Posts

pgmpy: Discrete Bayesian Network Trial-Inference
Bayesian network package ~ Pebl tutorial execution ~