Algorithms for Inference

Variable Elimination

class pgmpy.inference.ExactInference.VariableElimination(*args, **kwargs)[source]
induced_graph(elimination_order)[source]

Returns the induced graph formed by running Variable Elimination on the network.

elimination_order: list, array like
List of variables in the order in which they are to be eliminated.
>>> import numpy as np
>>> import pandas as pd
>>> from pgmpy.models import BayesianModel
>>> from pgmpy.inference import VariableElimination
>>> values = pd.DataFrame(np.random.randint(low=0, high=2, size=(1000, 5)),
...                       columns=['A', 'B', 'C', 'D', 'E'])
>>> model = BayesianModel([('A', 'B'), ('C', 'B'), ('C', 'D'), ('B', 'E')])
>>> model.fit(values)
>>> inference = VariableElimination(model)
>>> inference.induced_graph(['C', 'D', 'A', 'B', 'E'])
<networkx.classes.graph.Graph at 0x7f34ac8c5160>
induced_width(elimination_order)[source]

Returns the width (integer) of the induced graph formed by running Variable Elimination on the network. The width is the defined as the number of nodes in the largest clique in the graph minus 1.

elimination_order: list, array like
List of variables in the order in which they are to be eliminated.
>>> import numpy as np
>>> import pandas as pd
>>> from pgmpy.models import BayesianModel
>>> from pgmpy.inference import VariableElimination
>>> values = pd.DataFrame(np.random.randint(low=0, high=2, size=(1000, 5)),
...                       columns=['A', 'B', 'C', 'D', 'E'])
>>> model = BayesianModel([('A', 'B'), ('C', 'B'), ('C', 'D'), ('B', 'E')])
>>> model.fit(values)
>>> inference = VariableElimination(model)
>>> inference.induced_width(['C', 'D', 'A', 'B', 'E'])
3
max_marginal(variables=None, evidence=None, elimination_order=None)[source]

Computes the max-marginal over the variables given the evidence.

variables: list
list of variables over which we want to compute the max-marginal.
evidence: dict
a dict key, value pair as {var: state_of_var_observed} None if no evidence
elimination_order: list
order of variable eliminations (if nothing is provided) order is computed automatically
>>> import numpy as np
>>> import pandas as pd
>>> from pgmpy.models import BayesianModel
>>> from pgmpy.inference import VariableElimination
>>> values = pd.DataFrame(np.random.randint(low=0, high=2, size=(1000, 5)),
...                       columns=['A', 'B', 'C', 'D', 'E'])
>>> model = BayesianModel([('A', 'B'), ('C', 'B'), ('C', 'D'), ('B', 'E')])
>>> model.fit(values)
>>> inference = VariableElimination(model)
>>> phi_query = inference.max_marginal(['A', 'B'])
query(variables, evidence=None, elimination_order=None)[source]
variables: list
list of variables for which you want to compute the probability
evidence: dict
a dict key, value pair as {var: state_of_var_observed} None if no evidence
elimination_order: list
order of variable eliminations (if nothing is provided) order is computed automatically
>>> from pgmpy.inference import VariableElimination
>>> from pgmpy.models import BayesianModel
>>> import numpy as np
>>> import pandas as pd
>>> values = pd.DataFrame(np.random.randint(low=0, high=2, size=(1000, 5)),
...                       columns=['A', 'B', 'C', 'D', 'E'])
>>> model = BayesianModel([('A', 'B'), ('C', 'B'), ('C', 'D'), ('B', 'E')])
>>> model.fit(values)
>>> inference = VariableElimination(model)
>>> phi_query = inference.query(['A', 'B'])

Belief Propagation

class pgmpy.inference.ExactInference.BeliefPropagation(model)[source]

Class for performing inference using Belief Propagation method.

Creates a Junction Tree or Clique Tree (JunctionTree class) for the input probabilistic graphical model and performs calibration of the junction tree so formed using belief propagation.

model: BayesianModel, MarkovModel, FactorGraph, JunctionTree
model for which inference is to performed
calibrate()[source]

Calibration using belief propagation in junction tree or clique tree.

