As an example, I'll artificially divide the dipper data into 3 groups: "Y","J","A"
- Code: Select all
data(dipper)
# assign some artificial initial ages
dipper$ageclass="Y"
dipper$ageclass[100:150]="J"
dipper$ageclass[151:250]="A"
dipper$ageclass=factor(dipper$ageclass)
#process data and assign initial ages;note that the ageclass factor is ordered
#alphabetically so it is A,J,Y and the initial ages are then 2,1,0
dp=process.data(dipper,groups=c("sex","ageclass"),age.var=2,initial.age=c(2,1,0))
ddl=make.design.data(dp)
# look at some of the design data and see that the ages at
# each time are adjusted for each group according to their initial age and
# the time intervals
head(ddl$Phi[ddl$Phi$ageclass=="Y",])
head(ddl$Phi[ddl$Phi$ageclass=="J",])
head(ddl$Phi[ddl$Phi$ageclass=="A",])
Next we can fit a model and get predictions. Note this is nonsense data and I haven't made any attempt to do model selection or evaluate model fit. This is just an example. There are other examples like this in the documentation and workshop notes. Make sure you have the recent workshop notes as well.
- Code: Select all
model=mark(dp,ddl,model.parameters=list(Phi=list(formula=~age)))
Phiestimates=summary(model,se=TRUE)$reals$Phi
unique.estimates=unique(subset(Phiestimates,select=c("Age","estimate","lcl","ucl")))
with(unique.estimates[order(unique.estimates$Age),],
{
plot(Age,estimate,ylab="Survival",type="b",ylim=c(0,1))
lines(Age,lcl,lty=2)
lines(Age,ucl,lty=2)
})