La respuesta corta es que su δ está bien, pero su γ es malo. Con el fin de obtener el positivo estable de la distribución dada por la fórmula en R, es necesario establecer
γ=|1−itan(πα/2)|−1/α.
El ejemplo más temprano que pude encontrar de la fórmula que les dio fue que en (Feller, 1971), pero sólo he encontrado ese libro en forma física. Sin embargo (Hougaard, 1986) da la misma fórmula, junto con la transformada de Laplace
L(s)=E[exp(−sX)]=exp(−sα).
De la stabledist
manual (stabledist
se utiliza en fBasics
), el pm=1
parametrización es de (Samorodnitsky y Taqqu, 1994), otro recurso cuya reproducción ha eludido a mí. Sin embargo (Weron, 2001) da la función característica en Samorodnitsky y Taqqu la parametrización para α≠1
φ(t)=E[exp(itX)]=exp[δt−γα|t|α(1−iβsigno(t)\bronceadoπα2)].
He cambiado el nombre de algunos parámetros de Weron documento que coinside con la notación que estamos usando. Él usa μδσγ. En cualquier caso, conectar β=1δ=0, obtenemos
φ(t)=exp[−γα|t|α(1−isigno(t)tanπα2)].
Note that (1−itan(πα/2))/|1−itan(πα/2)|=exp(−iπα/2) for α\en(0,1) and that iα=exp(iπα/2). Formally, L(s)=φ(es), so by setting γ=|1−itan(πα/2)|−1/α in φ(t) obtenemos
φ(es)=exp(−sα)=L(s).
Un punto interesante a tener en cuenta es que el γ que corresponde a α=1/21/2, por lo que si usted fuera a tratar γ=α o γ=1−α, que en realidad no es una mala aproximación, se termina exactamente correcta para α=1/2.
He aquí un ejemplo en R para comprobar la corrección:
library(stabledist)
# Series representation of the density
PSf <- function(x, alpha, K) {
k <- 1:K
return(
-1 / (pi * x) * sum(
gamma(k * alpha + 1) / factorial(k) *
(-x ^ (-alpha)) ^ k * sin(alpha * k * pi)
)
)
}
# Derived expression for gamma
g <- function(a) {
iu <- complex(real=0, imaginary=1)
return(abs(1 - iu * tan(pi * a / 2)) ^ (-1 / a))
}
x=(1:100)/100
plot(0, xlim=c(0, 1), ylim=c(0, 2), pch='',
xlab='x', ylab='f(x)', main="Density Comparison")
legend('topright', legend=c('Series', 'gamma=g(alpha)'),
lty=c(1, 2), col=c('gray', 'black'),
lwd=c(5, 2))
text(x=c(0.1, 0.25, 0.7), y=c(1.4, 1.1, 0.7),
labels=c(expression(paste(alpha, " = 0.4")),
expression(paste(alpha, " = 0.5")),
expression(paste(alpha, " = 0.6"))))
for(a in seq(0.4, 0.6, by=0.1)) {
y <- vapply(x, PSf, FUN.VALUE=1, alpha=a, K=100)
lines(x, y, col="gray", lwd=5, lty=1)
lines(x, dstable(x, alpha=a, beta=1, gamma=g(a), delta=0, pm=1),
col="black", lwd=2, lty=2)
}
![Plot output]()
- Feller, W. (1971). Una Introducción a la Teoría de la Probabilidad y Sus Aplicaciones, 2, 2ª ed. Nueva York: Wiley.
- Hougaard, P. (1986). La supervivencia de los Modelos para Poblaciones Heterogéneas Derivadas de Distribuciones Estables, Biometrika 73, 387-396.
- Samorodnitsky, G., Taqqu, M. S. (1994). Estable No Aleatoria Gaussiana Procesos, Chapman & Hall, Nueva York, 1994.
- Weron, R. (2001). Levy-estable distribuciones revisited: cola índice > 2 no excluye la de Levy-régimen estable, Revista Internacional de la Física Moderna C, 2001, 12(2), 209-223.