>>> from pgmpy.models import BayesianModel
>>> from pgmpy.factors.discrete import TabularCPD
>>> from pgmpy.inference import BeliefPropagation
>>> G = BayesianModel([('diff', 'grade'), ('intel', 'grade'),
...                    ('intel', 'SAT'), ('grade', 'letter')])
>>> diff_cpd = TabularCPD('diff', 2, [[0.2], [0.8]])
>>> intel_cpd = TabularCPD('intel', 3, [[0.5], [0.3], [0.2]])
>>> grade_cpd = TabularCPD('grade', 3,
...                        [[0.1, 0.1, 0.1, 0.1, 0.1, 0.1],
...                         [0.1, 0.1, 0.1, 0.1, 0.1, 0.1],
...                         [0.8, 0.8, 0.8, 0.8, 0.8, 0.8]],
...                        evidence=['diff', 'intel'],
...                        evidence_card=[2, 3])
>>> sat_cpd = TabularCPD('SAT', 2,
...                      [[0.1, 0.2, 0.7],
...                       [0.9, 0.8, 0.3]],
...                      evidence=['intel'], evidence_card=[3])
>>> letter_cpd = TabularCPD('letter', 2,
...                         [[0.1, 0.4, 0.8],
...                          [0.9, 0.6, 0.2]],
...                         evidence=['grade'], evidence_card=[3])
>>> G.add_cpds(diff_cpd, intel_cpd, grade_cpd, sat_cpd, letter_cpd)
>>> bp = BeliefPropagation(G)
>>> bp.calibrate()
get_clique_beliefs()[source]

Returns clique beliefs. Should be called after the clique tree (or junction tree) is calibrated.

get_cliques()[source]

Returns cliques used for belief propagation.

get_sepset_beliefs()[source]

Returns sepset beliefs. Should be called after clique tree (or junction tree) is calibrated.

map_query(variables=None, evidence=None)[source]

MAP Query method using belief propagation.

variables: list
list of variables for which you want to compute the probability
evidence: dict
a dict key, value pair as {var: state_of_var_observed} None if no evidence
>>> from pgmpy.factors.discrete import TabularCPD
>>> from pgmpy.models import BayesianModel
>>> from pgmpy.inference import BeliefPropagation
>>> bayesian_model = BayesianModel([('A', 'J'), ('R', 'J'), ('J', 'Q'),
...                                 ('J', 'L'), ('G', 'L')])
>>> cpd_a = TabularCPD('A', 2, [[0.2], [0.8]])
>>> cpd_r = TabularCPD('R', 2, [[0.4], [0.6]])
>>> cpd_j = TabularCPD('J', 2,
...                    [[0.9, 0.6, 0.7, 0.1],
...                     [0.1, 0.4, 0.3, 0.9]],
...                    ['R', 'A'], [2, 2])
>>> cpd_q = TabularCPD('Q', 2,
...                    [[0.9, 0.2],
...                     [0.1, 0.8]],
...                    ['J'], [2])
>>> cpd_l = TabularCPD('L', 2,
...                    [[0.9, 0.45, 0.8, 0.1],
...                     [0.1, 0.55, 0.2, 0.9]],
...                    ['G', 'J'], [2, 2])
>>> cpd_g = TabularCPD('G', 2, [[0.6], [0.4]])
>>> belief_propagation = BeliefPropagation(bayesian_model)
>>> belief_propagation.map_query(variables=['J', 'Q'],
...                              evidence={'A': 0, 'R': 0, 'G': 0, 'L': 1})
max_calibrate()[source]

Max-calibration of the junction tree using belief propagation.

>>> from pgmpy.models import BayesianModel
>>> from pgmpy.factors.discrete import TabularCPD
>>> from pgmpy.inference import BeliefPropagation
>>> G = BayesianModel([('diff', 'grade'), ('intel', 'grade'),
...                    ('intel', 'SAT'), ('grade', 'letter')])
>>> diff_cpd = TabularCPD('diff', 2, [[0.2], [0.8]])
>>> intel_cpd = TabularCPD('intel', 3, [[0.5], [0.3], [0.2]])
>>> grade_cpd = TabularCPD('grade', 3,
...                        [[0.1, 0.1, 0.1, 0.1, 0.1, 0.1],
...                         [0.1, 0.1, 0.1, 0.1, 0.1, 0.1],
...                         [0.8, 0.8, 0.8, 0.8, 0.8, 0.8]],
...                        evidence=['diff', 'intel'],
...                        evidence_card=[2, 3])
>>> sat_cpd = TabularCPD('SAT', 2,
...                      [[0.1, 0.2, 0.7],
...                       [0.9, 0.8, 0.3]],
...                      evidence=['intel'], evidence_card=[3])
>>> letter_cpd = TabularCPD('letter', 2,
...                         [[0.1, 0.4, 0.8],
...                          [0.9, 0.6, 0.2]],
...                         evidence=['grade'], evidence_card=[3])
>>> G.add_cpds(diff_cpd, intel_cpd, grade_cpd, sat_cpd, letter_cpd)
>>> bp = BeliefPropagation(G)
>>> bp.max_calibrate()
query(variables, evidence=None)[source]

