FMStable (version 0.1-4)

Gstable: General Stable Distributions

Description

A procedure based on the R function integrate for computing the distribution function for stable distributions which may be skew but have standard location and scale parameters. This computation is not fast.

It is not designed to work for alpha near to 1.

Usage

tailsGstable(x, logabsx, alpha, oneminusalpha, twominusalpha, 
  beta, betaplus1, betaminus1, parametrization, lower.tail=TRUE)

Arguments

x

Value (scalar).

logabsx

Logarithm of absolute value of x. Must be used when x is outside the range over which numbers can be stored. (e.g. 1.e-5000)

alpha

Value of parameter of stable distribution.

oneminusalpha

Value of alpha. This should be specified when alpha is near to 1 so that the difference from 1 is specified accurately.

twominusalpha

Value of 2 - alpha. This should be specified when alpha is near to 2 so that the difference from 2 is specified accurately.

beta

Value of parameter of stable distribution.

betaplus1

Value of beta + 1. This should be specified when beta is near to -1 so that the difference from -1 is specified accurately.

betaminus1

Value of beta - 1. This should be specified when beta is near to 1 so that the difference from 1 is specified accurately.

parametrization

Parametrization: 0 for Zolotarev's M = Nolan S0, 1 for Zolotarev's A = Nolan S1 and 2 for Zolotarev's C = Chambers, Mallows and Stuck.

lower.tail

Logical: Whether the lower tail of the distribution is of primary interest. This parameter affects whether numerical integration is used for the lower or upper tail. The other tail is computed by subtraction.

Value

Returns a list with the following components:

left.tail.prob

The probability distribution function. I.e. the probability of being less than or equal to x.

right.tail.prob

The probability of being larger than x.

est.error

An estimate of the computational error in the previous two numbers.

message

A message produced by R's standard integrate routine.

References

Chambers, J.M., Mallows, C.L. and Stuck, B.W. (1976) A method for simulating stable random variables. Journal of the American Statistical Association 71, 340--344.

Examples

# NOT RUN {
# Check relationship between maximally skew and other stable distributions
# in paper by J.M. Chambers, C.L. Mallows and B.W. Stuck
alpha <- 1.9
beta <- -.5
k <- 1- abs(1-alpha)
denom <- sin(pi*k)
p <- (sin(.5*pi*k * (1+beta))/denom)^(1/alpha)
q <- (sin(.5*pi*k * (1-beta))/denom)^(1/alpha)
# Probability that p S1 - q S2 < x
S1 <- setParam(alpha=1.9, location=0, logscale =log(p), pm="C")
S2 <- setParam(alpha=1.9, location=0, logscale =log(q), pm="C")
S3 <- setParam(alpha=1.9, location=0, logscale =0, pm="C")
xgiven <- 1
f <- function(x) dEstable(x, S1) * pEstable(xgiven + x, S2)
print(integrate(f, lower=-Inf, upper=Inf, rel.tol=1.e-12)$value, digits=16)
f <- function(x) dEstable(x, S3) * pEstable((xgiven + p*x)/q, S3)
print(integrate(f, lower=-Inf, upper=Inf, rel.tol=1.e-8)$value, digits=16)
direct <- tailsGstable(x=xgiven, logabsx=log(xgiven),alpha=alpha,
  beta=beta, parametrization=2)
print(direct$left.tail.prob, digits=16)

# Compare Estable and Gstable
# List fractional discrepancies
disc <- function(tol){
  for(pm in pms) for (a in alphas) for(x in xs) {
    lx <- log(abs(x))
    beta <- if(pm==2 && a > 1) -1 else 1
    if(x > 0 || a > 1){
      a1 <- pEstable(x, setParam(alpha=a, location=0, logscale=0, pm=pm))
      a2 <- tailsGstable(x=x, logabsx=lx, alpha=a, beta=beta,
        parametrization=pm)$left.tail.prob
      print(paste("parametrization=", pm, "alpha=", a,"x=", x,
        "Frac disc=", a1/a2-1), quote=FALSE)
    }
  }
}

alphas <- c(.3, .8, 1.1, 1.5, 1.9)
pms <- 0:2
xs <- c(-2, .01, 4.3)
disc()
# }