### R code from vignette source 'tabulate-sol.rnw'

###################################################
### code chunk number 1: tabulate-sol.rnw:2-3
###################################################
library( Epi )


###################################################
### code chunk number 2: tabulate-sol.rnw:21-23
###################################################
lung <- read.table( "../data/lung5-M.txt", header=T )
head(lung) 


###################################################
### code chunk number 3: tabulate-sol.rnw:25-26
###################################################
attach( lung )


###################################################
### code chunk number 4: tabulate-sol.rnw:40-44
###################################################
D_table_nice <- stat.table(index=list(A,P), sum(D), data=lung, margin=T )
print( D_table_nice, digits=c(sum=0) ) 
Y_table_nice<-stat.table( index=list(A,P), sum(Y), data=lung, margin=T )  
print( Y_table_nice, digits=c(sum=2) )


###################################################
### code chunk number 5: tabulate-sol.rnw:53-56
###################################################
R_table_nice <- stat.table( index=list(A,P), list(Rate=ratio(D,Y,100000)),
                            data=lung, margin=T )  
print( R_table_nice, digits=c(sum=2) )


###################################################
### code chunk number 6: tabulate-sol.rnw:60-64
###################################################
cat( "tabulate-sol" )
D_table <- with( lung, tapply( D, list(A,P), sum ) )
Y_table <- with( lung, tapply( Y, list(A,P), sum ) )
R_table <- D_table/Y_table*(10^5)


###################################################
### code chunk number 7: tabulate-sol.rnw:73-74
###################################################
rateplot( R_table, which=c("AP"), ann=TRUE )


###################################################
### code chunk number 8: tabulate-sol.rnw:78-79
###################################################
rateplot( R_table, which=c("AC"), ann=TRUE )


###################################################
### code chunk number 9: tabulate-sol.rnw:83-84
###################################################
rateplot( R_table, which=c("PA"), ann=TRUE )


###################################################
### code chunk number 10: tabulate-sol.rnw:88-89
###################################################
rateplot( R_table, which=c("CA"), ann=TRUE )


###################################################
### code chunk number 11: RATEPLOT1
###################################################
par( mfrow=c(2,2) )
rateplot( R_table ,ann=TRUE)


###################################################
### code chunk number 12: tabulate-sol.rnw:121-126
###################################################
# age-specific rates remain still the same as in period 1943
R_table_no_change <- matrix( R_table[,1], dim(R_table)[1], dim(R_table)[2] )
colnames( R_table_no_change ) <- colnames( R_table )
rownames( R_table_no_change ) <- rownames( R_table )
R_table_no_change


###################################################
### code chunk number 13: RATEPLOT2_no_change
###################################################
par( mfrow=c(2,2) ) 
rateplot( R_table_no_change, log.ax="" )


###################################################
### code chunk number 14: tabulate-sol.rnw:149-157
###################################################
#age-specific rates are only influence by period
step <- 2
change_p <- matrix( rep(seq(1,11*step,step),10),10,11, byrow=T )
change_p
R_table_p <- R_table_no_change + change_p
colnames( R_table_p ) <- colnames( R_table )
rownames( R_table_p ) <- rownames( R_table )
R_table_p


###################################################
### code chunk number 15: RATEPLOT2_p_change
###################################################
par( mfrow=c(2,2) )
rateplot( R_table_p, log.ax="" ) 


###################################################
### code chunk number 16: tabulate-sol.rnw:178-191
###################################################
#age-specific rates are only influence by cohort
nr <- nrow( R_table )
nc <- 10
p <- c( rep(NA,nc ), R_table[,1] )
np <- length( p )
R_table_c <- cbind(p[(np-nr+1):np],p[(np-nr):(np-1)],p[(np-nr-1):(np-2)],
p[(np-nr-2):(np-3)],p[(np-nr-3):(np-4)],p[(np-nr-4):(np-5)],
p[(np-nr-5):(np-6)],p[(np-nr-6):(np-7)],p[(np-nr-7):(np-8)],
p[(np-nr-8):(np-9)],p[(np-nr-9):(np-10)]
)
colnames( R_table_c ) <- colnames( R_table )
rownames( R_table_c ) <- rownames( R_table )
R_table_c


###################################################
### code chunk number 17: RATEPLOT2_c_change
###################################################
par( mfrow=c(2,2) )
rateplot( R_table_c, log.ax="" )