Query method using belief propagation.

variables: list
list of variables for which you want to compute the probability
evidence: dict
a dict key, value pair as {var: state_of_var_observed} None if no evidence
>>> from pgmpy.factors.discrete import TabularCPD
>>> from pgmpy.models import BayesianModel
>>> from pgmpy.inference import BeliefPropagation
>>> bayesian_model = BayesianModel([('A', 'J'), ('R', 'J'), ('J', 'Q'),
...                                 ('J', 'L'), ('G', 'L')])
>>> cpd_a = TabularCPD('A', 2, [[0.2], [0.8]])
>>> cpd_r = TabularCPD('R', 2, [[0.4], [0.6]])
>>> cpd_j = TabularCPD('J', 2,
...                    [[0.9, 0.6, 0.7, 0.1],
...                     [0.1, 0.4, 0.3, 0.9]],
...                    ['R', 'A'], [2, 2])
>>> cpd_q = TabularCPD('Q', 2,
...                    [[0.9, 0.2],
...                     [0.1, 0.8]],
...                    ['J'], [2])
>>> cpd_l = TabularCPD('L', 2,
...                    [[0.9, 0.45, 0.8, 0.1],
...                     [0.1, 0.55, 0.2, 0.9]],
...                    ['G', 'J'], [2, 2])
>>> cpd_g = TabularCPD('G', 2, [[0.6], [0.4]])
>>> belief_propagation = BeliefPropagation(bayesian_model)
>>> belief_propagation.query(variables=['J', 'Q'],
...                          evidence={'A': 0, 'R': 0, 'G': 0, 'L': 1})

MPLP

class pgmpy.inference.mplp.Mplp(model)[source]

Class for performing approximate inference using Max-Product Linear Programming method.

We derive message passing updates that result in monotone decrease of the dual of the MAP LP Relaxation.

model: MarkovModel for which inference is to be performed. Examples ——– >>> import numpy as np >>> from pgmpy.models import MarkovModel >>> from pgmpy.factors.discrete import DiscreteFactor >>> student = MarkovModel() >>> student.add_edges_from([(‘A’, ‘B’), (‘B’, ‘C’), (‘C’, ‘D’), (‘E’, ‘F’)]) >>> factor_a = DiscreteFactor([‘A’], cardinality=[2], value=np.array([0.54577, 1.8323])) >>> factor_b = DiscreteFactor([‘B’], cardinality=[2], value=np.array([0.93894, 1.065])) >>> factor_c = DiscreteFactor([‘C’], cardinality=[2], value=np.array([0.89205, 1.121])) >>> factor_d = DiscreteFactor([‘D’], cardinality=[2], value=np.array([0.56292, 1.7765])) >>> factor_e = DiscreteFactor([‘E’], cardinality=[2], value=np.array([0.47117, 2.1224])) >>> factor_f = DiscreteFactor([‘F’], cardinality=[2], value=np.array([1.5093, 0.66257])) >>> factor_a_b = DiscreteFactor([‘A’, ‘B’], cardinality=[2, 2], value=np.array([1.3207, 0.75717, 0.75717, 1.3207])) >>> factor_b_c = DiscreteFactor([‘B’, ‘C’], cardinality=[2, 2], value=np.array([0.00024189, 4134.2, 4134.2, 0.00024189])) >>> factor_c_d = DiscreteFactor([‘C’, ‘D’], cardinality=[2, 2], value=np.array([0.0043227, 231.34, 231.34, 0.0043227])) >>> factor_d_e = DiscreteFactor([‘E’, ‘F’], cardinality=[2, 2], value=np.array([31.228, 0.032023, 0.032023, 31.228])) >>> student.add_factors(factor_a, factor_b, factor_c, factor_d, factor_e, factor_f, factor_a_b, ... factor_b_c, factor_c_d, factor_d_e) >>> mplp = Mplp(student)

