How can one represent distinct non-numeric symbols in sympy?

Viewed 77

I am experimenting with the representation of a trivial statistics problem in Sympy:

For a sample space S, there are 6 possible distinct outcomes a,b,c,d,e,f. We can define event A as having occurred if any of a,b,c have, and event B as having ocurred if any of b,c,d have. Given a probability mass function pmf defined over S, what is the probability of event A?

When attempting to implement this symbolically, I receive the following error:

~/project/.envs/dev/lib/python3.6/site-packages/sympy/stats/frv.py in _test(self, elem)
    164         elif val.is_Equality:
    165             return val.lhs == val.rhs
--> 166         raise ValueError("Undecidable if %s" % str(val))
    167 
    168     def __contains__(self, other):

ValueError: Undecidable if Eq(d, a) | Eq(d, b) | Eq(d, c)

The problem is implemented as follows with comments on the failing lines of code:

from sympy import Eq, Function, symbols
from sympy.logic import Or
from sympy.sets import FiniteSet, Union
from sympy.stats import FiniteRV, P

# 1. Define a sample space S with outcomes: a,b,c,d,e,f; Define events A, B
A = FiniteSet(*symbols('a b c'))
B = FiniteSet(*symbols('b c d'))
S = Union(A, B, FiniteSet(*symbols('e f')))
display("Sample Space", S)

pmfFunc = Function("pmf")
pmfDict = {v: pmfFunc(v) for v in S}
X = FiniteRV('X', pmfDict)

a,b = symbols('a b')
# 2. P(X = a) = pmf(a)
display(P(Eq(X,a)))

# 3. A.as_relational(X) yields `(X=a) \lor (X=b) \lor (X=c)`
display(A.as_relational(X))

# 4. P(X = a \lor X = b) = pmf(a) + pmf(b)
#    - Actual Output: ValueError: Undecidable if Eq(c, a) | Eq(c, b)
display(P(Or(Eq(X,a), Eq(X,b)))) # [FAILS]

# 5. P(A) = pmf(a) + pmf(b) + pmf(c)
#    - Actual Output: ValueError: Undecidable if Eq(d, a) | Eq(d, b) | Eq(d, c)
display(P(A.as_relational(X)))  # [FAILS]

I obtain expected output up to display(A.as_relational(X)):

enter image description here

Interpreting the failure message suggests that Sympy is unable to tell that the set members are distinct. Replacing the symbols with integers resolves the error and I get output similar to what I desire.

A = FiniteSet(1, 2, 3)
B = FiniteSet(2, 3, 4)
S = Union(A, B, FiniteSet(5, 6))

enter image description here

If I am not misunderstanding the error or the fundamental use of the library, is there a way to tell Sympy that a collection of symbols is entirely distinct? I have attempted to replace the symbols with Dummy instances without success, and I have also attempted to leverage the assumptions module without success:

facts = [Eq(a,b) if a is b else Not(Eq(a,b)) for a, b in itertools.product(S, S)]
with assuming(*facts):

I would like to avoid confusing mappings between integers and symbolic forms, as user error may not be apparent when the results are printed as latex. I am willing to adopt some burden in a workaround (e.g., as it would have been maintaining a collection of Dummy instances), but I have yet to find an acceptable workaround.

1 Answers

Interesting question. Maybe it can be done with with assuming(Ne(a,b), ...): context but I take a more pragmatic approach: replace symbols with cos(non-zero integer) which SymPy can easily distinguish as equal or not:

>>> reps = dict(zip(var('a:f'),(cos(i) for i in range(1,7))))
>>> ireps = {v:k for k,v in reps.items()}
>>> a,b,c,d,e,f = [reps[i] for i in var('a:f')]

Then remove your a, b = symbols... line and replace display(x) with display(x.subs(ireps) to get

('Sample Space', FiniteSet(a, b, c, d, e, f))
(pmf(a),)
(Eq(X, a) | Eq(X, b) | Eq(X, c),)
(pmf(a) + pmf(b),)

(I use cos(int) instead of int because I am not sure whether any computation would result in addition of two elements and I want to make sure they stay distinct.)

Another approach would be to define a constant class that derives from Symbol:

class con(Symbol):
    def __hash__(self):
        return id(self)
    def __eq__(a,b):
        if isinstance(b, con):
            return a.name == b.name
    _eval_Eq = __eq__

a,b,c,d,e,f=map(con,'abcdef')
display=lambda*x:print(x)

from sympy import Eq, Function, symbols
from sympy.logic import Or
from sympy.sets import FiniteSet, Union
from sympy.stats import FiniteRV, P

A = FiniteSet(a,b,c)
B = FiniteSet(b,c,d)
S = Union(A, B, FiniteSet(e,f))
pmfFunc = Function("pmf")
pmfDict = {v: pmfFunc(v) for v in S}
X = FiniteRV('X', pmfDict)
display("Sample Space", S)
display(P(Eq(X,a)))
display(A.as_relational(X))
display(P(Or(Eq(X,a), Eq(X,b))))
display(P(A.as_relational(X)))

gives

('Sample Space', FiniteSet(a, b, c, d, e, f))
(pmf(a),)
(Eq(X, a) | Eq(X, b) | Eq(X, c),)
(pmf(a) + pmf(b),)
(pmf(a) + pmf(b) + pmf(c),)
Related