Allison Gonzales — Oct 17, 2013, 7:52 PM
##Allison Gonzales
##PH 251D hw#3
##21 October 2013
##Prob 3.1.
status <- c("Dead", "Survived", "Dead", "Survived","Dead", "Survived","Dead", "Survived")
treatment <- c("Tolbutamide", "Tolbutamide", "Placebo", "Placebo", "Tolbutamide", "Tolbutamide", "Placebo", "Placebo")
agegrp <- c("<55", "<55", "<55", "<55", "55+", "55+", "55+", "55+")
freq <- c(8, 98, 5, 115, 22, 76, 16, 69)
df <- data.frame(status, treatment, agegrp, freq)
df
status treatment agegrp freq
1 Dead Tolbutamide <55 8
2 Survived Tolbutamide <55 98
3 Dead Placebo <55 5
4 Survived Placebo <55 115
5 Dead Tolbutamide 55+ 22
6 Survived Tolbutamide 55+ 76
7 Dead Placebo 55+ 16
8 Survived Placebo 55+ 69
##Prob 3.2.
##Prob 3.3.
aidsdata <- read.table("http://www.medepi.net/data/aids.txt", header=TRUE)
plot(aidsdata$year, aidsdata$cases, xlab = "Year", lwd = 2, ylab = "Cases", main = "AIDS cases per year")
##Prob 3.4.
measdata <- read.table("http://www.medepi.net/data/measles.txt", header=TRUE)
plot(measdata$year, measdata$cases, lwd = 2, xlab = "Year", ylab="Cases", main = "Measles cases per year")
plot(measdata$year, measdata$cases, log="y", lwd = 2, xlab = "Year", ylab="Cases", main = "Measles cases per year")
##Prob 3.5.
hepbdata <- read.table("http://www.medepi.net/data/hepb.txt", header=TRUE)
matplot(hepbdata$year, cbind(hepbdata$cases,aidsdata$cases), type = "l", lwd = 2, xlab = "Year", ylab = "Cases", main = "Reported cases of Hepatitis B and AIDS, United States, 1980-2003")
legend(1980, 100000, legend = c("Hepatitis B", "AIDS"), lwd = 2, lty = 1:2, col = 1:2)
##Prob 3.6.
evansdata <- read.table("http://www.medepi.net/data/evans.txt", header=TRUE)
## Recode binary to string labels
evansdata$chd2 <- factor(evansdata$chd, levels = 0:1, labels = c("No", "Yes"))
evansdata$cat2 <- factor(evansdata$cat, levels = 0:1, labels = c("Normal", "High"))
evansdata$smk2 <- factor(evansdata$smk, levels = 0:1, labels = c("Never Smoked", "Ever Smoked"))
evansdata$ecg2 <- factor(evansdata$ecg, levels = 0:1, labels = c("No Abnormality", "Abnormality"))
evansdata$hpt2 <- factor(evansdata$hpt, levels = 0:1, labels = c("No", "Yes"))
##Discretize age into levels
quantile(evansdata$age)
0% 25% 50% 75% 100%
40 46 52 60 76
evansdata$agecat <- cut(evansdata$age, quantile(evansdata$age), include.lowest=TRUE)
table(evansdata$agecat)
[40,46] (46,52] (52,60] (60,76]
169 148 140 152
##Create new hypertension categorical variable
hptcat <- rep(NA, nrow(evansdata))
normbp <- evansdata$sbp<120 & evansdata$dbp<80
hptcat[normbp] <- 1
prehpt <- (evansdata$sbp>=120 & evansdata$sbp<140) | (evansdata$dbp>=80 & evansdata$dbp<90)
hptcat[prehpt] <- 2
hpt1 <- (evansdata$sbp>=140 & evansdata$sbp<160) | (evansdata$dbp>=90 & evansdata$dbp<100)
hptcat[hpt1] <- 3
hpt2 <- evansdata$sbp>=160 & evansdata$dbp>=100
hptcat[hpt2] <- 4
##Table comparing old and new hypertension variables
table (evansdata$hpt, hptcat)
hptcat
1 2 3 4
0 56 165 133 0
1 0 18 124 108
##Prob 3.7.
wnvdata <- read.table("http://www.medepi.net/data/wnv/wnv2004raw.txt", header=TRUE, sep=",", na.strings= c(".", "Unknown"))
wnvdata$date.onset2 <- as.Date(wnvdata$date.onset, format="%m/%d/%Y")
wnvdata$date.tested2 <- as.Date(wnvdata$date.tested, format="%m/%d/%Y")
write.table(wnvdata, "c:/users/allie/documents/wnvdat.txt", sep=",", row.names=FALSE)
##Prob 3.8.
##Plot histogram of case by onset of illness
odata <- read.table("http://www.medepi.net/data/oswego.txt", header=TRUE, sep = "", na.strings=".")
dtm <- paste("4/18/1940", odata$meal.time)
dt.meal <- strptime(mdt, "%m/%d/%Y %I:%M %p")
Error: object 'mdt' not found
dto <- paste(paste(odat$onset.date,"/1940",sep = ""), odat$onset.time)
Error: object 'odat' not found
dt.onset <- strptime(dto, "%m/%d/%Y %I:%M %p")
Error: object 'dto' not found
hist(dt.onset, breaks=30, freq=TRUE)
Error: object 'dt.onset' not found