## Basic usage:
d <- density(faithful$waiting)
plot(d)

## Slide 4
x <- c(1,3,5,6,6.2,6.6,7,8,10)
n <- length(x)
xx <- seq(-2,max(x),len=101)
hist(x,border="white",main="",freq=FALSE)
rug(x,col="red",lwd=2)
for (i in 1:n)
  {
    lines(xx,dnorm(xx,x[i],0.75)/10,col="blue")
  }
d <- density(x,bw=0.75)
lines(d)

## Slide 6
faithful <- read.delim("../../data/faithful.txt")
par(mfrow=c(1,3))
plot(density(faithful$waiting,bw=10),main=expression(lambda==10),xlab="Waiting time (min)")
plot(density(faithful$waiting,bw=2.66),main=expression(lambda==2.66),xlab="Waiting time (min)")
plot(density(faithful$waiting,bw=1),main=expression(lambda==1),xlab="Waiting time (min)")

## Slide 10
col <- c("black","blue","red","green")
d <- vector("list",4)
plot(d[[1]] <- density(faithful$waiting,bw="nrd"),main="",xlab="Waiting time (min)",ylim=c(0,0.04),col=col[1],lwd=2)
lines(d[[2]] <- density(faithful$waiting,bw="sJ"),col=col[2],lwd=2)
lines(d[[3]] <- density(faithful$waiting,bw="ucv"),col=col[3],lwd=2)
lines(d[[4]] <- density(faithful$waiting,bw="bcv"),col=col[4],lwd=2)
legend("topleft",col=col,legend=paste(c("Normal reference","Sheather-Jones","Cross-validation","Biased CV"),formatC(unlist(lapply(d,function(x){x$bw})),2,format='f')),lty=1,lwd=2)

## Slider
fun <- function(x) plot(density(faithful$waiting,adjust=x))
slider(fun,0.1,3,0.01)
bmd <- read.delim("../../data/bmd.txt")
fun <- function(x){plot.loess(loess(spnbmd~age,bmd,span=x,deg=1),col="red",lwd=3,type="l")}
slider(fun,0.1,1,0.01)

## Slides 16-18
plot(faithful$eruptions,faithful$waiting,pch=19,xlab="Duration (min)",ylab="Waiting (min)",col="blue")
require(KernSmooth)
fit <- bkde2D(faithful, bandwidth=c(0.7, 7))
contour(fit$x1, fit$x2, fit$fhat,xlab="Duration (min)",ylab="Waiting (min)")
persp(fit$x1,fit$x2,fit$fhat,shade=0.25,col="lightblue",theta=30,phi=20,ticktype="detailed",xlab="duration",ylab="waiting",zlab="Density")

## Slide 22
heart <- read.delim("../../data/heart.txt")
trellis.par.set(superpose.line=list(col=c("gray","red"),lwd=2),superpose.symbol=list(col=c(rgb(.5,.5,.5,alpha=.4),rgb(1,0,0,alpha=.4)),pch=16,cex=.5))
densityplot(~typea,heart,group=chd,auto.key=list(columns=2,text=c("No CHD","CHD")),xlab="Type A (score)",bw="SJ")

## Slide 23
d0 <- with(heart,density(typea[chd==0],bw="SJ"))
d1 <- with(heart,density(typea[chd==1],bw="SJ"))
x <- seq(min(heart$typea),max(heart$typea),len=101)
f0 <- approx(d0$x,d0$y,xout=x)
f1 <- approx(d1$x,d1$y,xout=x)
p1 <- mean(heart$chd)
plot(x,p1*f1$y/(p1*f1$y+(1-p1)*f0$y),type="l",ylim=c(0,1),xlab="Type A (score)",ylab="Pr(CHD)",lwd=3)
abline(h=.5,lty=2)
