# RProject8_3_bootstrap_location.r
#   Problem 10.26 of Rice
#   Hampson and Walker data on heats of sublimation of
#       platinum,iridium, and rhodium

# To install these packages, uncomment the next two lines
#install.packages("MASS")
#install.packages("boot")

library(MASS)
## Warning: package 'MASS' was built under R version 3.1.3
library(boot)
## Warning: package 'boot' was built under R version 3.1.3
x.platinum=scan(file="Rice 3e Datasets/ASCII Comma/Chapter 10/platinum.txt")
x.iridium=scan(file="Rice 3e Datasets/ASCII Comma/Chapter 10/iridium.txt")
x.rhodium=scan(file="Rice 3e Datasets/ASCII Comma/Chapter 10/rhodium.txt")

# Parts (a)-(d)
x=x.platinum
hist(x)

stem(x)
##
##   The decimal point is at the |
##
##   132 | 7
##   134 | 134578889900224488
##   136 | 36
##   138 |
##   140 | 2
##   142 | 3
##   144 |
##   146 | 58
##   148 | 8
boxplot(x)

plot(x)

# (e)
# Measurements 8,9,10 are all very high as are measurements 14 and 15
# The data do not appear indepenent

# (f) Measures of location
mean(x); mean(x,trim=.1);mean(x,trim=.2);median(x); huber(x,k=1.5)[[1]]
## [1] 137.0462
## [1] 136.3091
## [1] 135.2875
## [1] 135.1
## [1] 135.3841
# (g) Standard error of the sample mean and approximate 90 percent conf. interval
x.stdev=sqrt(var(x))
x.mean.sterr=x.stdev/sqrt(length(x))

print(x.mean.sterr)
## [1] 0.8724542
alpha=0.10
z.upperalphahalf=qnorm(1-alpha/2)
x.mean=mean(x)
x.mean.ci<-c(-1,1)*z.upperalphahalf* x.mean.sterr + x.mean
print(x.mean); print(x.mean.ci)
## [1] 137.0462
## [1] 135.6111 138.4812
# parts (i), (j), (k).
# Bootstrap confidence intervals of different location measures

fcn.median<- function(x, d) {
return(median(x[d]))
}
fcn.mean<- function(x, d) {
return(mean(x[d]))
}

fcn.trimmedmean <- function(x, d, trim=0) {
return(mean(x[d], trim/length(x)))
}

fcn.huber<-function(x,d){
x.huber=huber(x[d],k=1.5)
return(x.huber[[1]])
}

# Bootstrap analysis of sample mean:
x.boot.mean= boot(x, fcn.mean, R=1000)
plot(x.boot.mean)
print(x.boot.mean$t0) ## [1] 137.0462 print(sqrt(var(x.boot.mean$t)))
##           [,1]
## [1,] 0.8406306
boot.ci(x.boot.mean,conf=.95, type="basic")
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = x.boot.mean, conf = 0.95, type = "basic")
##
## Intervals :
## Level      Basic
## 95%   (135.2, 138.5 )
## Calculations and Intervals on Original Scale
title(paste("Sample Mean:  StDev=",as.character(round(sqrt(var(x.boot.mean$t)),digits=4)),sep=""), adj=1) # Bootstrap analysis of sample median: x.boot.median= boot(x, fcn.median, R=1000) plot(x.boot.median) title(paste("Median: StDev=",as.character(round(sqrt(var(x.boot.median$t)),digits=4)),sep=""),
adj=1)

print(x.boot.median$t0) ## [1] 135.1 print(sqrt(var(x.boot.median$t)))
##          [,1]
## [1,] 0.238632
boot.ci(x.boot.median,conf=.95, type="basic")
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = x.boot.median, conf = 0.95, type = "basic")
##
## Intervals :
## Level      Basic
## 95%   (134.4, 135.4 )
## Calculations and Intervals on Original Scale
# Bootstrap analysis of sample trimmed mean:
x.boot.trimmedmean= boot(x, fcn.trimmedmean,trim=.2, R=1000)
plot(x.boot.trimmedmean)
title(paste("Trimmed Mean:  StDev=",as.character(round(sqrt(var(x.boot.trimmedmean$t)),digits=4)),sep=""), adj=1) print(x.boot.trimmedmean$t0)
## [1] 137.0462
print(sqrt(var(x.boot.trimmedmean$t))) ## [,1] ## [1,] 0.833879 boot.ci(x.boot.trimmedmean,conf=.95, type="basic") ## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS ## Based on 1000 bootstrap replicates ## ## CALL : ## boot.ci(boot.out = x.boot.trimmedmean, conf = 0.95, type = "basic") ## ## Intervals : ## Level Basic ## 95% (135.3, 138.6 ) ## Calculations and Intervals on Original Scale # Bootstrap analysis of Huber M estimate: x.boot.huber= boot(x, fcn.huber, R=1000) plot(x.boot.huber) title(paste("Huber M Est: StDev=",as.character(round(sqrt(var(x.boot.huber$t)),digits=4)),sep=""),
adj=1)

print(x.boot.huber$t0) ## [1] 135.3841 print(sqrt(var(x.boot.huber$t)))
##           [,1]
## [1,] 0.3827326
boot.ci(x.boot.huber,conf=.95, type="basic")
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = x.boot.huber, conf = 0.95, type = "basic")
##
## Intervals :
## Level      Basic
## 95%   (134.4, 135.9 )
## Calculations and Intervals on Original Scale