car_armando.jl

using Omega

# Beta distributed random variables in [0, 1]
p1 = beta(0.5)
pgood = beta(0.5)
pbad = beta(0.5)

# Function that simulates the car
function simulate_car(ω)
  car_is_good = Bernoulli(ω, p1(ω), Bool)
  cuts_you_off = ifelse(prob_car_is_good, Bernoulli(ω, pgood(ω), Bool), Bernoulli(ω, pbad, Bool))
  drive_sim = some_function(car_is_good, cuts_you_off)
  return drive_sim
end

# n data points
N = 10

# Random variable over vector of drive_sim
simulations  = rv(ω -> [simulate_car_(ω) for i = 1:n])

# Data
data = ...

# Posterior of pgood, pbad, p1
rand(cond((pgood, pbad, p1), data == simulations))

# Is it better (does variance decrease)
var(cond(pbad, data == simulations)) - var(cond(pbad, (data == simulations) & pgood > pbad))