Beta Distribution – Understanding Intuition for Beta Distribution with Alpha and/or Beta Less Than 1

bayesianbeta distributionintuitionjeffreys-prior

I am curious for myself, but also trying to explain this to others.

The beta distribution is often used as a Bayesian conjugate prior for a binomial likelihood. It is often explained with the example that $\left(\alpha-1\right)$ is analogous to the number of successes and $\left(\beta-1\right)$ is like the number of fails.

As expected, a beta distribution with $\alpha=\beta=1$ is equivalent to a uniform distribution.

But the beta distribution can have values less than 1 (any non-negative number). At the extreme case, $\alpha=\beta=0$ yields a bimodal PDF (probability density function) with values at only 0 and 1. I can still intuit this: it represents a case like flipping a coin – not the probability of heads or tails, but rather the outcomes: there are only 2 possibilities, 0 or 1 (or heads or tails).

But any $\alpha$ or $\beta$ value between 0 and 1 I cannot find a good way to explain or think about. I can calculate it, but not really grok it.

Bonus points for anyone who can help explain the difference between a conjugate prior using what to me seems it should provide no information, which would be a beta distribution with $\alpha=\beta=1$, and what is actually used as a prior with no information, the Jeffrey's Prior, which uses $\alpha=\beta=0.5$.

Addendum

Looks like I need to be clearer. I am looking to understand, conceptually what natural phenomenon might be represented by a beta distribution with $\alpha=\beta=\frac{1}{2}$.

For instance,

Binomial distribution with n=10 and k=4 "means": some phenomenon with a bimodal response experienced 4 "successes" in 10 attempts.
Poisson distribution with k=2 and $\lambda=4.5$ means: some phenomenon that "typically" happens 4.5 times per hour (or whatever unit of time) only happened twice in the interval.

Or even with positive integer beta distributions, I can say:

Beta distribution with $\alpha=4$ and $\beta=7$ means: some phenomenon with a bimodal response had 3 successes and 6 fails in 9 attempts.
- (I know this one is a bit inaccurate, since beta distributions are continuous and provide a probability density instead of mass, but this is often how it is conceptually viewed or explained, and why it is used as a conjugate prior.)

What sort of similar construct or meaning could I create for the beta distribution with $\alpha=\beta=\frac{1}{2}$?

I am not looking for a plot. As I said earlier, I know how to work with a beta distribution mathematically (plot it, calculate it, etc.) I am just trying to get some natural intuition.

Best Answer

Here is a frivolous example that may have some intuitive value.

