Thank you very much for your help, Greg!
Here's my ludicrous vision/attempt/whatever :o)
### Idea for making a density() with confidence interval ###
xx <- faithful$eruptions
xx.hist <- hist(xx, breaks="FD", freq=F)
# plot(xx.hist$mids, xx.hist$density)
# gives a rough "plot(density(xx))" output
# xx.hist$counts is later used to "repeat" observations
# for loess(), making the CI band narrower where there
# are more observations (or data in each histogram "bin")
x <- c() # preparing variables...
y <- c()
for (i in 1:length(xx.hist$mids)) {
# going through each xx.hist$mids/xx.hist$density pair
for (j in 1:xx.hist$counts[i]) {
# ...and repeating observations according to xx.hist$counts
x <- append(x, xx.hist$mids[i])
y <- append(y, xx.hist$density[i])
}
}
xx.loess <- loess(y ~ x, span=1)
xx.predict <- predict(xx.loess, se=T)
plot(x,y)
lines(x, xx.predict$fit, col="red")
lines(x, xx.predict$fit + 1.96 * xx.predict$s, col="red", lty=5)
lines(x, xx.predict$fit - 1.96 * xx.predict$s, col="red", lty=5)
# Kind regards, David
> Here is a simple approach that uses bootstrapping (this could
> probably be improved by using better bootstrap estimates and not
> ignoring the dependence between points):
>
> xx <- faithful$eruptions
>
> fit1 <- density(xx)
>
> fit2 <- replicate(10000, { x <- sample(xx, replace=TRUE); density(x,
> from=min(fit1$x), to=max(fit1$x))$y } )
>
> fit3 <- apply(fit2, 1, quantile, c(0.025,0.975) )
>
> plot(fit1, ylim=range(fit3)) polygon( c(fit1$x, rev(fit1$x)),
> c(fit3[1,], rev(fit3[2,])), col='grey', border=F) lines(fit1)
>
>
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
