42  Beta Prime Distribution

The Beta Prime distribution — also called the Inverted Beta or Beta Distribution of the Second Kind — extends the Beta distribution to the positive half-line. While the Beta models a proportion bounded to \([0,1]\), the Beta Prime models a positive ratio with no upper bound, arising naturally in Bayesian hierarchical models and as a distribution for odds ratios.

Formally, the random variate \(X\) defined for the range \(X > 0\), is said to have a Beta Prime Distribution (i.e. \(X \sim \text{BetaPrime}(\alpha, \beta, \theta)\)) with shape parameters \(\alpha > 0\) and \(\beta > 0\), and scale parameter \(\theta > 0\). If \(Y \sim \text{Beta}(\alpha, \beta)\) then \(\theta Y/(1-Y) \sim \text{BetaPrime}(\alpha, \beta, \theta)\).

42.1 Probability Density Function

\[ f(x) = \frac{(x/\theta)^{\alpha-1}(1+x/\theta)^{-\alpha-\beta}}{\theta\,\text{B}(\alpha,\beta)}, \quad x > 0 \]

where \(\text{B}(\alpha, \beta) = \Gamma(\alpha)\Gamma(\beta)/\Gamma(\alpha+\beta)\).

The figure below shows examples of the Beta Prime Probability Density Function for different parameter combinations with \(\theta = 1\).

Code 
dbetaprime <- function(x, alpha, beta, theta = 1) {
  ifelse(x > 0,
         (x/theta)^(alpha-1) * (1 + x/theta)^(-alpha-beta) /
           (theta * beta(alpha, beta)),
         0)
}

par(mfrow = c(2, 2))
x <- seq(0.001, 6, length = 500)

plot(x, dbetaprime(x, 2, 4, 1), type = "l", lwd = 2, col = "blue",
     xlab = "x", ylab = "f(x)", main = expression(paste(alpha==2, ", ", beta==4, ", ", theta==1)))

plot(x, dbetaprime(x, 4, 4, 1), type = "l", lwd = 2, col = "blue",
     xlab = "x", ylab = "f(x)", main = expression(paste(alpha==4, ", ", beta==4, ", ", theta==1)))

plot(x, dbetaprime(x, 2, 2, 1), type = "l", lwd = 2, col = "blue",
     xlab = "x", ylab = "f(x)", main = expression(paste(alpha==2, ", ", beta==2, ", ", theta==1)))

plot(x, dbetaprime(x, 5, 3, 1), type = "l", lwd = 2, col = "blue",
     xlab = "x", ylab = "f(x)", main = expression(paste(alpha==5, ", ", beta==3, ", ", theta==1)))

par(mfrow = c(1, 1))
Figure 42.1: Beta Prime Probability Density Function for various parameter combinations (scale = 1)

42.2 Purpose

The Beta Prime distribution models positive, right-skewed quantities for which the support extends to infinity but is bounded below by zero. It generalizes the Beta distribution to the full positive half-line via an odds-ratio transformation and arises naturally in Bayesian inference as a marginal distribution. Common applications include:

  • Odds ratios in clinical trials and epidemiology (positive, unbounded, right-skewed)
  • Bayesian hierarchical models: scale parameter priors with heavy right tails
  • Variance-ratio statistics: the F distribution is a scaled Beta Prime
  • Financial loss ratios: claim amount relative to policy limit
  • Reliability and survival: ratio of component strengths or lifetimes

Relation to the discrete setting. The Beta Prime arises in the same way as the Negative Binomial for count data: just as Poisson-Gamma mixture yields the Negative Binomial, a Beta-Prime-type distribution emerges as the marginal of Poisson-Gamma hierarchies for rates. Both generalize their simpler counterparts to handle overdispersion.

42.3 Distribution Function

\[ F(x) = I_{x/(x+\theta)}(\alpha, \beta) \]

where \(I_u(\alpha, \beta)\) is the regularized incomplete beta function. In R: pbeta(x/(x+theta), shape1 = alpha, shape2 = beta).

The figure below shows the Beta Prime Distribution Function for \(\alpha = 2\), \(\beta = 4\), \(\theta = 1\).

