Copulas20262/16_simular.jl at main · aerdely/Copulas20262 · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
### Simular cópula arquimediana
### Autor: Arturo Erdely
### Fecha: 2026-03-18

#=
    Cópula arquimediana con parámetro θ ≥ 1, definida mediante
    el generador φ(t) = (1-t)^θ, t ∈ [0,1] (familia 4.2.2 de Nelsen).
=#

## Paquetes instalados previamente

using Plots, LaTeXStrings, Distributions


# curva frontera de región cero
g(u,θ) = 1 - (1 - (1-u)^θ)^(1/θ)

# ψ(v|u) = ∂C(u,v)/∂u como función de v, para u y θ fijos
function ∂C∂u(u, v, θ)
    if g(u,θ) < v ≤ 1.0
        return ((1-u)^(θ-1)) * ((1-u)^θ + (1-v)^θ)^(1/θ-1)
    elseif 0.0 ≤ v < g(u,θ)
        return 0.0
    else
        return NaN
    end
end

# ψ⁻¹(t|u) como función de t, para u y θ fijos
function ψi(t, u, θ)
    if t > (1-u)^(θ-1)
        return 1 - (1-u)*(1/(t^(θ/(θ-1))) - 1)^(1/θ)
    else
        return g(u, θ)
    end
end

# Gráfica de ψ(v|u) para u y θ fijos
begin
    θ = 2.3   # θ ≥ 1
    u = 0.4   # 0 < u < 1
    vv = collect(range(0, 1, length = 1_000))
    ψ = [∂C∂u(u, v, θ) for v ∈ vv]
    scatter(vv, ψ, xlabel = L"v", ylabel = L"\partial C_{\theta}(u,v) / \partial u",
         title = "u = $u, θ = $θ", size = (500,400), label = L"v\mapsto \partial C(u,v) / \partial u",
         ms = 2.0, mcolor = :black)
    gu = g(u,θ)
    vg = (1-u)^(θ-1)
    vline!([gu], ls = :dash, color = :red, label = "g(u) = $(round(gu, digits=4))")
    hline!([vg], ls = :dash, color = :gray, label = "(1-u)^(θ-1) = $(round(vg, digits=4))")
end

# Agregar gráfica de ψ⁻¹(t|u) para u y θ fijos (para comprobar)
begin
    tt = collect(range(0, 1, length = 1_000))
    plot!([ψi(t, u, θ) for t ∈ tt], tt, color = :cyan,
               label = L"\psi^{(-1)}_{θ}(t|u)")
end


# Función para simular n observaciones de la cópula arquimediana con parámetro θ
function simular_copula(n, θ)
    U = rand(n)
    V = similar(U)
    for i ∈ 1:n
        u = U[i]
        w = rand()
        V[i] = ψi(w, u, θ)
    end
    return U, V
end

begin
    θ = 2.3
    U, V = simular_copula(1000, θ)
    scatter(U, V, xlabel = L"U", ylabel = L"V", size = (500,470),
            ms = 2.0, mcolor = :blue, label = "")
end
# agregar curva frontera de región cero
plot!(sort(U), g.(sort(U), θ), color = :red, lw = 2.0, label = L"g_{θ}(u)")


## Simular (X,Y) con cópula arquimediana y marginales dadas

function simular_XY(n, θ, distX, distY)
    U, V = simular_copula(n, θ)
    X = quantile(distX, U)
    Y = quantile(distY, V)
    return X, Y
end

begin
    θ = 2.3
    distX = Beta(0.8,0.5)
    distY = Gamma(2,3)
    X, Y = simular_XY(1000, θ, distX, distY)
    scatter(X, Y, xlabel = L"X", ylabel = L"Y", size = (500,470),
            ms = 2.0, mcolor = :green, label = "")
end

# comprobando marginales

histogram(X, bins = 30, xlabel = L"X", ylabel = "densidad",
          title = "Histograma de X ~ Beta(0.8,0.5)", size = (500,400), label = "empírica",
          mcolor = :orange, linecolor = :black, normalize = true)
plot!(x -> pdf(distX, x), 0:0.01:1, color = :red, lw = 2.0, label = "teórica")

histogram(Y, bins = 30, xlabel = L"Y", ylabel = "denidad",
          title = "Histograma de Y ~ Gamma(2,3)", size = (500,400), label = "empírica",
          mcolor = :purple, linecolor = :black, normalize = true)
plot!(x -> pdf(distY, x), 0:0.1:20, color = :red, lw = 2.0, label = "teórica")