load("C:/Users/Msi/Downloads/Datos.RData")

Simular al proceso Poisson compuesto Entregar una función (o funciones) en R que permita simular al procesos Poisson compuesto con innovaciones: exponenciales, normales y otra distribución (a elección del alumno)

#' Simulación de un Proceso Poisson Compuesto
#' @param T_max Tiempo total de la simulación
#' @param lambda Tasa de llegada (intensidad del proceso)
#' @param dist_func Función de distribución para la magnitud Y (ej: rexp, rnorm, rgamma)
#' @param params Lista con los parámetros de la distribución (ej: list(rate=0.5))
sim_pp_compuesto <- function(T_max, lambda, dist_func, params) {
  
  # 1. Generar los tiempos de salto (Proceso Poisson)
  # Usamos un número alto para asegurar que cubrimos T_max
  tiempos_entre <- rexp(n = ceiling(lambda * T_max * 5), rate = lambda)
  tiempos_llegada <- cumsum(tiempos_entre)
  
  # Filtrar solo los que ocurren dentro del horizonte T_max
  W <- tiempos_llegada[tiempos_llegada <= T_max]
  n_eventos <- length(W)
  
  # 2. Generar las magnitudes (los "pesos" o "costos" de cada evento)
  # Usamos do.call para aplicar la función elegida con sus parámetros
  params$n <- n_eventos # Aseguramos que genere n_eventos valores
  Y <- do.call(dist_func, params)
  
  # 3. Calcular la suma acumulada de las magnitudes (S_t)
  saltos_acumulados <- cumsum(Y)
  
  # 4. Preparar datos para la gráfica (escalonada)
  # Añadimos el inicio en (0,0)
  x_plot <- c(0, W)
  y_plot <- c(0, saltos_acumulados)
  
  # 5. Graficar el Proceso Poisson Compuesto
  plot(x_plot, y_plot, type = "s", col = "blue", lwd = 2,
       main = "Proceso Poisson Compuesto",
       xlab = "Tiempo (t)", ylab = "Suma Acumulada (St)")
  
  # Devolvemos los datos por si necesitas hacer cálculos después
  return(data.frame(tiempo = x_plot, acumulado = y_plot))
}
# Escenario A: Magnitudes Exponenciales (alpha = 0.5)
sim_pp_compuesto(T_max = 20, lambda = 3, dist_func = rexp, params = list(rate = 0.5))

##        tiempo   acumulado
## 1   0.0000000   0.0000000
## 2   0.5021782   0.5064709
## 3   0.5936316   3.7690461
## 4   0.6494565   4.3223622
## 5   0.6783799   8.5331389
## 6   0.7976330  11.3736624
## 7   0.8221237  12.1736700
## 8   0.8760109  18.5172924
## 9   1.2583551  22.2581499
## 10  1.5159967  26.0436961
## 11  1.9569392  26.1791477
## 12  2.0902872  29.7145178
## 13  3.0767095  29.7826311
## 14  3.2250586  30.3201211
## 15  3.7337128  30.9932450
## 16  4.1547430  31.9107763
## 17  4.5931514  36.9698069
## 18  4.7072266  37.3943406
## 19  4.8538821  38.5416304
## 20  5.1160350  40.6674628
## 21  5.4163876  41.6639838
## 22  5.6659583  42.2698557
## 23  5.8587673  44.6643462
## 24  5.9085276  45.6689167
## 25  6.3068918  47.0001534
## 26  6.3611420  47.8415288
## 27  6.7489819  50.3844582
## 28  6.9504202  51.4822677
## 29  7.1541608  58.6029362
## 30  7.1898250  59.6625810
## 31  7.5596872  59.9216406
## 32  8.0669303  61.9280801
## 33  8.2879320  61.9766898
## 34  8.3315058  63.7718055
## 35  8.6176111  64.9868031
## 36  8.6241604  68.2764132
## 37  9.5847119  72.1557206
## 38  9.9023330  73.3494273
## 39 10.1954459  73.3892047
## 40 10.3556580  73.8260531
## 41 10.4379393  74.5279404
## 42 10.8479805  80.2214416
## 43 11.0796695  81.6350440
## 44 11.2648022  82.3385369
## 45 11.7981404  87.6864838
## 46 12.0965090  89.1532503
## 47 12.1707122  94.5519206
## 48 12.2188372  97.5799377
## 49 12.3707336  98.1227486
## 50 12.7261755 100.5351990
## 51 13.0336107 108.8477449
## 52 13.0856318 109.4542519
## 53 13.6307115 109.4654613
## 54 14.0405825 109.5815674
## 55 14.1867451 113.3335928
## 56 14.1922238 117.3266616
## 57 14.1928494 118.4374399
## 58 14.4974292 120.9013365
## 59 14.5887078 122.5274983
## 60 15.1459805 131.3681602
## 61 15.2455457 135.8106679
## 62 15.6310783 136.7611345
## 63 16.0083133 138.5042589
## 64 17.0915647 142.1573591
## 65 17.3420705 147.5286060
## 66 17.8470093 150.2604097
## 67 17.8750278 152.7485592
## 68 17.9793925 152.9960987
## 69 18.2115590 153.7557755
## 70 18.8836468 155.4427981
## 71 19.2721364 157.2614077
## 72 19.2944783 159.2548702
## 73 19.3562134 161.0836876
## 74 19.4125624 161.4924412
## 75 19.5702455 162.0959947
## 76 19.8098408 163.6309557
# Escenario B: Magnitudes Normales (mu = 10, sd = 2)
sim_pp_compuesto(T_max = 20, lambda = 3, dist_func = rnorm, params = list(mean = 10, sd = 2))

##         tiempo  acumulado
## 1   0.00000000   0.000000
## 2   0.04375888   5.882389
## 3   0.39040304  16.253001
## 4   0.51223194  26.126894
## 5   0.62770897  35.768883
## 6   1.15212683  44.413558
## 7   1.69918128  53.707662
## 8   2.19828632  62.954474
## 9   3.11156905  72.229737
## 10  3.21796967  79.298497
## 11  3.69914846  88.792901
## 12  4.23651744  98.173204
## 13  4.29113431 108.803698
## 14  4.30947310 118.328477
## 15  4.34789362 129.415805
## 16  4.78919077 140.524738
## 17  4.86384086 146.887959
## 18  4.93746167 155.334259
## 19  5.14184157 164.573646
## 20  5.43263994 172.894213
## 21  5.80708530 184.214275
## 22  5.97122547 194.672405
## 23  5.97655565 204.769203
## 24  6.27440417 216.872041
## 25  6.32808977 227.173545
## 26  6.85930057 234.530389
## 27  7.07165817 242.093290
## 28  7.14966628 255.335918
## 29  7.74547044 265.847282
## 30  7.83524350 274.104277
## 31  7.84132248 280.304867
## 32  7.99140805 291.263313
## 33  8.46212781 301.513396
## 34  9.00748541 308.943933
## 35  9.40934791 319.264248
## 36 10.29314383 330.316480
## 37 10.30453618 340.745916
## 38 10.90108244 350.343588
## 39 11.19790886 359.498231
## 40 11.60627730 368.516494
## 41 12.31926531 377.379148
## 42 12.60998854 387.808801
## 43 12.67368735 398.838110
## 44 12.76097005 408.059663
## 45 13.07006330 416.815416
## 46 13.10857679 426.636280
## 47 13.13523474 437.285522
## 48 13.18525664 443.668581
## 49 13.36585356 451.445431
## 50 13.43110092 464.009003
## 51 13.64344361 476.669947
## 52 14.08599797 486.291839
## 53 14.20416883 493.721068
## 54 14.22839172 503.744589
## 55 14.23542830 513.446089
## 56 14.33876322 526.077018
## 57 14.50425215 535.710770
## 58 14.89469887 544.884602
## 59 15.09674280 554.874355
## 60 15.25168585 567.068828
## 61 15.52469480 576.333214
## 62 15.66380801 585.001389
## 63 15.90465551 595.664357
## 64 16.24634798 605.040218
## 65 17.24713553 614.309123
## 66 17.33861525 621.891786
## 67 17.36814729 630.666748
## 68 17.38853298 637.605767
## 69 18.44639671 648.592083
## 70 18.70626731 656.415214
## 71 19.03680655 663.403992
## 72 19.10017901 675.346432
## 73 19.59433582 684.590346
## 74 19.70294826 695.274281
# Escenario C: Magnitudes Gamma (shape = 2, rate = 0.5)
sim_pp_compuesto(T_max = 20, lambda = 3, dist_func = rgamma, params = list(shape = 2, rate = 0.5))

