# [R] Re: conf.int. for ecdfs {was "Two questions"}

Martin Maechler maechler at stat.math.ethz.ch
Mon Oct 22 19:15:35 CEST 2001

```>>>>> "Christian" == Christian Blouin <bongo at havana.biochem.dal.ca> writes:

Christian> I have two questions that I could not answer from the
Christian> documentation.

Christian> A - ecdf and confidence intervals : Is there a (simple) way
Christian> to generate confidence intervals (95%) for a ecdf?

About two weeks ago, Kjetil Halvorsen (CC'ed above)
sent me code to just do this, and I have improved on the code a bit
with the idea to incorporate into the stepfun package (part of "standard R").

However, it's content is also linked with the (hidden) pkstwo() function
in ctest's  kolomogorov-smirnov testing.
Also, I think part of this might become an option to ecdf() or plot.ecdf()
rather than the current setup.
Unfortunately, it's not finished for integration into R,
but here is the current code (without help pages) that works quite nicely:

## Fixme: In the following, computing and plotting should be separated

###--> ./ecdf.R plot.ecdf() should get conf.type and conf.int argument!!

ecdf.ksCI <- function(x, main = NULL, sub = NULL,
xlab = deparse(substitute(x)), ...)
{
xlab
if(is.null(main))
main <- paste("ecdf(",deparse(substitute(x)),") + 95% K.S. bands",
sep="")
n <- length(x)
if(is.null(sub))
sub <- paste("n = ", n)
ec <- ecdf(x)
xx <- get("x", envir=environment(ec))# = sort(x)
yy <- get("y", envir=environment(ec))
D <- approx.ksD(n)
yyu <- pmin(yy+D, 1)
yyl <- pmax(yy-D, 0)
ecu <- stepfun(xx, c(yyu, 1) )
ecl <- stepfun(xx, c(yyl, yyl[n]) )

## Plots -- all calling  plot.stepfun

plot(ec, main = main, sub = sub, xlab = xlab)## no "..." here ???  __hmm..__
col.hor="red" , col.vert="red", ...)
col.hor="red", col.vert="red", ...)
}

approx.ksD <- function(n)
{
## approximations for the critical level for Kolmogorov-Smirnov statistic D,
## for confidence level 0.95. Taken from Bickel & Doksum, table IX, p.483
## and Lienert G.A.(1975) who attributes to Miller,L.H.(1956), JASA
ifelse(n > 80,
1.358 /( sqrt(n) + .12 + .11/sqrt(n)),##Bickel&Doksum, table IX,p.483

splinefun(c(1:9, 10, 15, 10 * 2:8),# from Lienert
c(.975,   .84189, .70760, .62394, .56328,# 1:5
.51926, .48342, .45427, .43001, .40925,# 6:10
.33760, .29408, .24170, .21012,# 15,20,30,40
.18841, .17231, .15975, .14960)) (n))
}

###-- Example from Kjetil :

ecdf.ksCI( rchisq(50,3) )

--
Martin Maechler <maechler at stat.math.ethz.ch>	http://stat.ethz.ch/~maechler/
Seminar fuer Statistik, ETH-Zentrum  LEO D10	Leonhardstr. 27
ETH (Federal Inst. Technology)	8092 Zurich	SWITZERLAND
phone: x-41-1-632-3408		fax: ...-1228			<><
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

```