### R code from vignette source 'C:/Bendix/teach/NSCE/2022/pracs/Lexis-sol.rnw'

###################################################
### code chunk number 1: Lexis-sol.rnw:23-26
###################################################
library( Epi )
occup <- read.table("http://BendixCarstensen.com/NSCE/R/occup.txt", header=TRUE)
occup


###################################################
### code chunk number 2: Lexis-sol.rnw:30-36
###################################################
oL <- Lexis( entry = list( age=AoE, per=DoE ),
              exit = list(          per=DoX ),
      entry.status = factor( rep("W",nrow(occup)) ),
       exit.status = factor( Xst ),
              data = occup )
summary( oL )


###################################################
### code chunk number 3: Lexis-sol.rnw:42-45
###################################################
oL <- splitLexis( oL, time="age", breaks=seq(0,100,5) )
oL <- splitLexis( oL, time="per", breaks=seq(0,110,5)+1900 )
oL[order(oL$lex.id,oL$age),]


###################################################
### code chunk number 4: Lexis-sol.rnw:49-50
###################################################
table( timeBand(oL,"age","left"), timeBand(oL,"per","left"))


###################################################
### code chunk number 5: Lexis-sol.rnw:56-60
###################################################
FU <- with( oL, cbind(lex.Xst=="D",lex.dur) )
colnames(FU) <- c("D","Y")
Age    <- timeBand(oL,"age","left")
Period <- timeBand(oL,"per","left")


###################################################
### code chunk number 6: Lexis-sol.rnw:64-66
###################################################
FUtab <- xtabs( FU ~ Age + Period )
ftable(FUtab,col.vars=2:3)


###################################################
### code chunk number 7: Lexis-sol.rnw:72-78
###################################################
BC <- subset(oL,per-age>1952 & per-age<1962)
FU <- with( BC, cbind(lex.Xst=="D",lex.dur) )
colnames(FU) <- c("D","Y")
Age   <- timeBand(BC,"age","left")
FUctab <- xtabs( FU ~ Age )
FUctab


###################################################
### code chunk number 8: Lexis-sol.rnw:89-92
###################################################
CumRate <- sum(FUctab[,1]/FUctab[,2]*5)
CumRisk <- 1 - exp( - CumRate )
round( c(CumRate, CumRisk), 3)


###################################################
### code chunk number 9: Lexis-sol.rnw:112-116
###################################################
p.rates <- data.frame( rate=c(100,200,400), Age=c(40,45,50) )
oL$Age <- timeBand(oL,"age","left")
oL <- merge(oL,p.rates)
oL


###################################################
### code chunk number 10: Lexis-sol.rnw:119-122
###################################################
Obs <- with( oL, sum( lex.Xst=="D" ) )
Exp <- with( oL, sum( lex.dur*rate/10^5 ) )
round(c( Observed=Obs, Expected=Exp, SIR=Obs/Exp ), 3)


