### R code from vignette source 'C:/Bendix/teach/NSCE/2022/pracs/prev-trial1-sol.rnw'

###################################################
### code chunk number 1: prev-trial1-sol.rnw:10-12
###################################################
cases <- c(474, 402)
rates <- c(56.3, 47.5) # per 10000 years


###################################################
### code chunk number 2: prev-trial1-sol.rnw:15-16
###################################################
names( rates ) <- c("BetaCarotene","Placebo")


###################################################
### code chunk number 3: prev-trial1-sol.rnw:27-29
###################################################
pyears <- (cases/rates)*10000
pyears


###################################################
### code chunk number 4: prev-trial1-sol.rnw:38-43
###################################################
ratio <- rates[1]/rates[2]
SE.logr <- sqrt(sum(1/cases))
ratio.95low <- ratio/exp(1.96*SE.logr)
ratio.95up  <- ratio*exp(1.96*SE.logr)
cbind(ratio, SE.logr, ratio.95low, ratio.95up)


###################################################
### code chunk number 5: prev-trial1-sol.rnw:67-72
###################################################
diff <- rates[1] - rates[2]
SE.diff <- sqrt(sum(rates^2/cases))
diff.95low <- diff - 1.96*SE.diff
diff.95up <- diff + 1.96*SE.diff
cbind(diff, SE.diff, diff.95low, diff.95up)


###################################################
### code chunk number 6: prev-trial1-sol.rnw:81-85
###################################################
Z <- diff/SE.diff
P <- 1 - pchisq( Z^2, 1 )
test.diff <- cbind(Z, P)
test.diff


###################################################
### code chunk number 7: prev-trial1-sol.rnw:90-93
###################################################
Z <- log(ratio)/SE.logr
P <- 1 - pchisq( Z^2, 1 )
( test.ratio <- cbind(Z, P) )


###################################################
### code chunk number 8: prev-trial1-sol.rnw:97-100
###################################################
tt <- rbind( test.diff, test.ratio )
rownames( tt ) <- c("diff","ratio")
tt


###################################################
### code chunk number 9: prev-trial1-sol.rnw:135-141
###################################################
library( Epi )
D <- cases
Y <- pyears/10000
G <- factor(c("Beta","Plc"))
G <- Relevel( G, 2:1 )
data.frame( D, Y, G )


###################################################
### code chunk number 10: prev-trial1-sol.rnw:144-146
###################################################
mm <- glm( cbind(D,Y) ~ G, family=poisreg )
round( ci.exp( mm, pval=TRUE ), 3 )


###################################################
### code chunk number 11: prev-trial1-sol.rnw:149-151
###################################################
ma <- glm( cbind(D,Y) ~ G, family=poisreg(link=identity) )
round( ci.exp( ma, Exp=FALSE, pval=TRUE ), 3 )