Epi (version 2.34)

mod.Lexis: Fit intensity models to follow-up data in Lexis objects

Description

Modeling intensities based on Lexis objects, exploiting the structure of the Lexis objects where the events and risk time have predefined representations. This allows a simpler syntax than the traditional explicit modeling using glm, gam and coxph. But it is only a set of wrappers fro glm, gam and coxph.

Usage

glm.Lexis( Lx,         # Lexis object	
           from = preceding(Lx,to), # 'from' states
             to = absorbing(Lx)   , # 'to' states
        formula,         # ~ model	
         paired = FALSE, # only the pairwise
           link = "log", # link function
          scale = 1,     # scaling of PY
        verbose = TRUE,  # report what is done?
          … )        # further arguments to glm
  gam.Lexis( Lx,         # Lexis object	
           from = preceding(Lx,to), # 'from' states
             to = absorbing(Lx)   , # 'to' states	
        formula,         # ~ model	
         paired = FALSE, # only the pairwise
           link = "log", # link function
          scale = 1,     # scaling of PY
        verbose = TRUE,  # report what is done?
          … )        # further arguments to glm
coxph.Lexis( Lx,         # Lexis object	
           from = preceding(Lx,to), # 'from' states
             to = absorbing(Lx)   , # 'to' states	
        formula,         # timescale ~ model	
         paired = FALSE, # only the pairwise
        verbose = TRUE,  # report what is done?
          … )        # further arguments to glm

Arguments

Lx

A Lexis object representing cohort follow-up.

from

Character vector of states from which transitions are considered. May also be an integer vector in which case the reference will be to the position of levels of lex.Cst. Defaults to the collection of transient states immediately preceding the absorbing states.

to

Character vector of states to which a transition is considered an event. May also be an integer vector in which case the reference will be to the position of levels of lex.Xst. Defaults to the set of absorbing states.

formula

Model formula describing the model for the intensity(-ies). For glm and gam, the formula should be one-sided; for coxph the formula should be two-sided and have the name of the time-scale used for baseline as the l.h.s.

paired

Logical. Should the states mentioned in to, rep. from be taken as pairs, indicating the only transitions modeled. If FALSE all transitions from any of the states in from to any states in to are modeled.

link

Character; name of the link function used, allowed values are 'log' (the default), 'identity' and 'sqrt', see the family poisreg.

scale

Scalar. lex.dur is divided by this number before analysis, so that you can get resulting rates on a scale of your wish.

verbose

Print information on the states modeled?

…

Further arguments passed on to glm, glm or coxph

Value

glm.Lexis returns a glm object, which is also of class glm.lex, gam.Lexis returns a gam object, which is also of class gam.lex, and coxph.Lexis returns a coxph object, which is also of class coxph.lex. These extra class attributes are meant to facilitate the (still pending) implementation of an update function.

The returned objects all have an extra attribute, Lexis which is a list with names from, to and data; character vectors of state names from which and to which transitions are modeled, respectively, and the name of the Lexis object modeled. Note that it is not the object, only the name of it, which may not be portable.

The glm and gam objects also has a scale element in the list; a scalar indicating the scaling of lex.dur before modeling. The coxph object has a timescale element in the list, a character scalar with the name of the underlying timescale variable.

Details

The glm and gam models are fitted using the family poisreg which is a bit faster than the traditional poisson family. The response variable for this family is a two-column vector of events and person-time respectively, so the predictions, for example using ci.pred does not require lex.dur (and would ignore this) as variable in the newdata. ci.pred returns the estimated rates in units of the lex.dur in the Lexis object, scaled by scale, which has a default value of 1.

The default is to model all transitions into any absorbing state by the same model (how wise is that??). If only from is given, to is set to all states reachable from from, which may be a really goofy model and if so a warning is issued. If only to is given, from is set to the collection of states from which to can be reached directly --- see preceding and its cousins.

Occasionally you only want to model a subset of the possible transitions from states in from to states in to, in which case you specify from and to as character vectors of the same length and then only transitions from[i] to to[i], i=1,2,... will be modeled.

There is no working update functions for these objects (yet).

Strictly speaking, it is a bit counter-intuitive to have the time-scale on the l.h.s. of the formula for the coxph since the time scale is also a predictor of the occurrence rate. On the other hand, calling coxph directly would also entail having the name of the time scale in the Surv object on the l.h.s. of the formula. So the inconsistency is merely carried over from coxph.

See Also

Lexis, cutLexis, mcutLexis, addCov.Lexis, absorbing, transient

Examples

# NOT RUN {
library( Epi )
library( survival )
data( DMlate )

# Lexis object of total follow-up
mL <- Lexis( entry = list(age=dodm-dobth,per=dodm),
              exit = list(per=dox),
       exit.status = factor(!is.na(dodth),labels=c("Alive","Dead")),
              data = DMlate )

# Cut follow-up at start of insulin use
cL <- cutLexis( mL, cut = mL$doins,
              timescale = "per",
              new.state = "Ins",
       precursor.states = "Alive" )

# Split follow-up on age-axis
system.time( sL <- splitLexis( cL, breaks=0:25*4, time.scale="age") )
summary( sL )

# glm models for rates based on the time-split dataset by insulin and sex

# Proportional hazards model with insulin as time-dependent variable
# - uses the defaul of modeling all transitions from both transient
# states ("Alive" and "Ins") to the absorbing state ("Dead"). 
mt <- glm.Lexis( sL, formula= ~ sex + lex.Cst + Ns(age,knots=c(15,3:8*10)) )

# prediction of mortality rates from "Alive" with and without PH assumption
nA <- data.frame( age=40:70, sex="M", lex.Cst="Alive" )
nI <- data.frame( age=40:70, sex="M", lex.Cst="Ins" )
matshade( nA$age, cbind( ci.pred(mt,nA),
                         ci.pred(mt,nI) )*1000, plot=TRUE,
          lwd=3, lty=1, log="y", col=c("black","blue","red"),
          xlab="Age", ylab="Mortality per 1000 PY" )
 
# gam models takes quite some time so we leave it out
# }
# NOT RUN {
mt.gam <- gam.Lexis( sL, "Dead", ~ sex + lex.Cst + s(age), scale=1000 )
        
# }
# NOT RUN {
# Fit a Cox model with age as baseline time scale and insulin as time-dependent
mt.cox <- coxph.Lexis( sL, c("Alive","Ins"), "Dead", age ~ sex + lex.Cst )

# Pretty much the same results for regression paramters as the glm:
  ci.exp( mt    , subset="ex" )
# ci.exp( mt.gam, subset="ex" )
  ci.exp( mt.cox, subset="ex" )
# }