source("http://www.uvm.edu/~rsingle/stat221/data/scripts-221.R")

#Pulse data 
# Exercise: 1='little'; 2='moderate'; 3='high'
# Rest:     resting pulse
# Active:   pulse after brisk walk up & down stairs 3 times
# Gender:   0='male', 1='female'

 dat <- otherdata("pulse.txt")
 dat$Sex <- NA
 dat$Sex[dat$Gender==0] <- "M"
 dat$Sex[dat$Gender==1] <- "F"

#simple correlation between Active & Resting pulse
 plot(dat$Active,dat$Rest)
 cor(dat$Active,dat$Rest)
 cor(Active~Rest, data=dat)  #NOTE: this syntax does not work
 with(dat, cor(Active,Rest)) #      But, with() is an option
 
#test Ho: rho=0 
 cor.test(dat$Active, dat$Rest)  
 mod1 <- lm(Active~Rest, data=dat)
 summary(mod1)                 #Note: T-test for beta_1 gives same value
 
 
#test Ho: rho=.5 (using Fisher's Z-transformation)
#Suppose, based on previous data, our Ho is that the correlation is .5
 r   <- cor(dat$Active,dat$Rest) #sample estimate
 rho <- .5                       #value under Ho
 n <- length(dat$Active)
 z = ( 0.5*log((1+r)/(1-r)) - 0.5*log((1+rho)/(1-rho)) )/(1/sqrt(n-3))
 
 
#test Ho: rho=0 stratified by Sex
 cor.test(dat$Active[dat$Sex=="M"],dat$Rest[dat$Sex=="M"])
 cor.test(dat$Active[dat$Sex=="F"],dat$Rest[dat$Sex=="F"]) 
 
#test Ho: rho_male = rho_female 
 n1 <- length(dat$Active[dat$Sex=="M"])
 n2 <- length(dat$Active[dat$Sex=="F"]) 
 r1 <- cor(dat$Active[dat$Sex=="M"],dat$Rest[dat$Sex=="M"])
 r2 <- cor(dat$Active[dat$Sex=="F"],dat$Rest[dat$Sex=="F"])
 r1; n1
 r2; n2
 
 W1 <- 0.5 * log( (1+r1)/(1-r1) )
 W2 <- 0.5 * log( (1+r2)/(1-r2) )
 
 diff <- W1 - W2
 SE.diff <- sqrt( 1 / (n1 - 3) + 1 / (n2 - 3))
 z <- diff / SE.diff
 
 
 
# Smoking Status
 r1 <- cor(dat$Active[dat$Smoke==0],dat$Rest[dat$Smoke==0])
 r2 <- cor(dat$Active[dat$Smoke==1],dat$Rest[dat$Smoke==1])
 n1 <- length(dat$Active[dat$Smoke==0])
 n2 <- length(dat$Active[dat$Smoke==1])
 r1; n1
 r2; n2