# Condition Turing model on a variable created within the model

**URL:** <https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994>\
**Category:** Probabilistic Programming\
**Tags:** turing\
**Created:** [July 13, 2024, 10:53am UTC](https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994 "2024-07-13T10:53:32Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![DominiqueMakowski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dominiquemakowski/32/51410_2.png) [@DominiqueMakowski](https://discourse.julialang.org/u/DominiqueMakowski)\
**Post date:** [July 13, 2024, 10:53am UTC](https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994/1 "2024-07-13T10:53:32Z")

</div>

Assuming I have the following variable [0, 1]

```julia
using Downloads, CSV, DataFrames, Random
using Turing, Distributions
using GLMakie
using StatsFuns: logistic

y = vcat(rand(Beta(2, 10), 200),
    rand(Uniform(0, 1), 100),
    rand(Beta(10, 2), 200))
hist(y)

```

![image](https://global.discourse-cdn.com/julialang/original/3X/4/0/403d4b10178b1e13dfcdc377316b374421b80516.png)

I would like to know if it was possible, within a Turing model, to transform this variable (dichotomize it to 0 and 1 if y \> 0.5), and then fit a logistic model to this dichotomized variable.

I tried the following:

```julia
@model function model_choice(y)
    # Prior on probability
    choice_intercept ~ Normal(0, 1)

    for i in 1:length(y)
        # Logistic function
        v = logistic(choice_intercept)

        # Transform to binary choice
        choice = ifelse(y[i] > 0.5, 1, 0)

        # Likelihood    
        choice ~ Bernoulli(v)
    end
end

fit = model_choice(y)
chains = sample(fit, NUTS(), 400)

```

While this successfully samples:

```julia
Summary Statistics
        parameters mean std mcse ess_bulk ess_tail rhat ess_per_sec 
            Symbol Float64 Float64 Float64 Float64 Float64 Float64 Float64

  choice_intercept -0.5054 0.0000 0.0000 1.4125 20.8078 2.0369 0.1158
            choice 1.0000 0.0000 NaN NaN NaN NaN NaN

```

I am not sure I am doing this correctly as 1) `choice` shows up in the results as a parameter and it probably shouldn’t be there; and 2) I can’t use `predict()` on that model:

```julia
pred = predict(model_choice([(missing) for i in 1:length(y)]), chains)

```

```julia
ERROR: MethodError: no method matching ifelse(::Missing, ::Int64, ::Int64)

Closest candidates are:
  ifelse(::Bool, ::Any, ::Any)
   @ Base essentials.jl:647

```

I am assuming I am not doing this correctly. What is the best way to do that? Thanks for any pointers!

---

<div class="post-metadata">

