Let’s assume that I have a variable Z that is a sum of a standard normal variable X + Bernoulli variable Y.

\begin{equation} Z = X + Y, \text{where } X\sim N(0,1) \text{ and } Y\sim B(p) \end{equation}

I have a observations of Z and I want to estimate p.

Is it possible to define a Pymc estimation problem *without having to define the distributional properties* of Z, other than that it is a sum X+Y?