7
$\begingroup$

I'm trying to fit a weibull distribution to this but am having problems. Not sure why. What causes the NaNs?

temp <- dput(temp)
c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22, 759.46, 
1142.33, 134, 1232.23, 389.81, 7811.65, 992.11, 1152.4, 3139.01, 
2636.78, 3294.75, 2266.95, 32.12, 7356.84, 1448.54, 3606.82, 
465.39, 950.5, 3721.49, 522.01, 1548.62, 2196.3, 256.8, 2959.72, 
214.4, 134, 2307.79, 2112.74)


fitdist(temp, distr = "weibull", method = "mle")
Fitting of the distribution ' weibull ' by maximum likelihood 
Parameters:
          estimate  Std. Error
shape    0.8949019   0.1205351
scale 1803.8816283 357.9042207
Warning messages:
1: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
2: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
3: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
4: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
5: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
6: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
7: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
8: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
9: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
10: In dweibull(c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22,  :
  NaNs produced
$\endgroup$
  • 4
    $\begingroup$ Both parameters of the Weibull distribution are positive. If you type fitdist(temp, distr = "weibull", method = "mle", lower = c(0, 0)) no errors will be produced. Otherwise, optim tries negatives values which in turn produce the error. $\endgroup$ – COOLSerdash Jun 22 '15 at 19:43
  • $\begingroup$ Thanks! Why was it trying parameters below as a default -- I would assume that the Weibull would only fit positive shape and scale parameters? Similarly, what should the "lower" be for gamma and how do you know? It's not documented in ?fitdist. $\endgroup$ – user1357015 Jun 22 '15 at 20:11
  • $\begingroup$ Although this is a question about understanding an R error message (ie presumably off topic programming Q), it seems the underlying issue might be understanding the nature of the Weibull distribution. Thus, it may be on topic here. $\endgroup$ – gung Jun 22 '15 at 20:11
  • $\begingroup$ What package is the function fitdist coming from? I see fitdistr in MASS. CRAN can be quite lax with package evaluation. $\endgroup$ – AdamO Jun 22 '15 at 20:11
  • $\begingroup$ @COOLSerdash, since that seems to resolve this, why not turn your comment into an official answer (& maybe expand it a bit)? $\endgroup$ – gung Jun 22 '15 at 20:12
17
$\begingroup$

The Weibull distribution has two parameters, the scale $\lambda$ and shape $k$ (I'm following Wikipedias notation). Both parameters are positive real numbers.

The function fitdist from the fitdistrplus package uses the optim function to find the maximum likelihood estimations of the parameters. By default, optim imposes no limits to the parameters and tries out negative numbers as well. But negative values for the scale or shape produce NaNs for the Weibull distribution. By using the options lower and upper, you can impose limits on the parameter search space for optim.

The gamma distribution also has two parameters and as with the Weibull distribution, both are positive. So the same limits lower = c(0, 0) can be used for the gamma distribution.

Edit

Here is a small comparison of the Weibull and gamma fit for the posted data. The errors for the gamma distribution arise because of bad starting values. I provide them manually and then it works fine without errors.

library(fitdistrplus)

temp <- c(477.25, 2615.56, 1279.98, 581.57, 13.55, 80.4, 6640.22, 759.46, 
          1142.33, 134, 1232.23, 389.81, 7811.65, 992.11, 1152.4, 3139.01, 
          2636.78, 3294.75, 2266.95, 32.12, 7356.84, 1448.54, 3606.82, 
          465.39, 950.5, 3721.49, 522.01, 1548.62, 2196.3, 256.8, 2959.72, 
          214.4, 134, 2307.79, 2112.74)

fit.weibull <- fitdist(temp, distr = "weibull", method = "mle", lower = c(0, 0))
fit.gamma <- fitdist(temp, distr = "gamma", method = "mle", lower = c(0, 0), start = list(scale = 1, shape = 1))

Plot the fit for the Weibull:

plot(fit.weibull)

Weibull fit

And for the gamma distribution:

plot(fit.gamma)

Gamma fit

They are practically indistinguishable. The AICs are virtually the same for both fits:

gofstat(list(fit.weibull, fit.gamma))

Goodness-of-fit statistics
                             1-mle-weibull 2-mle-gamma
Kolmogorov-Smirnov statistic    0.07288424  0.07970184
Cramer-von Mises statistic      0.02532353  0.02361358
Anderson-Darling statistic      0.20489012  0.17609146

Goodness-of-fit criteria
                               1-mle-weibull 2-mle-gamma
Aikake's Information Criterion      601.7909    601.5659
Bayesian Information Criterion      604.9016    604.6766
$\endgroup$
  • $\begingroup$ lower = c(0,0) causes errors however... $\endgroup$ – user1357015 Jun 22 '15 at 20:38
  • $\begingroup$ @user1357015 When did the error occurr and what was it? $\endgroup$ – COOLSerdash Jun 22 '15 at 20:44
  • $\begingroup$ @COOLSerdah fitdist(temp, distr = "gamma", method = "mle", lower = c(0,0)) causes a "Error in fitdist(splitData[[1]], distr = "gamma", method = "mle", lower = c(0, : the function mle failed to estimate the parameters, with the error code 100 However, if I set c(0,1) it works fine though. $\endgroup$ – user1357015 Jun 22 '15 at 20:50
  • 1
    $\begingroup$ @user1357015 I have updated my answer accordingly. The error arises because of bad starting values. $\endgroup$ – COOLSerdash Jun 22 '15 at 21:00
  • $\begingroup$ Ok, so that works but now bootstrap also throws warnings, namely `bootdist(temp)'. I think some iterations fail to converge which mathematically seems fine since 993, out of 1001 iterations did converge -- but thought I'd double check :-) $\endgroup$ – user1357015 Jun 22 '15 at 23:37

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.