Difficulté à adapter la distribution gamma avec R

Je cherche à estimer les paramètres pour un gamma de distribution de l'adapter à la densité écologique (c'est à dire de la biomasse par unité de surface) de données. J'ai été en utilisant le fitdistr() la commande à partir de la MASSE paquet dans R (version 3.0.0 : x86_64-w64-mingw32/x64 (64-bit)). C'est une estimation du maximum de vraisemblance de commande pour la répartition des paramètres.

Les vecteurs de données sont assez grandes, mais des statistiques sommaires sont comme suit:

Min. = 0; 1st Qu. = 87.67; Médiane = 199.5; Moyenne = 1255; Variance =
2.79 E+07; 3rd Qu. = 385.6; Max. = 33880

Le code j'utilise pour exécuter le MLE procédure est la suivante:

gdist <- fitdistr(data, dgamma, 
                  start=list(shape=1, scale=1/(mean(data))),lower=c(1,0.1))

R me donne l'erreur suivante:

Erreur dans optim(x = c(6.46791148085828, 4060.54750836902,
99.6201565968665, : non-finis, différences finies valeur [1]

D'autres qui ont vécu ce type de problème et se sont tournés vers stackoverflow pour obtenir de l'aide semblent avoir trouvé la solution en ajoutant le "lower=" argument de leur code, et/ou de supprimer les zéros. Je trouve que le R sera de fournir des paramètres pour un ajustement si je supprimer le zéro observations, mais j'étais sous l'impression que le gamma distributions couverts de la gamme 0 <= x > inf (Forbes et al. 2011. Les Distributions Statistiques)?

Ai-je eu l'impression, à tort, sur la gamme de la distribution gamma? Ou est-il un autre problème je suis manquantes concernant la MLE (dans lequel je suis un novice).

source d'informationauteur niafall

  1. 22

    Obtenir une estimation par la méthode des moments (correspondant à la moyenne=forme*à l'échelle et de la variance=forme*à l'échelle^2) nous avons

    mean <- 1255
var <- 2.79e7
shape = mean^2/var   ## 0.056
scale = var/mean     ## 22231

    Maintenant produire des données à partir de cette distribution:

    set.seed(101)
x = rgamma(1e4,shape,scale=scale)
summary(x)
##     Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
##     0.00      0.00      0.06   1258.00     98.26 110600.00 

MASS::fitdistr(x,"gamma")  ## error

    Je soupçonne fortement que le problème est que le sous-jacent optim appel assume les paramètres de forme et d'échelle, ou de la forme et de taux) sont approximativement de la même importance, qu'ils ne le sont pas. Vous pouvez contourner ce problème par la mise à l'échelle de vos données:

    (m <- MASS::fitdistr(x/2e4,"gamma"))  ## works fine
##      shape           rate    
##  0.0570282411   0.9067274280 
## (0.0005855527) (0.0390939393)

    fitdistr donne un taux de paramètre plutôt qu'un paramètre d'échelle: pour revenir au paramètre de forme que vous souhaitez, retourner et re-scale ...

    1/coef(m)["rate"]*2e4  ## 22057

    Par ailleurs, le fait que les quantiles de la simulation, les données ne correspondent pas à vos données très bien (par exemple, la médiane de x=0.06 vs une moyenne de 199 dans vos données) suggèrent que vos données peuvent ne pas convenir à un Gamma de tout ce que le bien -- par exemple, il pourrait y avoir quelques valeurs aberrantes affectant la moyenne et la variance de plus que les quantiles?

    PS ci-dessus j'ai utilisé le built-in 'gamma' estimateur en fitdistr plutôt que d'utiliser dgamma: avec les valeurs de départ basé sur la méthode des moments, et de mise à l'échelle les données par 2e4, cela fonctionne (bien qu'il donne un avertissement au sujet de NaNs produced à moins que nous préciser lower)

     m2 <- MASS::fitdistr(x/2e4,dgamma,
        start=list(shape=0.05,scale=1), lower=c(0.01,0.01))