##        tiempo  acumulado
## 1   0.0000000   0.000000
## 2   0.3943965   6.531906
## 3   0.9734487  12.087706
## 4   1.2341116  15.539120
## 5   1.2842248  17.432943
## 6   1.5437230  19.263391
## 7   1.7172444  25.729575
## 8   1.9400745  35.425073
## 9   1.9552358  37.039119
## 10  2.1323187  40.611799
## 11  2.6760518  44.469242
## 12  3.3453383  49.147412
## 13  3.6588789  51.092353
## 14  4.1122722  53.169536
## 15  4.8904432  59.435302
## 16  5.2273925  66.466132
## 17  5.5993189  70.808675
## 18  5.7089761  74.741688
## 19  6.1114539  78.767084
## 20  6.7749505  88.510462
## 21  6.7794973  89.636404
## 22  8.5729546  93.197666
## 23  9.6580601  95.016961
## 24  9.6712507  98.256122
## 25 10.5587991  99.244087
## 26 10.5951253 103.480331
## 27 10.9345034 104.007858
## 28 11.0865666 105.580520
## 29 12.2615468 111.476853
## 30 12.3475436 121.723685
## 31 13.3083522 124.455779
## 32 13.5715381 125.205177
## 33 13.6762911 130.710041
## 34 14.0290225 139.197583
## 35 14.6954809 141.083330
## 36 14.8256659 145.539457
## 37 15.1441801 146.427757
## 38 15.4468638 150.198278
## 39 15.9231178 152.292786
## 40 16.0224291 155.406406
## 41 16.0929669 156.916742
## 42 16.5828305 161.014511
## 43 16.6133239 163.627229
## 44 17.5306554 164.907880
## 45 18.3279503 166.401550
## 46 18.4404524 174.860213
## 47 19.6130675 182.139743

Simular el modelo de Cramer - Lundberg Entregar una función (o funciones) en R que permita simular el modelo de Cramer - Lundberg utilizando distribución exponencial y log-normal para las reclamaciones.

#' Simulador de Modelo de Cramér-Lundberg
#' @param u Capital inicial
#' @param c Tasa de cobro de primas
#' @param lambda Frecuencia de accidentes (ritmo Poisson)
#' @param T_max Horizonte temporal
#' @param tipo_dist Tipo de distribución: "exp" o "logn"
#' @param params Parámetros de la distribución (list)
simular_CL <- function(u, c, lambda, T_max, tipo_dist = "exp", params = list()) {
  
  # 1. Generar los tiempos de ocurrencia (Proceso Poisson)
  tiempos_entre <- rexp(n = ceiling(lambda * T_max * 2), rate = lambda)
  tiempos_llegada <- cumsum(tiempos_entre)
  tiempos_llegada <- tiempos_llegada[tiempos_llegada <= T_max]
  n_eventos <- length(tiempos_llegada)
  
  # 2. Generar las magnitudes de los reclamos (Y) según la distribución elegida
  if (tipo_dist == "exp") {
    Y <- rexp(n_eventos, rate = params$rate)
  } else if (tipo_dist == "logn") {
    Y <- rlnorm(n_eventos, meanlog = params$meanlog, sdlog = params$sdlog)
  }
  
  # 3. Construir la trayectoria del capital
  # Capital(t) = u + (c * t) - Suma_Reclamos_hasta_t
  trayectoria_capital <- u + (c * tiempos_llegada) - cumsum(Y)
  
  # Crear data frame para graficar
  df <- data.frame(tiempo = c(0, tiempos_llegada), 
                   capital = c(u, trayectoria_capital))
  
  # 4. Graficar
  plot(df$tiempo, df$capital, type = "s", col = "blue", lwd = 2,
       main = paste("Modelo de Cramér-Lundberg (Dist:", tipo_dist, ")"),
       xlab = "Tiempo", ylab = "Capital (u + ct - St)")
  abline(h = 0, col = "red", lty = 2) # Línea de ruina en 0
  
  return(df)
}
# Escenario 1: Distribución Exponencial
simular_CL(u = 10, c = 2, lambda = 1, T_max = 20, 
           tipo_dist = "exp", params = list(rate = 0.5))

##        tiempo   capital
## 1   0.0000000 10.000000
## 2   0.8261201  4.626591
## 3   2.6127783  6.503786
## 4   4.4043318  7.979823
## 5   5.4344095  9.199588
## 6   6.1449655  9.082012
## 7   6.3470129  6.414435
## 8   8.2065683  8.829179
## 9   9.7463932  6.984187
## 10 10.3174224  8.003865
## 11 10.7826974  8.264508
## 12 11.9949234 10.172970
## 13 13.5593075 12.108249
## 14 19.8177477 24.024744
# Escenario 2: Distribución Log-Normal (más volatilidad)
simular_CL(u = 10, c = 2, lambda = 1, T_max = 20, 
           tipo_dist = "logn", params = list(meanlog = 0, sdlog = 0.5))

