нижний предел гамма-распределения в R

Я хотел бы подогнать гамма-распределение к набору данных, состоящему из 338 элементов, с фиксированным низким порогом (я использую R). Чтобы представить нижний предел, я решил использовать гамму с тремя параметрами, добавив местоположение. Вот мой код:

library(FAdist)
library(fitdistrplus)
fit <- fitdist(mydata,"gamma3",start=list(1,1,mythreshold))

Каждый раз, когда я запускаю код, я получаю ту же ошибку:

<simpleError in optim(par = vstart, fn = fnobj, fix.arg = fix.arg, obs = data,     gr = gradient, ddistnam = ddistname, hessian = TRUE, method = meth,     lower = lower, upper = upper, ...): non-finite finite-difference value [3]>
Error in fitdist(Hs, "gamma3", start = list(1, 1, 3)) : 
  the function mle failed to estimate the parameters, 
                with the error code 100

Поэтому я попытался изменить начальные значения, используя те, которые были получены путем подгонки гаммы с двумя параметрами, что дало такой результат:

fit <- fitdist(mydata,"gamma")
fit
Fitting of the distribution ' gamma ' by maximum likelihood 
Parameters:
       estimate Std. Error
shape 21.417503  1.6348313
rate   5.352422  0.4133735

Но все равно не работает.. Я бы понял, если бы даже гамма с двумя параметрами не работала, но это не так и я не могу себе объяснить. Более того, график QQ и ecdf для гаммы с двумя параметрами не так хороши... но если я подгоню распределение к исходному набору данных, масштабированному относительно нижнего порога, подгонка выглядит идеально:

fit <- fitdist(mydata-mythreshold,"gamma")
fit
Fitting of the distribution ' gamma ' by maximum likelihood 
Parameters:
      estimate Std. Error
shape 1.059540 0.07212832
rate  1.058007 0.09117620

Но тогда я не знаю, правильно ли это делать... параметры очень разные, и они мне нужны для расчета периодов возврата, связанных с моими данными! вот почему я подумал о гамме с параметром местоположения.

p.s. я не приложил данные, потому что их слишком много, но я могу сообщить их сводку:

 summary(mydata)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  3.003   3.282   3.753   4.001   4.444   8.087 

person F. De Leo    schedule 26.12.2016    source источник


Ответы (1)


Настройте воспроизводимый пример;

library(FAdist); library(fitdistrplus)
set.seed(101)
x <- rgamma3(1000,shape=1,scale=1,thres=1)

Для значений порогового параметра, превышающих минимальное значение в наборе данных, вероятность бесконечна (поскольку эти значения предполагаются невозможными/имеют нулевое значение при усеченном распределении). Мы можем заставить все работать, установив верхнюю границу для этого параметра:

fitdist(x,dgamma3,start=list(shape=1,scale=1,thres=0.5),
        upper=c(Inf,Inf,min(x)))
## Fitting of the distribution ' gamma3 ' by maximum likelihood 
## Parameters:
##        estimate Std. Error
## shape 0.9321949         NA
## scale 1.0586079         NA
## thres 1.0000012         NA

Обратите внимание: (1) при максимальной вероятности этот тип порогового параметра всегда будет иметь наибольшее значение, которое он может принимать (т. Е. Минимальное значение в наборе данных); (2) стандартные ошибки равны NA, потому что стандартные ошибки Вальда не могут быть вычислены, когда параметры достигают границы.

В качестве альтернативы вы можете исправить пороговый параметр, определив функцию-оболочку:

dgconstr <- function(x,shape,scale,...) {
    dgamma3(x,shape,scale,thres=0.5,...)
}
pgconstr <- function(...) {
    pgamma3(...,thres=0.5)
}

fitdist(x,dgconstr,start=list(shape=1,scale=1))
person Ben Bolker    schedule 26.12.2016
comment
Большое спасибо Бен! Я получил это, и это окончательно сработало. - person F. De Leo; 29.12.2016