# logConCI: Compute pointwise confidence interval for a density assuming... In logcondens: Estimate a Log-Concave Probability Density from Iid Observations

## Description

Compute approximate confidence interval for the true log-concave density, on a grid of points. Two main approaches are implemented: In the first, the confidence interval at a fixed point is based on the pointwise asymptotic theory for the log-concave maximum likelihood estimator (MLE) developed in Balabdaoui, Rufibach, and Wellner (2009). In the second, the confidence interval is estimated via the boostrap.

## Usage

 1 2 3 logConCI(res, xx0, conf.level = c(0.8, 0.9, 0.95, 0.99), type = c("DR", "ks", "nrd", "ECDFboot", "NPMLboot"), htype = c("hscv", "hlscv", "hpi", "hns"), BB = 500) 

## Arguments

 res An object of class dlc, usually a result of a call to logConDens. xx0 Vector of grid points at which to calculate the confidence interval. conf.level Confidence level for the confidence interval(s). The default is 95%. type Vector of strings indicating type of confidence interval to compute. When type = ks is chosen, then htype should also be specified. The default is type = ks. htype Vector of strings indicating bandwidth selection method if type = ks. The default is htype = hns. BB number of iterations in the bootstrap if type = NPMLboot or type = ECDFboot. The default is BB = 500.

## Details

In Balabdaoui et al. (2009) it is shown that (if the true density is strictly log-concave) the limiting distribution of the MLE of a log-concave density \widehat f_n at a point x is

n^{2/5}(\widehat f_n(x)-f(x)) \to c_2(x) \bar{C}(0).

The nuisance parameter c_2(x) depends on the true density f and the second derivative of its logarithm. The limiting process \bar{C}(0) is found as the second derivative at zero of a particular operator (called the "envelope") of an integrated Brownian motion plus t^4.

Three of the confidence intervals are based on inverting the above limit using estimated quantiles of \bar{C}(0), and estimating the nuisance parameter c_2(x). The options for the function logConCI provide different ways to estimate this nuisance parameter. If type = "DR", c_2(x) is estimated using derivatives of the smoothed MLE as calculated by the function logConDens (this method does not perform well in simulations and is therefore not recommended). If type="ks", c_2(x) is estimated using kernel density estimates of the true density and its first and second derivatives. This is done using the R package ks, and, with this option, a bandwidth selection method htype must also be chosen. The choices in htype correspond to the various options for bandwidth selection available in ks. If type = "nrd", the second derivative of the logarithm of the true density in c_2(x) is estimated assuming a normal reference distribution.

Two of the confidence intervals are based on the bootstrap. For type = "ECDFboot" confidence intervals based on re-sampling from the empirical cumulative distribution function are computed. For type = "NPMLboot" confidence intervals based on re-sampling from the nonparametric maximum likelihood estimate of log-concave density are computed. Bootstrap confidence intervals take a few minutes to compute!

The default option is type = "ks" with htype = "hns". Currently available confidence levels are 80%, 90%, 95% and 99%, with a default of 95%.

Azadbakhsh et al. (2014) provides an empirical study of the relative performance of the various approaches available in this function.

## Value

The function returns a list containing the following elements:

 fhat MLE evaluated at grid points. up_DR Upper confidence interval limit when type = DR. lo_DR Lower confidence interval limit when type = DR. up_ks_hscv Upper confidence interval limit when type = ks and htype = hscv. lo_ks_hscv Lower confidence interval limit when type = ks and htype = hscv. up_ks_hlscv Upper confidence interval limit when type = ks and htype = hlscv. lo_ks_hlscv Lower confidence interval limit when type = ks and htype = hlscv. up_ks_hpi Upper confidence interval limit when type = ks and htype = hpi. lo_ks_hpi Lower confidence interval limit when type = ks and htype = hpi. up_ks_hns Upper confidence interval limit when type = ks and htype = hns. lo_ks_hns Lower confidence interval limit when type = ks and htype = hns. up_nrd Upper confidence interval limit when type = nrd. lo_nrd Lower confidence interval limit when type = nrd. up_npml Upper confidence interval limit when type = NPMLboot. lo_npml Lower confidence interval limit when boot = NPMLboot. up_ecdf Upper confidence interval limit when boot = ECDFboot. lo_ecdf Lower confidence interval limit when boot = ECDFboot.

## Author(s)

Mahdis Azadbakhsh

Hanna Jankowski, hkj@mathstat.yorku.ca,
http://www.math.yorku.ca/~hkj/

## References

Azadbakhsh, M., Jankowski, H. and Gao, X. (2014). Computing confidence intervals for log-concave densities. Comput. Statist. Data Anal., to appear.

Baladbaoui, F., Rufibach, K. and Wellner, J. (2009) Limit distribution theory for maximum likelihood estimation of a log-concave density. Ann. Statist., 37(3), 1299–1331.

Tarn Duong (2012). ks: Kernel smoothing. R package version 1.8.10. http://CRAN.R-project.org/package=ks

## 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 ## Not run: ## =================================================== ## Confidence intervals at a fixed point for the density ## =================================================== data(reliability) x.rel <- sort(reliability) # calculate 95 grid <- seq(min(x.rel), max(x.rel), length.out = 200) res <- logConDens(x.rel) ci <- logConCI(res, grid, type = c("nrd", "ECDFboot")) par(las = 1, mar = c(2.5, 3.5, 0.5, 0.5)) hist(x.rel, n = 25, col = gray(0.9), main = "", freq = FALSE, xlab = "", ylab = "", ylim = c(0, 0.0065), border = gray(0.5)) lines(grid, ci$fhat, col = "black", lwd = 2) lines(grid, ci$lo_nrd, col = "red", lwd = 2, lty = 2) lines(grid, ci$up_nrd, col = "red", lwd = 2, lty = 2) lines(grid, ci$lo_ecdf, col = "blue", lwd = 2, lty = 2) lines(grid, ci\$up_ecdf, col = "blue", lwd = 2, lty = 2) legend("topleft", col = c("black", "blue", "red"), lwd = 2, lty = c(1, 2, 2), legend = c("log-concave NPMLE", "CI for type = nrd", "CI for type = ECDFboot"), bty = "n") ## End(Not run) 

logcondens documentation built on May 2, 2019, 6:11 a.m.