diff --git a/Project.toml b/Project.toml index 0fe1c4e0..5f9b79c2 100644 --- a/Project.toml +++ b/Project.toml @@ -14,6 +14,7 @@ Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" TaylorSeries = "6aa5eb33-94cf-58f4-a9d0-e4b2c4fc25ea" +WilliamsonTransforms = "48feb556-9bdd-43a2-8e10-96100ec25e22" [compat] Cubature = "1.5" @@ -21,12 +22,13 @@ Distributions = "0.25" ForwardDiff = "0.10" GSL = "1" MvNormalCDF = "0.2, 0.3" +Random = "1.6" Roots = "1, 2" SpecialFunctions = "2" StatsBase = "0.33, 0.34" TaylorSeries = "0.12, 0.13, 0.14, 0.15" +WilliamsonTransforms = "0.0" julia = "1.6" -Random = "1.6" [extras] InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240" diff --git a/src/ArchimedeanCopulas/WilliamsonCopula.jl b/src/ArchimedeanCopulas/WilliamsonCopula.jl new file mode 100644 index 00000000..23c3b029 --- /dev/null +++ b/src/ArchimedeanCopulas/WilliamsonCopula.jl @@ -0,0 +1,28 @@ +struct WilliamsonCopula{d,Tϕ,TX} <: ArchimedeanCopula{d} + ϕ::Tϕ + X::TX +end +function WilliamsonCopula(X::Distributions.UnivariateDistribution, d) + ϕ = WilliamsonTransforms.𝒲(X,d) + return WilliamsonCopula{d,typeof(ϕ),typeof(X)}(ϕ,X) +end +function WilliamsonCopula(ϕ::Function, d) + X = WilliamsonTransforms.𝒲₋₁(ϕ,d) + return WilliamsonCopula{d,typeof(ϕ),typeof(X)}(ϕ,X) +end +function WilliamsonCopula(ϕ::Function, X::Distributions.UnivariateDistribution, d) + return WilliamsonCopula{d,typeof(ϕ),typeof(X)}(ϕ,X) +end +function Distributions._rand!(rng::Distributions.AbstractRNG, C::CT, x::AbstractVector{T}) where {T<:Real, CT<:WilliamsonCopula} + r = rand(rng,C.X) + Random.rand!(rng,x) + for i in 1:length(C) + x[i] = -log(x[i]) + end + sx = sum(x) + for i in 1:length(C) + x[i] = ϕ(C,r * x[i]/sx) + end + return x +end +ϕ( C::WilliamsonCopula, t) = C.ϕ(t) \ No newline at end of file diff --git a/src/Copulas.jl b/src/Copulas.jl index cdcff2ac..99e852c6 100644 --- a/src/Copulas.jl +++ b/src/Copulas.jl @@ -18,6 +18,7 @@ end Copulas import ForwardDiff import Cubature import MvNormalCDF + import WilliamsonTransforms # Standard copulas and stuff. include("utils.jl") @@ -59,12 +60,13 @@ end Copulas include("ArchimedeanCopulas/GumbelCopula.jl") include("ArchimedeanCopulas/FrankCopula.jl") include("ArchimedeanCopulas/AMHCopula.jl") + include("ArchimedeanCopulas/WilliamsonCopula.jl") export IndependentCopula, ClaytonCopula, JoeCopula, GumbelCopula, FrankCopula, - AMHCopula - + AMHCopula, + WilliamsonCopula end diff --git a/test/williamson_test.jl b/test/williamson_test.jl new file mode 100644 index 00000000..f1ef321e --- /dev/null +++ b/test/williamson_test.jl @@ -0,0 +1,18 @@ + +@testitem "williamson test" begin + using Distributions, Random + Random.seed!(123) + taus = [0.0, 0.1, 0.5, 0.9, 1.0] + + ϕ_clayton(x, θ) = max((1 + θ * x),zero(x))^(-1/θ) + + Cops = ( + WilliamsonCopula(Dirac(1),10), + WilliamsonCopula(x -> exp(-x),10), + WilliamsonCopula(x -> ϕ_clayton(x,2),2), + WilliamsonCopula(x -> ϕ_clayton(x,-0.3),2) + ) + for C in Cops + x = rand(C,1000) + end +end