#R/Splus program to plot empirical distribution functions and pointwise (1-alpha)% confidence band of edf

#

#given two files:

#	data file contains data points separated by spaces

#	ci file contains the pointwise (1-alpha)% confidence bands (L,H) for the order statistics in a table

#

#this program will plot the empirical distribution function of the data along with the (1-alpha)% confidence band



#read in data and ci bands

data <- scan("data.dat")  		#should contain pathway to data file

table <- read.table("ci.dat",header=T)  	#should contain pathway to ci file



data <- sort(data)		#orders data

a <- table[,1]			#a is lower confidence band for order stats

b <- table[,2]			#b is upper confidence band for order stats

n <- length(data)		#n is number of data points



#edf associated with each data point

s <- 1:n				

edf <- s/n



prec <- sqrt(var(data))/100   		#precision (prec) tells how fine the plot should be.  this default is sd(data)/100

x <- seq(data[1]-5*prec, data[n]+5*prec, prec)

y <- rep(NA,length(x))				#y is edf(x)

h <- rep(NA,length(x))				#h is lower confidence band for edf

l <- rep(NA,length(x))				#l is upper confidence band for edf

for(i in 1:length(x)) {
	j <- sum(data <= x[i])
	if (j == 0) {
		y[i] <- 0
		h[i] <- b[1]
		l[i] <- 0
	}
	if (j > 0 & j < n) {
		y[i] <- edf[j]
		h[i] <- b[j+1]
		l[i] <- a[j]
	}
	if (j == n) {
	y[i] <- 1
	h[i] <- 1
	l[i] <- a[n]
	
}
}



plot(data,edf, type="n", xlab="Data",ylab="", xlim=c(data[1]-5*prec, data[n]+5*prec),ylim=c(0,1)) 		#sizes plot

lines(x,y,lwd=3)										#plots edf

lines(x,h,lwd=1,lty=2)										#plots upper band

lines(x,l,lwd=1,lty=2)										#plots lower band