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
|
# Demo for prediction and confidence intervals
par(ask = TRUE)
# Classical Gaussian Version
sub <- "Classical Gaussian Version"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
ss = seq(0,27,len=100)
flm <- lm(dist ~ speed + I(speed^2),data = cars)
cflm <- predict(flm,list(speed=ss), data = cars, interval = "confidence")
pflm <- predict(flm,list(speed=ss), data = cars, interval = "prediction")
matlines(ss, pflm,lty = c(1,2,2),col=c("black","red","red"))
matlines(ss, cflm,lty = c(1,2,2),col=c("black","blue","blue"))
# Now try several different ways to compute alternative bands based on QR
fl=rq(dist~speed+I(speed^2),data=cars,tau= 0.25)
fu=rq(dist~speed+I(speed^2),data=cars,tau= 0.75)
# 1. Covariance matrix estimation using "nid"
sub <- "Covariance matrix estimation using nid"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
pl <- predict(fl,list(speed=ss),interval = "confidence", se = "nid")
pu <- predict(fu,list(speed=ss),interval = "confidence", se = "nid")
matlines(ss,pl,lty = c(1,2,0), col=c("black","red","red"))
matlines(ss,pu,lty = c(1,0,2), col=c("black","red","red"))
# 2. Covariance matrix estimation using xy-bootstrap
sub <- "Covariance matrix estimation using xy-bootstrap"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
pl <- predict(fl,list(speed=ss),interval = "confidence", se = "boot",bsmethod = "xy")
pu <- predict(fu,list(speed=ss),interval = "confidence", se = "boot",bsmethod = "xy")
matlines(ss,pl,lty = c(1,2,0), col=c("black","red","red"))
matlines(ss,pu,lty = c(1,0,2), col=c("black","red","red"))
# 3. Covariance matrix estimation using weighted xy-bootstrap
sub <- "Covariance matrix estimation using weighted xy-bootstrap"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
pl <- predict(fl,list(speed=ss),interval = "confidence", se = "boot",bsmethod = "wxy")
pu <- predict(fu,list(speed=ss),interval = "confidence", se = "boot",bsmethod = "wxy")
matlines(ss,pl,lty = c(1,2,0), col=c("black","red","red"))
matlines(ss,pu,lty = c(1,0,2), col=c("black","red","red"))
# 4. Percentile method using xy-bootstrap
sub <- "Percentile method using xy-bootstrap"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
pl <- predict(fl,list(speed=ss),interval = "confidence",
type="percentile", se = "boot", bsmethod = "xy")
pu <- predict(fu,list(speed=ss),interval = "confidence",
type="percentile", se = "boot", bsmethod = "xy")
matlines(ss,pl,lty = c(1,2,0), col=c("black","red","red"))
matlines(ss,pu,lty = c(1,0,2), col=c("black","red","red"))
# 5. Direct method using "iid" covariance matrix
sub <- "Direct method using iid covariance matrix"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
pl <- predict(fl,list(speed=ss),interval = "confidence", type="direct", se = "iid")
pu <- predict(fu,list(speed=ss),interval = "confidence", type="direct", se = "iid")
matlines(ss,pl,lty = c(1,2,0), col=c("black","red","red"))
matlines(ss,pu,lty = c(1,0,2), col=c("black","red","red"))
# 6. Direct method using "nid" covariance matrix
sub <- "Direct method using nid covariance matrix"
plot(cars$speed,cars$dist,xlim=c(0,27),ylim=c(-10,130), main = "Ezekiel Data",
sub = sub, xlab = "speed (mph)",ylab = "stopping distance (ft)")
pl <- predict(fl,list(speed=ss),interval = "confidence", type="direct", se = "nid")
pu <- predict(fu,list(speed=ss),interval = "confidence", type="direct", se = "nid")
matlines(ss,pl,lty = c(1,2,0), col=c("black","red","red"))
matlines(ss,pu,lty = c(1,0,2), col=c("black","red","red"))
par(ask = FALSE)
|