In US Major League Baseball each team plays 162 games per season. Suppose a team is equally likely to win or lose each of its games. What proportion of the time will such a team have more wins than losses? (In order to have symmetry, if a team's wins and losses are tied at any point, we say it is ahead if it was ahead just before the tie occurred, otherwise behind.)

Suppose we look at a team's win-loss record as the season progresses. For our team with wins and losses are as if determined by tosses of a fair coin, you might think a team would most likely be ahead about half the time throughout a season. Actually, half the time is the least likely proportion of time for being ahead.

The "bathtub shaped" histogram below shows the approximate distribution of the proportion of time during a season that such a team is ahead. The curve is the PDF of $\mathsf{Beta}(.5,.5).$ The histogram is based on 20,000 simulated 162-game seasons for a team where wins and losses are like independent tosses of a fair coin, simulated in R as follows:

set.seed(1212);  m = 20000;  n = 162;  prop.ahead = numeric(m)
for (i in 1:m)
 {
 x = sample(c(-1,1), n, repl=T);  cum = cumsum(x)
 ahead = (c(0, cum) + c(cum,0))[1:n]  # Adjustment for ties
 prop.ahead[i] = mean(ahead >= 0)
 }

cut=seq(0, 1, by=.1); hdr="Proportion of 162-Game Season when Team Leads"  
hist(prop.ahead, breaks=cut, prob=T, col="skyblue2", xlab="Proportion", main=hdr)
curve(dbeta(x, .5, .5), add=T, col="blue", lwd=2)

Note: Feller (Vol. 1) discusses such a process. The CDF of $\mathsf{Beta}(.5,.5)$ is a constant multiple of an arcsine function, so Feller calls it an 'Arcsine Law'.

Related Solutions

Solved – Approximation for Beta distribution when alpha is less than 10

Johnson and Kotz point out that

The $B$ (Beta) and $F$ distributions are closely related: When $k$ has a $B(\nu_1, \nu_2)$ distribution, then $f$ has an $F(2\nu_1, 2\nu_2)$ distribution where

$$f = \frac{k \nu_2}{(1-k)\nu_1}.$$
The Wilson-Hilferty approximation to $\chi^2$ (or, equivalently, a $\Gamma$ distribution)--which has "remarkable accuracy"--can be applied to $F$, which itself is a ratio of $\chi^2$ distributions.

The result is that

$$z = \frac{\sqrt[3]{f} \left(1-\frac{2}{9 \nu_2}\right)+\left(\frac{2}{9 \nu_1}-1 \right)}{\sqrt{\frac{2 f^{2/3}}{9 \nu_2}+\frac{2}{9 \nu_1}}}$$

has approximately a standard Normal distribution. This provides relatively fast calculations of the CDF of either the $B$ or the $F$ distribution and is not difficult to invert.

This approximation works fairly well for $\nu_2 \ge 1000$ and $\nu_1=1$ and progressively better as $\nu_1$ increases. Here is a set of plots of the difference between the approximate CDF and the correct CDF for $\nu_2=1000$ and $\nu_1 = 1,2,3,6,10$, shown with identical scales on their axes:

CDF error plots

The error is extremely low in the right tail, where interest often focuses, but even in the left tail it rarely exceeds $0.002$ once $\nu_1 \ge 2$. As $\nu_2$ increases, the errors remain about the same.

Given the level of accuracy of this approximation, it makes little sense to compute the standard Normal CDF (or its inverse) to double-precision accuracy. This will allow for even faster computation.

One approach is to create a table of values and interpolate within it. Because the Normal CDF is so curved near the tails, it would seem a largish table might be needed for high accuracy. However, by considering Mills Ratio one might be led instead to create a (very) small table of pairs $(z, \sqrt{-2\log(\Phi(z))}$ for $z$ in the lower tail, say for $-8 \le z \le 0$, where $\Phi$ is the standard Normal CDF. Cubic spline interpolation does a great job even when these pairs are evaluated only for integral $z$ (nine values total in the entire table). To further simplify the programming (if a cubic spline interpolation routine is not readily available) one could implement a quadratic spline. Even better--and barely more difficult--is to form the weighted average of a linear interpolator (which tends to underestimate the CDF) and a quadratic interpolator (which overestimates the CDF). A weight of $0.238$ applied to the linear interpolator works well:

Comparison with normal CDF

This plot shows the differences between three interplators and the standard Normal CDF: the red one is for the linear interpolator, the gold for the quadratic interpolator, and the blue (in between) for their weighted average.

(The same technique of averaging in a linear interpolator also decreases the error by an order of magnitude when using a cubic spline.)

The weighted average is another quadratic interpolator but it is no longer a spline (it fails to have continuous first derivatives at the negative integers)--but so what? In effect, it computes the CDF with absolute accuracy better than $0.00005$ for $z\le -1$ and, by symmetry, for $z\ge 1$. Quadratic or cubic interpolation between a small number of pairs of $(z, \Phi(z))$ should work fine for $-1 \lt z \lt 1$ if such calculations are needed. (Typically only very low accuracy, say to about $0.01$, is needed in this range anyway.) The upshot is that the CDF (or, using similar methods, its inverse) can be computed to sufficient accuracy using about a dozen simple arithmetic operations and one exponentiation, at a cost of storing around a dozen precomputed values: even on an interpreted platform like Java that should go very quickly, taking a small fraction of a microsecond.

Examples

Mathematica code to produce the first figure follows.

fDist[f_, {ν1_, ν2_}] := CDF[NormalDistribution[]][((1-2/(9 ν2)) f^(1/3) - (1-2/(9 ν1))) /
    Sqrt[2 f^(2/3)/(9 ν2) + 2/(9 ν1)]];
f[k_, {ν1_, ν2_}] := (k ν2)/((1 - k) ν1);
With[{ν1 = 1, ν2 = 1000},
 Plot[fDist[f[k,{ν1,ν2}], {2 ν1,2 ν2}] - CDF[FRatioDistribution[2 ν1,2 ν2]][f[k,{ν1,ν2}]], 
    {k, 0, 10 ν1/(ν1 + ν2)}, PlotRange -> Full]
 ]

Code to perform the interpolation for $|z| \ge 1$ and plot its error:

d = Table[{z, Sqrt[-2 Log[CDF[NormalDistribution[], z]]]}, {z, -8, -0, 1}];
f1 = Interpolation[d, InterpolationOrder -> 1];
f2 = Interpolation[d, InterpolationOrder -> 2];
f[z_] := 0.238 f1[z] + (1 - 0.238) f2[z];
Plot[{Exp[-f[z]^2/2 ] - CDF[NormalDistribution[], z]}, {z, -4, -1}, PlotRange -> Full]

References

Norman L. Johnson and Samuel Kotz, Continuous Univariate Distributions - 2. J. Wiley & Sons, New York, 1970.

Solved – Conjugate priors for Gamma distribution of unknown $\alpha$ and $\beta$

If you take a Gamma Ga$(\alpha,\beta)$ likelihood $$\dfrac{\alpha^{n\beta}}{\Gamma(\beta)^n}\exp\left\{\beta\sum_{i=1}^n\log x_i -\alpha \sum_{i=1}^n x_i\right\}$$ the distribution with density $$\pi(\alpha,\beta) \propto \dfrac{\alpha^{\lambda\beta}}{\Gamma(\beta)^\lambda}\exp\left\{\beta\xi -\alpha \mu\right\}$$ is conjugate, if non-standard.

Furthermore, the conditional of $\alpha$ given $\beta$ is a Gamma, hence can be integrated out. A posteriori and a priori. The symmetric property for $\beta$ given $\alpha$ does not hold though.