Compute log-concave density estimator and related quantities
Compute the log-concave and smoothed log-concave density estimator.
logConDens(x, xgrid = NULL, smoothed = TRUE, print = FALSE,
gam = NULL, xs = NULL)x |
Vector of independent and identically distributed numbers, not necessarily unique. |
xgrid |
Governs the generation of weights for observations. See |
smoothed |
If |
print |
|
gam |
Only necessary if |
xs |
Only necessary if |
See activeSetLogCon for details on the computations.
logConDens returns an object of class "dlc", a list containing the
following components:
xn, x, w, phi, IsKnot, L, Fhat, H,
n, m, knots, mode, and sig as generated
by activeSetLogCon. If smoothed = TRUE, then the returned object additionally contains
f.smoothed, F.smoothed, gam, and xs as generated by evaluateLogConDens. Finally, the
entry smoothed of type "logical" returnes the value of smoothed.
The methods summary.dlc and plot.dlc are used to obtain a summary and generate plots of the estimated
density.
Kaspar Rufibach, kaspar.rufibach@gmail.com,
http://www.kasparrufibach.ch
Duembgen, L, Huesler, A. and Rufibach, K. (2010). Active set and EM algorithms for log-concave densities based on complete and censored data. Technical report 61, IMSV, Univ. of Bern, available at http://arxiv.org/abs/0707.4643.
Duembgen, L. and Rufibach, K. (2009). Maximum likelihood estimation of a log–concave density and its distribution function: basic properties and uniform consistency. Bernoulli, 15(1), 40–68.
Duembgen, L. and Rufibach, K. (2011). logcondens: Computations Related to Univariate Log-Concave Density Estimation. Journal of Statistical Software, 39(6), 1–28. http://www.jstatsoft.org/v39/i06
## ===================================================
## Illustrate on simulated data
## ===================================================
## Set parameters
n <- 50
x <- rnorm(n)
res <- logConDens(x, smoothed = TRUE, print = FALSE, gam = NULL,
xs = NULL)
summary(res)
plot(res, which = "density", legend.pos = "topright")
plot(res, which = "log-density")
plot(res, which = "CDF")
## Compute slopes and intercepts of the linear functions that
## compose phi
slopes <- diff(res$phi) / diff(res$x)
intercepts <- -slopes * res$x[-n] + res$phi[-n]
## ===================================================
## Illustrate method on reliability data
## Reproduce Fig. 2 in Duembgen & Rufibach (2009)
## ===================================================
## Set parameters
data(reliability)
x <- reliability
n <- length(x)
res <- logConDens(x, smooth = TRUE, print = TRUE)
phi <- res$phi
f <- exp(phi)
## smoothed log-concave PDF
f.smoothed <- res$f.smoothed
xs <- res$xs
## compute kernel density
sig <- sd(x)
h <- sig / sqrt(n)
f.kernel <- rep(NA, length(xs))
for (i in 1:length(xs)){
xi <- xs[i]
f.kernel[i] <- mean(dnorm(xi, mean = x, sd = h))
}
## compute normal density
mu <- mean(x)
f.normal <- dnorm(xs, mean = mu, sd = sig)
## ===================================================
## Plot resulting densities, i.e. reproduce Fig. 2
## in Duembgen and Rufibach (2009)
## ===================================================
plot(0, 0, type = 'n', xlim = range(xs), ylim = c(0, 6.5 * 10^-3))
rug(res$x)
lines(res$x, f, col = 2)
lines(xs, f.normal, col = 3)
lines(xs, f.kernel, col = 4)
lines(xs, f.smoothed, lwd = 3, col = 5)
legend("topleft", c("log-concave", "normal", "kernel",
"log-concave smoothed"), lty = 1, col = 2:5, bty = "n")
## ===================================================
## Plot log-densities
## ===================================================
plot(0, 0, type = 'n', xlim = range(xs), ylim = c(-20, -5))
legend("bottomright", c("log-concave", "normal", "kernel",
"log-concave smoothed"), lty = 1, col = 2:5, bty = "n")
rug(res$x)
lines(res$x, phi, col = 2)
lines(xs, log(f.normal), col = 3)
lines(xs, log(f.kernel), col = 4)
lines(xs, log(f.smoothed), lwd = 3, col = 5)
## ===================================================
## Confidence intervals at a fixed point for the density
## see help file for logConCI()
## ===================================================Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.