class Cluster(intersection_set_variables, cluster_potential)[source]

Inner class for representing a cluster. A cluster is a subset of variables.

set_of_variables: tuple
This is the set of variables that form the cluster.
intersection_set_variables: set containing frozensets.
collection of intersection of all pairs of cluster variables.

For eg: {{C_1 cap C_2}, {C_2 cap C_3}, {C_3 cap C_1} } for clusters C_1, C_2 & C_3.

cluster_potential: DiscreteFactor
Each cluster has a initial probability distribution provided beforehand.
Mplp.find_triangles()[source]

Finds all the triangles present in the given model

>>> from pgmpy.models import MarkovModel
>>> from pgmpy.factors.discrete import DiscreteFactor
>>> from pgmpy.inference import Mplp
>>> mm = MarkovModel()
>>> mm.add_nodes_from(['x1', 'x2', 'x3', 'x4', 'x5', 'x6', 'x7'])
>>> mm.add_edges_from([('x1', 'x3'), ('x1', 'x4'), ('x2', 'x4'),
...                    ('x2', 'x5'), ('x3', 'x6'), ('x4', 'x6'),
...                    ('x4', 'x7'), ('x5', 'x7')])
>>> phi = [DiscreteFactor(edge, [2, 2], np.random.rand(4)) for edge in mm.edges()]
>>> mm.add_factors(*phi)
>>> mplp = Mplp(mm)
>>> mplp.find_triangles()
Mplp.get_integrality_gap()[source]
Returns the integrality gap of the current state of the Mplp algorithm. The lesser it is, the closer we are
towards the exact solution.
>>> from pgmpy.models import MarkovModel
>>> from pgmpy.factors.discrete import DiscreteFactor
>>> from pgmpy.inference import Mplp
>>> mm = MarkovModel()
>>> mm.add_nodes_from(['x1', 'x2', 'x3', 'x4', 'x5', 'x6', 'x7'])
>>> mm.add_edges_from([('x1', 'x3'), ('x1', 'x4'), ('x2', 'x4'),
...                    ('x2', 'x5'), ('x3', 'x6'), ('x4', 'x6'),
...                    ('x4', 'x7'), ('x5', 'x7')])
>>> phi = [DiscreteFactor(edge, [2, 2], np.random.rand(4)) for edge in mm.edges()]
>>> mm.add_factors(*phi)
>>> mplp = Mplp(mm)
>>> mplp.map_query()
>>> int_gap = mplp.get_integrality_gap()
Mplp.map_query(init_iter=1000, later_iter=20, dual_threshold=0.0002, integrality_gap_threshold=0.0002, tighten_triplet=True, max_triplets=5, max_iterations=100, prolong=False)[source]

MAP query method using Max Product LP method. This returns the best assignment of the nodes in the form of a dictionary.

init_iter: integer
Number of maximum iterations that we want MPLP to run for the first time.
later_iter: integer
Number of maximum iterations that we want MPLP to run for later iterations
dual_threshold: double
This sets the minimum width between the dual objective decrements. If the decrement is lesser than the threshold, then that means we have stuck on a local minima.
integrality_gap_threshold: double
This sets the threshold for the integrality gap below which we say that the solution is satisfactory.
tighten_triplet: bool
set whether to use triplets as clusters or not.
max_triplets: integer
Set the maximum number of triplets that can be added at once.
max_iterations: integer
Maximum number of times we tighten the relaxation. Used only when tighten_triplet is set True.
prolong: bool

If set False: The moment we exhaust of all the triplets the tightening stops. If set True: The tightening will be performed max_iterations number of times irrespective of the

triplets.

Reference: Section 3.3: The Dual Algorithm; Tightening LP Relaxation for MAP using Message Passing (2008) By Sontag Et al.

