0
$\begingroup$

I'm trying to find a distribution that fits my data in order that I can then predict a 5th percentile but none of the standard distributions seem to fit.

I'll explain my approach so far with examples below. I have then been trying to fit a Burr distribution in fitdistrplus but cannot find any suitable initial values. But I am unsure it this is because this is an incorrect distribution for the data or if its me. I've also experimented with prefit function but am getting a similar problem of not being able to choose feasible starting values. Maybe there is a more appropriate distribution I haven't tried?

My data are:

0.0001900 0.0002100 0.0002200 0.0003000 0.0007800 0.0008400 0.0011000 0.0011300 0.0012000 0.0016000 0.0016000 0.0020000 0.0020000 0.0031000 0.0056500 0.0059000 0.0082449 0.0130000 0.0180000 0.0191000 0.0510000

The Cullen and Frey graph using the following code is as follows:

descdist(Data, boot = 500)

enter image description here

from this I thought the Beta distribution may work best but it isn't quite right. I can only post 2 images so have only included the QQ plot here:

fitln <- fitdist(SSD2$NOEC,"lnorm")
fitW <- fitdist(SSD2$NOEC, "weibull")
fitg <- fitdist(SSD2$NOEC, "gamma")
fitn <- fitdist(SSD2$NOEC, "norm")
fitexp <- fitdist(SSD2$NOEC,"exp")
fitB <- fitdist(SSD2$NOEC,"beta")
fitP <- fitdist(SSD2$NOEC,"pareto")

cdfcomp(list(fitW, fitg, fitln, fitn, fitexp, fitB, fitP), 
    legendtext=c("Weibull", "gamma", "lognormal", "norm", "exp", "Beta", "Pareto"))
denscomp(list(fitW, fitg, fitln, fitn, fitexp, fitB, fitP), 
     legendtext=c("Weibull", "gamma", "lognormal", "norm", "exp", "Beta", "Pareto"))
qqcomp(list(fitW, fitg, fitln, fitn, fitexp, fitB, fitP), 
   legendtext=c("Weibull", "gamma", "lognormal", "norm", "exp", "Beta", "Pareto"))

enter image description here

ppcomp(list(fitW, fitg, fitln, fitn, fitexp, fitB, fitP), 
   legendtext=c("Weibull", "gamma", "lognormal", "norm", "exp", "Beta", "Pareto"))
gofstat(list(fitW, fitg, fitln, fitn, fitexp, fitB, fitP))

Goodness-of-fit statistics
                             1-mle-weibull 2-mle-gamma 3-mle-lnorm 4-mle-norm 5-mle-exp 6-mle-beta
Kolmogorov-Smirnov statistic    0.17899327   0.2201135  0.13002018  0.2888904 0.3552864  0.2208923
Cramer-von Mises statistic      0.08976409   0.1567659  0.04560304  0.5880214 0.6158410  0.1589649
Anderson-Darling statistic      0.53379880   0.8341498  0.30557002  3.1553792 3.5184364  0.8457680
                             7-mle-pareto
Kolmogorov-Smirnov statistic   0.12410532
Cramer-von Mises statistic     0.04296412
Anderson-Darling statistic     0.29429372

Goodness-of-fit criteria
                               1-mle-weibull 2-mle-gamma 3-mle-lnorm 4-mle-norm 5-mle-exp 6-mle-beta
Akaike's Information Criterion     -173.8742   -171.8001   -177.3758  -124.3444 -167.3059  -171.7056
Bayesian Information Criterion     -171.7851   -169.7111   -175.2868  -122.2554 -166.2614  -169.6166
                               7-mle-pareto
Akaike's Information Criterion    -176.0548
Bayesian Information Criterion    -173.9658

The code I have used so far to try burr distribution is as follows but am struggling with the initial values - I have tried a few variations of shape 1 and 2:

fitBurr <- fitdist(SSD2$NOEC,"burr", start = list(shape1 = 1, shape2 = 3, rate = 1))
prefit(SSD2$NOEC,"burr", method = "mle", start = list(shape1 = 1, shape2 = 3))

Any help in getting a suitable distribution for this data would be greatly appreciated.

Improve this question
$\endgroup$
3
  • 1
    $\begingroup$ What's the parameterization of this Burr? Are "shape1" and "shape2" the $c$ and $k$ (respectively) in the Wikipedia page on the Burr type Xii? Is the rate the inverse of the scale ($\lambda$) there? Looking at it, I'd expect that a best fitting $c$ would be less than 1 and $k$ a bit greater than one, if anything. $\endgroup$
    – Glen_b
    Commented Aug 8, 2017 at 3:10
  • $\begingroup$ It may be that most "standard" type distributions won't fit the lower tail all that well; the bulk of the distribution is close to lognormal but the left tail is a bit shorter than it would be for a lognormal or even a shifted lognormal; this would suggest a gamma probably wouldn't suit either (since its lower tail will tend to be longer, not shorter, than a lognormal), this would then rule out the less general exponential as well. On similar grounds you could rule out a number of the others. $\endgroup$
    – Glen_b
    Commented Aug 8, 2017 at 3:11
  • $\begingroup$ Yes - I agree that its the lower tail that is is the cause of the deviation from the distributions. I believe that the shapes 1 and 2 are c and k respectively as you say but im still failing to find feasible starting values. $\endgroup$
    – Gareth Le Page
    Commented Aug 9, 2017 at 10:13

1 Answer 1

Reset to default
0
$\begingroup$

I tried beta and got what I think was a good result:

Parameters: a = 1.6589757232166452E-01 b = 1.0014191317056027E+00 location = 1.8999999999999998E-04 scale = 5.2536884151069704E-02

the scipy documentation for beta is here for the distribution details: http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.beta.html

I fit your data to the 80+ continuous statistical distributions in scipy using an open source statistical distribution fitter I had written years ago. The beta distribution was indeed at the top of the results list.

beta

Improve this answer
$\endgroup$
3
  • $\begingroup$ Beta doesn't fit the lower tail at all. Look at the density of the log of this beta random variable compared to the log-data. I believe this will actually be worse at estimating the low quantiles (like the 5th percentile) than the lognormal $\endgroup$
    – Glen_b
    Commented Aug 9, 2017 at 10:16
  • $\begingroup$ My result for your suggestion of lognormal shows that it is indeed quite a good fit to this data: Parameters: s = 2.3544583746739800E+00 location = 1.8727097624690790E-04 scale = 1.2329959394857546E-03 with scipy documentation for the distribution at docs.scipy.org/doc/scipy/reference/generated/… $\endgroup$
    – James Phillips
    Commented Aug 9, 2017 at 11:20
  • $\begingroup$ Yes, I know it's a good fit overall (I tried it myself on your data the other day), but the problem is all these possibilities still get the low percentiles wrong - right where you need the fit to be good; the three parameter lognormal isn't terrible but it's not great either. You may need a truncated distribution or a mixture. Perhaps a four-parameter beta might do pretty well but I've been assuming you have a two parameter beta there. $\endgroup$
    – Glen_b
    Commented Aug 9, 2017 at 12:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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