# Read in data Diarrhea
DiaData <- read.delim("http://myweb.uiowa.edu/pbreheny/data/diarrhea.txt")


# Boxplot
boxplot(DiaData$Stool ~ DiaData$Group, col = c("burlywood4", "darkgoldenrod3"), main = "Effect of bismuth salicylate", horizontal = T)


# Summary Statistics
#Find n in each group
table(DiaData$Group)

#Find mean of each group
by(DiaData$Stool, DiaData$Group, mean)

#Find standard deviation of each group 
by(DiaData$Stool, DiaData$Group, sd)



# 
x1 <- 260.2976 # Mean of control group
x2 <- 181.8706 # Mean of treatment group

sd1 <- 253.6135 # SD of control group
sd2 <- 197.3636 # SD of treatment group

n1 <- 84 # n for control group
n2 <- 85 # n for treatment group

sd_pooled <- sqrt(((n1 - 1)*sd1^2 + (n2 - 1)*sd2^2) / (n1 + n2 - 2))
SE <- sd_pooled * sqrt((1 / n1) + (1 / n2))
df <- n1 + n2 - 2

t <- (x1 - x2) / (SE)
2*pt(t, df = df, lower.tail = F)


# 
(x1 - x2) - qt(0.975, df = df)*SE
(x1 - x2) + qt(0.975, df = df)*SE


# Equal Variance t test
t.test(DiaData$Stool~DiaData$Group, var.equal = TRUE)


# 
by(DiaData$Stool, DiaData$Group, sd)


# Unequal Varaince t test
t.test(DiaData$Stool~DiaData$Group, var.equal = FALSE)


# read in data heart
heart <- read.delim("http://myweb.uiowa.edu/pbreheny/data/infant-heart.txt")