Distribution and Quantile Functions for Unimodal Distributions

Share:

Description

Given the density function of a unimodal distribution specified by the root of the density function name, returns the distribution function and quantile function of the specified distribution.

Usage

1
2
3
4
5
6
7
pDist(densFn = "norm", q, param = NULL, subdivisions = 100,
      lower.tail = TRUE, intTol = .Machine$double.eps^0.25,
      valueOnly = TRUE, ...)
qDist(densFn = "norm", p, param = NULL,
      lower.tail = TRUE, method = "spline", nInterpol = 501,
      uniTol = .Machine$double.eps^0.25,
      subdivisions = 100, intTol = uniTol, ...)

Arguments

densFn

Character. The name of the density function for which the distribution function or quantile function is required.

q

Vector of quantiles.

p

Vector of probabilities.

param

Numeric. A vector giving the parameter values for the distribution specified by densFn. If no param values are specified, then the default parameter values of each distribution are used instead.

method

Character. If "spline" quantiles are found from a spline approximation to the distribution function. If "integrate", the distribution function used is always obtained by integration.

lower.tail

Logical. If lower.tail = TRUE, the cumulative density is taken from the lower tail.

subdivisions

The maximum number of subdivisions used to integrate the density and determine the accuracy of the distribution function calculation.

intTol

Value of rel.tol and hence abs.tol in calls to integrate. See integrate.

valueOnly

Logical. If valueOnly = TRUE calls to pDist only return the value obtained for the integral. If valueOnly = FALSE an estimate of the accuracy of the numerical integration is also returned.

nInterpol

Number of points used in qDist for cubic spline interpolation of the distribution function.

uniTol

Value of tol in calls to uniroot. See uniroot.

...

Passes additional arguments to integrate, distMode or distCalcRange. In particular, the parameters of the distribution.

Details

The name of the unimodal density function must be supplied as the characters of the root for that density (e.g. norm, ghyp).

pDist uses the function integrate to numerically integrate the density function specified. The integration is from -Inf to x if x is to the left of the mode, and from x to Inf if x is to the right of the mode. The probability calculated this way is subtracted from 1 if required. Integration in this manner appears to make calculation of the quantile function more stable in extreme cases.

qDist provides two methods to calculate quantiles both of which use uniroot to find the value of x for which a given q is equal to F(x) where F(.) denotes the distribution function. The difference is in how the numerical approximation to F is obtained. The more accurate method, which is specified as "integrate", is to calculate the value of F(x) whenever it is required using a call to pDist. It is clear that the time required for this approach is roughly linear in the number of quantiles being calculated. The alternative (and default) method is that for the major part of the distribution a spline approximation to F(x) is calculated and quantiles found using uniroot with this approximation. For extreme values of some heavy-tailed distributions (where the tail probability is less than 10^(-7)), the integration method is still used even when the method specified as "spline".

If accurate probabilities or quantiles are required, tolerances (intTol and uniTol) should be set to small values, i.e 10^(-10) or 10^(-12) with method = "integrate". Generally then accuracy might be expected to be at least 10^(-9). If the default values of the functions are used, accuracy can only be expected to be around 10^(-4). Note that on 32-bit systems .Machine$double.eps^0.25 = 0.0001220703 is a typical value.

Value

pDist gives the distribution function, qDist gives the quantile function.

An estimate of the accuracy of the approximation to the distribution function can be found by setting valueOnly = FALSE in the call to pDist which returns a list with components value and error.

Author(s)

David Scott d.scott@auckland.ac.nz Joyce Li xli053@aucklanduni.ac.nz

Examples

1
2
3
4
5
6
7
pDist("norm", q = 2, mean = 1, sd = 1)
pDist("t", q = 0.5, df = 4)
require(GeneralizedHyperbolic)
pDist("ghyp", q = 0.1)
require(SkewHyperbolic)
qDist("skewhyp", p = 0.4, param = c(0, 1, 0, 10))
qDist("t", p = 0.2, df = 4)