Skip to content

Commit 02b8a7b

Browse files
authored
[Internals] Simplify \rho for Clayton and InvGaussian (#331)
1 parent 939969e commit 02b8a7b

File tree

2 files changed

+4
-30
lines changed

2 files changed

+4
-30
lines changed

src/Generator/ClaytonGenerator.jl

Lines changed: 2 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -81,25 +81,8 @@ function Distributions._logpdf(C::ClaytonCopula{d,TG}, u) where {d,TG<:ClaytonGe
8181
S1==d-1 && return eltype(u)(-Inf)
8282
return log+ 1) * (d - 1) -+ 1) * S2 + (-1 / θ - d) * log(S1 - d + 1)
8383
end
84-
### only for test...
85-
@inline function _C_clayton(u::Float64, v::Float64, θ::Real)
86-
s = u^(-θ) + v^(-θ) - 1
87-
if θ < 0
88-
return (s <= 0) ? 0.0 : s^(-1/θ) # soporte recortado para θ<0
89-
else
90-
return s^(-1/θ) # para θ>0 siempre s≥1
91-
end
92-
end
93-
# Spearman (vía CDF) — with _safett Integral
94-
function ρ(G::ClaytonGenerator; rtol=1e-8, atol=1e-10)
95-
θ = float(G.θ)
96-
θ -1 && throw(ArgumentError("Para Clayton: θ > -1."))
97-
iszero(θ) && return 0.0
98-
I = HCubature.hcubature(x -> _C_clayton(x[1], x[2], θ),
99-
[0.0,0.0], [1.0,1.0];
100-
rtol=rtol, atol=atol)[1]
101-
return 12I - 3
102-
end
84+
85+
ρ(G::ClaytonGenerator) = @invoke ρ(ArchimedeanCopula(2, G)::Copula)
10386

10487
# Inverse ρ → θ for Clayton (without trimming to [0,1])
10588
function ρ⁻¹(::Type{<:ClaytonGenerator}, ρ̂; atol=1e-10)

src/Generator/InvGaussianGenerator.jl

Lines changed: 2 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -141,21 +141,12 @@ function τ⁻¹(::Type{<:InvGaussianGenerator}, τ)
141141
end
142142
frailty(G::InvGaussianGenerator) = Distributions.InverseGaussian(G.θ,1)
143143

144-
145-
function _rho_invgaussian(θ; rtol=1e-7, atol=1e-9, maxevals=10^6)
146-
θeff = clamp(θ, 1e-12, Inf)
147-
= Copulas.ArchimedeanCopula(2, InvGaussianGenerator(θeff))
148-
f(x) = _cdf(Cθ, (x[1], x[2])) # <- tu _cdf
149-
I = HCubature.hcubature(f, (0.0,0.0), (1.0,1.0); rtol=rtol, atol=atol, maxevals=maxevals)[1]
150-
return 12I - 3
151-
end
152-
153-
ρ(G::InvGaussianGenerator) = _rho_invgaussian(G.θ)
144+
ρ(G::InvGaussianGenerator) = @invoke ρ(ArchimedeanCopula(2, G)::Copula)
154145
function ρ⁻¹(::Type{<:InvGaussianGenerator}, rho)
155146
# Numerically inverts _rho_invgaussian using Brent's method.
156147
# Spearman's rho for InvGaussian: [0, 1/2)
157148
rho 0 && return zero(rho)
158149
rho log(2) && return Inf * rho
159-
xhat = Roots.find_zero(x -> _rho_invgaussian(-log(x)) - rho, (0, 1))
150+
xhat = Roots.find_zero(x -> ρ(InvGaussianGenerator(-log(x))) - rho, (0, 1))
160151
return -log(xhat)
161152
end

0 commit comments

Comments
 (0)