>>> from pgmpy.models import MarkovModel
>>> from pgmpy.factors.discrete import DiscreteFactor
>>> import numpy as np
>>> student = MarkovModel()
>>> student.add_edges_from([('A', 'B'), ('B', 'C'), ('C', 'D'), ('E', 'F')])
>>> factor_a = DiscreteFactor(['A'], cardinality=[2], value=np.array([0.54577, 1.8323]))
>>> factor_b = DiscreteFactor(['B'], cardinality=[2], value=np.array([0.93894, 1.065]))
>>> factor_c = DiscreteFactor(['C'], cardinality=[2], value=np.array([0.89205, 1.121]))
>>> factor_d = DiscreteFactor(['D'], cardinality=[2], value=np.array([0.56292, 1.7765]))
>>> factor_e = DiscreteFactor(['E'], cardinality=[2], value=np.array([0.47117, 2.1224]))
>>> factor_f = DiscreteFactor(['F'], cardinality=[2], value=np.array([1.5093, 0.66257]))
>>> factor_a_b = DiscreteFactor(['A', 'B'], cardinality=[2, 2], value=np.array([1.3207, 0.75717, 0.75717, 1.3207]))
>>> factor_b_c = DiscreteFactor(['B', 'C'], cardinality=[2, 2], value=np.array([0.00024189, 4134.2, 4134.2, 0.0002418]))
>>> factor_c_d = DiscreteFactor(['C', 'D'], cardinality=[2, 2], value=np.array([0.0043227, 231.34, 231.34, 0.0043227]))
>>> factor_d_e = DiscreteFactor(['E', 'F'], cardinality=[2, 2], value=np.array([31.228, 0.032023, 0.032023, 31.228]))
>>> student.add_factors(factor_a, factor_b, factor_c, factor_d, factor_e, factor_f, factor_a_b,
...    factor_b_c, factor_c_d, factor_d_e)
>>> mplp = Mplp(student)
>>> result = mplp.map_query()
Return: {'B': 0.93894, 'C': 1.121, 'A': 1.8323, 'F': 1.5093, 'D': 1.7765, 'E': 2.12239}

Dynamic Bayesian Network Inference

Elimination Ordering

class pgmpy.inference.EliminationOrder.BaseEliminationOrder(model)[source]

Base class for finding elimination orders.

cost(node)[source]

The cost function to compute the cost of elimination of each node. This method is just a dummy and returns 0 for all the nodes. Actual cost functions are implemented in the classes inheriting BaseEliminationOrder.

node: string, any hashable python object.
The node whose cost is to be computed.
fill_in_edges(node)[source]

Return edges needed to be added to the graph if a node is removed.

node: string (any hashable python object)
Node to be removed from the graph.
get_elimination_order(nodes=None)[source]

Returns the optimal elimination order based on the cost function. The node having the least cost is removed first.

nodes: list, tuple, set (array-like)
The variables which are to be eliminated.
>>> import numpy as np
>>> from pgmpy.models import BayesianModel
>>> from pgmpy.factors.discrete import TabularCPD
>>> from pgmpy.inference.EliminationOrder import MinFill
>>> model = BayesianModel([('c', 'd'), ('d', 'g'), ('i', 'g'),
...                        ('i', 's'), ('s', 'j'), ('g', 'l'),
...                        ('l', 'j'), ('j', 'h'), ('g', 'h')])
>>> cpd_c = TabularCPD('c', 2, np.random.rand(2, 1))
>>> cpd_d = TabularCPD('d', 2, np.random.rand(2, 2),
...                   ['c'], [2])
>>> cpd_g = TabularCPD('g', 3, np.random.rand(3, 4),
...                   ['d', 'i'], [2, 2])
>>> cpd_i = TabularCPD('i', 2, np.random.rand(2, 1))
>>> cpd_s = TabularCPD('s', 2, np.random.rand(2, 2),
...                   ['i'], [2])
>>> cpd_j = TabularCPD('j', 2, np.random.rand(2, 4),
...                   ['l', 's'], [2, 2])
>>> cpd_l = TabularCPD('l', 2, np.random.rand(2, 3),
...                   ['g'], [3])
>>> cpd_h = TabularCPD('h', 2, np.random.rand(2, 6),
...                   ['g', 'j'], [3, 2])
>>> model.add_cpds(cpd_c, cpd_d, cpd_g, cpd_i, cpd_s, cpd_j,
...                cpd_l, cpd_h)
>>> WeightedMinFill(model).get_elimination_order(['c', 'd', 'g', 'l', 's'])
['c', 's', 'l', 'd', 'g']
>>> WeightedMinFill(model).get_elimination_order(['c', 'd', 'g', 'l', 's'])
['c', 's', 'l', 'd', 'g']
>>> WeightedMinFill(model).get_elimination_order(['c', 'd', 'g', 'l', 's'])
['c', 's', 'l', 'd', 'g']