Migrate from Cubature.jl to HCubature.jl (#236)

Author: DilumAluthge

Project.toml

@@ -1,13 +1,13 @@
 name = "Copulas"
 uuid = "ae264745-0b69-425e-9d9d-cf662c5eec93"
 authors = ["Oskar Laverny"]
-version = "0.1.25"
+version = "0.1.26"
 
 [deps]
 Combinatorics = "861a8166-3701-5b0c-9a16-15d98fcdc6aa"
-Cubature = "667455a9-e2ce-5579-9412-b964f529a492"
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
+HCubature = "19dc6840-f33b-545b-b366-655c7e3ffd49"
 InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 LogExpFunctions = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
@@ -25,9 +25,9 @@ WilliamsonTransforms = "48feb556-9bdd-43a2-8e10-96100ec25e22"
 [compat]
 Aqua = "v0.8"
 Combinatorics = "1"
-Cubature = "1.5"
 Distributions = "0.25"
 ForwardDiff = "0.10"
+HCubature = "1"
 HypothesisTests = "v0.11"
 InteractiveUtils = "1.6"
 LinearAlgebra = "1.6"

src/Copula.jl

@@ -1,14 +1,14 @@
 abstract type Copula{d} <: Distributions.ContinuousMultivariateDistribution end
 Base.length(::Copula{d}) where d = d
 
-# The potential functions to code: 
+# The potential functions to code:
 # Distributions._logpdf
 # Distributions.cdf
 # Distributions.fit(::Type{CT},u) where CT<:Mycopula
 # Distributions._rand!
 # Base.eltype
 # τ, τ⁻¹
-# Base.eltype 
+# Base.eltype
 function Distributions.cdf(C::Copula{d},u::VT) where {d,VT<:AbstractVector}
     length(u) != d && throw(ArgumentError("Dimension mismatch between copula and input vector"))
     if any(iszero,u)
@@ -25,13 +25,13 @@ end
 function _cdf(C::CT,u) where {CT<:Copula}
     f(x) = Distributions.pdf(C,x)
     z = zeros(eltype(u),length(C))
-    return Cubature.hcubature(f,z,u,reltol=sqrt(eps()))[1]
+    return HCubature.hcubature(f,z,u,rtol=sqrt(eps()))[1]
 end
 function ρ(C::Copula{d}) where d
     F(x) = Distributions.cdf(C,x)
     z = zeros(d)
     i = ones(d)
-    r = Cubature.hcubature(F,z,i,reltol=sqrt(eps()))[1]
+    r = HCubature.hcubature(F,z,i,rtol=sqrt(eps()))[1]
     return 12*r-3
 end
 function τ(C::Copula)
@@ -40,7 +40,7 @@ function τ(C::Copula)
     return 4*r-1
 end
 function StatsBase.corkendall(C::Copula{d}) where d
-    # returns the matrix of bivariate kendall taus. 
+    # returns the matrix of bivariate kendall taus.
     K = ones(d,d)
     for i in 1:d
         for j in i+1:d
@@ -51,7 +51,7 @@ function StatsBase.corkendall(C::Copula{d}) where d
     return K
 end
 function StatsBase.corspearman(C::Copula{d}) where d
-    # returns the matrix of bivariate spearman rhos. 
+    # returns the matrix of bivariate spearman rhos.
     K = ones(d,d)
     for i in 1:d
         for j in i+1:d
@@ -71,11 +71,11 @@ function measure(C::CT, u,v) where {CT<:Copula}
     # We use a gray code according to the proposal at https://discourse.julialang.org/t/looping-through-binary-numbers/90597/6
 
     eval_pt = copy(u)
-    d = length(C) # is known at compile time 
+    d = length(C) # is known at compile time
     r = zero(eltype(u))
     graycode = 0    # use a gray code to flip one element at a time
     which = fill(false, d) # false/true to use u/v for each component (so false here)
-    r += Distributions.cdf(C,eval_pt) # the sign is always 0. 
+    r += Distributions.cdf(C,eval_pt) # the sign is always 0.
     for s = 1:(1<<d)-1
         graycode′ = s ⊻ (s >> 1)
         graycomp = trailing_zeros(graycode ⊻ graycode′) + 1

src/Copulas.jl

@@ -10,22 +10,22 @@ module Copulas
     import StatsFuns
     import TaylorSeries
     import ForwardDiff
-    import Cubature
+    import HCubature
     import MvNormalCDF
     import WilliamsonTransforms
     import Combinatorics
     import LogExpFunctions
     import QuadGK
     import LinearAlgebra
 
-    # Standard copulas and stuff. 
+    # Standard copulas and stuff.
     include("utils.jl")
     include("Copula.jl")
     include("SklarDist.jl")
     export pseudos,
            SklarDist
 
-    # Others. 
+    # Others.
     include("MiscellaneousCopulas/SurvivalCopula.jl")
     include("MiscellaneousCopulas/PlackettCopula.jl")
     include("MiscellaneousCopulas/EmpiricalCopula.jl")
@@ -41,10 +41,10 @@ module Copulas
     include("EllipticalCopula.jl")
     include("EllipticalCopulas/GaussianCopula.jl")
     include("EllipticalCopulas/TCopula.jl")
-    export GaussianCopula, 
+    export GaussianCopula,
            TCopula
 
-    # These three distributions might be merged in Distrbutions.jl one day. 
+    # These three distributions might be merged in Distrbutions.jl one day.
     include("UnivariateDistribution/Sibuya.jl")
     include("UnivariateDistribution/Logarithmic.jl")
     include("UnivariateDistribution/AlphaStable.jl")
@@ -67,7 +67,7 @@ module Copulas
     include("Generator/UnivariateGenerator/GumbelGenerator.jl")
     include("Generator/UnivariateGenerator/InvGaussianGenerator.jl")
     include("Generator/UnivariateGenerator/JoeGenerator.jl")
-    
+
     # Archimedean copulas
     include("ArchimedeanCopula.jl")
     include("BivariateArchimedeanMethods.jl")
@@ -94,9 +94,9 @@ module Copulas
     include("ExtremeValueCopulas/HuslerReissCopula.jl")
     include("ExtremeValueCopulas/LogCopula.jl")
     include("ExtremeValueCopulas/MixedCopula.jl")
-    include("ExtremeValueCopulas/MOCopula.jl")    
+    include("ExtremeValueCopulas/MOCopula.jl")
     include("ExtremeValueCopulas/tEVCopula.jl")
-    
+
     export AsymGalambosCopula,
            AsymLogCopula,
            AsymMixedCopula,
@@ -106,12 +106,12 @@ module Copulas
            HuslerReissCopula,
            LogCopula,
            MixedCopula,
-           MOCopula,   
+           MOCopula,
            tEVCopula
 
 
     # Subsetting
-    include("SubsetCopula.jl") # not exported yet. 
+    include("SubsetCopula.jl") # not exported yet.
 
     using PrecompileTools
     @setup_workload begin
@@ -143,8 +143,8 @@ module Copulas
                 # WCopula(2),            ################ <<<<<<<<<-------------- Does not work and I cannot explain why !
                 # RafteryCopula(2, 0.2), ################ <<<<<<<<<<------------- BUGGY
                 # RafteryCopula(3, 0.5), ################ <<<<<<<<<<------------- BUGGY
-                # We should probably add others to speed up again. 
-            ) 
+                # We should probably add others to speed up again.
+            )
                 u1 = rand(C)
                 u = rand(C,2)
                 if applicable(Distributions.pdf,C,u1) && !(typeof(C)<:EmpiricalCopula)