##       tiempo  capital
## 1   0.000000 10.00000
## 2   2.150242 13.17957
## 3   3.076975 14.58545
## 4   4.891565 17.53394
## 5   5.544502 17.26532
## 6   7.039786 19.10215
## 7   7.858436 20.17375
## 8  10.022099 23.68567
## 9  11.085228 24.71519
## 10 12.364948 26.63528
## 11 12.874780 25.13881
## 12 13.030337 24.80347
## 13 13.167782 24.42516
## 14 13.306633 23.85359
## 15 14.210776 24.97031
## 16 14.717796 25.18575
## 17 15.702528 26.07970

Estimar probabilidad de ruina Entregar una función en R que permita realizar machas simulaciones de el modelo de Cramer - Lundberg y calcule el número de simulaciones que terminaron en ruina para estimar la probabilidad de ruina como se describió en clase.

#' Estimar Probabilidad de Ruina mediante simulación de Monte Carlo
#' @param n_sim Número de simulaciones a realizar
#' @param u Capital inicial (reserva)
#' @param c Tasa de prima por periodo
#' @param lambda Frecuencia de eventos
#' @param T_max Horizonte de tiempo
#' @param tipo_dist Tipo de distribución ("exp" o "lnorm")
#' @param params Lista de parámetros de la distribución
estimar_probabilidad_ruina <- function(n_sim, u, c, lambda, T_max, tipo_dist, params) {
  
  # Usamos replicate para repetir el escenario n_sim veces
  resultados_ruina <- replicate(n_sim, {
    
    # 1. Generar tiempos de llegada (Poisson)
    tiempos_entre <- rexp(n = ceiling(lambda * T_max * 2), rate = lambda)
    tiempos_llegada <- cumsum(tiempos_entre)
    tiempos_llegada <- tiempos_llegada[tiempos_llegada <= T_max]
    n_eventos <- length(tiempos_llegada)
    
    # 2. Generar montos de reclamaciones (Severidad)
    if (tipo_dist == "exp") {
      Y <- rexp(n_eventos, rate = params$rate)
    } else if (tipo_dist == "lnorm") {
      Y <- rlnorm(n_eventos, meanlog = params$meanlog, sdlog = params$sdlog)
    }
    
    # 3. Calcular la trayectoria del capital
    # Capital(t) = u + (c * t) - Suma(Y)
    trayectoria_capital <- u + (c * tiempos_llegada) - cumsum(Y)
    
    # 4. Verificar si hubo ruina (el capital bajó de 0 en algún momento)
    # min() nos da el punto más bajo alcanzado en toda la trayectoria
    ifelse(min(trayectoria_capital) < 0, 1, 0)
  })
  
  # Calculamos la proporción de veces que ocurrió la ruina
  prob_ruina <- sum(resultados_ruina) / n_sim
  return(prob_ruina)
}

Estimación de parámetros Usando los datos que les dejo en esta misma publicación deberán estimar los parámetros lambda (con las fechas de ocurrencia) y mu (con el valor del reclamo) para hacer una aproximación de la política de cobro por primas por periodo (constante c del modelo C-L).

tiempos_entre <- diff(Fechas)
lambda_est <- 1 / mean(tiempos_entre)
mu_est <- mean(Y)
list(lambda = lambda_est, mu = mu_est)
## $lambda
## [1] 1.888491
## 
## $mu
## [1] 26.00487
theta <- 0.2
prima_sugerida <- (1 + theta) * lambda_est * mu_est

• Análisis de sensibilidad - Usando los mismos datos, realicen un ajuste a una de las distribuciones ya trabajadas (exponencial y log-normal) a los montos de reclamaciones. - Con las funciones hechas anteriormente calculen las probabilidades de ruina usando los parámetros estimados y modificando los valores de u (capital inicial) y de c (cobro de primas por periodo) para obtener gráficas que muestren el comportamiento de estas probabilidades como función de esos parámetros.

library(fitdistrplus)
## Warning: package 'fitdistrplus' was built under R version 4.5.3
## Cargando paquete requerido: MASS
## Cargando paquete requerido: survival
# 1. Ajustar la distribución Log-Normal a tus datos
ajuste_lognorm <- fitdist(Y, "lnorm")

