2.2. Modelos de distribucion de perdidas: cuantia de los siniestros.

Distribución log-normal

Es una distribución conocida como distribución de colas anchas, lo que quiere decir, que su densidad refleja una caracteristica que tambien se observa en su varianza: hay una probabilidad alta de que la variable tome valores alejados de la media. Nos servirá para modelizar siniestros con grandes costes o cuantias por reclamación.

Se suele definir la variable aleatoria \(X\) que sigue una distribucion log-normal como aquella variable cuyo logaritmo neperiano \(Y=\ln X\) sigue una distribución normal \(N(\mu, \sigma^2)\). Se obtiene su función de densidad a traves de la transformación de variables:

\[f(x)=\frac{1}{x}\left(\frac{1}{\sqrt{2\pi}\sigma} e^{\frac{(-\log x-\mu)^2}{2\sigma^2}}\right) \] La esperanza y la varianza son:

\(E(X)=e^{\mu+\frac{\sigma^2}{2}}, var(X)=e^{2\mu+2\sigma^2}-e^{2\mu+\sigma^2}\)

Estimadores de la distribución lognormal

\[ \hat{\mu}=\sum_{i=1}^{n}log(x_i)/n, \hat{\sigma}^2=\sum_{i=1}^{n}(\log x_i-\hat{\mu})^2/n \]

Ejercicio: Se estudia la proporción de rentistas por encima de 18.000 anuales para un sector económico cuya distribución salarial medida en miles de dolares sigue un modelo logarítmico normal con parámetros \(\mu=2\) y \(\sigma=2\).

p=1-plnorm(18,2,1.2)
p

## [1] 0.2290508

Un siniestro tiene para su cuantía monetaria una función de distribución log-normal con \(\mu=7\) y \(\sigma=1.5\). ¿cuál es la probabilidad de tener un siniestro de cuantía inferior a 200? ¿cuál es la probabilidad de tener un siniestro de cuantía superior a 1000?

p=plnorm(200,7,1.5)

p=1-plnorm(1000,7,1.5)

Distribución de Pareto

Al igual que la distribución anterior, es una ley de colas anchas por lo que la utilizaremos para modelizar cuantias elevadas. La distribución de Pareto surge al considerar que la probabilidad de que una v.a. tome un valor superior a un \(x\) determinado, tiene la siguiente forma funcional:

\[ P(X>x)=\left(\frac{\lambda}{x}\right)^\alpha \]

si \(x\geq \lambda\) y \(P(X>x)=1\), en otro caso. A partir de aqui se puede obtener facilmente la función de de densidad:

\[ f(x)=\frac{\alpha\lambda^\alpha}{x^{\alpha+1}}, x>0 \]

La función de cola viene dada por:

