cf2DistGP: Evaluating CDF/PDF/QF from CF of a continous distribution F...

Description Usage Arguments Value Note See Also Examples

Description

cf2DistGP(cf, x, prob, options) evaluates the CDF/PDF/QF (Cumulative Distribution Function, Probability Density Function, Quantile Function) from the Characteristic Function CF by using the Gil-Pelaez inversion formulae:
cdf(x) = 1/2 + (1/pi) * Integral_0^inf imag(exp(-1i*t*x)*cf(t)/t)*dt.
pdf(x) = (1/pi) * Integral_0^inf real(exp(-1i*t*x)*cf(t))*dt.

Usage

1
2
cf2DistGP(cf, x, prob, option, isCompound, isCircular, N, SixSigmaRule, xMin,
  xMax, isPlot)

Arguments

cf

function handle for the characteristic function CF

x

vector of values from domain of the distribution F, if missing, cf2DistFFT automatically selects vector x from the domain

prob

vector of values from [0,1] for which the quantiles will be estimated

option

list with the following (default) parameters:

  • option$isCompound = FALSE treat the compound distributions

  • option$isCircular = FALSE treat the circular distributions

  • option$isInterp = FALSE create and use the interpolant functions for PDF/CDF/RND

  • option$N = 2^10 set N points used by FFT

  • option$xMin = -Inf set the lower limit of X

  • option$xMax = Inf set the upper limit of X

  • option$xMean = NULL set the MEAN value of X

  • option$xStd = NULL set the STD value of X

  • option$dt = NULL set grid step dt = 2*pi/xRange

  • option$T = NULL set upper limit of (0,T), T = N*dt

  • option$SixSigmaRule = 6 set the rule for computing domain

  • option$tolDiff = 1e-4 tol for numerical differentiation

  • option$isPlot = TRUE plot the graphs of PDF/CDF

  • option$DIST list with information for future evaluations. option$DIST is created automatically after first call:

    • option$DIST$xMean

    • option$DIST$xMin

    • option$DIST$xMax

    • option$DIST$xRange

    • option$DIST$cft

    • option$DIST$N

    • option$DIST$k

    • option$DIST$t

    • option$DIST$dt

isCompound

treat the compound distributions, default FALSE

isCircular

treat the circular distributions, default FALSE

N

N points used by GP, default 2^10 (2^14 for compound distribution)

SixSigmaRule

set the rule for computing domain, default 6

xMin

set the lower limit of X

xMax

set the upper limit of X

isPlot

plot the graphs of PDF/CDF, default TRUE

Value

result - list with the following results values:

result$x

<- x

result$cdf

<- cdf (Cumulative Distribution Function)

result$pdf

<- pdf (Probability Density Function)

result$prob

<- prob

result$qf

<- qf (Quantile Function)

result$SixSigmaRule

<- option$SixSigmaRule

result$N

<- option$N

result$dt

<- dt

result$T

<- t[(length(t))]

result$PrecisionCrit

<- PrecisionCrit

result$myPrecisionCrit

<- options$crit

result$isPrecisionOK

<- isPrecisionOK

result$xMean

<- xMean

result$xRange

<- xRange

result$xSTD

<- xStd

result$xMin

<- xMin

result$xMAx

<- xMax

result$cf

<- cf

result$const

<- const

result$isCompound

<- option$isCompound

result$details$count

<- count

result$details$correction

<- correction

result$option

<- option

result$tictoc

<- toc

Note

If options.DIST is provided, then cf2DistGP evaluates CDF/PDF based on this information only (it is useful, e.g., for subsequent evaluation of the quantiles). options.DIST is created automatically after first call. This is supposed to be much faster, bacause there is no need for further evaluations of the characteristic function. In fact, in such case the function handle of the CF is not required, i.e. in such case set cf = []

The required integrals are evaluated approximately by using the simple trapezoidal rule on the interval(0,T), where T = N * dt is a sufficienly large integration upper limit in the frequency domain.

If the optimum values of N and T are unknown, we suggest, as a simple rule of thumb, to start with the application of the six-sigma-rule for determining the value of dt = (2*pi)/(xMax-xMin), where xMax = xMean + 6*xStd, and xMin = xMean - 6*xStd, see [1].

Please note that THIS (TRAPEZOIDAL) METHOD IS AN APPROXIMATE METHOD: Frequently, with relatively low numerical precision of the results of the calculated PDF/CDF/QF, but frequently more efficient and more precise than comparable Monte Carlo methods.

However, the numerical error (truncation error and/or the integration error) could be and should be properly controled!

CONTROLING THE PRECISION: Simple criterion for controling numerical precision is as follows: Set N and T = N*dt such that the value of the integrand function Imag(e^(-1i*t*x) * cf(t)/t) is sufficiently small for all t > T, i.e. PrecisionCrit = abs(cf(t)/t) <= tol, for pre-selected small tolerance value, say tol = 10^-8. If this criterion is not valid, the numerical precission of the result is violated, and the method should be improved (e.g. by selecting larger N or considering other more sofisticated algorithm - not considered here).

See Also

For more details see: https://arxiv.org/pdf/1701.08299.pdf

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
## EXAMPLE1 (Calculate CDF/PDF of N(0,1) by inverting its CF)
cf <- function(t)
  exp(-t ^ 2 / 2)
result <- cf2DistGP(cf)

## EXAMPLE2 (PDF/CDF of the compound Binomial-Exponential distribution)
n <- 25
p <- 0.3
lambda <- 5
cfX <- function(t)
  cfX_Exponential(t, lambda)
cf <- function(t)
  cfN_Binomial(t, n, p, cfX)
x <- seq(0, 5, length.out = 101)
prob <- c(0.9, 0.95, 0.99)
result <- cf2DistGP(cf, x, prob, isCompound = TRUE)

## EXAMPLE3 (PDF/CDF of the compound Poisson-Exponential distribution)
lambda1 <- 10
lambda2 <- 5
cfX <- function(t)
  cfX_Exponential(t, lambda2)
cf <- function(t)
  cfN_Poisson(t, lambda1, cfX)
x <- seq(0, 8, length.out = 101)
prob <- c(0.9, 0.95, 0.99)
result <- cf2DistGP(cf, x, prob, isCompound = TRUE)

Simkova/CharFun documentation built on May 9, 2019, 1:30 p.m.