From 101e8004784359189c98a2421f6f5b21589b4507 Mon Sep 17 00:00:00 2001 From: Oskar Laverny Date: Tue, 24 Oct 2023 16:46:36 +0200 Subject: [PATCH] Adding special routes for the independent copulas sampling cdf and pdf, fix a few typos, add a test --- src/ArchimedeanCopulas/IndependentCopula.jl | 19 ++++-- src/Copulas.jl | 4 +- src/MiscellaneousCopulas/MCopula.jl | 2 +- test/sampling_test | 71 +++++++++++++++++++++ 4 files changed, 88 insertions(+), 8 deletions(-) create mode 100644 test/sampling_test diff --git a/src/ArchimedeanCopulas/IndependentCopula.jl b/src/ArchimedeanCopulas/IndependentCopula.jl index 014c5970..cbf237cc 100644 --- a/src/ArchimedeanCopulas/IndependentCopula.jl +++ b/src/ArchimedeanCopulas/IndependentCopula.jl @@ -19,13 +19,22 @@ It happends to be an Archimedean Copula, with generator : ``` """ struct IndependentCopula{d} <: ArchimedeanCopula{d} end -IndependentCopula(d) = IndependentCopula{d} -function Distributions._logpdf(::IndependentCopula{d},u) where d - return all(0 .<= u .<= 1) ? 1 : 0 +IndependentCopula(d) = IndependentCopula{d}() +function Distributions._logpdf(C::IndependentCopula{d},u) where d + return all(0 .<= u .<= 1) ? zero(eltype(u)) : -Inf end -function Distributions.cdf(::IndependentCopula{d},u) where d - return all(0 .<= u .<= 1) ? prod(u) : 0 +function Distributions.cdf(C::IndependentCopula{d},u) where d + return prod(u) end ϕ(::IndependentCopula,t) = exp(-t) ϕ⁻¹(::IndependentCopula,t) = -log(t) τ(::IndependentCopula) = 0 + +# Exceptionally we overload the functions as we dont want to take the slow route of archemedean copulas for the independant copula. + +function Distributions._rand!(rng::Distributions.AbstractRNG, C::IndependentCopula{d}, x::AbstractVector{T}) where {T<:Real,d} + Random.rand!(rng,x) +end +function Base.rand(rng::Distributions.AbstractRNG,C::IndependentCopula{d}) where {d} + rand(rng,length(C)) +end \ No newline at end of file diff --git a/src/Copulas.jl b/src/Copulas.jl index 6b16b457..eca3ed81 100644 --- a/src/Copulas.jl +++ b/src/Copulas.jl @@ -61,8 +61,8 @@ end Copulas include("MiscellaneousCopulas/WCopula.jl") include("MiscellaneousCopulas/SurvivalCopula.jl") export MCopula, - Wcopula, - SurvivalCopula + WCopula, + SurvivalCopula # PrecompileTools stuff @setup_workload begin diff --git a/src/MiscellaneousCopulas/MCopula.jl b/src/MiscellaneousCopulas/MCopula.jl index d7b50d66..1de15040 100644 --- a/src/MiscellaneousCopulas/MCopula.jl +++ b/src/MiscellaneousCopulas/MCopula.jl @@ -1,6 +1,6 @@ struct MCopula{d} <: Copula{d} end MCopula(d) = MCopula{d}() -Distributions.cdf(::MCopula{d},u) where {d} = min(u) +Distributions.cdf(::MCopula{d},u) where {d} = minimum(u) function Distributions._rand!(rng::Distributions.AbstractRNG, ::MCopula{d}, x::AbstractVector{T}) where {d,T<:Real} x .= rand(rng) end diff --git a/test/sampling_test b/test/sampling_test new file mode 100644 index 00000000..3985b49b --- /dev/null +++ b/test/sampling_test @@ -0,0 +1,71 @@ +# @testitem "test default parametrisation of every concrete copula types" begin +# for T in subtypes(Copulas) +# if isconcretetype(T) +# rand(T(),100) +# @test true +# end +# end +# end + +# @testitem "test sampling of every concrete copula type" begin +# for T in subtypes(Copulas) +# if isconcretetype(T) +# rand(T(),100) +# @test true +# end +# end +# end + + +@testitem "small sampling test" begin + # test constructed from https://github.com/lrnv/Copulas.jl/issues/35 + using Copulas, Distributions + rand(GaussianCopula([1.0 0.5; 0.5 1.0]),1) + rand(IndependentCopula(2),1) + rand(MCopula(2),1) + rand(WCopula(2),1) + + @test true +end + + +@testitem "test cdf for three copulas." begin + # test constructed from https://github.com/lrnv/Copulas.jl/issues/35 + using Copulas, Distributions + + u = range(0, stop=1, length=100) + + for G in ( + GaussianCopula([1.0 0.5; 0.5 1.0]), + IndependentCopula(2), + MCopula(2), + WCopula(2) + ) + for uᵢ in u + for vᵢ in u + cdf(G,[uᵢ,vᵢ]) + end + end + end + @test true +end + +@testitem "test pdf for three copulas." begin + # test constructed from https://github.com/lrnv/Copulas.jl/issues/35 + using Copulas, Distributions + + u = range(0, stop=1, length=100) + + for G in ( + GaussianCopula([1.0 0.5; 0.5 1.0]), + IndependentCopula(2), + # MCopula(2) + ) + for uᵢ in u + for vᵢ in u + pdf(G,[uᵢ,vᵢ]) + end + end + end + @test true +end \ No newline at end of file