treat <- c(0,0,1,1)
tr <- tcctomat(treat)
dose <-
matrix(c(9,13,16,7,12,6,9,10,11,9,10,10,7,9,9,9,8,10,15,4),
ncol=5,byrow=TRUE)
dd <- tvctomat(dose)
y <- restovec(structure(c(6, 4, 0, 0, 3, 6, 1, 1, 1, 5, 0, 0, 0, 4, 0, 1, 0,
13, 0, 3), dim = 4:5))
reps <- rmna(y, ccov=tr, tvcov=dd)
kalcount(y, intensity="log normal", dep="independence",
preg=0.3, pshape=0.1)
kalcount(y, intensity="log normal", dep="frailty", pinitial=0.1,
preg=1, psh=0.1)
kalcount(y, intensity="log normal", dep="serial", pinitial=0.1,
preg=1, pdep=0.75, psh=0.1)
# random effect and autoregression (commented out: AR difficult to estimate)
#kalcount(y, intensity="log normal", dep="frailty", pinitial=0.1,
# pdep=0.5, preg=1, psh=0.1)
# add time-constant variable
kalcount(y, intensity="log normal", pinitial=1, psh=1,
preg=c(0.8,0.11), ccov=treat)
# or
kalcount(y, intensity="log normal", mu=~b0+b1*treat,
pinitial=1, psh=.1, preg=c(0.4,-0.04), envir=reps)
# add time-varying variable
kalcount(y, intensity="log normal", pinitial=1, psh=1,
preg=c(-1,2), ccov=treat, ptvc=0, tvc=dose)
# or equivalently, from the environment
dosev <- as.vector(t(dose))
kalcount(y, intensity="log normal", mu=~b0+b1*rep(treat,rep(5,4))+b2*dosev,
pinitial=1, psh=1, ptvc=c(-1,2,0))
# or from the reps data object
kalcount(y, intensity="log normal", mu=~b0+b1*treat+b2*dose,
pinitial=1, psh=1, ptvc=c(-1,2,0), envir=reps)
# try power variance family
kalcount(y, intensity="log normal", mu=~b0+b1*treat+b2*dose,
pinitial=1, psh=1, ptvc=c(-1,2,0.1), envir=reps,
pfamily=0.8)
Run the code above in your browser using DataLab