Code 
pbetaprime <- function(x, alpha, beta, theta = 1) {
  ifelse(x > 0, pbeta(x / (x + theta), shape1 = alpha, shape2 = beta), 0)
}

x <- seq(0, 6, length = 500)
plot(x, pbetaprime(x, 2, 4, 1), type = "l", lwd = 2, col = "blue",
     xlab = "x", ylab = "F(x)", main = "Beta Prime Distribution Function",
     sub = expression(paste(alpha==2, ", ", beta==4, ", ", theta==1)))
Figure 42.2: Beta Prime Distribution Function (alpha = 2, beta = 4, scale = 1)

42.4 Moment Generating Function

The moment generating function does not exist for \(t > 0\) due to the heavy right tail.

42.5 1st Uncentered Moment

\[ \mu_1' = \frac{\theta\,\alpha}{\beta - 1}, \quad \beta > 1 \]

42.6 2nd Uncentered Moment

\[ \mu_2' = \frac{\theta^2\,\alpha(\alpha+1)}{(\beta-1)(\beta-2)}, \quad \beta > 2 \]

42.7 3rd Uncentered Moment

\[ \mu_3' = \frac{\theta^3\,\alpha(\alpha+1)(\alpha+2)}{(\beta-1)(\beta-2)(\beta-3)}, \quad \beta > 3 \]

42.8 4th Uncentered Moment

\[ \mu_4' = \frac{\theta^4\,\alpha(\alpha+1)(\alpha+2)(\alpha+3)}{(\beta-1)(\beta-2)(\beta-3)(\beta-4)}, \quad \beta > 4 \]

In general: \(\mu_n' = \theta^n\,\text{B}(\alpha+n,\,\beta-n)/\text{B}(\alpha,\beta)\) for \(n < \beta\).

42.9 2nd Centered Moment

\[ \mu_2 = \frac{\theta^2\,\alpha(\alpha+\beta-1)}{(\beta-1)^2(\beta-2)}, \quad \beta > 2 \]

42.10 3rd Centered Moment

Obtained by expanding from raw moments; requires \(\beta > 3\).

42.11 4th Centered Moment

Obtained by expanding from raw moments; requires \(\beta > 4\).

42.12 Expected Value

\[ \text{E}(X) = \frac{\theta\,\alpha}{\beta - 1}, \quad \beta > 1 \]

42.13 Variance

\[ \text{V}(X) = \frac{\theta^2\,\alpha(\alpha+\beta-1)}{(\beta-1)^2(\beta-2)}, \quad \beta > 2 \]

42.14 Median

The median has no closed form and is computed numerically via qbeta.

42.15 Mode

\[ \text{Mo}(X) = \begin{cases} \theta\,\dfrac{\alpha-1}{\beta+1} & \alpha > 1 \\[4pt] 0 & \alpha \leq 1 \end{cases} \]

42.16 Coefficient of Skewness

\[ g_1 = \frac{2(2\alpha+\beta-1)}{\beta-3}\sqrt{\frac{\beta-2}{\alpha(\alpha+\beta-1)}}, \quad \beta > 3 \]

42.17 Coefficient of Kurtosis

\[ g_2 = \frac{3(\beta-2)\bigl(\alpha^2\beta + 5\alpha^2 + \alpha\beta^2 + 4\alpha\beta - 5\alpha + 2\beta^2 - 4\beta + 2\bigr)}{\alpha(\beta-4)(\beta-3)(\alpha+\beta-1)}, \quad \beta > 4 \]

The kurtosis is always greater than 3 for valid parameter values.

42.18 Parameter Estimation

Maximum likelihood estimation is numerical. Method-of-moments provides starting values from the mean and variance:

set.seed(42)
alpha_true <- 2; beta_true <- 4; theta_true <- 1

# Simulate BetaPrime via Beta transformation
y <- rbeta(100, shape1 = alpha_true, shape2 = beta_true)
x_obs <- theta_true * y / (1 - y)

# Summary statistics
cat("Sample mean:", round(mean(x_obs), 4), "\n")
cat("Theoretical mean:", theta_true * alpha_true / (beta_true - 1), "\n")
cat("Sample var:", round(var(x_obs), 4), "\n")
cat("Theoretical var:",
    theta_true^2 * alpha_true * (alpha_true + beta_true - 1) /
    ((beta_true - 1)^2 * (beta_true - 2)), "\n")
Sample mean: 0.6892 
Theoretical mean: 0.6666667 
Sample var: 0.5095 
Theoretical var: 0.5555556

42.19 R Module

42.19.1 RFC

The Beta Prime Distribution module is available in RFC under the menu “Distributions / Beta Prime Distribution”.

42.19.3 R Code

The following code demonstrates Beta Prime probability calculations:

alpha <- 2; beta <- 4; theta <- 1

pbetaprime <- function(x, alpha, beta, theta = 1) {
  ifelse(x > 0, pbeta(x / (x + theta), shape1 = alpha, shape2 = beta), 0)
}

# P(X <= 1)
pbetaprime(1, alpha, beta, theta)

# Mean and mode
cat("Mean:", theta * alpha / (beta - 1), "\n")
cat("Mode:", theta * (alpha - 1) / (beta + 1), "\n")
[1] 0.8125
Mean: 0.6666667 
Mode: 0.2

42.20 Example

An odds ratio for a treatment effect is modeled as \(X \sim \text{BetaPrime}(\alpha = 2, \beta = 4, \theta = 1)\). The mean odds ratio is \(\alpha/(\beta-1) = 2/3\) and the mode is \((\alpha-1)/(\beta+1) = 1/5\).

alpha <- 2; beta <- 4; theta <- 1

pbetaprime <- function(x, alpha, beta, theta = 1) {
  ifelse(x > 0, pbeta(x / (x + theta), shape1 = alpha, shape2 = beta), 0)
}

# P(odds ratio <= 1)
cat("P(odds ratio <= 1):", round(pbetaprime(1, alpha, beta, theta), 4), "\n")

# Mean and mode
cat("Mean odds ratio:", theta * alpha / (beta - 1), "\n")
cat("Mode:", theta * (alpha - 1) / (beta + 1), "\n")
P(odds ratio <= 1): 0.8125 
Mean odds ratio: 0.6666667 
Mode: 0.2
Interactive Shiny app (click to load).
Open in new tab

42.21 Random Number Generator

Beta Prime random variates are generated via the Beta transformation:

\[ \text{If } Y \sim \text{Beta}(\alpha, \beta) \text{ then } X = \frac{\theta Y}{1-Y} \sim \text{BetaPrime}(\alpha, \beta, \theta) \]

set.seed(123)
n <- 1000; alpha <- 2; beta <- 4; theta <- 1

# Generate via Beta transformation
y <- rbeta(n, shape1 = alpha, shape2 = beta)
x_sim <- theta * y / (1 - y)

cat("Simulated mean:", round(mean(x_sim), 4), "\n")
cat("Theoretical mean:", theta * alpha / (beta - 1), "\n")
Simulated mean: 0.6751 
Theoretical mean: 0.6666667
Interactive Shiny app (click to load).
Open in new tab

42.22 Property 1: Beta Transformation

If \(Y \sim \text{Beta}(\alpha, \beta)\) then \(\theta Y/(1-Y) \sim \text{BetaPrime}(\alpha, \beta, \theta)\). This is the defining property and provides the simplest interpretation: the Beta Prime models odds-ratio-type quantities derived from bounded Beta variates (see Chapter 30).

42.23 Property 2: F Distribution Relationship

If \(X \sim F(2\alpha, 2\beta)\) (Fisher-Snedecor F), then \((\alpha/\beta) X \sim \text{BetaPrime}(\alpha, \beta, 1)\). This means the F distribution is a scaled Beta Prime, and tables of the F distribution implicitly cover the Beta Prime (see Chapter 26).

42.24 Property 3: Lomax (Shifted Pareto) Special Case

When \(\alpha = 1\): \(\text{BetaPrime}(1, \beta, \theta)\) is the Lomax distribution (also known as the shifted Pareto or Pareto Type II). See Chapter 32.