using Turing

@model function linear_regression(x)
    # Priors
    α ~ Normal(0, 1)
    β ~ Normal(0, 1)
    σ² ~ truncated(Normal(0, 1), lower=0)

    # Likelihood
    μ = α .+ β .* x
    y ~ MvNormal(μ, σ² * I)
end

x, y = rand(10), rand(10)
posterior = linear_regression(x) | (; y = y)
chain = sample(posterior, NUTS(), 1000)