Comment puis-je adapter un ensemble de données à une distribution de Pareto dans R?

Eh bien, si vous avez un échantillon d'une distribution de Pareto avec les paramètres et (où est le paramètre de limite inférieure et est le paramètre de forme), la probabilité logarithmique de cela l'échantillon est: $X_1, ..., X_n$ $m>0$ $\alpha>0$ $m$ $\alpha$

n \log (α) + n α \log (m) - (α + 1) \sum_{i = 1}^{n} \log (X_{i})

$n \log(\alpha) + n \alpha \log(m) - (\alpha+1) \sum_{i=1}^{n} \log(X_i)$

il s'agit d'une augmentation monotone de , de sorte que le maximiseur est la plus grande valeur cohérente avec les données observées. Puisque le paramètre définit la borne inférieure du support pour la distribution de Pareto, l'optimum est $m$ $m$

\hat{m} = min_{i} X_{i}

$\hat{m} = \min_{i} X_i$

qui ne dépend pas de . Ensuite, en utilisant des astuces de calcul ordinaires, le MLE pour doit satisfaire $\alpha$ $\alpha$

\frac{n}{α} + n \log (\hat{m}) - \sum_{i = 1}^{n} \log (X_{i}) = 0

$\frac{n}{\alpha} + n \log( \hat{m} ) - \sum_{i=1}^{n} \log(X_i) = 0$

une algèbre simple nous dit que le MLE de est $\alpha$

\hat{α} = \frac{n}{\sum_{i = 1}^{n} \log (X_{i} / \hat{m})}

$\hat{\alpha} = \frac{n}{\sum_{i=1}^{n} \log(X_i/\hat{m})}$

Dans de nombreux sens importants (par exemple, efficacité asymptotique optimale en ce qu'elle atteint la limite inférieure de Cramer-Rao), c'est la meilleure façon d'ajuster les données à une distribution de Pareto. Le code R ci - dessous calcule le MLE pour un ensemble de données, X.

pareto.MLE <- function(X)
{
   n <- length(X)
   m <- min(X)
   a <- n/sum(log(X)-log(m))
   return( c(m,a) ) 
}

# example. 
library(VGAM)
set.seed(1)
z = rpareto(1000, 1, 5) 
pareto.MLE(z)
[1] 1.000014 5.065213

Edit: Sur la base du commentaire de @cardinal et I ci-dessous, nous pouvons également noter que est l'inverse de la moyenne de l'échantillon des , ce qui arrive à ont une distribution exponentielle. Par conséquent, si nous avons accès à un logiciel qui peut s'adapter à une distribution exponentielle (ce qui est plus probable, car il semble se produire dans de nombreux problèmes statistiques), l'ajustement d'une distribution de Pareto peut être accompli en transformant l'ensemble de données de cette manière et en l'adaptant à une distribution exponentielle à l'échelle transformée. $\hat{\alpha}$ $\log(X_i /\hat{m})$

Macro
la source

(+1) Nous pouvons écrire les choses de façon un peu plus suggestive en notant que est distribué de façon exponentielle avec le taux . à partir de cela et de l'invariance des MLE en transformation, nous concluons immédiatement que , où nous remplaçons par dans cette dernière expression. Cela indique également comment nous pourrions utiliser un logiciel standard pour adapter un Pareto même si aucune option explicite n'est disponible.

Y_{i} = \log (X_{i} / m)

$Y_i = \log(X_i/m)$

α

$\alpha$

\hat{α} = 1 / \bar{Y}

$\hat\alpha = 1/\bar Y$

m

$m$

\hat{m}

$\hat m$

cardinal

@cardinal - Donc, est l'inverse de la moyenne de l'échantillon des , qui se trouvent avoir une distribution exponentielle. Comment cela nous aide-t-il?

\hat{α}

$\hat{\alpha}$

\log (X_{i} / \hat{m})

$\log(X_i/\hat{m})$

Macro

Salut, macro. Le point que j'essayais de faire était que le problème d'estimation des paramètres d'un Pareto peut être (essentiellement) réduit à celui de l'estimation du taux d'une exponentielle: via la transformation ci-dessus, nous pouvons convertir nos données et problème en un (peut-être) plus familier et extrayez immédiatement la réponse (en supposant que nous, ou notre logiciel, savons déjà quoi faire avec un échantillon d'exponentielles).

cardinal

Comment puis-je mesurer l'erreur de ce type d'ajustement?

emanuele

@emanuele, la variance approximative d'un MLE est l'inverse de la matrice d'informations du pêcheur, qui vous obligera à calculer au moins une dérivée de la log-vraisemblance. Ou, vous pouvez utiliser une sorte de rééchantillonnage bootstrap pour estimer l'erreur standard.

Macro

Comment puis-je adapter un ensemble de données à une distribution de Pareto dans R?

Réponses: