mydat <- read.csv(file="http://sites.williams.edu/bklingen/files/2012/02/smoking.csv")
mydat
rownames(mydat) <- mydat[,1]
mydat <- mydat[,-1]
mydat
## marginal table:
colSums(mydat)
## Preparations for CMH test:
mydat1 <- mydat[,c(1,3,2,4)]
my3dtable <- array(unlist(t(mydat1)), dim=c(2,2,4),
dimnames=list(Smoking=c("Yes", "No"),
              Survival=c("Yes","No"),
              Age=c("18-34", "35-54", "55-64", "65+")))
mantelhaen.test(my3dtable, correct=FALSE)
###############
table.long <- data.frame(ftable(my3dtable))
table.long
index <- table.long$Survival=="Yes"
table1 <- cbind(table.long[index,-2], table.long[!index,4])
names(table1)[3:4] <- c("S","F")
table1
logistic.fit <- glm(cbind(S,F)~ Age + (Smoking=="Yes"), family="binomial", data=table1)
summary(logistic.fit)
anova(logistic.fit, test="Chisq")
#### Estimating Common Odds Ratio
confint(logistic.fit)
exp(c(-0.7984353, -0.1087374))
exp(-0.45)
## Is homogeneous association plausible?
logistic.fit1 <- glm(cbind(S,F)~ Age + (Smoking=="Yes") + Age*(Smoking=="Yes"), family="binomial", data=table1)
summary(logistic.fit1)
anova(logistic.fit, test="Chisq")