#rejection sampling
#P(j|-a)
import random
samples = 0
PJA = 0
PnotB = 0
iterations = 2000000
for i in range(iterations): #not number of samples, number of samples will be sampled
    rEQ = random.random()
    rB = random.random()
    rA = random.random()

    Pb = rB <= 0.001
    Peq = rEQ <= 0.002

    if Pb and Peq:
        Pa = rA <= 0.95
    elif Pb and not Peq:
        Pa = rA <= 0.94
    elif not Pb and Peq:
        Pa = rA <= 0.29
    elif not Pb and not Peq:
        Pa = rA <= 0.001

    if not Pb:
        PnotB += 1

    if Pa:
        samples += 1
    else:
        continue    #rejection

    rJ = random.random()
    PJ = rJ <= 0.9

    if  PJ:
        PJA += 1

PnotB = PnotB/iterations
PJA = 0 if samples == 0 else PJA / samples
print("P(-b)", PnotB)
print("P(j|a) =", PJA)
print("P(j|a, -b) =", PnotB*PJA)