lmoms.cov {lmomco}

Distribution-Free Variance-Covariance Structure of Sample L-moments

Description

Compute the distribution-free, variance-covariance matrix (\widehat{\mathrm{var}}(\lambda)) of the sample L-moments (\hat\lambda_r) or alternatively the sample probability-weighted moments (\hat\beta_k, Elamir and Seheult, 2004, sec. 5). The \widehat{\mathrm{var}}(\lambda) is defined by the matrix product

\widehat{\mathrm{var}}(\lambda) = \mathbf{C}\,\mathbf{\hat\Theta}\,\mathbf{C}^{\mathrm{T}}\mbox{,}

where the r \times r matrix \mathbf{C} for number of moments r represents the coefficients of the linear combinations converting \beta_k to \lambda_r and the rth row in the matrix is defined as

\mathbf{C}[r,]_{k{=}0:(r-1)} = (-1)^{(r-1-k)} {r-1 \choose k} {r-1+k \choose k}\mbox{,}

where the row is padded from the right with zeros for k < r to form the required lower triangular structure. Elamir and Seheult (2004) list the \mathbf{C} matrix for r = 4.

Letting the falling factorial be defined (matching Elamir and Seheult's nomenclature) as

a^{(b)} = \Gamma(b+1) {a \choose b}\mbox{,}

and letting an entry in the \mathbf{\hat\Theta} matrix denoted as \hat\theta_{kl} be defined as

\hat\theta_{kl} = \hat\beta_k\hat\beta_l - \frac{A}{n^{(k+l+2)}}\mbox{,}

where \hat\beta_k are again the sample probability-weighted moments and are computed by pwm, and finally A is defined as

A = \sum_{i=1}^{n-1}\sum_{j=i+1}^{n} \bigl[ (i-1)^{(k)} (j-k-2)^{(l)} + (i-1)^{(l)} (i-l-2)^{(k)} \bigr] X_{i:n}X_{j:n}\mbox{,}

where X_{i:n} are the sample order statistics for a sample of size n.

Incidentally, the matrix \mathbf{\hat\Theta} is the variance-covariance structure (\widehat{\mathrm{var}}) of the \hat\beta, thus \widehat{\mathrm{var}}(\beta) = \mathbf{\hat\Theta}, which can be returned by a logical function argument (as.pwm=TRUE) instead of \widehat{\mathrm{var}}(\lambda). The last example in Examples provides a demonstration.

Usage

lmoms.cov(x, nmom=5, as.pwm=FALSE, showC=FALSE,
             se=c("NA", "lamse", "lmrse", "pwmse"), ...)

Arguments

x

A vector of data values.

nmom

The number of moments to compute. Default is 5.

as.pwm

A logical controlling whether the distribution-free, variance-covariance of sample probability-weighted moments (\mathbf{\hat\Theta}) is returned instead.

showC

A logical controlling whether the matrix \mathbf{C} is printed during function operation, and this matrix is not returned as a presumed safety feature.

se

Compute standard errors (SE) for the respective moments. The default of "NA" retains the return of either \widehat{\mathrm{var}}(\beta) or \widehat{\mathrm{var}}(\lambda) depending on setting of as.pwm. The "lamse" returns the square root of the diagonal of \widehat{\mathrm{var}}(\lambda), and notationally these are \lambda_r^{SE}. Similarly, "pwmse" returns the square root of the diagonal of \widehat{\mathrm{var}}(\beta) by internally setting as.pwm to TRUE, and notationally these are \beta_{r-1}^{SE}. (Remember that \beta_0 \equiv \lambda_1—the indexing of the former starts at 0 and at the later at 1). The "lmrse" returns the square root of the first two terms of the \widehat{\mathrm{var}}(\lambda) diagonal (\lambda_{1,2}^{SE}) but computes SE for the L-moment ratios (\tau_r^{SE}) for r \ge 3 using the Taylor-series-based approximation (see Note) shown by Elamir and Seheult (2004, p. 348). (Remember that L-moment ratios are \tau_r = \lambda_r/\lambda_2 for r \ge 3 and that \tau_2 = \lambda_2/\lambda_1 [coefficient of L-variation].)

...

Other arguments to pass should they be needed (none were at first implementation).

Value

An R matrix is returned. In small samples and substantially sized r, one or more \hat\theta_{kl} will be NaN starting from the lower right corner of the matrix. The function does not test for this nor reduce the number of moments declared in nmom itself. To reiterate, the square roots along the \widehat{\mathrm{var}}(\lambda) diagonal are SE for the respective L-moments.

Note

Function lmoms.cov was developed as a double check on the evidently separately developed r \le 4 (nmom) implementations of \widehat{\mathrm{var}}(\lambda) in packages Lmoments and nsRFA. Also the internal structure closely matches the symbolic mathematics by Elamir and Seheult (2004), but this practice comes at the expense of more than an order of magnitude slower execution times than say either of the functions Lmomcov() (package Lmoments) or varLmoments() (package nsRFA). For a high speed and recommended implementation, please use the Lmoments package by Karvanen (2016)—Karvanen extended this implementation to larger r for the lmomco package.

For se="lmrse", the Taylor-series-based approximation is suggested by Elamir and Seheult (2004, p. 348) to estimate the variance of an L-moment ratio (\tau_r for r \ge 3) is based on structure of the variance of the ratio of two uniform variables in which the numerator is the rth L-moment and the denominator is \lambda_2:

\mathrm{var}(\tau_r) \cong \biggl[ \frac{\mathrm{var}(\lambda_r)}{\mathrm{E}(\lambda_r)^2} + \frac{\mathrm{var}(\lambda_2)}{\mathrm{E}(\lambda_2)^2} - \frac{2\mathrm{cov}(\lambda_r,\lambda_2)}{\mathrm{E}(\lambda_r)\mathrm{E}(\lambda_2)} \biggr] \biggl[\frac{\mathrm{E}(\lambda_r)}{\mathrm{E}(\lambda_2)} \biggr]^2\mbox{,}

where \mathrm{var}(\cdots) are the along the diagonal of \widehat{\mathrm{var}}(\lambda) and \mathrm{cov}(\cdots) are the off-diagonal covariances. The expectations \mathrm{E}(\cdots) are replaced with the sample estimates. Only for se="lmrse" the SE of the coefficient of L-variation (\tau_2^{SE}) is computed but retained as an attribute (attr() function) of the returned vector and not housed within the vector—the \lambda_2^{SE} continues to be held in the 2nd position of the returned vector.

Author(s)

W.H. Asquith

References

Elamir, E.A.H., and Seheult, A.H., 2004, Exact variance structure of sample L-moments: Journal of Statistical Planning and Inference, v. 124, pp. 337–359.

Karvanen, Juha, 2016, Lmoments—L-moments and quantile mixtures: R package version 1.2-3, accessed February 22, 2016 at https://cran.r-project.org/web/packages/Lmoments/index.html

See Also

lmoms, pwm

Examples

## Not run: 
nsim <- 1000; n <- 10 # Let us compute variance of lambda_3
VL3sample <- mean(replicate(nsim, { zz <- lmoms.cov(rexp(n),nmom=3); zz[3,3] }))
falling.factorial <- function(a, b) gamma(b+1)*choose(a,b)
VL3exact  <- ((4*n^2 - 3*n - 2)/30)/falling.factorial (10, 3) # Exact variance is from
print(c(VL3sample, VL3exact)) # Elamir and Seheult (2004, table 1, line 8)
#[1] 0.01755058 0.01703704  # the values obviously are consistent
## End(Not run)
## Not run: 
# Data considered by Elamir and Seheult (2004, p. 348)
library(MASS); data(michelson); Light <- michelson$Speed
lmoms(Light, nmom=4)$lambdas # 852.4, 44.3, 0.83, 6.5 # matches those authors
lmoms.cov(Light) # [1, ] ==> 62.4267, 0.7116, 2.5912, -3.9847 # again matches
# The authors report standard error of L-kurtosis as 0.03695, which matches
lmoms.cov(Light, se="lmrse")[4] # 0.03695004 
## End(Not run)
## Not run: 
D <- rnorm(100) # Check results of Lmoments package.
lmoms.cov(D, rmax=5)[,5]
#        lam1         lam2         lam3         lam4         lam5
#3.662721e-04 3.118812e-05 5.769509e-05 6.574662e-05 1.603578e-04
Lmoments::Lmomcov(D, rmax=5)[,5]
#          L1           L2           L3           L4           L5
#3.662721e-04 3.118812e-05 5.769509e-05 6.574662e-05 1.603578e-04
## End(Not run)