I am teaching myself about Bayesian graphical networks. I'm attempting to use the python package pgmpy to generate the networks in python. This seems like a great resource.

For my first test, I generated a simple network depicted below (I set the known probabilities and conditional probabilities to infer the unconditional probabilities): enter image description here

Now, I entered the probabilities for $A$ and $B$, as well as the probabilities for $P(C|A,B)$ into a bayesian graphical model structure in pgmpy:

#These are based on the the Monte Hall example found at https://github.com/pgmpy/pgmpy/blob/dev/examples/Monte%20Hall%20Problem.ipynb

from pgmpy.models import BayesianModel
from pgmpy.factors import TabularCPD

# Defining the network structure
model = BayesianModel([('A', 'C'), ('B', 'C')])

# Defining the CPDs:
cpd_p = TabularCPD('A', 2, [[0.99, 0.01]])
cpd_a = TabularCPD('B', 2, [[0.9, 0.1]])
cpd_t = TabularCPD('C', 2, [[0.9, 0.5, 0.4, 0.1], 
                                  [0.1, 0.5, 0.6, 0.9]],
                  evidence=['A', 'B'], evidence_card=[2, 2])

# Associating the CPDs with the network structure.
model.add_cpds(cpd_p, cpd_a, cpd_t)

# Some other methods

[<TabularCPD representing P(A:2) at 0x10e24cfd0>,
 <TabularCPD representing P(B:2) at 0x10e24cf10>,
 <TabularCPD representing P(C:2 | A:2, B:2) at 0x10df9a750>]

However, when I calculate the probabilities in pgmpy I get:

# Infering the posterior probability 
from pgmpy.inference import VariableElimination

print 'P(B|A=1,C=1)'
infer = VariableElimination(model)
posterior_p = infer.query(['B'], evidence={'A': 1, 'C': 1})

print 'P(B|C=1)'
posterior_p = infer.query(['B'], evidence={'C': 1})

print 'probs'
posterior_p = infer.query(['B','C','A'])
for entry in posterior_p:
    print posterior_p[entry]

| B   |   phi(B) |
| B_0 |   0.8333 |
| B_1 |   0.1667 |
| B   |   phi(B) |
| B_0 |   0.6082 |
| B_1 |   0.3918 |
| A   |   phi(A) |
| A_0 |   0.9900 |
| A_1 |   0.0100 |
| C   |   phi(C) |
| C_0 |   0.8461 |
| C_1 |   0.1539 |
| B   |   phi(B) |
| B_0 |   0.9000 |
| B_1 |   0.1000 |

However, when I calculate the conditional probabilities (or total probabilities for variable $C$) shown above by hand (calculated from the unconditional probabilities), I get different answers:

$$P(B=1|A=1,C=1) = \frac{P(B=1|A=1,C=1)}{P(B=0|A=1,C=1)+P(B=1|A=1,C=1)} = \frac{0.0009}{0.0054+0.0009} = 0.143$$

$$P(B=1|C=1) = \frac{P(B=1,C=1)}{P(C=1)} = \frac{0.0495+0.0009}{0.0891+0.0495+0.0054+0.0009} = 0.348$$

$$P(C=1) = 0.0891 + 0.0495 + 0.0054 + 0.0009 = 0.145 $$

Obviously this is different than the probabilities output by pgmpy.

Does anyone have a clue about where I'm going wrong? Either with my "by-hand" calculations or with the coding.



Seems as thought I figured this out. I'm posting it here just in case it helps someone.

Well, it turns out that I wasn't making an error with my "by-hand" calculation, but it was indeed a problem with the pgmpy package. For some reason, it was inferring independence between $A$ and $B$ given $C$ (they should be independent without knowing C).

I updated the package to the developer edition and it gives the desired result.

  • $\begingroup$ Thank you SO much for answering this here, I've been pulling my freaking HAIR OUT trying to figure out why pgmpy is giving me the wrong conditional probabilities! I updated to dev like you said, and it works now. It's really unacceptable to have such a major issue in the production version of such a widely used package. It was causing some seriously egregious errors. See here: stackoverflow.com/questions/47087068/… $\endgroup$ – jon_simon Nov 3 '17 at 5:37

I tried installing the new package and the bug remains. So i looked at the Bayesian network code itself and found the bug.

Its in line 530 of BayesianModel.py when finding the descendents. Because neighbors is an iterator, it is empty after visit.extend and so descendants.extend does nothing and the result is empty. I fixed it in place by casting neighbors as a list.

    Since Bayesian Networks are acyclic, this is a very simple dfs
    which does not remember which nodes it has visited.
    descendents = []
    visit = [node]
    while visit:
        n = visit.pop()
        neighbors = self.neighbors(n)
    return descendents

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.