Given a model defined in Turing I would like to calculate the joint probability function over all the variables. Say I have the model

```
@model function model(y)
mu ~ Normal(0, 1)
y ~ Normal(mu, 1)
end
y_obs = 1
mu = 1
```

With the prob macro I can do

```
prob"mu = mu, y = y_obs | model = model"
```

and it gives me exactly what I want. However, now I would like to be able to not have to explicitly pass in the value of `y_obs`

and instead do

```
m = model(y_obs)
prob"mu = mu | model = m"
```

Is that somehow possible?

My end goal is to create a function which only takes as input an instantiated model and returns the joint probability function. So I want to have a function `get_joint_prob`

such that:

```
m = model(y_obs)
joint_prob = get_joint_prob(m)
joint_prob((mu=mu,)) == prob"mu = mu, y = y_obs | model = model"
```

(`joint_prob`

takes in a named tuple so that it can deal with multiple variables if necessary).

I have already read through the documentation on the Turing side and looked a bit at the implementation of the prob macro and the VarInfo stuff but I didn’t quite understand all of it well enough to get something to work. Any pointers would be really helpful!