dlomax: The Lomax distribution.

Description

Lomax density, distribution and quantile functions, random number generation, and expectation and variance.

Usage

dlomax(x, shape=1, scale=1, log=FALSE)
  plomax(q, shape=1, scale=1)
  qlomax(p, shape=1, scale=1)
  rlomax(n, shape=1, scale=1)
  elomax(shape=1, scale=1)
  vlomax(shape=1, scale=1)

Value

‘dlomax()’ gives the density function,

‘plomax()’ gives the cumulative distribution function (CDF),

‘qlomax()’ gives the quantile function (inverse CDF), and ‘rlomax()’ generates random deviates. The ‘elomax()’ and ‘vlomax()’

functions return the corresponding Lomax distribution's expectation and variance, respectively.

Arguments

x,q: quantile.
p: probability.
n: number of observations.
shape: shape parameter ($\alpha > 0$).
scale: scale parameter ($\lambda > 0$).
log: logical; if TRUE, logarithmic density will be returned.

Author

Christian Roever christian.roever@med.uni-goettingen.de

Details

The Lomax distribution is a heavy-tailed distribution that also is a special case of a Pareto distribution of the 2nd kind. The probability density function of a Lomax distributed variable with shape $\alpha>0$ and scale $\lambda>0$ is given by $$p(x) = (\alpha / \lambda) (1 + x / \lambda)^{-(\alpha+1)}.$$ The density function is monotonically decreasing in $x$. Its mean is $\lambda / (\alpha-1)$ (for $\alpha>1$) and its median is $\alpha(2^{1/\alpha}-1)$. Its variance is finite only for $\alpha > 2$ and equals $(\lambda^2 \alpha) / ((\alpha-1)^2 (\alpha-2))$. The cumulative distribution function (CDF) is given by $$P(x) = 1-(1+ x / \lambda)^{-\alpha}.$$

The Lomax distribution also arises as a gamma-exponential mixture. Suppose that $X$ is a draw from an exponential distribution whose rate $\theta$ again is drawn from a gamma distribution with shape $a$ and scale $s$ (so that $\mathrm{E}[\theta]=as$ and $\mathrm{Var}(\theta)=as^2$, or $\mathrm{E}[1/\theta]=\frac{1}{s(a+1)}$ and $\mathrm{Var}(1/\theta)=\frac{1}{s^2(a-1)^2(a-2)}$). Then the marginal distribution of $X$ is Lomax with scale $1/s$ and shape $a$. Consequently, if the moments of $\theta$ are given by $\mathrm{E}[\theta]=\mu$ and $\mathrm{Var}(\theta)=\sigma^2$, then $X$ is Lomax distributed with shape $\alpha=\left(\frac{\mu}{\sigma}\right)^2$ and scale $\lambda=\frac{\mu}{\sigma^2}=\frac{\alpha}{\mu}$. The gamma-exponential connection is also illustrated in an example below.

References

C. Roever, R. Bender, S. Dias, C.H. Schmid, H. Schmidli, S. Sturtz, S. Weber, T. Friede. On weakly informative prior distributions for the heterogeneity parameter in Bayesian random-effects meta-analysis. Research Synthesis Methods, 12(4):448-474, 2021. tools:::Rd_expr_doi("10.1002/jrsm.1475").

N.L. Johnson, S. Kotz, N. Balakrishnan. Continuous univariate distributions, volume 1. Wiley, New York, 2nd edition, 1994.

Examples

Run this code

#######################
# illustrate densities:
x <- seq(0,6,le=200)
plot(x, dexp(x, rate=1), type="l", col="cyan", ylim=c(0,1),
     xlab=expression(tau), ylab=expression("probability density "*f(tau)))
lines(x, dlomax(x), col="orange")
abline(h=0, v=0, col="grey")

# show log-densities (note the differing tail behaviour):
plot(x, dexp(x, rate=1), type="l", col="cyan", ylim=c(0.001,1), log="y",
     xlab=expression(tau), ylab=expression("probability density "*f(tau)))
lines(x, dlomax(x), col="orange")
abline(v=0, col="grey")

######################################################
# illustrate the gamma-exponential mixture connection;
# specify a number of samples:
N <- 10000
# specify some gamma shape and scale parameters
# (via mixing distribution's moments):
expectation <- 2.0
stdev       <- 1.0
gammashape <- (expectation / stdev)^2
gammascale <- stdev^2 / expectation
print(c("expectation"=expectation, "stdev"=stdev,
        "shape"=gammashape, "scale"=gammascale))
# generate gamma-distributed rates:
lambda <- rgamma(N, shape=gammashape, scale=gammascale)
# generate exponential draws according to gamma-rates:
y <- rexp(N, rate=lambda)
# determine Lomax quantiles accordingly parameterized:
x <- qlomax(ppoints(N), scale=1/gammascale, shape=gammashape)
# compare distributions in a Q-Q-plot:
plot(x, sort(y), log="xy", main="quantile-quantile plot",
     xlab="theoretical quantile", ylab="empirical quantile")
abline(0, 1, col="red")

Run the code above in your browser using DataLab