# Feedback request: \`rand\` and \`predict\` semantics

**URL:** <https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977>\
**Category:** Probabilistic Programming\
**Created:** [June 18, 2022, 5:26pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977 "2022-06-18T17:26:39Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)\
**Post date:** [June 18, 2022, 5:26pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/1 "2022-06-18T17:26:39Z")

</div>

Hi, I’m working on [Tilde.jl](https://github.com/cscherrer/Tilde.jl), and thinking some more about how `rand` and `predict` should be set up.

In current development, this works:

```julia
julia> m = @model n begin
           μ ~ Normal(σ=10)
           σ ~ Exponential()
           y ~ Normal(μ, σ) ^ n
       end;

julia> obs = (y = [3.0, 4.0],);

julia> post = m(2) | obs
ModelPosterior given
    arguments (:n,)
    observations (:y,)
@model n begin
        μ ~ Normal(σ = 10)
        σ ~ Exponential()
        y ~ Normal(μ, σ) ^ n
    end

julia> r = rand(post)
(μ = 4.99758715428609, σ = 1.582563929898804)

julia> predict(post, r)
(y = [3.6962760705758146, 3.661799255763297],)

```

I think this makes sense for this example. But what if you had, say `y = [3.0, missing]`? I think conditioning on `missing` should be the same as not conditioning at all. So for arrays, conditioning on `y` ought to really mean conditioning on _the known elements_ of `y`, so `rand` should sample the others.

So this would look something like

```julia
julia> rand(m(2) | (y = [3.0, missing]))
(μ = 4.99758715428609, σ = 1.582563929898804, y = [3.0, 6.346809940724782])

```

But it seems weird to have a completely different return type when `y` is fully observed. So maybe the fully-observed case should instead work like this:

```julia
julia> rand(post)
(μ = 4.99758715428609, σ = 1.582563929898804, y = [3.0, 4.0])

```

This is equivalent to calling `rand` on the model where `y ~ Normal(μ, σ) ^ n` is replaced with `y ~ Dirac([3.0, 4.0])`, or could also be seen as calling `rand` on the original (unconditioned) model, and updating that result with any observations.

## `logdensity` / `rand` interaction

This seems fine, but we need to be sure we have consistent semantics. Since `rand` and `logdensityof` are both defined for a `post`, we should be able to call

```julia
logdensityof(post, rand(post))

```

So… maybe the semantics are that `logdensityof` can be called without including `Dirac` components, and if they _are_ given they must match what the model specifies. Does this make sense?

## `predict` / `rand` interaction

Currently we can do

```julia
julia> predict(post, rand(post))
(y = [17.638620683052217, 17.891076168944927],)

```

What does this look like if `rand` also has a `y` slot? I think this could be something like

> Among all sampled (`~`) variables in a posterior model `post`, `predict(post, pars)` samples those that are either (1) not included in `pars` or (2) _are_ conditioned upon in `post`

## Thoughts?

I think it’s important to get this right, so I’d appreciate any feedback on this. Thanks 🙂

---

<div class="post-metadata">

**Author:** ![EvoArt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/evoart/32/25357_2.png) [@EvoArt](https://discourse.julialang.org/u/EvoArt)\
**Post date:** [June 18, 2022, 6:04pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/2 "2022-06-18T18:04:14Z")

</div>

"and if they are given they must match what the model specifies. "

Presumably this would throw an error if they don’t. It’s tempting to return -Inf in this case, but I can’t think of a real situation where people pass in the wrong y values and it’s not just an error.

---

<div class="post-metadata">

**Author:** ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)\
**Post date:** [June 18, 2022, 6:40pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/3 "2022-06-18T18:40:54Z")

</div>

I think rather than decide what happens if they don’t match, it’s helpful to ask what `rand` means on a posterior.

If it means “sample form the unconditioned model, except for observed values which are treated as Dirac” then the case where the input doesn’t match that needs to return a log-density of `-Inf`, to match Dirac semantics.

Maybe it should mean something different, but I don’t have any alternatives yet.

---

<div class="post-metadata">

**Author:** ![marius311](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/marius311/32/3953_2.png) [@marius311](https://discourse.julialang.org/u/marius311)\
**Post date:** [June 18, 2022, 6:50pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/4 "2022-06-18T18:50:15Z")

</div>

My (maybe bad) 2c as a casual onlooker but I was indeed curious what `rand(post)` meant, not knowing for Tilde and not being able to remember for other PPLs. It almost makes the most sense to me that it would mean a sample from the posterior, since I think of the conditioned posterior as some implicit distribution, and `rand` samples distributions. But ofc this requires sophisticated sampling. So one thing it could be is an error unless its trivially sampleable. Not sure if that makes any of these other downstream issues irrelevant.

---

<div class="post-metadata">

**Author:** ![EvoArt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/evoart/32/25357_2.png) [@EvoArt](https://discourse.julialang.org/u/EvoArt)\
**Post date:** [June 18, 2022, 7:05pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/5 "2022-06-18T19:05:05Z")

</div>

I would tend to agree that calling rand() on the posterior makes me think of sampling from the posterior.

---

<div class="post-metadata">

**Author:** ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)\
**Post date:** [June 18, 2022, 7:55pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/6 "2022-06-18T19:55:34Z")

</div>

This is a great point, makes me think we should just disallow `rand` on conditional measures. You can always get the prior and call `rand` on that if. The extra step is easy, and it makes things a lot clearer.

---

<div class="post-metadata">

**Author:** ![mschauer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mschauer/32/13946_2.png) [@mschauer](https://discourse.julialang.org/u/mschauer)\
**Post date:** [June 19, 2022, 8:28am UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/7 "2022-06-19T08:28:59Z")

</div>

Not so much disallow but only define it if sampling the conditional is day O(1) and then rand should give such samples

---

<div class="post-metadata">

**Author:** ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)\
**Post date:** [June 19, 2022, 1:12pm UTC](https://discourse.julialang.org/t/feedback-request-rand-and-predict-semantics/82977/8 "2022-06-19T13:12:03Z")

</div>

@Tamas_Papp made a similar point a while back, and I agree. To try to nail this down a little, I think a sampling procedure should only be callable using `rand` if it’s

1. “fast” in some sense ^\dagger
2. exact
3. iid

^\dagger though O(1) might be too restrictive, for example `MvNormal` requires matrix-vector multiply