qua2ci: Estimate Confidence Interval for a Single Quantile of a Parent Distribution

Description

This function estimates the lower and upper limits of a specified confidence interval for an aribitrary quantile value of a specified parent distribution [quantile function $Q(F,\theta)$ with parameters $\theta$] using Monte Carlo simulation. The quantile value, actually the nonexceedance probability ($F$ for $0 \le F < 1$) of the value, is specified by the user. The user also provides the parameters of the parent distribution (see lmom2par).

For nsim simulation runs (ideally a very large number), samples of size $n$ are drawn from $Q(F,\theta)$ using the R-function runif(n) to generate the simulated $F$ values and the lmomco function par2qua function to compute the quantile. The L-moments of the simulated sample are computed using lmoms and a distribution of the same type as the parent is fit to these simulated L-moments. The $F$-quantile of the just-fitted distribution is computed and placed into a vector. The process of simulating the sample, computing the L-moments, computing the parameters, and solving for the $F$-quantile is repeated for the specified number of simulation runs.

To estimate the confidence interval, the L-moments of the vector simulated quantiles are computed. Subsequently, the parameters of a user-specified distribution error distribution are computed. The two quantiles of this error distribution for the specified confidence interval are computed. These two quantiles represent the estimated lower and upper limits for the confidence interval of the parent distribution for samples of size $n$.

Finally, it is often useful to have vectors of lower and upper limits for confidence intervals for a vector of $F$ values. The function genci does just that and uses the qua2ci function as the computational engine.

Usage

qua2ci(f, para, n, ci=0.90, edist='nor', nsim=1000, 
verbose=FALSE, showpar=FALSE, maxlogdiff=6)

Arguments

Nonexceedance probability ($0 \le F \le 1$) of the quantile for which the confidence interval is needed.

para

The parameters from lmom2par or similar.

The sample size that the Monte Carlo simulation will use.

The confidence interval ($0.5 \le$ ci $< 1$). The interval is specified as the size of the interval. The default is 0.90 or the 90th percentile. The function will return the 5th (1-0.90)/2 and 95th (1-(1-0.90)/2) percentile cumulative probabi

edist

The model for the error distribution. Although the normal (the default) is typically assumed in error analyses, it need not be, as support for other distributions supported by the lmomco package is available. However, one should seriously consi

nsim

The number of simulations for the sample size n to perform. Large numbers produce more refined confidence limit estimates at the cost of CPU time. The default is anticipated to be large enough for evaluative-useage without too much computatio

verbose

The verbosity of the operation of the function.

showpar

The parameters of the edist for each simulation are printed.

maxlogdiff

The maximum permitted difference in log10 space between a simulated quantile and the true value. It is possible that a well fit simulated sample to the parent distribution type provides crazy quantile estimates in the far reaches of either tail. The defau

Value

An R list is returned.
lowerThe lower value of the confidence interval having nonexceedance probability equal to (1-ci)/2.
trueThe value returned by par2qua(f,para).
upperThe upper value of the confidence interval having nonexceedance probability equal to 1-(1-ci)/2.
elmomsThe L-moments of the distribution of simulated of quantiles.
eparaThe parameters of the error distribution fit using the elmoms.
ifailA diagnostic value. A value of zero means that successful exit was made.
ifailtextA descriptive message related to the ifail value.

Examples

Run this code

# It is well known that standard deviation (sigma) of the
  # sample mean is equal to sigma/sample_size.
  MEAN  <- 0
  SIGMA <- 100
  PAR   <- vec2par(c(MEAN,SIGMA),type='nor')
  CI    <- qua2ci(0.5,PAR,n=10,nsim=20) # F=0.5-->median=mean
  # Theoretrical sample mean sigma = 100/10 = 10
  # L-moment theory: L-scale*sqrt(pi) = sigma
  # Thus, it follows that
  CI$elmoms$lambdas[2]/sqrt(pi)
  # approaches 10 as nsim --> Inf.


  # Another example.
  D   <- c(123,34,4,654,37,78, 93, 95, 120) # fake sample
  lmr <- lmoms(D)    # compute the lmoments of the data
  WEI <- parwei(lmr) # estimate parameters of the weibull
  CI  <- qua2ci(0.75,WEI,20,nsim=20,ci=0.95)
  # CI contains the estimate 95-percent 
  # confidence interval for the 75th-percent of the parent
  # weibull distribution for size 20 samples from the parent

Run the code above in your browser using DataLab

Get 50% off unlimited learning