1

I want to simulate from a beta distribution which has the following properties:

  1. The mode (peak) of the distribution is closer to 0.6.
  2. The first quartile and the third quartile should be 0.3 and 0.9 respectively.

Is there is a way to calculate the parameters of beta distribution ($\alpha$,$\beta$) based on this 2 conditions using R? Since the CDF of beta is not in closed form, I couldn't think of a way to do this.

Any help would be highly appreciated.

student_R123
  • 751
  • 4
  • 13
  • Maybe using [ABC algorithm](https://en.wikipedia.org/wiki/Approximate_Bayesian_computation) with quantile and mode as summary function ? – Kledou Sep 16 '21 at 16:02

1 Answers1

1

Since the beta distribution has only two parameters, but you want to satisfy three conditions, you will probably not be able to get a perfect fit.

One possibility would be to define an objective function that takes two parameters $(\alpha,\beta)$ and calculates the deviation between the quantiles and the mode of the corresponding beta distribution on the one hand, and your target values on the other hand. Then run a straightforward optimizer over it. If you use the sum of the squared deviations, this is easy to numerically optimize.

foo <- optim(par=c(1.2,1.2),
    fn=function(par)(qbeta(0.25,par[1],par[2])-0.3)^2+ # first quartile
                    (qbeta(0.75,par[1],par[2])-0.9)^2+ # third quartile
                    ((par[1]-1)/(sum(par)-2)-0.6)^2,   # mode
    method="L-BFGS-B",lower=1)

Unfortunately, it turns out that the result fits your conditions on the first quartile and the mode rather well - but the third quartile quite badly. Essentially, making the beta distribution conform to the first two conditions makes it very hard to also satisfy the third one.

> qbeta(c(0.25,0.75),foo$par[1],foo$par[2])
[1] 0.2964684 0.7504854
> (foo$par[1]-1)/(sum(foo$par)-2)
[1] 0.6273578

Here is the density:

beta

xx <- seq(0,1,by=0.01)
plot(xx,dbeta(xx,foo$par[1],foo$par[2]),type="l",las=1,xlab="",ylab="")

You could start playing around with weights on your three conditions, but then of course, as you improve the third quartile, the other two conditions suffer. Or you could try looking for a three-parameter generalization of the beta, which you might be able to match exactly to your requirements.

Stephan Kolassa
  • 95,027
  • 13
  • 197
  • 357
  • 2
    The OP was wise to specify a *fuzzy* condition ("closer to"), implicitly recognizing your point about conditions and parameters. It therefore suffices to find the Beta distribution with the given quartiles and check whether its mode is sufficiently close to 0.6. – whuber Sep 16 '21 at 15:56
  • @whuber Can we say with certainty that only one beta distribution will have those particular quartiles? – Dave Sep 16 '21 at 15:58
  • @Dave: it appears that yes, judging from a skim of the first duplicate. (However, it the problem may be quite ill conditioned, so numerics might be problematic.) – Stephan Kolassa Sep 16 '21 at 16:00
  • 1
    @Dave Yes: see the duplicates. Stephan: see the duplicates for remarks on avoiding numerical instability and working code that performs surprisingly well even with severely ill-conditioned inputs (such as two very close extreme quantiles). – whuber Sep 16 '21 at 16:03