1

I'm fitting a beta distribution to my data using MATLAB.

To do so, there are several options like fitdist, which provide a set of a and b and the stats (mean and std) are accessible using betastat or simply mean and std commands.

But what can I do if I want my distribution have a given mean, and only the std be unknown.

There is an example. I have a variable-bin-size semi-histogram with defined central ratios and their probabilities as follow:

central_ratios = [.005 .1 .4 .8   1]; 
probabilities  = [.5   .3 .1 .08 .02];
vul = sum(central_ratios.*probabilities);

Now I want to fit a 'beta distribution' to the data, bounded between [0,1] with a mean equal to vul. all I need now is the best std or a and b pair.

parisa
  • 29
  • 6
  • Is this the Beta distribution bounded between [0,1] or the generalized Beta distribution which is bounded between [x1,x2 where x1,x2 are real numbers? The answer depends on this. In both cases, the Beta distribution also has two shape parameters, either denoted alpha1, alpha2, or in MATLAB's case, a & b. – SecretAgentMan May 07 '19 at 20:54
  • Also, see [How to Ask](https://stackoverflow.com/help/how-to-ask). If you [edit] your question to include a [Minimal, Complete, and Verifiable example](https://stackoverflow.com/help/mcve), you'll get more specific feedback and/or assistance. – SecretAgentMan May 07 '19 at 21:01
  • I did edit the question. hope it helps getting better answers.@SecretAgentMan – parisa May 08 '19 at 09:01

1 Answers1

0

fitdist for a beta distribution gets the distribution parameters from betafit, which sets up an appropriate likelihood function for the distribution given your data, some heuristic initial guesses for a and b and then optimises log(a) and log(b) to maximise likelihood using fminsearch.

Your constraint defining the distribution mean establishes an enforced relationship between a and b. From Wikipedia the mean mu is related to a and b thus:

mu = 1 / (1 + b/a)

This can be rearranged to give one distribution parameter given the other:

b = a * (1/mu - 1)

To examine the unconstrained implementation of beta distribution fitting that is available in MATLAB, and that you're seeking to constraint, you can view the betafit source code using:

edit betafit

At least in MATLAB R2018b you will find that the optimisation of log(a) and log(b) takes place at a point which declares:

phat = fminsearch(negloglike,pstart,opts);
phat = exp(phat);

Your constrained distribution fitting problem can be described in terms of the optimised objective function used here, which could allow you to reuse the other aspects of betafit's behaviour:

negloglike1 = @(loga) negloglike([loga log(exp(loga) * (1/mu - 1))]);

You could either create your own duplicate of betafit which makes this declaration before the call to fminsearch, or stop the built-in betafit on a breakpoint after the line has been called and declare the new likelihood function from the command line. Either way, you can then replace the parameters with constrained ones which maximise likelihood within this constraint:

loga = fminsearch(negloglike1,pstart(1),opts);
phat = exp(loga) * [1 (1/mu - 1)];

The resulting beta distribution parameters phat will be guaranteed to result in a distribution with mean mu, and locally maximise the likelihood function for your data given this constraint.

Will
  • 1,835
  • 10
  • 21
  • thanks but the steps are unclear.. and I believe it is "negloglik". I would appreciate it if you could edit the answer with a complete executed sample. – parisa May 07 '19 at 15:36
  • 1
    I've tried to clarify the answer a bit. `negloglike` is a function handle used internally in `betafit`, which you can inspect as described. This is distinct from `negloglik` which is a method of distribution objects, not used internally to `betafit`. You could, alternatively, make an objective function which makes a call to `makedist` using trial parameters and then calls `negloglik`, but you would be incurring a lot of overhead on each iteration in `makedist`. If you want answers with complete executed samples, your question needs a complete executed example from which to reference. – Will May 07 '19 at 16:22