I’m trying to create a third variable based on two prior distributions. Say that I’m trying to apply this function

```
def f(x1,x2):
p = np.zeros(3,3)
for i in range(3):
for j in range(3):
p[i,j] += x1+x2-i-j
return p
```

where x1,x2 both come from pm.Normal() with shape 5. I should be expecting the returning variable to be 5 3X3 matrices.

Is there a way of doing so? I’m aware of pymc.Deterministic() but it doesn’t seem to take in functions.