Description Usage Arguments Details Value Author(s) References Examples
Compute the logconcave and smoothed logconcave density estimator.
This function is as it is in the logcondens package except we've added
the 'prec' variable as an argument, and modified the values returned
as output, to be in line with the activeSetLogCon.mode
function.
1 2 
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 
prec 
Governs precision of various subfunctions, e.g. the NewtonRaphson procedure. 
See activeSetLogCon
for details on the
computations.
logConDens
returns an object of class
"dlc"
, a list containing the following components: xn
,
x
, w
, L
, IsKnot
, knots
, phi
,
fhat
,
Fhat
, H
, n
, m
, mode
, dlcMode
,
sig
, phi.f
, fhat.f
, Fhat.f
, E.f
,
phiPL
,
phiPR
, phiPL.f
, and phiPR.f
,
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
Lutz Duembgen, duembgen@stat.unibe.ch,
http://www.staff.unibe.ch/duembgen
Duembgen, L, Huesler, A. and Rufibach, K. (2010). Active set and EM algorithms for logconcave 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 LogConcave Density Estimation. Journal of Statistical Software, 39(6), 1–28. http://www.jstatsoft.org/v39/i06
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82  ## ===================================================
## 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 = "logdensity")
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 logconcave 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("logconcave", "normal", "kernel",
"logconcave smoothed"), lty = 1, col = 2:5, bty = "n")
## ===================================================
## Plot logdensities
## ===================================================
plot(0, 0, type = 'n', xlim = range(xs), ylim = c(20, 5))
legend("bottomright", c("logconcave", "normal", "kernel",
"logconcave 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()
## ===================================================

Loading required package: logcondens
Attaching package: 'logcondens.mode'
The following objects are masked from 'package:logcondens':
activeSetLogCon, intF, logConDens
The following object is masked from 'package:base':
dir.exists
Warning messages:
1: In rgl.init(initValue, onlyNULL) : RGL: unable to open X11 display
2: 'rgl_init' failed, running with rgl.useNULL = TRUE
3: .onUnload failed in unloadNamespace() for 'rgl', details:
call: fun(...)
error: object 'rgl_quit' not found
Estimation of a logconcave density from i.i.d. data
Number of initial observations: n = 50
Number of unique observations (or grid points): m = 50
loglikelihood: 2.31
Maximum likelihood estimate:
Mode: x[17] = 0.43
Value of logdensity at mode: 0.73
Value of density at mode: 0.48
Smoothed maximum likelihood estimate:
Mode: 0.32
Value of logdensity at mode: 0.86
Value of density at mode: 0.42
Number of knots of the MLE: 6
Knots of the MLE:
x[1, 14, 17, 39, 44, 50] =
1.71, 0.58, 0.43, 0.75, 1.25, 2.71
[1] "ASLC: Beginning"
[1] "iter1=2 / L=6.7569 / max(H)=2.163 / #knots = 3"
[1] "iter1=3 / L=6.7337 / max(H)=1.3981 / #knots = 4"
[1] "iter1=4 / L=6.7266 / max(H)=0.4054 / #knots = 5"
[1] "iter1=5 / L=6.7212 / max(H)=0.2512 / #knots = 6"
[1] "iter1=6 / L=6.719 / max(H)=0.1588 / #knots = 7"
[1] "iter1=7 / L=6.7177 / max(H)=0.0748 / #knots = 7"
[1] "iter1=7 / L=6.7177 / max(H)=0.0748 / #knots = 7"
[1] "iter1=8 / L=6.7176 / max(H)=0.0363 / #knots = 8"
[1] "iter1=9 / L=6.7168 / max(H)=0.0321 / #knots = 8"
[1] "iter1=9 / L=6.7168 / max(H)=0.0321 / #knots = 8"
[1] "iter1=10 / L=6.7166 / max(H)=0.0207 / #knots = 9"
[1] "iter1=11 / L=6.7165 / max(H)=0.0174 / #knots = 9"
[1] "iter1=11 / L=6.7165 / max(H)=0.0174 / #knots = 9"
[1] "iter1=12 / L=6.7164 / max(H)=0.0074 / #knots = 9"
[1] "iter1=12 / L=6.7164 / max(H)=0.0074 / #knots = 9"
[1] "iter1=13 / L=6.7163 / max(H)=0.005 / #knots = 9"
[1] "iter1=13 / L=6.7163 / max(H)=0.005 / #knots = 9"
[1] "iter1=14 / L=6.7163 / max(H)=2e04 / #knots = 10"
[1] "iter1=15 / L=6.7163 / max(H)=0 / #knots = 10"
[1] "iter1=15 / L=6.7163 / max(H)=0 / #knots = 10"
[1] "iter1=16 / L=6.7163 / max(H)=0 / #knots = 10"
[1] "iter1=16 / L=6.7163 / max(H)=0 / #knots = 10"
[1] "ASLC: returning"
[1] "10% of computation of smooth estimates done"
[1] "20% of computation of smooth estimates done"
[1] "30% of computation of smooth estimates done"
[1] "40% of computation of smooth estimates done"
[1] "50% of computation of smooth estimates done"
[1] "60% of computation of smooth estimates done"
[1] "70% of computation of smooth estimates done"
[1] "80% of computation of smooth estimates done"
[1] "90% of computation of smooth estimates done"
[1] "100% of computation of smooth estimates done"
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.