\[ \overline{F}=\left\{\begin{array}{ll} 1, & \hbox{x<0;}\\ \left(\frac{\lambda}{\lambda+x}\right)^\alpha, & \hbox{x>0.} \end{array}\right. \]

Propiedades

Su esperanza y varianza son:

\(E(X)=\frac{\lambda}{\alpha-1}, var(X)=\frac{\alpha\lambda^2}{(\alpha-1)^2(\alpha-2)}\)

Normalmente \(\lambda\) es un parámetro prefijado de antemano, a partir de los datos iniciales de la aplicación. El parámetro \(\alpha\) es un parámetro de ajuste cuyo valor se estima o determina a posteriori.

Si \(X\sim P( \alpha,\lambda)\) , para \(k>0\), entonces \(kX\sim P(\alpha,k\lambda)\). Se puede interpretar en términos de inflación o cambio de moneda.
Si \(X\sim P(\alpha,\lambda)\) y \(X>M\), siendo (\(M>0\)) una constante, entonces \(X-M \sim P(\alpha,\lambda+M)\) . Se puede interpretar en terminos de franquicias o reaseguros.
Si \(X \sim P(\alpha, \lambda)\) entonces \(X\equiv_{i.d.}\lambda(U^{-1/\alpha}-1)\) , con \(U\sim Uniforme(0,1)\). Puede ser interesante para simulaciones.

Ejercicio: Genera una muestra de tamaño 300 de una distribución de Pareto con \(\alpha = 3\) y \(\lambda = 800\). Visualiza la muestra y calcula su media y varianza.

sample<-800*((runif(n=300))^(-1/3)-1) 
fix(sample) 
mean(sample)

## [1] 425.5489

var(sample)

## [1] 446364.9

EMM de la distribución Pareto

Los estimadores del método de los momentos (EMM) se obtienen resolviendo las ecuaciones:

\[ \overline{x}=\frac{\lambda}{\alpha-1}; s^2=\frac{\alpha\lambda^2}{(\alpha-1)^2(\alpha-2)} \]

Los estimadores son:

\(\tilde{\alpha}=\frac{2s^2}{s^2-\overline{x}^2}\)

\(\tilde{\lambda}=(\tilde{\alpha}-1)\overline{x}\)

donde \(\overline{x}\) y \(s^2\) son, respectivamente, la media y la varianza de la muestra.

EMV de la distribución Pareto

Los estimadores de máxima verosimilitud (EMV) tienen mejores propiedades asintóticas. Para una muestra \(x_1, x_2, …, x_n\) de una distribución de Pareto \(P(\alpha,\lambda)\), la función de verosimilitud es:

\[ L(\alpha,\lambda)=\prod_{i=1}^n \frac{\alpha\lambda^\alpha}{(\lambda+x_i)^{\alpha+1}} \]

Los EMV verifican:

\[ \frac{\partial}{\partial\alpha}L(\alpha,\lambda)\Rightarrow\hat{\alpha}=\frac{n}{\sum_{i=1}^n\log(1+x_i/\hat{\lambda)}} \]

\[ \frac{\partial}{\partial\lambda}L(\alpha,\lambda)\Rightarrow\hat{\alpha}=\frac{\sum_{i=1}^n1/(\hat{\lambda}+x_i)}{\sum_{i=1}^nx_i/(\hat{\lambda}(\hat{\lambda}+x_i))} \]

Distribución Gamma

Una variable aleatoria continua \(X\) sigue una distribución gamma de parámetros \(\alpha>0\) (parámetro de forma) y \(\lambda>0\) (parámetro de escala) y se denota por \(X \sim G(\alpha,\lambda)\), si su función de densidad viene dada por:

\[ f(x)=\frac{\lambda^{\alpha}}{\Gamma(\alpha)}x^{\alpha-1}e^{-\lambda x}, x>0 \]

donde \(\Gamma(\alpha)=\int_{0}^{\infty}e^{-x}x^{\alpha-1}dx\)

Su esperanza y varianza son:

\[ E(X)=\frac{\alpha}{\lambda}, var(X)= \frac{\alpha}{\lambda^2} \]

Cuando \(\alpha=1\), se tiene la distribución exponencial
Cuando \(\alpha=\frac{r}{2}\) y \(\lambda=\frac{1}{2}\), \(X\) se distribuye como una \(\chi^2\) con \(r\) grados de libertad.

Representación gráfica

curve(dgamma(x, shape=2, scale=0.5), xlab = "x", ylab = "densidad gamma", col=4,
xlim=c(0,15), ylim=c(0,0.8), main = "Funciones de densidad Gamma") 
curve(dgamma(x, shape=2, scale=1), xlab = "x", ylab = "densidad gamma", col=3,
xlim=c(0,15), ylim=c(0,0.8), main = "Gráfico distribución Gamma", add=T)
curve(dgamma(x, shape=2, scale=2), xlab = "x", ylab = "densidad gamma", col=2,
xlim=c(0,15), ylim=c(0,0.8), add=T)
legend (x="topright", legend=c("G(2,1)", "G(2,0.5)","G(2,2)"), fill=c("3", "2", "4"))

Estimadores de la distribución Gamma

Los estimadores por el método de los momentos son:

\(\tilde{\alpha}=\frac{\overline{x}^2 }{S^2}, \tilde{\lambda}=\frac{\overline{x} }{S^2}\)

donde \(\overline{x}\) y \(S^2\) son, respectivamente, la media y la varianza de la muestra.
No hay formas cerradas para los EMV. Se obtienen, tras una reparametrización, por métodos numéricos.

Distribución exponencial

Es un caso particular de la distribución Gamma cuando \(\alpha=1\). Su función de densidad viene dada por:

\[ f(x)=\lambda e^{-\lambda x}, x \geq 0 \]

y \(f(x)=0\) , en otro caso.

Su esperanza y varianza son:

\[ E(X)=\frac{1}{\lambda}, var(X)=\frac{1}{\lambda^2} \]

Distribución de Weibull

La distribución de Weibull es de uso habitual en análisis de supervivencia, teoria de valores extremos, seguros, etc. Es otra generalización de la distribución exponencial, en la que intervienen dos parámetros, el parámetro \(c>0\) y \(\gamma>0\), que se denota como \(X\sim W(c,\gamma)\) , y su función de densidad viene dada por:

\[ f(x)=c\gamma x^{\gamma-1}e^{-cx^{\gamma}}, x>0 \]

Su función de cola es \(\overline{F}(x)=e^{-cx{\gamma}}, x>0\) y sus momentos centrales son:

\[ EX^k=\frac{1}{c^{\frac{k}{\alpha}}}\Gamma\left(1+\frac{k}{\gamma}\right) \]

Cuando \(\gamma=1\) se tiene una distribución exponencial de parámetro \(c\)
La distribución Weibull es de cola pesada si y sólo si \(\gamma<1\).

Estimación de los parámetros

Paquetes de R para la estimación

Usando el paquete fitdistrplus con la función fitdist. En este caso las distribuciones disponibles son: norm, lnorm, exp, pois, cauchy, gamma, logis, nbinom, geom, beta y weibull del paquete stats; invgamma, llogis, invweibull, pareto1 y pareto del paquete actuar. Se puede realizar mediante el método de los momentos (mme) o mediante el método de máxima verosimilitud (mle).

Contrastes de ajuste

Test de Kolmogorov-Smirnov

H\(_0\): Los datos proceden de una distribución F\(_0\)

Se trata de un test no parámetrico que se usa en distribuciones continuas
Compara la función de distribución teórica F\(_0\) con la empírica F\(_n\)
Se rechaza H\(_0\) si \(d_n=\sup |F_n(x)-F_0(x)|>k, -\infty<x<\infty\) donde dado un conjunto de datos \(x_1, x_2,...,x_n\): \[F_n(x)=\frac{1}{n}\sum_{i=1}^{n}I_{(-\infty,x)}(x_i)\]
El test K-S es invariante ante transformaciones.
Ventaja: La distribución del estadístico bajo la hipótesis nula no depende de F\(_0\).
Inconveniente: no detecta bien discrepacias en la cola. El test de Anderson-Darling es una alternativa que da más peso a la cola. Tiene como inconveniente que los valores criticos deben calcularse para cada distribución.

Selección de distribuciones

Una vez que se realiza los test de bondad de ajuste, se pueden usar técnicas alternativas para seleccionar la distribución con mejor ajuste.

Criterio de información de Akaike (AIC)

\[AIC =-2log(likelihood)+2k\]

Criterio de información bayesiano (BIC)

\[BIC=-2log(likelihood)+log(n)×k\]

Ambas tienen en cuenta el log-likelihood y penalizan el número de parámetros de la distribución (grados de libertad) de manera diferente.

Para todas ellas, cuanto menor sea el valor, mejor el ajuste. Es importante tener en cuenta que, ninguna de estas métricas, sirven para cuantificar la calidad del modelo en un sentido absoluto, sino para comparar la calidad relativa entre modelos/ajustes. Si todos los ajustes candidatos son malos, no proporcionan ningún aviso de ello.

Aplicación práctica

Se consideran 120 reclamaciones (en euros) por robo en cierto vecindario. ¿Qué distribución de probabilidad ajusta adecuadamente los datos? El interés se centra en modelar la cola derecha de la distribución: la forma y el tamaño de las reclamaciones grandes. Se pide:

Realizar un análisis exploratorio de los datos.
Estimar los parámetros para la distribución exponencial, Weibull, Pareto y Lognormal
Analizar en detalle el ajuste con la distribución exponencial.
Contrastar hipótesis referidas al ajuste de modelos propuestos.
Analisis exploratorio de los datos

datosrobos <- read_excel("datosrobos.xlsx")
Robos<-datosrobos$R
summary(Robos)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##     3.0   271.0   868.5  2020.3  1733.0 32043.0

skewness(Robos)

## [1] 5.033963

kurtosis(Robos)

## [1] 31.10723

mean=summary(Robos)[4]
n=length(Robos)

hist(Robos, main="Histograma de Reclamaciones por robo", xlab= "Cuantía
de reclamación", ylab = "frecuencia", col="red")

Estimación de los parámetros

library(fitdistrplus)
library(actuar)
fitexp<- fitdist(Robos, "exp", method = "mme")
fitexp

## Fitting of the distribution ' exp ' by matching moments 
## Parameters:
##         estimate
## rate 0.000494978

fitp<- fitdist(Robos, "pareto")
fitp

## Fitting of the distribution ' pareto ' by maximum likelihood 
## Parameters:
##          estimate  Std. Error
## shape    1.880085   0.4897099
## scale 1871.122666 699.1626464

fitw <- fitdist(Robos, "weibull")
fitw

## Fitting of the distribution ' weibull ' by maximum likelihood 
## Parameters:
##           estimate   Std. Error
## shape    0.7158071   0.04733433
## scale 1557.0525771 210.29278380

fitln <- fitdist(Robos, "lnorm")
fitln

## Fitting of the distribution ' lnorm ' by maximum likelihood 
## Parameters:
##         estimate Std. Error
## meanlog 6.624172 0.13795726
## sdlog   1.511246 0.09755032

Gráficos

E<-fitdist(Robos,"exp",method = "mme")
plot(E,histo=FALSE,demp=TRUE)

P<-fitdist(Robos,"pareto")
plot(P,histo=FALSE,demp=TRUE)

W<-fitdist(Robos,"weibull")
plot(W,histo=FALSE,demp=TRUE)

L<-fitdist(Robos,"lnorm")
plot(L,histo=FALSE,demp=TRUE)

3.Análisis ajuste distribución Exponencial

hist(Robos, main="Ajuste exponencial", xlab= "Cuantía de reclamación", ylab
= "frecuencia relativa", col="red",breaks="Sturges", freq=FALSE)
curve(dexp(x,rate=1/mean), col="black", from=0, to=35000, add=T)
legend (x="topright", legend=c("Frec.rel.","Exp(1/2020.29)"), fill=c("red","black"))

¿Es razonable asumir que los datos proceden de una distribució exponencial?

p8000<- sum(Robos > 8000)/n
p8000

## [1] 0.05

pexp(8000, rate=1/mean, lower.tail = F)

## [1] 0.01906646

p10000<- sum(Robos > 10000)/n
p10000

## [1] 0.025

pexp(10000, rate=1/mean, lower.tail = F)

## [1] 0.007084965

p20000<- sum(Robos > 20000)/n
p20000

## [1] 0.01666667

pexp(20000, rate=1/mean, lower.tail = F)

## [1] 5.019673e-05

Contrastes de ajuste

Representación de la función de distribución empírica.

ECDF=ecdf(Robos)
plot(ECDF, main="Cumulative distribution", xlab="cuantía de reclamación")

Test de Kolmogorov Smirnov

ks.test(Robos, "pexp", 1/mean(Robos))

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  Robos
## D = 0.20133, p-value = 0.0001192
## alternative hypothesis: two-sided

ks.test(Robos, "plnorm", 1/mean(Robos))

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  Robos
## D = 0.98341, p-value < 2.2e-16
## alternative hypothesis: two-sided

ks.test(Robos, "pweibull", 1/mean(Robos))

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  Robos
## D = 0.63232, p-value < 2.2e-16
## alternative hypothesis: two-sided

#Adaptación para realizar el contraste para la distribución Pareto
ks.test(1-(1872.13/(1872.13+Robos))**(1.88) ,"punif")

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  1 - (1872.13/(1872.13 + Robos))^(1.88)
## D = 0.056026, p-value = 0.8456
## alternative hypothesis: two-sided

Selección de distribución

library(fitdistrplus)
library(actuar)
fitexp<- fitdist(Robos, "exp", method = "mme")
fitp<- fitdist(Robos, "pareto")
fitw <- fitdist(Robos, "weibull")
fitln <- fitdist(Robos, "lnorm")
gofstat(list(fitexp,fitp, fitw, fitln))

## Goodness-of-fit statistics
##                              1-mme-exp 2-mle-pareto 3-mle-weibull 4-mle-lnorm
## Kolmogorov-Smirnov statistic 0.2013282   0.05619364     0.1006106   0.0866735
## Cramer-von Mises statistic   1.6440506   0.06790826     0.2126488   0.1358042
## Anderson-Darling statistic   8.6398126   0.36596900     1.2740786   0.6967935
## 
## Goodness-of-fit criteria
##                                1-mme-exp 2-mle-pareto 3-mle-weibull 4-mle-lnorm
## Akaike's Information Criterion  2068.639     2028.423      2038.858    2033.451
## Bayesian Information Criterion  2071.427     2033.998      2044.433    2039.026

Matemáticas aplicada a la Economía y la Empresa. Modelos para la variabilidad y administración del riesgo. Sesión 2

Catalina Garcia

20 de marzo de 2025

Previo: Paquetes necesarios

2.2. Modelos de distribucion de perdidas: cuantia de los siniestros.

Distribución log-normal

Estimadores de la distribución lognormal

Distribución de Pareto

Propiedades

EMM de la distribución Pareto

EMV de la distribución Pareto

Distribución Gamma

Representación gráfica

Estimadores de la distribución Gamma

Distribución exponencial

Distribución de Weibull

Estimación de los parámetros

Paquetes de R para la estimación

Contrastes de ajuste

Test de Kolmogorov-Smirnov

Selección de distribuciones