################
# Section 12
################

twosamp <- function(yvec,trtvec,alpha=0.05,header="") {
#################################################
# A function to compute a two-sample t-test and confidence 
# interval (equal-variance, independent samples).  yvec is 
# a numeric vector containing both samples' data.  trtvec 
# is a vector, same length as yvec, of treatment
# identifiers for the data in yvec.  A boxplot comparing
# the treatments' data is constructed.  Output is a one-row
# data frame reporting the results of the test and 
# confidence interval
##################################################
trtvec<-as.factor(trtvec)
boxplot(split(yvec,trtvec))
title(header)
ybar<-tapply(yvec,trtvec,mean)
varvec<-tapply(yvec,trtvec,var)
nvec<-table(trtvec)
error.df<-nvec[1]+nvec[2]-2
pooled.var<-((nvec[1]-1)*varvec[1]+(nvec[2]-1)*varvec[2])/error.df
diff12estimate<-ybar[1]-ybar[2]
stderr<-sqrt(pooled.var*((1/nvec[1])+(1/nvec[2])))
tratio<-diff12estimate/stderr
twosidedP<-2*(1-pt(abs(tratio),error.df))
tcrit<-qt(1-alpha/2,error.df)
lower<-diff12estimate-tcrit*stderr
upper<-diff12estimate+tcrit*stderr
out<-data.frame(diff12estimate,stderr,tratio,twosidedP,lower,upper,alpha)
out
}

twosamp

fix(twosamp)

my.yvec <- c(3.1, 4.2, 6.2, 2.7, 3.4, 3.1, 2.1, 7.5, 5.4, 2.2, 2.5, 4.9, 6.3, 2.4, 3.5)

my.trtvec <- c("med", "control", "med", "control", "med", "control", "med", "control",
"med", "med", "med", "med", "control", "control", "med")

twosamp(my.yvec, my.trtvec)

twosamp(my.yvec, my.trtvec, alpha=0.01, header="Results of Clinical Trial")

data(sleep)

?sleep

results<-twosamp(sleep$extra, sleep$group)

results

# Using the source command to read in code that is saved in a file:
# Save twosamp.txt then change file directory, and then:

source("twosamp.txt")


#Using t.test to analyze data; two different ways
plot(extra~group,data=sleep)
t.test(extra~group,data=sleep, var.equal=T)
with(sleep,t.test(extra[group==1],extra[group==2], var.equal=T))

# go through browser() example:

#Experiment with browser()
fix(twosamp) #Place browser() after ybar statement
twosamp(sleep$extra,sleep$group)
ybar
varvec
nvec
n         # this enters the "step-through debugger"
n
ybar
varvec
nvec
n
ybar
varvec
nvec
Q   #or use c for normal execution (c stands for "continue")