# 2. Ver el resumen (esto te mostrará los valores de estimación)
summary(ajuste_lognorm)
## Fitting of the distribution ' lnorm ' by maximum likelihood 
## Parameters : 
##          estimate Std. Error
## meanlog 2.8363102 0.03308130
## sdlog   0.8645555 0.02339187
## Loglikelihood:  -2806.931   AIC:  5617.862   BIC:  5626.915 
## Correlation matrix:
##         meanlog sdlog
## meanlog       1     0
## sdlog         0     1
# 3. Extraer los valores para tu script
# Esto guarda los valores directamente en la lista que pide tu script
params_real <- list(
  meanlog = as.numeric(ajuste_lognorm$estimate["meanlog"]),
  sdlog   = as.numeric(ajuste_lognorm$estimate["sdlog"])
)

# Verifica que se hayan guardado correctamente
print(params_real)
## $meanlog
## [1] 2.83631
## 
## $sdlog
## [1] 0.8645555
# 1. DEFINICIÓN DE PARÁMETROS REALES (MODIFICA ESTOS VALORES)
# ------------------------------------------------------------
# Sustituye lambda_real por el resultado de tu función de estimación (1 / promedio de tiempos)
lambda_real <- 1.88          

# Sustituye meanlog y sdlog por los resultados obtenidos con fitdist() o cálculo manual
# (Deben ser los parámetros de la distribución Log-Normal que mejor se ajustó a tus datos)
params_real  <- list(meanlog = 2.83631, sdlog = 0.8645555) 

# Horizonte de tiempo (ejemplo: si tus datos son mensuales, mantén 30 días)
T_max_real   <- 30            

# 2. DEFINICIÓN DE RANGOS DE SENSIBILIDAD (MODIFICA ESTOS VALORES)
# ------------------------------------------------------------
# Ancla tus rangos al costo puro: (lambda_real * media_de_tus_reclamos)
# Ejemplo: si el costo puro es 10, probamos rangos de 10 a 25 para primas
costo_puro <- lambda_real * 15 # <-- Sustituye 15 por la media de tus reclamos reales

vector_u <- seq(20, 100, length.out = 30) # Rango de capital inicial
vector_c <- seq(costo_puro, costo_puro * 2, length.out = 30) # Rango de primas (c)

# 3. CÁLCULO DE LA MATRIZ DE RIESGO
# ------------------------------------------------------------
matriz_ruina <- outer(vector_u, vector_c, Vectorize(function(u_val, c_val) {
  estimar_probabilidad_ruina(
    n_sim     = 1000, 
    u         = u_val, 
    c         = c_val, 
    lambda    = lambda_real, 
    T_max     = T_max_real, 
    tipo_dist = "lnorm", 
    params    = params_real
  )
}))

# 4. VISUALIZACIÓN
# ------------------------------------------------------------
filled.contour(vector_u, vector_c, matriz_ruina,
               color.palette = function(n) hcl.colors(n, "YlOrRd", rev = TRUE),
               main = "Mapa de Calor: Probabilidad de Ruina",
               xlab = "Capital Inicial (u)",
               ylab = "Prima por periodo (c)")

• Realicen las conclusiones que consideren pertinentes.

Para parametrizar el modelo, se utilizó la tasa de siniestralidad λ calculada a partir de los intervalos de llegada de los siniestros históricos y se ajustó una distribución Log-Normal a los montos de reclamación mediante máxima verosimilitud, resultando en los parámetros meanlog y sdlog empleados en la simulación.

Hemos contrastado el modelo bajo una distribución Normal, que asume una variabilidad simétrica, frente a la distribución Log-Normal, que captura mejor la asimetría positiva de los siniestros. Esta comparativa nos permite cuantificar cuánto subestima el riesgo un modelo que no considera adecuadamente las colas pesadas en la severidad de las reclamaciones.

Los rangos fueron definidos en función del costo puro esperado (Cp) y el monto promedio de los siniestros observados, asegurando que el análisis de sensibilidad cubra escenarios desde el riesgo extremo hasta la solvencia operativa”.

Respecto al mapa de calor, este nos permite identificar la frontera de solvencia. Las zonas con colores cálidos representan configuraciones de capital y primas donde el riesgo es inaceptable, mientras que las zonas con colores fríos (o claros) definen las estrategias que garantizan la sostenibilidad de la aseguradora.