Density, distribution function, quantile function and random generation for the Pearson type IV distribution.
dpearsonIV(x, m, nu, location, scale, params, log = FALSE)ppearsonIV(q, m, nu, location, scale, params, lower.tail = TRUE,
log.p = FALSE, tol = 1e-08, ...)
qpearsonIV(p, m, nu, location, scale, params, lower.tail = TRUE,
log.p = FALSE, tol = 1e-08, ...)
rpearsonIV(n, m, nu, location, scale, params)
vector of quantiles.
vector of probabilities.
number of observations.
first shape parameter of Pearson type IV distribution.
second shape parameter (skewness) of Pearson type IV distribution.
location parameter of Pearson type IV distribution.
scale parameter of Pearson type IV distribution.
vector/list of length 4 containing parameters m, nu,
location, scale for Pearson type IV distribution
(in this order!).
logical; if TRUE, probabilities p are given as log(p).
logical; if TRUE, probabilities are \(P[X\le x]\),
otherwise, \(P[X>x]\).
relative tolerance for evaluation of hypergeometric function 2F1
(ppearsonIV) or absolute target q-error for Newton method
(qpearsonIV).
further parameters for underlying hypergeometric function.
dpearsonIV gives the density, ppearsonIV gives the distribution
function, qpearsonIV gives the quantile function, and rpearsonIV
generates random deviates.
If at all possible, package gsl should be installed.
Otherwise, calculations for distributions with (more or less) extreme
parameters may slow down by factors of more than 1000 and/or may become
unstable.
The Pearson type IV distribution with location parameter
location\(=\lambda\),
scale parameter scale\(=a\), and shape parameters
\(m\) and \(\nu\) can
be obtained by its probability density function
$$p(x) = \frac{\left|\frac{\Gamma(m+\frac{\nu}{2}i)}{\Gamma(m)}\right|^2}
{a B(m-\frac{1}{2},\frac{1}{2})}
\left[1+\left(\frac{x-\lambda}{a}\right)^2\right]^{-m}
e^{-\nu \mathrm{arctan}\left(\frac{x-\lambda}{a} \right)}$$
for \(a>0\), \(m>\frac{1}{2}\), \(\nu\ne 0\)
(\(\nu=0\) corresponds to the Pearson type VII distribution family).
The normalizing constant, which involves the complex Gamma function, is
calculated with lngamma_complex (see Gamma) of
package gsl, if package gsl is installed.
Section 5.1 of [2] contains an algorithm (C code) for the calculation of the
normalizing constant, which is used otherwise, but this will be very slow for
large absolute values of \(\nu\).
The generation of random numbers (rpearsonIV) uses the C code from
section 7 of [2]. It is (thus) restricted to distributions with \(m>1\).
For the cumulative distribution function (ppearsonIV), numerical
integration of the density function is used, if package gsl is not
available.
If package gsl is installed, three different
methods are used, depending on the parameter constellation (the corresponding
parameter regions were obtained by comprehensive benchmarks):
numerical integration of the density function
cdf representation of Heinrich [2]
cdf representation of Willink [4]
The hypergeometric functions involved in the latter two representations are
approximated
via partial sums of the corresponding series (see [1], 15.1.1, p. 556).
Depending on the parameter constellation, transformation 15.3.5 of [1]
(p. 559) is applied for Heinrich's method.
The evaluation of the partial sums is first carried out in (ordinary)
double arithmetic. If cancellation reduces accuracy beyond tol, the
evaluation is redone in double-double arithmetics. If cancellation still
reduces accuracy beyond tol, the evaluation is again redone in
quad-double arithmetic. Code for double-double and quad-double arithmetics
is based on [3].
For Willink's representation, the hypergeometric function in the denominator
of \(R\) in equation (10) is evaluated via complex gamma
functions (see [1], 15.1.20, p. 556), which is fast and much more stable.
A warning is issued if the approximation of the hypergeometric function
seems to fail (which should not
happen, since numerical integration should be carried out for critical
parameter constellations).
The quantile function (qpearsonIV) is obtained via Newton's method.
The Newton iteration begins in the (single) mode of the distribution,
which is easily calculated (see [2], section 3). Since the mode of the
distribution is the only inflection point of the cumulative distribution
function, convergence is guaranteed.
The Newton iteration stops when the target q-error is reached
(or after a maximum of 30 iterations).
[1] Abramowitz, M. and I. A. Stegun (1972) Handbook of mathematical functions, National Bureau of Standards, Applied Mathematics Series - 55, Tenth Printing, Washington D.C.
[2] Heinrich, J. (2004) A Guide to the Pearson Type IV Distribution, Univ. Pennsylvania, Philadelphia, Tech. Rep. CDF/Memo/Statistics/Public/6820 http://www-cdf.fnal.gov/physics/statistics/notes/cdf6820_pearson4.pdf
[3] Hida, Y., X. S. Li and D. H. Bailey (2000) Algorithms for quad-double precision floating point arithmetic, Lawrence Berkeley National Laboratory. Paper LBNL-48597.
[4] Willink, R. (2008) A Closed-form Expression for the Pearson Type IV Distribution Function, Aust. N. Z. J. Stat. 50 (2), pp. 199-205
# NOT RUN {
## define Pearson type IV parameter set with m=5.1, nu=3, location=0.5, scale=2
pIVpars <- list(m=5.1, nu=3, location=0.5, scale=2)
## calculate probability density function
dpearsonIV(-2:2,params=pIVpars)
## calculate cumulative distribution function
ppearsonIV(-2:2,params=pIVpars)
## calculate quantile function
qpearsonIV(seq(0.1,0.9,by=0.2),params=pIVpars)
## generate random numbers
rpearsonIV(5,params=pIVpars)
# }
Run the code above in your browser using DataLab