Run a program file (filename1.R) using the 'source' command;
setwd("~/Dropbox/PH251D/final-project-1/")
source("filename1.R", echo = TRUE)
##
## > X <- matrix(rnorm(100))
##
## > Y <- matrix(rnorm(100))
##
## > Data <- data.frame(Y, X)
##
## > plot(Data)
##
## > fit <- lm(Y ~ X)
##
## > plot(fit)
require(MASS)
## Loading required package: MASS
require(lattice)
## Loading required package: lattice
data(anorexia)
data(leuk)
Demonstrate reading an ASCII data file (filename2.dat) to create a 'data frame';
Data2 <- read.csv("filename2.csv", header = TRUE)
Data2
## group strikes
## 1 a 3
## 2 b 5
## 3 c 2
## 4 b 4
## 5 a 4
## 6 a 3
## 7 a 4
## 8 c 2
## 9 a 1
## 10 b 5
## 11 d 3
## 12 c 2
## 13 d 7
## 14 b 5
## 15 a 4
## 16 a 4
## 17 d 3
## 18 b 7
## 19 c 3
## 20 c 4
## 21 b 3
## 22 a 2
## 23 b 1
## 24 a 4
## 25 a 5
## 26 b 7
## 27 d 3
## 28 c 5
## 29 a 4
## 30 b 7
## 31 d 5
## 32 a 6
summary(Data2)
## group strikes
## a:12 Min. :1.00
## b: 9 1st Qu.:3.00
## c: 6 Median :4.00
## d: 5 Mean :3.97
## 3rd Qu.:5.00
## Max. :7.00
table(Data2)
## strikes
## group 1 2 3 4 5 6 7
## a 1 1 2 6 1 1 0
## b 1 0 1 1 3 0 3
## c 0 3 1 1 1 0 0
## d 0 0 3 0 1 0 1
hist((Data2$strikes))
class(Data2)
## [1] "data.frame"
Demonstrate simple data manipulation (e.g., variable transformation, recoding, etc.);
require(epitools)
## Loading required package: epitools
data(anorexia)
head(anorexia)
## Treat Prewt Postwt
## 1 Cont 80.7 80.2
## 2 Cont 89.4 80.1
## 3 Cont 91.8 86.4
## 4 Cont 74.0 86.3
## 5 Cont 78.1 76.1
## 6 Cont 88.3 78.1
anorexia$wtgain <- anorexia$Postwt - anorexia$Prewt
fivenum(anorexia$wtgain)
## [1] -12.20 -2.45 1.65 9.20 21.50
anorexia$wtgainCat <- cut(anorexia$wtgain, breaks = c(-13, -3, 2, 9, 10, 20),
right = FALSE)
table(anorexia$Treat, anorexia$wtgainCat)
##
## [-13,-3) [-3,2) [2,9) [9,10) [10,20)
## CBT 5 13 5 0 5
## Cont 10 7 5 0 4
## FT 2 2 4 2 6
Demonstrate the use of calendar and Julian dates;
data(leuk)
attach(leuk)
head(leuk)
## wbc ag time
## 1 2300 present 65
## 2 750 present 156
## 3 4300 present 100
## 4 2600 present 134
## 5 6000 present 16
## 6 10500 present 108
class(time) #survival time in weeks
## [1] "integer"
leuk$days <- time * 7
leuk$deathDate <- as.Date(days, origin = "1965-01-01") #arbitrarily chosen date- reference from 1965; assuming treatment was given all the same origin date and within 7 days.
## Error: object 'days' not found
head(leuk)
## wbc ag time days
## 1 2300 present 65 455
## 2 750 present 156 1092
## 3 4300 present 100 700
## 4 2600 present 134 938
## 5 6000 present 16 112
## 6 10500 present 108 756
detach(leuk)
Conduct a simple analysis using existing functions (from R, colleagues, etc.);
require(epitools)
data3 <- matrix(c(53, 274, 121, 640), 2, 2)
rownames(data3) <- c("maternal smoke", "no smoke")
colnames(data3) <- c("premature", "full term")
data3
## premature full term
## maternal smoke 53 121
## no smoke 274 640
oddsratio(data3)
## $data
## premature full term Total
## maternal smoke 53 121 174
## no smoke 274 640 914
## Total 327 761 1088
##
## $measure
## NA
## odds ratio with 95% C.I. estimate lower upper
## maternal smoke 1.000 NA NA
## no smoke 1.025 0.7153 1.451
##
## $p.value
## NA
## two-sided midp.exact fisher.exact chi.square
## maternal smoke NA NA NA
## no smoke 0.8931 0.9282 0.8989
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"
Conduct a simple analysis demonstrating simple programming (e.g., a 'for' loop);
x <- 1:7
for (i in 1:5) {
sum <- sum(i)
sumplus <- sum(i + 2)
result <- x * sum(i + 1)
}
sum
## [1] 5
sumplus
## [1] 7
result
## [1] 6 12 18 24 30 36 42
Conduct a simple analysis demonstrating an original function created by student;
values <- function(x) {
matrix((c(log(x), exp(x), (x)^2)), 4, 3)
}
table <- values(1:4)
table
## [,1] [,2] [,3]
## [1,] 0.0000 2.718 1
## [2,] 0.6931 7.389 4
## [3,] 1.0986 20.086 9
## [4,] 1.3863 54.598 16
Create a simple graph with title, axes labels and legend, and output to file;
matplot(1:4, table, type = "l", xlab = "x", ylab = "values(x)", col = (2:4))
legend("top", c("log(x)", "exp(x)", "(x)^2"), lty = 1, col = (2:4))
pdf("plot.pdf")
matplot(1:4, table, type = "l", xlab = "x", ylab = "values(x)", col = (2:4))
legend("top", c("log(x)", "exp(x)", "(x)^2"), lty = 1, col = (2:4))
dev.off()
## pdf
## 2
Demonstrate the use of regular expressions;
grep("C", anorexia$Treat)
## [1] 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
## [24] 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
## [47] 47 48 49 50 51 52 53 54 55
anorexia$Treat[grep("C", anorexia$Treat)]
## [1] Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont
## [15] Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont Cont CBT CBT
## [29] CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT
## [43] CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT CBT
## Levels: CBT Cont FT
grep("CBT", anorexia$Treat)
## [1] 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49
## [24] 50 51 52 53 54 55
Demonstrate the use of the 'sink' function to generate an output file;
sink("sinkfile")
values(4:7)
## [,1] [,2] [,3]
## [1,] 1.386 54.6 16
## [2,] 1.609 148.4 25
## [3,] 1.792 403.4 36
## [4,] 1.946 1096.6 49
sink()