12

Given a dataset: 

x <- c(4.9958942,5.9730174,9.8642732,11.5609671,10.1178216,6.6279774,9.2441754,9.9419299,13.4710469,6.0601435,8.2095239,7.9456672,12.7039825,7.4197810,9.5928275,8.2267352,2.8314614,11.5653497,6.0828073,11.3926117,10.5403929,14.9751607,11.7647580,8.2867261,10.0291522,7.7132033,6.3337642,14.6066222,11.3436587,11.2717791,10.8818323,8.0320657,6.7354041,9.1871676,13.4381778,7.4353197,8.9210043,10.2010750,11.9442048,11.0081195,4.3369520,13.2562675,15.9945674,8.7528248,14.4948086,14.3577443,6.7438382,9.1434984,15.4599419,13.1424011,7.0481925,7.4823108,10.5743730,6.4166006,11.8225244,8.9388744,10.3698150,10.3965596,13.5226492,16.0069239,6.1139247,11.0838351,9.1659242,7.9896031,10.7282936,14.2666492,13.6478802,10.6248561,15.3834373,11.5096033,14.5806570,10.7648690,5.3407430,7.7535042,7.1942866,9.8867927,12.7413156,10.8127809,8.1726772,8.3965665)

.. I would like to determine the most fitting probability distribution (gamma, beta, normal, exponential, poisson, chi-square, etc) with an estimation of the parameters. I am already aware of the question on the following link, where a solution is provided using R: https://stackoverflow.com/questions/2661402/given-a-set-of-random-numbers-drawn-from-a-continuous-univariate-distribution-f the best proposed solution is the following:

> library(MASS)
> fitdistr(x, 't')$loglik                                                              #$
> fitdistr(x, 'normal')$loglik                                                         #$
> fitdistr(x, 'logistic')$loglik                                                       #$
> fitdistr(x, 'weibull')$loglik                                                        #$
> fitdistr(x, 'gamma')$loglik                                                          #$
> fitdistr(x, 'lognormal')$loglik                                                      #$
> fitdistr(x, 'exponential')$loglik                                                    #$

And the distribution with the smallest loglik value is selected. However, other distrubtions such as beta distribution require the specification of some addition parameters in the fitdistr() function:

   fitdistr(x, 'beta', list(shape1 = some value, shape2= some value)).

Given that i am trying to determine the best distribution without any prior information, i don't know what the value of the parameters can possibly be for each distribution. Is there another solution that takes this requirement into account? it does not have to be in R.

Community
  • 1

2 Answers2

19

What do you do about the infinity of distributions that aren't in the list?

What do you do when none of the ones in your list fit adequately? e.g. if your distribution is strongly bimodal

How are you going to deal with the fact that the exponential is just a special case of the gamma, and as such, the gamma must always fit any set of data better, since it has an additional parameter, and hence must have a better likelihood?

How do you deal with the fact that the likelihood is only defined up to a multiplicative constant and that the likelihood for different distributions might not automatically be comparable unless defined consistently?

It's not that these are necessarily insoluble, but doing this stuff in a sensible way is nontrivial; certainly more thought is required than just bunging everything through the calculation of a MLE and comparison of likelihoods.

Glen_b
  • 257,508
  • 32
  • 553
  • 939
  • 3
    I only care for the distributions on the list, and if none of the distributions fit, then i ll deal with that problem next. but for now reaching that point is good enough for me. As for the last question regarding the gamma distribution, yes it can fit better with some parameter, my question is exactly about this, is there an algorithm that allows me to loop through different parameter values for the different distributions on the list? and return the most fitting distribution with the appropriate paramaters? –  Jun 14 '12 at 10:09
  • 1
    Well, yes and no. You can come up with a "figure of merit" to calculate which parameter set gives you the best fit, and write a "hill-climb" loop to optimize the value of the figure of merit. One example of a FOM is the R-value for regression fits. – Carl Witthoft Jun 14 '12 at 12:01
  • i subscribe to this +1. –  Jun 14 '12 at 15:05
  • 4
    @shachem You missed the point about the gamma. Distributions with additional parameters will always have a better likelihood, *even when the data come from the distribution with fewer parameters*. You need to consider this. Some measures of fit adjust for this effect. To be honest, I think the thing you're trying to do is likely 'answering the wrong question', somewhat like asking 'how do I figure out which of these hammers is the best one for pounding in this screw' – Glen_b Jun 15 '12 at 00:54
  • 1
    If "the likelihood is only defined up to a multiplicative constant," Glen, then how can one make sense of your statement that "distributions with additional parameters will always have a better likelihood"? Indeed, how could one possibly compare likelihoods that are so ill-defined? I suspect some key idea has been left unstated... – whuber Jun 15 '12 at 12:57
  • @whuber +1 I just noticed I never responded to that comment for some reason. You're right, of course. The difference is that in the case of distributions which are special cases of some other family, if you treat both distributions as members of the broader family (write it the same way), the likelihoods are comparable if you scale by the same constants. So it can make sense to compare two gammas or even an exponential nested in the gamma family, in a way that at least some people (Ripley is one example, IIRC) argue may not as readily make sense when comparing, say, gamma and lognormal. – Glen_b Mar 20 '14 at 20:10
3

I have found a function that answers my question using matlab. It can be found on this link: http://www.mathworks.com/matlabcentral/fileexchange/34943

I takes a data vector as input

   allfitdist(data)

and returns the following information for the best fitting distribution:

   DistName- the name of the distribution
   NLogL - Negative of the log likelihood
   BIC - Bayesian information criterion (default)
   AIC - Akaike information criterion
   AICc - AIC with a correction for finite sample sizes 
   ParamNames
   ParamDescription
   Params
   etc.
  • 1
    OK, now all you need to do is port the m-file to an R-file. :-) – Carl Witthoft Jun 14 '12 at 15:25
  • Not even! As i had mentioned it does not need to be an R-file, so the matlab function completely solves my problem :-) –  Jun 14 '12 at 15:40
  • AIC, BIC, AICc, etc. may be one way of deciding, but whether that makes sense really depends on what you are trying to do. E.g. what do you do if one distribution is really close to the others in terms of log L? If what you are trying to do is to predict, then usually ignoring models that are close in terms of your criterion should not just be discarded. E.g. model-averaging is one way of taking them into account. – Björn Nov 17 '15 at 11:56