Hi all,
I am a big fan of both libraries, but when I try to include copulas in probabilistic (Turing) models, I often get errors related to quantities being out of bounds. I have included a simple example below.
The below Turing model runs fine 
using Distributions, Random, Copulas
μ₁_prior = Normal(3.0, 1/2)
σ₁_prior = Normal(0, 1/2) |> dist -> truncated(dist, lower = 0)
μ₂_prior = Normal(1.15, 1/3)
σ₂_prior = Normal(0, 1/3) |> dist -> truncated(dist, lower = 0)
ρ_prior = Normal(-1/3, 1/4) |> dist -> truncated(dist, lower = -1, upper = 1)
@model function multi_var_meas_error(; x_meas, ϵ::Float64)
# starting points (priors)
ρ ~ ρ_prior
μ₁ ~ μ₁_prior; σ₁ ~ σ₁_prior
μ₂ ~ μ₂_prior; σ₂ ~ σ₂_prior
# Likelihood
# Defining marginal distributions and copula
Σ = [
σ₁^2 ρ * σ₁ * σ₂
ρ * σ₁ * σ₂ σ₂^2
]
# For numerical stability, and a positive-definite covariance matrix:
Σ = [i == j ? Σ[i, j] + eps() : Σ[i, j] for i ∈ 1:2, j ∈ 1:2]
# The input x_meas includes the measurement error ϵ
# Estimate the true levels without the measurement error
x₁ ~ Normal(x_meas[1], ϵ) |> d -> truncated(d, lower = 0)
x₂ ~ Normal(x_meas[2], ϵ) |> d -> truncated(d, lower = 0)
Turing.@addlogprob! logpdf(
MvNormal([μ₁, μ₂], Σ),
[x₁, x₂])
# Turing.@addlogprob! logpdf(
# SklarDist(
# GaussianCopula(Σ),
# (Normal(μ₁, σ₁) |> d -> truncated(d, lower = 0),
# Normal(μ₂, σ₂) |> d -> truncated(d, lower = 0))
# ),
# [x₁, x₂])
end
n_draws = 1_000; n_warmup = 10_000; n_chains = 4; target_acc = 0.65;
mcmc_sampler = Turing.NUTS(n_warmup, target_acc)
posterior = multi_var_meas_error(x_meas = [2.5, 1.0], ϵ = 0.2) |>
x -> sample(x, mcmc_sampler, MCMCThreads(),
n_draws, n_chains, seed = 231123)
…however, when I comment out the the MvNormal, and uncomment the Sklar distribution with the Gaussian copula (using the same covariance matrix)…
...
# Turing.@addlogprob! logpdf(
# MvNormal([μ₁, μ₂], Σ),
# [x₁, x₂])
Turing.@addlogprob! logpdf(
SklarDist(
GaussianCopula(Σ),
(Normal(μ₁, σ₁) |> d -> truncated(d, lower = 0),
Normal(μ₂, σ₂) |> d -> truncated(d, lower = 0))
),
[x₁, x₂])
end
…then I get the following error:
ERROR: TaskFailedException
nested task error: TaskFailedException
nested task error: DomainError with 1.0000000000000058:
`abs(x)` cannot be greater than 1.
and your example did trigger that.
, let us check that it works and then I’ll merge & push onto the registry.