**Author:** ![p-gw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/p-gw/32/210518_2.png) [@p-gw](https://discourse.julialang.org/u/p-gw)\
**Post date:** [July 15, 2024, 1:50pm UTC](https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994/2 "2024-07-15T13:50:09Z")

</div>

This fails because `missing` is passed to `ifelse` which is not defined:

```julia
ifelse(missing, 0, 1)

```

```julia
ERROR: MethodError: no method matching ifelse(::Missing, ::Int64, ::Int64)

```

Since your model does not depend on `y` but only `choice` it’s probably easiest to pass it directly.

```julia
using Distributions, Turing

y = vcat(rand(Beta(2, 10), 200),
    rand(Uniform(0, 1), 100),
    rand(Beta(10, 2), 200))

choice = y .> 0.5

@model function choice_model(choice)
    choice_intercept ~ Normal(0, 1)
    
    for i in eachindex(choice)
        choice[i] ~ BernoulliLogit(choice_intercept)
    end
end

fit = choice_model(choice)
chains = sample(fit, NUTS(), 400)

```

```julia
Summary Statistics
        parameters mean std mcse ess_bulk ess_tail rhat ess_per_sec 
            Symbol Float64 Float64 Float64 Float64 Float64 Float64 Float64 

  choice_intercept -0.0185 0.0905 0.0063 211.6200 194.0690 1.0137 25.6292

Quantiles
        parameters 2.5% 25.0% 50.0% 75.0% 97.5% 
            Symbol Float64 Float64 Float64 Float64 Float64 

  choice_intercept -0.1860 -0.0812 -0.0138 0.0509 0.1486

```

```julia
pred = predict(choice_model([missing for _ in 1:3]), chains)

```

```julia
Summary Statistics
  parameters mean std mcse ess_bulk ess_tail rhat ess_per_sec 
      Symbol Float64 Float64 Float64 Float64 Float64 Float64 Missing 

   choice[1] 0.5000 0.5006 0.0244 420.1271 NaN NaN missing
   choice[2] 0.5200 0.5002 0.0272 337.7291 NaN 1.0011 missing
   choice[3] 0.5175 0.5003 0.0251 397.3654 NaN 1.0017 missing

Quantiles
  parameters 2.5% 25.0% 50.0% 75.0% 97.5% 
      Symbol Float64 Float64 Float64 Float64 Float64 

   choice[1] 0.0000 0.0000 0.5000 1.0000 1.0000
   choice[2] 0.0000 0.0000 1.0000 1.0000 1.0000
   choice[3] 0.0000 0.0000 1.0000 1.0000 1.0000

```

---

<div class="post-metadata">

**Author:** ![DominiqueMakowski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dominiquemakowski/32/51410_2.png) [@DominiqueMakowski](https://discourse.julialang.org/u/DominiqueMakowski)\
**Post date:** [July 15, 2024, 2:54pm UTC](https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994/3 "2024-07-15T14:54:26Z")

</div>

**TLDR; How to, in a model predicting `y`, also estimate another parameter based on a derivative (a transformation done in the model) of `y`**

Thanks! Unfortunately this is just a minimal working example, but in reality the model is slightly more complex, and creating a variable to condition (the “choice” here) inside the model would be necessary.

Regardless of the `predict()` method not working (as this might be caused by the model not returning anything) - which is a somewhat secondary issue, how to correctly define the variable inside, as the version below is probably wrong, as evidenced by the presence of the `choice` parameter in the output and very high and nonsensical parameter values for `choice_intercept`

In the following model, I also condition some other parameters on the outcome variables.

```julia
using Downloads, CSV, DataFrames, Random
using Turing, Distributions
using GLMakie
using StatsFuns: logistic

y = vcat(rand(Beta(2, 10), 200),
    rand(Uniform(0, 1), 100),
    rand(Beta(10, 2), 200))

@model function model_choice(y)
    choice_intercept ~ Normal(0, 1)
    μ ~ Normal(0.5, 0.2)
    σ ~ truncated(Normal(0, 1); lower=0)

    for i in 1:length(y)
        # Transform to binary choice
        choice = ifelse(y[i] > 0.5, 1, 0)
        # Likelihood    
        choice ~ BernoulliLogit(choice_intercept)
        y[i] ~ Normal(μ, σ)
    end
end

fit = model_choice(y)
chains = sample(fit, NUTS(), 400)

```

But the convergence is bad…

---

<div class="post-metadata">

**Author:** ![p-gw](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/p-gw/32/210518_2.png) [@p-gw](https://discourse.julialang.org/u/p-gw)\
**Post date:** [July 15, 2024, 4:10pm UTC](https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994/4 "2024-07-15T16:10:35Z")

</div>

The proper way to use transformed data is to either transform it outside the model and pass it as an argument or for simple transformations like this code the transformation as the default value for the function argument:

```julia
@model function model_choice(y, choice = y .> 0.5)
    choice_intercept ~ Normal(0, 1)
    μ ~ Normal(0.5, 0.2)
    σ ~ truncated(Normal(0, 1); lower=0)

    for i in eachindex(y)
        choice[i] ~ BernoulliLogit(choice_intercept)
        y[i] ~ Normal(μ, σ)
    end
end

fit = model_choice(y)
chains = sample(fit, NUTS(), 400)

```

Whatever your approach, `choice` must be a part of the model function definition. Otherwise Turing.jl will treat it as a parameter and it will show up in your `chains` object.

---

<div class="post-metadata">

**Author:** ![DominiqueMakowski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dominiquemakowski/32/51410_2.png) [@DominiqueMakowski](https://discourse.julialang.org/u/DominiqueMakowski)\
**Post date:** [July 15, 2024, 4:19pm UTC](https://discourse.julialang.org/t/condition-turing-model-on-a-variable-created-within-the-model/116994/5 "2024-07-15T16:19:17Z")

</div>

Got it! Thanks 🙂