2013-04-12 13 views
8

mam problem z fitdistr {} function MASS w R. Mam ten wektor:Błędy podczas próby dopasować rozkład gamma z R fitdistr {MASS}

a <- c(26,73,84,115,123,132,159,207,240,241,254,268,272,282,300,302,329,346,359,367,375,378, 384,452,475,495,503,531,543,563,594,609,671,687,691,716,757,821,829,885,893,968,1053,1081,1083,1150,1205,1262,1270,1351,1385,1498,1546,1565,1635,1671,1706,1820,1829,1855,1873,1914,2030,2066,2240,2413,2421,2521,2586,2727,2797,2850,2989,3110,3166,3383,3443,3512,3515,3531,4068,4527,5006,5065,5481,6046,7003,7245,7477,8738,9197,16370,17605,25318,58524) 

i chcę dopasować rozkład gamma dane z polecenia:

fitted.gamma <- fitdistr(a, "gamma") 

ale mam taki błąd:

Error in optim(x = c(26, 73, 84, 115, 123, 132, 159, 207, 240, 241, 254, : 
non-finite finite-difference value [1] 
In addition: Warning messages: 
1: In densfun(x, parm[1], parm[2], ...) : NaNs produced 
2: In densfun(x, parm[1], parm[2], ...) : NaNs produced 
3: In densfun(x, parm[1], parm[2], ...) : NaNs produced 
4: In densfun(x, parm[1], parm[2], ...) : NaNs produced 

Więc próbowałem z inicjalizacji parametry:

(fitted.gamma <- fitdistr(a, "gamma", start=list(1,1))) 

Przedmiotem fitted.gamma jest tworzony, ale po wydrukowaniu, tworzy błąd:

Error in dn[[2L]] : subscript out of bounds 

Czy wiesz, co się dzieje, czy może znasz jakieś inne funkcje R, aby dopasować jednowymiarowe rozkłady autorem MLE?

Z góry dziękuję za pomoc lub odpowiedź.

Kuba

Odpowiedz

9

zawsze wykreślić swoje rzeczy po pierwsze, skalowanie jest daleko offfffffff.

library(MASS) 
a <- c(26,73,84,115,123,132,159,207,240,241,254,268,272,282,300,302,329,346,359,367,375,378, 384,452,475,495,503,531,543,563,594,609,671,687,691,716,757,821,829,885,893,968,1053,1081,1083,1150,1205,1262,1270,1351,1385,1498,1546,1565,1635,1671,1706,1820,1829,1855,1873,1914,2030,2066,2240,2413,2421,2521,2586,2727,2797,2850,2989,3110,3166,3383,3443,3512,3515,3531,4068,4527,5006,5065,5481,6046,7003,7245,7477,8738,9197,16370,17605,25318,58524) 
## Ooops, rater wide 
plot(hist(a)) 
fitdistr(a/10000,"gamma") # gives warnings 
# No warnings 
fitted.gamma <- fitdistr(a/10000, dgamma, start=list(shape = 1, rate = 0.1),lower=0.001) 

Teraz można zdecydować, co zrobić z skalowanie

+0

Dzięki za odpowiedź. Widzę, że dodanie "niższego" argumentu ze skalowaniem sprawiło, że trik został wykonany. Czy to oznacza, że ​​optymalizując parametry gamma, przyjmujemy w pewnym momencie wartości ujemne? Jeśli chodzi o skalowanie, dlaczego konieczne jest skalowanie wartości (parametr szybkości jest niski)? Kuba – kuba

+0

Tak, podczas optymalizacji gradientu łatwo wpadamy w obszary o słabym gradiencie dla niektórych próbek. Gamma może nie być właściwą dystrybucją, po prostu spróbuj wykreślić kilka przykładów. Jednak log (a) wygląda prawie normalnie ... –

+0

Bardzo dziękuję za pomoc :) – kuba

1

przypadku danych, które wyraźnie pasuje rozkład gamma, ale jest na niewłaściwej skali (czyli jak gdyby zostały pomnożone/dzielone przez duże numer), tutaj jest alternatywne podejście do montażu rozkład gamma:

fitgamma <- function(x) { 
    # Equivalent to `MASS::fitdistr(x, densfun = "gamma")`, where x are first rescaled to 
    # the appropriate scale for a gamma distribution. Useful for fitting the gamma distribution to 
    # data which, when multiplied by a constant, follows this distribution 
    if (!requireNamespace("MASS")) stop("Requires MASS package.") 

    fit <- glm(formula = x ~ 1, family = Gamma) 
    out <- MASS::fitdistr(x * coef(fit), "gamma") 
    out$scaling_multiplier <- unname(coef(fit)) 
    out 
} 

Zastosowanie:

set.seed(40) 
test <- rgamma(n = 100, shape = 2, rate = 2)*50000 
fitdistr(test, "gamma") # fails 
dens_fit <- fitgamma(test) # successs 
curve(dgamma(x, 2, 2), to = 5) # true distribution 
curve(dgamma(x, dens_fit$estimate['shape'], dens_fit$estimate['rate']), add=TRUE, col=2) # best guess 
lines(density(test * dens_fit$scaling_multiplier), col = 3) 

plot of density

Powiązane problemy