Multi-state models for biomedical research : new contributions in statistical modelling, software development, and applications

UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Departamento de Estatística e Investigación Operativa MULTI-STATE MODELS FOR BIOMEDICAL RESEARCH: NEW CONTRIBUTIONS IN STATISTICAL MODELLING, SOFTWARE DEVELOPMENT, AND APPLICATIONS Luís Filipe Meira Machado July 2005 uárez, Profesora Titular del Departamento de Estatística e vestigación Operativa de la Universidade de Santiago de Compostela, y D. Jacobo de ha sido alizada bajo su dirección por Luís Filipe Meira Machado, estimando que el interesado encuentra en condiciones de optar al grado de doctor en Ciencias Matemáticas, por lo que solicitan que ésta sea admitida a trámi Santiago de Compostela, 12 de Julio de 2005. D. Carmen Cadarso S In Uña Álvarez, Profesor Titular del Departamento de Estatística e Investigación Operativa de la Universidade de Vigo, INFORMAN que la memoria titulada “Multi-state models for biomedical research: new contributions in statistical modelling, software development, and applications” re se te para su lectura y defensa pública. D. Luís Filipe Meira Machado, profesor Assistente de la Universidade do Minho (Portugal), DECLARA que ha realizado la Tesis Doctoral titulada “Multi-state models for biomedical research: new contributions in statistical modelling, software development, nd applications” y solicita su admisión a tramite en el Departamento de Estatística e vestigación Operativa de la Universidade de Santiago de Compostela. a In To my family ix Summary This thesis contains some co ntributions for statistical models studying the disease d, focusing on possible advantages and disadvantages for each method. Softwa progression. Methods developed in this thesis are largely motivated by the applications to the medical sciences. Disease progression can be well described using multi-state models. These models may be considered as a generalization of the survival process where several (intermediate) events occur successively over time. Multi-state models offer a better understanding of the process of the illness, leading to a better knowledge of the evolution of the disease over time. Issues of interest include the estimation of progression rates, assessing the effects of individual risk factors, survival rates or prognostic forecasting. The influence of these intermediate events in survival is often analyzed using the Cox time-dependent regression model. This thesis contains a comprehensive review of several multi-state models for studying disease progression. Differences between these approaches and the time-dependent Cox regression model are discusse re implementation of these models has been developed in the form of an R library. Traditionally, statistical methods for analyzing such models depend on the Markov assumption, for which future evolution only depends on the current state. By ignoring disease history behaviour, these models may carry severe limitations which can make the model inappropriate. One alternative approach is to use a semi-Markov assumption, in which the future of the process does not depend on the current time, but only on the duration in the current state. Finally, we developed a new non-Markov approach for which the future of the process depends on the current time, but also on the xvi 3.5 Comparison of probability survival curves………………………………... Multi-state model for Psoriatic Arthritis data.…………………………….. 71 3.6 74 4.1 Curves obtained for setting 3 with 25 s = , 200 = N and 32% of censored observations for Aalen-Johansen estimator (solid line) and non- 120 4.2 s-death model for PROVA trial……………………………………… 124 4.3 Estimated transition probabilities for Aalen-Johansen estimator (solid 5.1 Extended illness-death model…..………………………………………….. 141 5.2 raphical output for the Cox Markov model. Stanford Heart 148 5.3 ulti-state model for Stomach Cancer data….……………………………. 149 5.4 Graphical output for the homogeneous Markov model. Stomach Cancer data……………………………………….….…………………………….. 151 Markovian estimator (solid bold line)……………………………………... Illnes line) and non-Markov estimator (bold solid line)…………………………. 125 4.4 The four-state progressive model………………………………………….. 134 G Transplantation data………..…..………………………………………….. M xvii List of Tables n ssion odel, according to sample size, over 500 replicates……………………... 57 3.2 ox time-dependent model for a single sample with size 3.1 Mean estimate of the covariate effect for the homogeneous Markov model, Cox Markov model and for the Cox time-depe dent regre m Covariate effect in homogeneous Markov model, Cox Markov model and C500 = n.............. 58 3.3 orresponding standard errors. Stanford Heart Transplantation data……... 61 3.4 tanford Heart Transplantation data………………………………………. 63 3.5 ransplantation data ………………………………………………………. 64 3.6 66 3.7 omparison of observed and predicted percentages of patients in each of ransplantation data………………………………………………………. 72 3.8 on-homogeneous model. Estimated transition rates and hazard rates on 73 3.9 bserved transitions in the Psoriatic Arthritis data……….……………….. 74 4.1 Time-dependent Cox regression models. Effect estimates and c Estimated effects on the mortality intensity in transplanted patients. S Final Cox Markov model for all transitions. Stanford Heart T Multi-state homogeneous Markov model. Estimated transition rates and hazard rates. Stanford Heart Transplantation data………………………... C the three states for Cox Markov Model (CMM), Homogeneous Markov Model (HMM) and Non-Homogeneous Model (NHM). Stanford Heart T N transition intensities. Stanford Heart Transplantation data………………... O Summary statistics measuring integrated bias, integrated variance and integrated mean square error……................................................................. 118 4.2 Estimates of integrated absolute bias, integrated variance and integrated mean square error for setting 1, according to fixed value s, censoring and sample size……………………...…………………………………………. 121 4.3 Estimates of integrated absolute bias, integrated variance and integrated mean square error for setting 2, according to fixed value s, censoring and sample size………………………………………………............................ 122 xviii 4.4 Estimates of integrated absolute bias, integrated variance and integrated mean square error for setting 3, according to fixed value s, censoring and 123 5.1 put data file for the Stanford Heart Transplantation data……………….. 142 5.2 tanford Heart Transplantation data……………………...……………….. 145 5.3 ransplantation data………...……………………...……………………… 146 5.4 ransplantation data…….………...……………………...……………….. 147 5.5 ata………………...…….………...……………………...……………….. 150 5.6 Output for the homogeneous Markov model. Stomach Cancer data………………...…….………...……………………...……………….. 151 sample size………………………………………………………………… In Sample of the output for the time-dependent Cox regression model. S Sample of the output for the homogeneous Markov model. Stanford Heart T Sample of the output for the Cox Markov model. Stanford Heart T Sample of the output for the Cox Markov model. Stomach Cancer d xix Chapter 1 Introduction Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 2 1.1 Survival analysis This chapter is concerned with studying the time, T, between a well-defined time origin and a subsequent event. In biomedical applications, this is known as survival analysis, and the times may represent the survival time of a living organism or the time until a disease is cured. Although we are mainly interested with the application of these methodologies to data from epidemiological and clinical studies, these methods can also be applied to data from different areas, such as the social sciences, economics and engineering. It is common in survival analysis for some of the data to be censored. That is, the time elapsed before the occurrence of the event is not known; it is only known that this time exceeds some value (right-censoring) or that it is inferior to some value (leftcensoring). Some observations in survival analysis may not be observed for the full time to the event. For example, in a clinical study, some patients will have survived to the end of the study, while for others the survival status at the end of study is unknown because there is either loss to follow-up or they have died due to other causes unrelated to the study. In these studies, complete survival information will be available for some patients, while for others we do not know the exact survival time but rather the time elapsed from the entry in the study until last known survival time. In such cases, we say that the observation is right-censored. Sometimes the survival time is less than that which is observed. This might occur if subjects are observed at intermittent visits, and only then it is observed that the event had occurred some time before. A typical example is the study concerning the time to recurrence of a tumour. When a patient is found to have a recurrence, the actual Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 3 of time of recurrence (since the operation) is less than the observed, and the observation is said to be left-censored. Another form of censoring is interval-censoring, which is commonly present in clinical studies where the observations are made at intermittent visits. In these cases, the exact survival time is not observed; it is only known that the event occurred within an interval of time. The observation is said to be interval-censored. Censored observations cannot be ignored since they carry important information about the survival, but because of censoring, special methods may be required to apply some model and to analyze the data. Generic survival data is in the form ( ) , δ Y ) C , where YT and () ( min ,= δ =≤I TC for which C are the censoring time associated with T . With random censoring, we will assume that the time of censoring and the survival time T are independent. This assumption specifies that the probability of an event occurring in the erval [[ ,tt dt+, given that no event occur until time t, is the same for individuals in general and those whose censoring time is s greater than t; that is, small int ( ) ( ) ,tTtdtTt t≤TtdtTtCt≤<+ ≥= <+ ≥ ≥PP . In addition to censoring, observations may be incomplete because of lefttruncation, that is, the observation is only considered if some condition bearing on subject history is met. For example, subjects may not be followed from time 0, but only from a later entry time, conditional on being alive at this entry time. Basic concepts Let T be a random non-negative variable representing the individual survival time from a homogeneous population. Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 4 Ft l function, , is the probability that the survival time is greater than som Ft Assume T is an absolutely continuous variable with underlying density function . The distribution function of T is given by () ft () ( ) () =≤= t T t fudu . [1.1] 0 The surviva ∫ P () St e value t, () ( ) 1=>=−P St T t . [1.2] The hazard function, () t () α , is the probability that an individual “dies” at some time t, conditional that he survived until that time. Thus, the hazard function represents the instantaneous probability that the event will occur at a given time t, and can be written as () 0 α → ≤<+ ≥ ⎡⎤ ⎣⎦ =P| lim dt tT tdt tdt or , representing the probability that an event will T t , [1.3] () () α ≅≤<+ ≥P| tdt tTtdtTt occur in the small time interval [ [ , tt dt + , given that no event occur until time t. The hazard function is also called risk, the failure rate, the mortality rate, or intensity. The integrated or cumulative hazard function is defined as Wse functions: () () t At udu . [1.4] 0 α =∫ e can now obtain some useful relationships between each one of the () () () ft tSt α =, [1.5] , [1.6] () () 0 exp t St uduα ⎛⎞ =− ⎜⎟ ⎝⎠ ∫ Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 5 () () ( ) =−exp St At , [1.7] and () () () 0 αα ⎛⎞ =× Because of relations [1.1]-[1.8], the distribution of T can be univocally specified through any one of the follo − ⎜⎟ ⎝⎠ ∫ exp t ft t udu . [1.8] wing functions: ( ) ( ) ( ) ( )() and α ⋅ ⋅⋅,,, fFS ⋅ ⋅ A . A basic task in survival analysis is the estimation of survival in the presence of red observations in a sample of dimension n with r for survival is the The Kaplan-Meier estimator censoring. If there are no censo observed survival times, , the most natural estimato 12 , ,..., n tt t empirical estimator, given by  () () 1= ni i n which is unbiased and is a consistent estimate of 1 n =≥ ∑I St t t , ( ) St . Whenever observations are censored, alternative estimators are needed. Suppose th w Kaplan and Meier (1958), obtained a nonparametric estimate of the survival function, called product-limit, which is the generalization of the empirical estimator for censored data. ere are n subjects with observed survival times, 12 , ,..., n tt t . Some of these observations may be right-censored or there may be tied observations. If we now assume that e have k events () kn ≤ and nk − censored observations, then we may write () ( ) ( ) 12 ... k tt t <<< for the t tim s (order statistics ed bjects at risk (alive and uncensored) just prior to k even e ) arranged in increas order. Let i n denote the number of su Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 6 nt at that e time. the event time () i t , and i d the number of subjects who experienced the eve sam The Kaplan-Meier estimate of ( ) St , is defined by,  () 1= ⎛⎞ ⎝⎠ j i ii d n ( ) 1 j tt + ≤< and 12, ,..., 1 =− St , [1.9] for j t ⎜⎟ ∏ () j k =. The Kaplan-Meier function is a consistent estimator of () St which, under general conditions, can be considered a nonparametric maximum likelihood estimator. This estimator is a step function which steps d n at each event time (only), with ow  () 1= t for () 1 tt <. T S his estimator can be used in the context of survival data subjected to inde presence of lefttruncation, the estimation [1.9] is pe subjects who have entered the study before time The variance of the Kaplan-Meier estimator is estimated by pendent right-censoring and for left-truncation as well. In the rformed as earlier, but now the i n is the number of () i t and are still in the study just prior to () i t (Borgan, 1998b).  ()  () () 22i d tSt nn d σ =− ∑, [1.10] 1 j ii i i= for jj ttt ≤< and 12, ,..., j k = () ( ) 1+. This estimator is known as Greenwood’s formula The cumulative hazard may be estimated by the Nelson estimator (Nelson, 1969), originally introduced to check graphically the fit of parametric models. This estimator is also known as Nelson-A (Greenwood, 1926). alen estimator and takes the form,  () 1 j i i i d n = At =∑, Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 13 corresponding to the event of interest. No further observations are nec s are illustra in survival is often important in the patient’s outcome and can be handled using the Cox time-dependent regression model. In such cases, the event of interest is considered as the main event, while intermediate events are often included as time-dependent covariates in a proportional hazards model. The process ends when the patient enters the (absorbing, “dead”) state essary after this event occurs. Thus, in survival analysis only two states are considered, and the event of interest is the passage from one state to another. As mentioned before, the proportional hazards model has the advantage of being easily interpreted and available in the majority of statistical packages. Multi-state models may be considered as a generalization of the basic framework dealing with survival data where survival is the ultimate outcome of interest but where intermediate (transient) states are identified. In contrast to survival data, in these models, a sequence of events is observed, leading to more than one observation per individual. In medicine, the states might be based on clinical symptoms (e.g. bleeding episodes), biological markers (e.g. CD4 T-lymphocyte cell counts; serum immunoglobulin levels), some scale of the disease (e.g. stages of cancer or HIV infection) or a non-fatal complication in the course of the illness (e.g. transplantation). A change of state is called a transition, or an event. States can be transient or absorbing, if no transitions emerge from the state (for example, death). The state structure identifies the states and the transitions allowed between states. There is no unique structure for the series of states. Some of the commonly-used state structure ted in Chapter 2. Multi-state models have several advantages over the Cox proportional hazards model with time-dependent covariates. First, they provide a comprehensive view of the disease process, putting the incomplete information to more efficient use when portions Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 14 regression model for example. In fact, it is very unlikely that the ris interest in statistical methods for studying disease progression. Some diseases that h ve been studied using multi-state models include: HIV infection (Lagakos et al., 1988 Longini et al., 1989 and 1991 ; Gentleman et al., 1994; Satten and Longini, 1996; Aalen et al., 1997), breast cancer (Duffy and Chen, 1995; Chen et al., 1996; Pérez-ócon et al., 2001), cirrhosis of the liver (Kay, 1986; Andersen et al., 2000), leukaemia (Klein et al., 1984 and 1994; Andersen et al., 1999; Keiding et al., 2001; Chevret et al., 2000), asthma (Saint-Pierre et al., 2003), alzheimer’s (Commenges et al., 2004), transplantations (Hansen et al., 1994; Klotz and Sharples, 1994; Jackson and Sharples, 2002), diabetes (Andersen, 1988), diabetic retinopathy (Andersen, 1991; Marshall and Jones, 1995), malaria (Gottschau and Hogh, 1995) and multiple sclerosis (Esbjerg et al., 1999). Multi-state models have been also of the history of an individual’s illness are known (efficiently handling heavily censored data). Secondly, in this framework, the so-called transition intensities provide the instantaneous hazard for movement out of one state into another. These intensity functions can be used to determine the mean sojourn time in a given state of illness, the number of individuals in different states at a certain moment, and survival proportions in each state. Finally, covariates in the transition intensities can also explain differences in the course of the illness among the population. In this way, multi-state models can reveal that different covariates affect different transitions, which would not be possible with other models, the Cox k of death in patients that received different treatments would be the same. Furthermore, the prognostic factors associated with the risk of death can be different in these groups of patients. In conclusion, multi-state models dynamically evaluate the patient’s progress of disease, depending on the occurrence of intermediate events. Recent years have seen an increasing a ; Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 15 used in other fields, such as signal processing (Juang and Rabiner, 1991) and reliability ngly used for modelling survival data in the case where interest is to study how covariates affect survival. However, in spite f the advantages mentioned previously, there are still few applications in the literature sing multi-state models. analysis (Su et al., 2000). For example, a search on PUBMED using the terms “multistate model” and/or "multi-state model" takes us to 86 published articles. Figure 1.1 shows the evolution of the use of multi-state models up to September 2005. This figure shows that multi-state models have been increasi several events occur successively over time. Figure 1.1. Number of published articles in recent years using multi-state models. The Kaplan-Meier estimate of the survival function is today well-known and 45 14 24 10 0 10 15 25 1987 Number of Arti es 4 15 5 20 Before 1988-1990 1991-1993 1994-1996 1997-1999 2000-2002 2003-2005 Period cl extensively used in the medical sciences; the proportional hazards model is also widely used whenever the o u Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 16 plications of multi-state models have been done with and Lawless, 1985; Kay, 1986). While the researc t of them are self-contained. Moreover, 1.3 Applications ov (when we The earlier ap homogeneous Markov models (Kalbfleisch h on such models has been increasing in recent years, there are still few applications on non-homogeneous models or non-Markov models. Furthermore, when using multi-state models, different patterns of censoring and the different state structure have to be considered. For these purposes, flexible estimation methods have been proposed along with software. Even so, we feel that the available software for these models is still scarce, difficult to use, and mos although the use of such models has been increasingly suggested by several authors, we feel that the applications of such models are still insufficient in the literature. In conclusion, and even with the limitations mentioned above, it is our belief that because of the development of more flexible estimation methods, and with the increasing power of numerical methods and computers, that multi-state models will be more and more applied in the medical sciences. A wide literature on these models includes books by Andersen et al. (1993) and Hougaard (2000). Throughout this thesis we will use some databases to illustrate all presented methodologies. The Stanford Heart Transplant dataset is one of the most widely used examples in survival analysis. An additional advantage of this data is its wide availability to researchers. Because Stanford Heart Transplant data is Mark Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 17 herapy versus that of no treatment on risk of bleeding and survival. This database shall use our own databases. A sample of 419 su Portugal. ansplant study began in October 1967. The available data in Crowley and Hu’s article (Crowley and Hu, 1977) covers the period up to April 1, 1974. Some patients died befor found. Of the 103 patients, 69 consider time since entry in study as the time scale), we may use this database to compare some reviewed methods in Chapter 3. One of the major aims of this thesis is to study non-Markov multi-state models. To illustrate proposed methodologies, we use data from a clinical trial of cirrhosis of the liver. PROVA clinical trial was conceived to evaluate the effect of propanolol and/or sclerot was previously used by Andersen et al., (2000), who concluded (among other things) that the Markov assumption is not valid in this case. In addition to these two databases we bjects suffering from Psoriatic Arthritis from Galicia (Spain), and a database provided by the Portuguese Oncology Institute of Oporto to study the stomach cancer treatment in the north of Stanford Heart Transplantation data The Stanford Heart Tr e an appropriate heart was received a heart transplant. The total number of deaths was 75; the remaining 28 patients contributed with censored survival times. For each individual were recorded the following: an indicator of his final vital status (censored or not), the survival times from the patient’s entry into the study (in days), and a covariate vector including age at acceptance, year of acceptance, previous surgery (coded as yes1 =; no0 =), and transplant (coded as yes1 = ;no0 =). The covariate transplant is the only time-dependent covariate, while the other covariates included are fixed. Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 18 state models to investigate e effect of the transplantation on survival. For the Stanford data, it i the hazard of death before and ter th l conceived to evaluate the effect s entered the study. Since their entry to the study, some patients had bleeding episodes, while others did not. The occurrence of these intermedi patient outcome and can be include umber of deaths was 75. In the context of multi-state modelling, we may consider the covariate ‘transplant’ as an associated state of risk, and then use multith s interesting to compare af e transplantation. We may also explore the potential fixed covariate effects in each of the transitions. Among others, Turnbull, Brown and Hu (1974), Mantel and Byar (1974) and Crowley and Hu (1977) have studied the Stanford Heart data to evaluate the influence of transplantation on hazard. PROVA clinical trial PROVA was a Danish multi-centre clinical tria of propanolol and/or sclerotherapy versus that of no treatment on risk of bleeding episodes and survival in patients with cirrhosis of the liver. This data cover the period between 1985 and 1989. In this period, 286 patient in ate events may affect the d as a transient state in a multi-state model. Before 1 January 1990, only 50 of the 286 patients developed bleeding episodes and from these, 29 died. The total n Psoriatic Arthritis Due to large number of people affected by Psoriasis and Psoriatic Arthritis (PsA), there is much demand for information on these chronic diseases. Here, we consider a sample of 419 subjects from Galicia (Spain). All these patients were Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 19 esearch goals will determine the time scale to use. For example, in Psoriatic Arthrit We consider a sample of 314 subjects drawn from the database provided by the ortuguese Oncology Institute of Oporto. This study is concerned with stomach cancer eatment in the north of Portugal. val times are recorded from the surgical removal of e tumour. Along with the survival time, the vital status and a covariate vector ing age, sex, the initial state (surgical intervention time) through one of two mutually exclusive states (‘m diagnosed as suffering from PsA, and then some treatment was considered. This data was obtained retrospectively and conceived to investigate the age of developing PsA. Here, we are interested in investigating the evolution of the disease from birth up to the current state. For this purpose, a series of stages was considered: Psoriatic, Arthritis (inflammatory arthritis), PsA, Remission and Active. The r is data, since our interest focuses on the age of developing of Psoriatic Arthritis, then age is the appropriate choice for a time scale. Stomach cancer study P tr In this study, patient’s survi th includ and clinical staging were recorded for each patient. Following removal of the tumour, 41 patients have a recurrence of stomach cancer while 53 patients fall ill as a result of metastasis to other solid organs. In this context, patients may pass from etastases’ and ‘recurrence’) to an absorbing state (death). In this study, our main goals include: (a) to obtain some informative events such as the mean sojourn time, transition rates and survival rates for each state; (b) to compare the mortality intensities; (c) to explore the potential fixed covariate effects in each of the transitions. Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 20 we present some standard multi-state models, we conside 1986); and (iii) Markov models with piecewise constant intensities (Pérez-Ocón et al., transition intensities, regression coefficients, goodness-of-fit, etc). These checks are also It should be notice that patients only enter our study if the surgical procedure was considered curative, that is to say, if after the surgery there is no evidence of illness. During the 5-year follow-up period, the number of deaths was 68; 38 of which occurred from state ‘Metastases’ while the remaining 30 were from the recurrent state. 1.4 Outline of the thesis. Chapter 2 provides a deeper review of multi-state models. We introduce these models as stochastic processes; r some censoring patterns, and we discuss the most common simplifying assumptions. In Chapter 3, we review some possible approaches following the methodology of the multi-state Markov models. Specifically, we review the following models: (i) Cox Markov models (Andersen et al., 2000), where the transition intensities are modeled using different Cox models separately; (ii) time-homogeneous Markov models (Kay, 2001). The reviewed methods are illustrated with the data of Stanford heart transplant study, providing some guidance about the use of these methodologies for studying the evolution of the disease. We also use this database to discuss hypothesis testing procedures and methods for model checking, such as the time-homogeneous assumption or the Markov assumption. Furthermore, when examining the data, one is often interested in testing several hypotheses about the model (e.g. hypothesis about the Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 21 , focusing on possible advant order to illustrate possible benefits of usin illustrated (discussed) using Stanford database. Differences between these approaches and the time-dependent Cox regression model are discussed ages and disadvantages for each method. Through this illustration, we show that multi-state models can provide new insights while confirming the results obtained when using the Cox regression model. Special attention must be paid to evaluating the covariate effect on the hazard. As mentioned before, the Cox model provides only constant estimates of the covariate effect over the whole study period. To avoid this problem, we used spline function methods (P-splines) to obtain a dynamic Cox model. Furthermore, we introduce these spline methods in multi-state models to find out possible non-linear covariate effect on the transition intensities. The use of these methods in the scope of multi-state models is new. While Multi-state Markov models provide non-biased estimates of the importance and ability of covariates to predict the course of the illness, these issues cannot be fully explored using Cox models alone. In g multi-state models, simulation studies were undertaken. Through these simulation studies we provide a useful description of why the Cox model can be misleading (with difficult interpretation) when used in the presence of time-dependent covariates. In this way, we show that a multi-state approach can provide a more conclusive analysis of the covariate effect than the Cox model. Because, in some cases, Markov models do not fit satisfactorily, alternative models must be considered. These models are considered in Chapter 4. One alternative approach is to use a semi-Markov assumption. In addition to the Cox semi-Markov models, we propose a new approach for illness-death models based on less restrictive assumptions than those based on the Markov property. To date, we believe this to be the Chapter 1. Introduction _____________________________________________________________________________________________________________________________________________ 22 me clinical studies, a model with th model is preferable. In several cases, a homogeneous model will be satisfactory, while in others not. Furthermore, possible comparisons between different multi-state models are rather difficult because each of the current programs requires its own input data structure. In addition, most of the available programs only provide the regression parame rs estimates and do not supply graphical output for the survival estimates nor for transition probabilities estimates. To overcome these difficulties, we have developed a user-friendly R library, called tdc.surv, which can be used to fit almost all the proposed models (the non-Markov model in Chapter 4 is the exception). Advantages of this software include the same data input for fitting the different models while providing the corresponding numerical and graphical outputs obtained: parameter estimates with standar errors for the covariates; transition rates; survival estimates; transition probabilities estimates; and flexible p-spline hazard ratio estimates for continuous covariates (Eilers and Marx, 1996). Moreover, users are able to include any number of first non-parametric modelling approach completely free from the Markov assumption. Along with theoretical properties, we have presented estimated transition probabilities and cumulative incidence functions for such models. A simulation study is performed to compare both Markovian and non-Markov approaches under a variety of scenarios. For illustration purposes, we have applied our methodology to data from the PROVA clinical trial, described above. In this chapter, we also show that these methods can be extended to more complex multi-state models. As a result of this work, a manuscript is being written and will be submitted for publication. As previously mentioned, one main difficulty in the implementation of multistate models is the lack of available software for these models; most of the current ones present some difficulties and limitations in practice. In so e Markov assumption may be appropriate, while for others the semi-Markov te d Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 29 transient state, “alive”, and one possible transition from this state to an absorbing state, “dead” nsition, , representing some terminal event of interest. Patients enter the study at the “alive” state and then transfer to the “dead” state by some coefficient of tra ( ) t α , at time ss as transition intensity, is a generalization of the hazard function of survival analysis, as they provide the instantaneous hazard for movement out of one state into another. Note that for the mortality model, the only relevant information in t. This coefficient of transition, which we shall expre t − F is the individual status (alive or dead) along with the covariate history. Figure 2.1: Mortality model for survival analysis. These models have been discussed in more detail in Chapter 1. Note that the transition (occupation) probabilities are expressed by 0, exp () ( ) ( ) () 111 t 0 p tp t uduSt α ⎧⎫ ⎪⎪ ==− ⎨⎬ ⎪⎪ ⎩⎭ ∫ , [2.3] = and () ( )() 212 0, 1 p tp t St==−. [2.4] Splitting the “Alive” state from the simple mortality model for survival data into two transient states, we therefore obtain the simplest progressive three-state model, illustrated in Figure 2.2. It has three states and the only possible transitions are 12→ and 23→. For this model, the entry time into state 2 is a relevant information in t − F. Note that for the progressive three-state model we assume that the transition intensity 1. Alive 2. Dead ( ) t α 1. Alive 2. Dead ( ) t α Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 30 : Progressive three-state model from state 2 into state 3 might depend, in some way, on the entry time in state 2, denoted by 12 t. ( ) 12 is the entry time into state 2tFigure 2.2 . When we refer to pr to all processes for which each state has at most one transition into it, and none into the initial state. ogressive models, we are referring The state occupation probabilities are now given by, t () ( ) ( ) 111 12 0 0, exp p tp t udu α ==− ⎨⎬ ⎪ ⎧⎫ ⎪⎪ ⎪ ⎩⎭ ∫ 212 12 12 00 0, exp e tu t u () ( ) ( ) () ( ) 23 xp , p tp t u vdv αα ⎡⎤ ⎧⎫ ⎪⎪ ⎢⎥ == − ⎨⎬ ⎢⎥ ⎪⎪⎪ ⎪ ⎩⎭⎩ ⎭ ⎣⎦ ∫∫ vudvdu α ⎧ ⎫ ⎪ ⎪ − ⎨ ⎬ ∫ () () () 312 1 p tptpt=− − The mortality model and the progressive three-state model are particular cases of the k-state progressive model. This model is illustrated in Figure 2.3. Figure 2.3: k-state progressive model. State 1 State 2 State 3State 1 State 2 State 3 () 12 t α ( ) 23 12 ,tt α State 1 State 2 State 3State 1 State 2 State 3 () 12 t α ( ) 23 12 ,tt α State 1 State 2 State 3 State KState 1 State 2 State 3 State K Another possible and extensively used model in the literature to describe the disease progression is the illness-death model (also known as disability model). This model is fully characterized by three transition intensities, each one describing the Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 31 the incidence of the disease as well as death. One important roblem here is to evaluate whether previously diseased subjects have the same risk of death as those who have been healthy all their lives. This model is applied to the Heart Transplant data in Chapter 3, and to the PROVA data in Chapter 4. instantaneous risk of moving from one state to another, namely: the disease intensity (incidence), the death intensity without disease, and the death intensity after the occurrence of the disease. These models are widely used in the medical literature and can be used to study p Figure 2.4: Illness-death model. 1. Healthy 2. Diseased 3. Dead ( ) 12 t α () 13 t α ( ) 23 12 ,tt α 1. Healthy 2. Diseased 3. Dead 1. Healthy 2. Diseased 3. Dead ( ) 12 t α () 13 t α ( ) 23 12 ,tt α Again, we assume that the transition intensity from state 2 into state 3 might depend on the entry time in state 2. For such a model, the probability of staying in the healthy state is given by t () ( ) () () () 111 12 13 0 0, exp p tp t u udu αα ⎧⎫ ⎪ ==− + ⎨⎬ ⎪⎪ ⎩⎭ ∫ The probability of being in the diseased state is vdu ⎫ ⎪ ⎬ ⎪ ⎭ and finally, the probability of having died is 2 . ⎪ () ( ) () () () () () 212 12 12 13 23 00 0, exp exp , tu t u pt p t uvvdvvud ααα α = ⎡⎤ ⎧⎫⎧ ⎪⎪⎪ ⎢⎥ =−+ − ⎨⎬⎨ ⎢⎥ ⎪⎪⎪ ⎩⎭⎩ ⎣⎦ ∫∫ ∫ () () () 31 1pt pt pt=− − Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 32 model can become more transparent (Hougaard, 2000). For example, when sing an illness-death model, it is not possible to know the previously visited states for a patient found in the ‘dead’ state. We can, however, extend these models to become progressive by adding an extra state, representing death after disease. This model is depicted in Figure 2.5. In this mo the healthy state, the current state includes the information on the states previously visited, and the order in which they have been visited. ay be is will m individuals to move back and forth between the two states. Such a model is referred to as the illness-death model with recovery (see Figure 2.6). Note that the illness-death model of Figure 2.4 is not progressive since there are two transitions into the ‘death’ state. One disadvantage with non-progressive models is the impossibility to know the states visited before and the visiting order just by knowing the current occupied state. Even so, in some cases the state structure can be changed so that the u 1. Healthy 3. Dead following disease 2. Diseased 4. Dead prior to disease 1. Healthy 3. Dead following disease 2. Diseased 4. Dead prior to disease Figure 2.5: Progressive illness-death model. del, if we assume that all patients enter the study in Reversibility between the “healthy” state and the “diseased” state m considered although th ake the analysis more difficult. These models allow Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 33 Figure 2.6: Illness-death model with recovery The competing risks model (Andersen et al., 2002) is another multi-state model which extends the simple mortality model for survival data in which each individual may ‘die’ due to any of several causes. As shown in Figure 2.7, there is one transient “alive” state and several ab ing to different causes of eath. Figure 2.7: Competing risks model. For the competing risks model, the hazard for the total mortality is given by . Thus, the survival function is given by sorbing “death” states correspond d () () 1 2 N j j tt αα = =∑ () ( ) 0 exp t St udu α ⎧ ⎫ ⎪ ⎪ =− ⎨ ⎬ ⎪ ⎪ ⎩ ∫⎭ 1. “Alive” N. “Dead, cause N-1” 3. “Dead, cause 2” 2. “Dead, cause 1” 1. “Alive” N. “Dead, cause N-1” 3. “Dead, cause 2” 2. “Dead, cause 1” ( ) 12 t α ( ) 13 t α ( ) 1Nt α 1. “Alive” N. “Dead, cause N-1” 3. “Dead, cause 2” 2. “Dead, cause 1” 1. “Alive” N. “Dead, cause N-1” 3. “Dead, cause 2” 2. “Dead, cause 1” ( ) 12 t α ( ) 13 t α ( ) 1Nt α 1. Healthy 2. Diseased 3. Dead 1. Healthy 2. Diseased 3. Dead Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 34 hen, the occupation probabilities are given by T () ( )() 0,pt p t St== () ( ) ( ) () 11 1 0, exp jj j l p t p t u v dv du αα 111 and tu N l= ⎧ ⎫ ⎪ ⎪ == − ⎨ ⎬ 2 00 ⎪ ⎪ ⎩⎭ model (Hougaard, 2000), depicted in Figure 2.8, is the multi-state odel for bivariate parallel data, with states ‘both alive’, ‘individual 1 dead’, ndividual 2 dead’ and ‘both dead’. ∑ ∫∫ The bivariate m ‘i Figure 2.8: The bivariate model. Although model depicted in Figure 2.8 is non-progressive, we can however, make this model progressive using a similar procedure as those shown above for the illness-death model. For the bivariate model, the occupation probabilities are given by () ( ) () () () 0 0, exp t pt p t u u du αα ⎧⎫ ⎪⎪ 2. Individual 1 dead 111 12 13 == ⎪ − + ⎨⎬ ⎪ ⎩⎭ 2 pt= ∫, () ( ) () () () () () 12 12 12 13 24 00 0, exp exp , , tu t u p t uvvdvvudvdu ααα α ⎡⎤ ⎧⎫⎧⎫ ⎪⎪⎪⎪ ⎢⎥ ⎢⎥ ⎪⎪⎪⎪ =−+ − ⎨⎬⎨⎬ ⎩⎭⎩⎭ ⎣⎦ ∫∫ ∫ 1. Both alive 3. Individual 2 dead 4. Both dead 1. Both alive 2. Individual 1 dead 3. Individual 2 dead 4. Both dead () t α 12 () 13 t α ( ) 24 12 ,tt α ( ) ,tt α 34 13 1. Both alive 2. Individual 1 dead 3. Individual 2 dead 4. Both dead 1. Both alive 2. Individual 1 dead 3. Individual 2 dead 4. Both dead () t α 12 () 13 t α ( ) 24 12 ,tt α ( ) ,tt α 34 13 Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 35 () ( ) () () () () () 313 00 0, exp exp , tu t u pt p t uvvdvvudvdu ααα α = ⎡⎤ ⎧⎫⎧⎫ ⎪⎪⎪⎪ 13 12 13 34 ⎢⎥ =−+ − ⎨⎬⎨⎬ ⎢⎥ ⎪⎪⎪⎪ ⎩⎭⎩⎭ ⎣⎦ 3 . The choice of the appropriate state structure depends on the available data as well as the research goal. Goals for the data analysis will be determinant for the choice of the events and the time scale. For example, in the Psoriatic Arthritis data we are interested in studying the evolution of the disease through age. Thus, age is the appropriate choice for the time scale. An overview of different time scale can be found in Keiding (1991). s was not observed, d from its origin, leading to left-censored observations. Usually, atients are observed at intermittent follow-up visits, at which time individual data and ovariate information are collected, but the information from the periods between visits ∫∫ ∫ and () () () () 412 1pt pt pt pt=− − − 2.3. Sampling times. In medical applications, measurements of a disease are often made at several times, providing incomplete observations in some way. Because the process cannot be observed over an infinite time period, at times it will be ended before an absorbing state is reached. In these situations the whole trajectory of the proces providing in this way right-censored observation times. It can also happen that the process is not observe p c Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 36 not available. In such cases the transition times of movement between states are not visits are not known; these e process has not been observed at Figure 2.9: Evolution of an illness-death multi-state model. Also frequent is the presence of left-truncation when dealing with multi-state models. For example, a common occurrence in the illness-death model is that for one patient to be selected, he must be in the “Healthy” state. Note that censoring and truncation are two different concepts. Censoring indicates that the information about the process is only partially observed, whereas truncation represents an exclusion of a part 0 3 6 Time is exactly observed and the states occupied in between observations are said to be interval-censored. Figure 2.9, illustrates a possible evolution of a patient in the illness-death model. The process of Figure 2.9 is observed at four occasions during the time period of 9 months, being the states occupied at times 0, 3, 6 and 9 months, the only information available. Following this figure we do not know the exact time of transition between the “Healthy” state and the “Diseased” state. More, if th time 6, for some reason, the transition between the “Healthy” state and the “Diseased” state was not observed at all. 0 3 6 9 Time 0 3 6 Time 0 3 6 9 Time Death Diseased Healthy Death Diseased Healthy Death Diseased Healthy Death Diseased Healthy Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 37 of the information because it is not observed. In most cases, the truncation and censoring mechanisms are assumed independent from the process. Ideal data collection should record every event or transition along with transition times and covariates. 2.4. C . arious possible models for the transition rates between states can be accomm ost common characterized by one of the following assumptions: 1. Time-Homogeneity: the intensities are constant over time, that is, transition intensities are independent of . Therefore we have ommon simplifying assumptions V odated in expression [2.1]. Transition intensities can depend on the states previously visited, the time since the last event, covariates, etc. Furthermore, they may be constant over time or not. The m models (see Figure 2.10) are t ( ) ( ) αα = hj hj tZ Z . on the previous history of the individual. Th 2. The Markov assumption: future evolution only depends on the current state and not at is, transition intensities do not depend on other information prior to time t. 3. The semi-Markov assumption: future evolution not only depends on the time t since origin, but also on the time spent in the current state h, that is, tt hj − , where t is the transition time from h to j. If we assume, in addition, that the transitions do not depend hj directly on t, we will have intensity functions of the general form ( ) α − hj hj ttZ . Chapter 2. Generalities on multi-state processes _____________________________________________________________________________________________________________________________________________ 38 Figure 2.10: Three-state progressive models. (a) Homogeneous Markov model; (b) Non-homogeneous Markov model; (c) Semi-Markov model ( is the transition time from state 1 to state 2). Because of their simplicity, Markov models are the most used in the literature. In ct, if the state structure has a large number of states, non-Markovian models can pidly become complicated. Note that 12 α 1 2 3 12 t fa ra ( ) { } , 0Xt t≥ is a Markov process, if for any with , and , st 0st≤< () { } , , 1,2,...,hjxu n∈ with hj ≤ , we have () () ( ) ( ) ( )() ,,0 X tjXshXuxu us XtjXsh ⎡⎤⎡ ===≤<= == ⎣⎦⎣ PP ⎤ ⎦ . e s depends only on the state occupied at time s. Traditionally, semi-Markov models are considered when the Markov assumption odels the transition intensities depend on Thus, the future of the process after tim is violated. For these m − t F only through the ate currently occupied, the elapsed time since the transition and the covariates. In hapter 4 we consider a different approach to non-Markov modelling, for which, in ddition to the current time, the transition intensities are allowed to depend also on the ntry time on the current state. By choosing the best state structure for the data, the assumption about the model an be more transparent, for example, making the process a Markov process. For st C a e c 1 2 3 1 2 3 (a) 23 α (b) 12 α ( ) t 12 α ( ) t 23 α (c) ( ) t 12 α ( ) 1223 tt − α 1 2 3 1 2 3 1 2 3 1 2 3 1 2 3 1 2 3 (a) 23 α (b) ( ) t 12 α ( ) t 23 α (c) ( ) t 12 α ( ) 1223 tt − α Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 45 , () t etp 23 22 α − = () () ( ) ( ) 231312 12 13 231312 ααα ααα −+ −−+− tee and . plex and require using special software. to relate the transition intensities 1312 1 α αα +−= +− t t etp , () t etp 23 1 23 α − −= Equations for more than three states are more com For homogeneous processes, expression [2.1] reduces to () ( ) () 1 α =∼; , ,..., ihj Xt MSM Zhj N , and then we may use Cox proportional hazards models of type () () T ααβ =exp , [3.4] hj hj hj hj ZZ α with covariates Z. Under model [3.4], inference is conducted by a general likelihood, which is erived as follows. Assuming that the stochastic process General Likelihood d ( ) i X ⋅ is observed at times ,1,0, ... <<< , where are the indexed individuals, let us consider a secutively observed i miii tt ni ,...,1=t general multi-state model, with a pair of states con ( ) ( ) ( ) 1,, ,+riirii tXtX . The general likelihood is then the product of all the terms over all the individuals and all the transitions (Kay, 1986), , ∏∏ = − = =n i m r ri i lL 1 1 0 , Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 46 where )()( ( ) riri tt r,1, 1, − + + is the contribution t h individual for the pair of states tXtXri pl iirii ,, , =o the likelihood for the it ( ) ( ) ( ) 1, +ri t observed. ,,irii XtX ogeneous Markov models. In some applications, the hypothesis of homogeneity may be unrealistic since the illness ogenous model is recommended. Several non-parametric approaches to non-homogeneous processes have ature (Frydma ; Joly et al., 2002). An alternative intervals and then ecewise constant intensities model (Pérez-Ocón et al., 2001; Saint-Pierre et al., 2003), leading to transition intensity functions as step functions. Such a model will be described in detail in the next section. 3.3.1. Piecewise homogeneous Markov model. Here we consider a piecewise process where the transition intensities are defined by stepwise constant functions of type: 3.3. Non-hom tends to evolve over time. In these occasions, a non-hom then been proposed in the liter n, 1995 (parametric) procedure consists of partitioning the whole study period in two or more fitting a pi () () () 1 1 1 , 0 , hj hj l hj l l Z t θ tZ Z t θ α ααθ − ⎧ ≤ ≤ ⎪ =⎨ < ≤ ⎪ ⎩ , kl ,...,3 ,2 = Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 47 whe 1,..., re 1− =k () θ θ θ is the vector of points cut-off ( ) ∞ =<< < < −kk θ θ θ θ 121 ... . Transition intensities () l hj Z α , kl ,...,3 ,2 = are often expressed according to the Cox model [3.4]. The likelihood function is built following parametric methodological procedures. Let us assume k-1 cut points: 012 0 , , ,..., k θ θθ θ = =∞ and let us define one piecewise intensities matrix: ( () ) () 101 1 , , lll Q QQ θ θ θ θ − ≤< ⎧ =⎨≤< ⎩ tZ Zt Zt kl ,...,3 ,2 = Assuming that the process ( ) i X ⋅ is observed at times i miii ttt ,1,0, ... << < , the likelihood is expressed as r ri lL 1 1 0 ,, where ∏∏ = − = =n i mi ()( ) ( ) riritXtXri ttpl riirii ,1,,, 1,, − = + + th individual for the pair of states is the contribution to the likelihood for the i ( ) ( ) ( ) 1,, ,+riirii tXtX observed. Let us define the following intervals: 1 θθ + ⎡ ⎡ = ⎣ ⎣ , qqq B and ] ] 1 θθ + =, sss C with . Then, each contribution is constructed as follows: 1. and , then 1,...,2,1, −= ksq ri l, if qri Bt ∈ ,qri Ct ∈ +1, ()( ) () ( ) ,,1 , , ir ir ir Xt Xt tt +− ,,1 qi iir iir QZ lp + =; 2. if qri Bt ∈ , and 11, ++ ∈qri Ct , then ( ()() () ) ()( ) () ( ) 1 , qi q i QZ Q Z ir tt lp p + −− 1, ,1 1 ,, ,, ,,1 qir ir q Xt Xt Xt Xt iir iir iir iir θθ +++ + =× ; Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 48 2≥ 3. if qri Bt ∈ , and lri Ct ∈ +1, with − ql , then ()() ()() () () () ()( ) 1 , 1, ,,, ,, ,1 1 , ,,1 2 uu Xt Xt qir Xt Xt iir iir iir iir ir l Xt Xt iir iir l tt + ++− + − () () () 11 , qi ui li QZ QZ ir uq QZ lp pp θθ θ θ +− − − = − =× ⎛⎞ ⎜⎟ Assuming an illness-death model, the survival function, × ⎜⎟ ⎝⎠ ∏. ( ) StZ is expressed as: 1. If [ [ 01 ,t θ θ ∈, ( ) ( ) ( ) ( ) ( ) 11 11 |StZ p t=12 || QZ QZ Z p tZ+; 2. If [ [ 12 ,t θ θ ∈, () () () () () () () { } () () () () 12 2 12 11 1 11 1 12 1 12 1 22 1 || | | | | ; QZ QZ Q Z QZ QZ StZ p Z p t Z p t Z pZptZ θθθ θθ =×−+−+ +×− 1 3. If θ θ − ∈⎡⎡ ⎣⎣ , qq t with , 3≥q () () () () () () () {} () () () () () 1 12 1 22 1 22 1 | | | | | { q ui u Q uu q QZ QZ pZp ZptZ pZpZ θθθ θ θθ θθ − = −− +× −×−+ +−×−× ⎝⎠ ⎜⎟ ∏ ∑∏ ) () () 1 22 1 1 | .}q q QZ q ui ZptZ θ − − =+ ×− The assumption of time-homogeneity may be assessed using a piecewise model (Kay, 1986). Likelihood ratio tests can be used to compare the piecewise model with the time ho has approximately a () 11 1 11 1 12 1 1 2 11 11 1 12 1 21 |||| qq u u q QZ QZ QZ uu q q q Z QZ QZ u qi uu ii iu StZ p Z p t Z p t Z θθ θ θ −−− − = −− −− == =−×−+−+ ⎛⎞ ⎜⎟ ⎛⎞ ⎝⎠ ∏ () ( 22 1 | u QZ uu p θθ − ×− ∏ () () () () 1 1 mogeneous model. Under the null hypothesis of homogeneity, the test statistic 2 kr χ − distribution, wher der . Alternatively, the local score test can be used. A complete description of this method can be found in Kalbfleisch and Lawless (1985). e r is the number of parameters under 0 and k the number of parameters un H 1 H Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 49 Cox Markov models. For simplicity, assume the illness-death model of Figure 2.4. The transition intensities, 3.3.2. () α hj tZ , are modelled using Cox-like models of the form () () () T 0 ααβ =exp hji i hj hj i tZ t Z , [3.5] assuming the process to be Markovian. the hazard of ‘death’ without disease, For ( ) 13 α tZ , survival times from diseased patients are taken as censored in disease time. Patients who are healthy also contribute with censored survival times. For the disease intensity, () 12 α tZ , the final the beginning of the disease. Survival times of patients who did not become diseased are taken as censo whether they have died without having been diseased. Finally, to model point is the time of red, whether they are alive or ( ) 23 the occurrence of the disease, we only ente α tZ , the death intensity after r the survival times truncated on disease time, censored or not, of the individuals that experienced the disease. Note that patients are at risk only after entering state 2. Let now  ()  ()  ( ) T β =exp AtZ A t Z be the estimate of the cumulative 0 hj hj hj intensity function with  () 0hj A⋅ the Breslow estimator for du. transition probabilities () ( ) 00 0 t hj hj tu Αα =∫ The estimation of the ) () () () , ( , hj p stZ X t jX s hZ===P st ≤ and jh ≤ for a given covariate vector Z , are expressed in the following way:  ()  () 3 11 1 2 = <≤ ⎛⎞ ⎜⎟ ⎝⎠ j j sut 1 =− ⎜⎟ ∑ ∏ , p stZ dA uZ , Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 50 that is  ()  ()  ( )  ( ) ( ) 11 11 12 13 1=−− −,, p stZ p st Z dA tZ dA tZ ,  ()  ()  ( )  ( ) 12 11 12 22 ≤ =− + ∑ , | , | , | ut p st Z p su Z d uZ p u t Z , where A  ()  ( ) () 22 23 1 <≤ sut =− ∏ ,| p st Z dA uZ , therefore,  ()  ()  ( )  ( )  () ( ) 12 11 12 12 23 1=− +−−,| ,| ,| p st Z p st Z dA tZ p st Z dA tZ ,  ()  ()  ( ) 13 11 12 1=− −, | , | , | p st Z p st Z p st Z , and  ()  () 23 22 1=−, | , | p st Z p st Z . Note, however, that the probability of having died without ever being diseased can be estimated by  ()  ()  () 13 11 13 =− ∑ ,| , | nd <≤ sut p st Z p su ZdA uZ . ( ) StZ The conditional survival probability which is defined as () ()() 11 12 0, 0,StZ p tZ p tZ=+, can now be estimated by  ()  ( )  ( ) 0, 0,StZ p tZ p tZ=+, hazards. For example, for the illness-death model, one approach that is often used is to the bas zards for transition 13→ and for the 23→ transition to be proportional. In such 11 12 with changes in its values at observed failure times. Note that in some cases we may assume some conditions about the baseline assume eline ha cases, the model is given by Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 51 () () () T 13 130 13 αα β =exp tZ t Z and () () () T 23 130 23 α αβ =+exp tZ t Z These models can be fitted using most of the statistical packages as long as we use a counting process notation, representing each patient with several observations. More details on these models, such as the variance of the estimated transition probabilities, were discussed in the book of Andersen et al. (1993). In this section we will present the simple case of a nonparametric Markov model without covariates. This approach can be thought of as the generalization of the Kaplan- Meier estimate of the simple mortality model. For survival data, the transition probability from the ‘alive’ state in Kaplan-Meier estimator (see [1.9] and [2.4]). In this section, we propose a generalization of this approach states. Such a generalization was considered by Aalen and Johansen (1978), and is denoted as the Aalen-Johansen estimator. Again, for simplicity, we start by assuming the illness-death model, and later extend such an for the illness-death model, it is enough to consider the transition probabilities ,pst, process is Markov, δ . 3.3.3. Non-parametric Markov models. to the ‘dead’ state may be estimated using the to general multi-state models with a finite number of approach to general models. Note that () 11 () 12 ,pst and () 22 ,pst; all the others can be obtained from these ones. Because the Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 52 () () () 11 12 13 ,exp t s pst u udu αα ⎧⎫ ⎪⎪ =− + ⎨⎬ ⎪⎪ ⎩⎭ ∫ () () 22 23 ,exp t s pst udu α ⎧⎫ ⎪⎪ =− ⎨⎬ ⎪⎪ ⎩⎭ ∫ () ( ) () () 12 11 12 22 ,, , t s p st p su u p ut du α =∫ Assume that we have a sample of n subjects with observe e and death times, 12 , ,..., n tt t . If we now assume th we have events and nk −censored observations, then we may write () ( ) 1 k t t d diseas at k 2 ... t ( ) < << ber of diseas subjec ely, just prior to the event time . Further, let subjects who become di numbers of healthy and diseased subjects who die at that same time. Then the above transition probabilities may be estimated by for the k event times arranged in increased order. Let now 1 i n and 2 i n denote the num healthy and e () i t , while 13 i d and 23 i d denote, respectively, the ts, respectiv () i t 12 i d be the number of seased at time  () () i st t <≤ ⎝ 12 13 11 1 ,1 ii i dd pst n ⎛⎞ + =− ⎜⎟ ⎠ ∏, [3.6]  () () 23 2 ,1 i i i st t d pst n <≤ ⎛⎞ 22 =− ⎜⎟ nd ⎝⎠ ∏, [3.7] a  ()  () () ()  () ( ) 12 12 11 22 1 1 ,, i i ii i st t d pst pst ptt n − <≤ =∑, . [3.8] Because [3.6] and [3.7] are Kaplan-Meier estimators, we may use Greenwood’s rmula [1.10] to achieve a variance estimator for such transition probabilities, whereas, fo Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 53 () ()  () ()  () () ()   () ()   () ()  () () ()  ()  () () () () 2 21 11 12 12 13 1 2 1 12 22 23 13 1 , , , i i i i iii i st t i st t i ii n 12 11 var , p st p = ∑ 1 11 12 13 13 22 2 1 , , 1 , , i i i ii i i st t st ptt ptt d n − <≤ <≤ − − ⎡⎤ − ⎣⎦ ⎣⎦ ate j at time , and let be the num er of subjects in state h just prior to time . For each n pst ptt d n n pst ptt d n − <≤ − ⎡⎤ +⎣⎦ − ⎡⎤ + ∑ ∑ For general cases, explicit expressions for the transition probabilities like those mentioned above cannot be given. For such cases, consider hji d the number of transitions from state h into st () i thi nb () i t { } 1,...,ik∈ consider the intensity matrix with i Q entry () ,hj given by hji α . Consider now , with entry  i Q ( ) ,hj given by  hji hji hi d n α = for hji hj≠; and  jh hhi α ≠ =− hi d n ∑ for () ,hh . We may then express the transition probability rms of the in trix. The Aalen-Johansen estimator takes the form  () ,Pst I Q matrix P in te tensity ma () () i st t i π <≤ where I is the identity matrix. An indispensable assumption here is that censoring is independent, so that censoring times do not carry information nsition between states. =+, on the hazard of tra Further details can be found in Andersen et al. (1993). Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 54 3.4. S e est n of an arise is whe events through a multi-state model provides a better knowledge about how the covariate affects the hazard. In this section, we intend to illustrate how the covariate effect can be wrongly interpreted when using simple regression models such as the Cox timedependent regression model. For simplicity, we have chosen to generate a homogeneous Markov process. We then compare the homogeneous Markov model and Cox Markov model with the Cox An individu imulation study. In situations where th imatio the covariate effect is of interest, one question that c ther the inclusion of information on the intermediate time-dependent regression model. We consider the illness-death model and the population vector as () 231312 ,, TTT where, hj T is the potential sojourn time in state h prior to transition to state j. al in state 1 is then exposed to two mutually exclusive events, for which only the first (small) of these events is observed. Therefore, patient histories can be decoupled in two groups, 321 →→ ( ) 1312 if TT ≤ or 31 → 12 T<e assume also that individual’s times are at of being right-censored by some censoring variable C, denoting the potential censoring time. Thus, individuals times can be censored in state 1 () 13 if T. W risk ( ) ( ) 1312 ,min if TTC < or in state 2 ) if ,<C C. e illness-d odel leads to the tim which is observed and 0 otherwise. () and min+ ≤T T T T Note that th eath m e-dependent Cox model, for the time-dependent covariate is coded as 1, if a transition from state 2 into state 3 ( 12 23 12 13 . Time-dependent Cox regression models. Effect estimates and corresponding standard errors. Stanford Heart Transplantation data. Model Estimate (SE) transplant Age year surgery Transplant by age Transplant by year AIC Table 3.3 I  () Pvalue β − SE 0.127 (0.301) 0.67 598.07 II  () Pvalue β − SE -0.004 (0.312) 0.99 0.031 (0.015) 0.03 595.07 III  () Pvalue β − SE 0.123 (0.303) 0.68 -0.191 (0.070) 0.01 592.55 IV  () Pvalue β − SE 0.158 (0.297) 0.59 -0.749 (0.360) 0.04 594.88 V  () Pvalue β − SE -0.031 (0.318) 0.92 0.027 (0.014) 0.06 -0.179 (0.070) 0.01 590.58 VI  () Pvalue β − SE 0.016 (0.309) 0.96 0.031 (0.014) 0.03 -0.773 (0.360) 0.03 591.52 VII  () Pvalue β − SE 0.076 (0.321) 0.81 0.012 (0.018) 0.51 0.041 (0.028) 0.15 594.97 VIII  () Pvalue β − SE -0.282 (0.514) 0.58 -0.265 (0.105) 0.01 0.137 (0.141) 0.33 593.60 IX  () Pvalue β − SE -0.010 (0.314) 0.97 0.027 (0.014) 0.05 -0.146 (0.070) 0.04 -0.637 (0.367) 0.08 589.13 X  () Pvalue β − SE -0.606 (0.540) 0.26 0.029 (0.014) 0.04 -0.279 (0.106) 0.01 0.186 (0.143) 0.19 590.84  () XI Pvalue β − 0.24 0.03 0.02 0.07 0.16 SE -0.621 (0.531) 0.030 (0.014) -0.253 (0.105) -0.664 (0.368) 0.197 (0.140) 589.10 SE=Standard error; AIC=Akaike Information Criterion. u .2 ) ard ratio estimation es for year o ce c it % ntw co azard o m n penalized splines fo ta (with 95% p n an alo w atic fit ( a line). n H splantation data. d ratio estimation es for year o ce c it % ntw co azard o m n penalized splines fo ta (with 95% p n an alo w atic fit ( a line). n H splantation data. f ac ptan e (w h 95 poi ise nfidence bands). (b) H ointwise confide ce b ds), ng ith quadr d shed f ac ptan e (w h 95 poi ise nfidence bands). (b) H ointwise confide ce b ds), ng ith quadr d shed -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 Fig re 3 : (a Haz rati esti atio with Sta ford eart Tran rati esti atio with Sta ford eart Tran with r ag with r ag pe e a pe e a nali t ac nali t ac zed cep zed cep splin nce splin nce -2 -1 0 1 2 3 year of acceptance log hazard 0123456 age log hazard 18 38 (b) -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 28 48 54 62 (a) -2 -1 0 1 2 3 year of acceptance log hazard 0123456 age log hazard 18 38 (b) 28 48 54 62 (a) Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 63 umpti Validity of Markov’s assumption The models studied in this chapter rely on the Markov assumption. Before the analyses of the Stanford data, we will check the validity of this ass on, verifying that the transition rate, () t 23 α , is not affected by the time since transplant. We, therefore, consider 1"age"Z= and 2"time since transplant"Z = as the covariates, and then we fit the following two hazards: ( ) tZ Z αβββ ++ for patients with 01123 exp( ) 2 transplant, and ( ) 011 for patients without. The results obtained (see Table 3.4) show that effect the of time since the transplant is not significant (p-value=0.93), indicating that the transplant do clude that Markov’s model is d Hear tation data. ramete Es te SE exp( )tZ αβ es not lead to an increase in mortality in the period immediately after surgery. This allows us to con satisfactory for the Stanford data (with time since entry in study as the time scale). Further details about this issue are discussed in Chapter 6. Table 3.4. Estimated effects on the mortality intensity in transplanted patients. Stanfor Pa t Transplan r tima P-value (age) .0 5 0.031 0 15 <0.0 1 β (indicator of transplant .3 ) -0.015 0 36 0.96 2 β (Tim t) 0 0 e since transplan .0001 0. 016 0.93 3 β =Standard x v el The Cox Markov model [3.5] allows us to observe how the covariate effects behave when transition intensities are modelled separately. Through this multi-state methodology we pretend to show that these methods can provide new biological insights while confirming some of the results obtained using the time-dependent Cox SE error. Co Marko Mod Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 64 model. Results obtained f Table 3.5. It can be seen rom fitting this model are given in that year of acceptance, which revealed a strong effect on survival in the Cox model, under the Cox Markov model only shows a significant effect on () 13 t α (HR:0.753; 95%CI: 0.606-0.936), although this covariate ceases to be a significant predictor when the patients admitted in the first year were excluded from the study. Age showed a linear effect on each transition, and is the best predictor for the mortality transition () t 23 α of transplanted patients ( ) HR:1.050;95%CI:1.008 1.100−, as well as for transition () , corresponding to receive a new heart 12 t α () HR:1.032;95%CI:1.004-1.060 mortality intensity in patients without transplant . Interestingly, age has no significant effect on the ( ) HR:1.020;95%CI:0.984 1.057−, suggesting that the effect of age enters through the transplant incidence. As shown in Table 3.5, there is no effect of a previous surgery in any transition. Table 3.5. Final Cox Markov model for all transitions. Stanford Heart Transplantation data. Transitions Covariate  β SE HR 95% CI for HR p-value →12 Age Year Surgery (yes=1,no=0) 0.031 0.001 0.047 0.014 0.070 0.315 1.032 1.001 1.048 1.004 – 1.060 0.873 – 1.150 0.565 – 1.940 0.03 0.99 0.88 Age Surgery (yes=1,no=0) 0.020 -0.229 0.018 0.636 1.020 0.796 0.984 – 1.057 0.229 – 2.768 0.27 0.72 Year -0.023 0.097 0.977 08 – 1.100 0.808 – 1.180 0.02 0.81 13→ Year -0.283 0.111 0.753 0.606 – 0.936 0.01 Age 0.050 0.021 1.050 1.0 23→ Surgery (yes=1,no=0) -0.817 0.455 0.442 0.181 – 1.080 0.07 S E = Standard Error, HR = Hazard Ratio; CI=Confidence Interval. Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 65 model offer ts obtained from the fitted model (see Table 3.6) are in good agreement with those obtained when using the Cox Markov mode in any of the previously studied m transplanted patients . Results ind te that age is the only covariate showing a significant linear effect in all transitions. This occurs even for the mo in patients without a transplant, something that did not occur when using the Cox Markov modelling approach. We have also observed that the acceptance time in the study (year) is a significant predictor, though only for mortality intensity in patients without transplant Homogeneous Markov Model In comparison with the preceding studied models, the homogeneous Markov s a detailed description of the survival process, making use of all the available information to estimate the transition probabilities and intensity rates. By applying this modelling approach, we refit the Stanford data including the potential effects of age, year and surgery in all transitions. Resul l. Some exceptions are the surgery effect, which was not a significant factor odels but now reveals a significant effect on survival for () HR:0.306;95%CI:0.128 0.730−ica rtality intensity ( ) HR:0.739;95%CI:0.595 0.919−. We have also observed that the fitted homogeneous Markov model leads to very similar effects of age in the three transitions (see Table 3.6). To assume that these effects are equal, we use the Wald test statistic, yielding a value of 1.163, revealing non significant differences between them (Marshall and Jones, 1995). eous Ma ko surgery in the mortality transition of transplanted patients. Notice that this new model e over-fitting the process, lik In view of the results obtained, we then consider a simplified homogen r v model, including the (same) potential effect of age in all transitions, the effect of year only in the mortality transition of patients without transplant, and the effect of prevents th of data with redundant parameters. In elihood Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 66 ratio tests (LRT) w from zero (Kay,1986) . This statistic has an approximately ere used to test if the regression parameters are statistically different 1 2 χ distribution under H 0 hj 0: β =. Estimated parameters for the simplified model (results not shown) did not differ substantially from those obtained by the initial model (shown in Table 3.6). Therefore, hypothesis testing can be carried out using either LRT or Wald tests although opinions regarding which is optimal lack consistency. Table 3.6. Multi-state homogeneous Markov model. Estimated transition rates and hazard rates. Stanford Heart Transplantation data. TR (SE)  12 α ( ) 0017.00137.0  13 α ( ) 0011.00054.0  23 α ( ) 0003.00018.0 HR (95%CI) Age 21 → ( ) 1.068 1.039 1.098− 31→ ( ) 1.056 1.020 1.093− 32 → ( ) 125.1030.1076.1 − Year 21 → ( ) 0.975 0.852 1.116− 31 → ( ) 0.739 0.595 0.919− 32 → ( ) 325.1928.0109.1 − Surgery 21 → ( ) 539.2737.0368.1 − 31 → ( ) 0.959 0.277 3.315− 32 → ( ) 730.0128.0306.0 − TR=T ror; HR=Hazard ratio; CI=Confidence Interval. given by ransition rate; SE=Standard er Further, we have used Wald’s test to verify whether or not a relation between transplant and survival exists (Kay, 1986). Formally, the hypothesis of no relation is 23130   ( ) 2 αα 13 23 11 Wv − =, being : α α = H, and then Wald’s test reduces to Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 67   () rv11 13 23 va α α −. With our data, under the null hypothesis the W statistic (which survival probabilities, 13 = follows a ) yields a value of 18.5, suggesting that the transplant is significantly associat inishing in mortality risk. In Figure 3.3 we also compare fitted 2 1 χ ed to a dim  p and 23  p , confirming the beneficial effect of the transplant. 0 500 1000 1500 0.0 0. 0.4 0.6 0.8 02 1. Time ted sur ival probabiliFit v ty From state 1 From state 2 Figure 3.3. Fitted survival probability for the mortality intensity from Heart Transplantation data. Note however that likelihood ratio tests (fitting unrestricted and restricted models) can also be used for constructing a test of H against the general alternative a multi-state homogeneous continuous-time Markov model. Stanford (under the null hypothesis the test statistic has an approximately distribution) (Gentleman et al., 1994). 0 2 1 χ Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 68 The goodness-of-fit of a multi-state model can be assessed by comparing the observed and predicted number of patients undergoing each transition. Table 3.7 (page 72) reports the observed percentages of patients in states “own heart”, “new heart” and “dead”, together with the corresponding expected percentages obtained from the fitted homogeneous Markov model. For comparison purposes, we also included the expected percentages obtained from all fitted multi-state models. In this table, we can observe that, for lower survival times, the mortality is underestimated from the fitted homogeneous Markov model. In many cases, these discrepancies can be explained by the failure of the Markov assumption. Another possibility is that the transition rates vary heart’ to state ‘new heart’ (approximately 73% from the total) occur up to 51 days of survival. Figure 3.4 compares the resulting transition probabilities with time, so that the model is non-homogeneous. This is the case for the Stanford data; it is seen that most of the transitions from state ‘own ( ) 12 0, ptz from the fitted Cox Markov model and homogeneous Markov model, respectively. This figure seems to point out that the Cox Markov model explains in a better way how transition rates vary with time, indicating a rapid increase of the transition probability up to 51 days, then decreasing quickly afterwards. Taken as a whole, these results suggest that a homogeneous model may be inappropriate. To assess the assumption of time homogeneity, Kay (1986) suggests the use of a piecewise model. Again, likelihood ratio tests can be used to compare the piecewise model with the homogeneous model. For the Stanford Heart Transplantation data the test statistic (which follows a 2 4 χ ) suggests the use of a nonres from time-homogeneity is via a local score test (K homogeneous model. Another method for verifying departu albfleisch and Lawless, 1985). Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 69 ated transition probabilities from state 1 (‘own heart’) to state 2 (‘new heart’), obtained from the Cox Markov model (solid line), and the homogeneous Markov model (dash Figure 3.4. Estim ed line). Stanford Heart Transplantation data. Non-homogeneous model In this section, we construct a piecewise constant intensities model with one cutoff point, specified from the Stanford data covariates that showed a significant effect when fitting the homogeneous model. After examining the likelihood for several cut-off points θ , a value of 90 θ = days was selected, and two intervals (days 90 time ≤,time >90 days ) were then considered. Focusing mainly on the shortterm survival period of 90 days, results in Table 3.7 indicate that, for the nonhomogeneous model, the agreement between predicted and observed percentages of patient ach transition is globally satisfactory, being clearly better than that which was obtained previously by using an homogeneous model. Estimated mortality rates and covariate effects produced from fitting the nonhomogeneous model are presented in Table 3.8. It is seen that in both intervals the 0 0, 0,2 0,3 0,4 0 200 400 600 P babilit 1 0,5 0,6 800 1000 1200 1400 Time (days) ro y s in e Chapter 3. Multi-state Markov models ____________________________________________________________________________________________________________________________________________ 70 resulting estimates for the mortality intensity were lower in transp nts, though only in the second interval (time >90 days ) a significant difference was found. When examining the fixed covariate effects, we see that, for time ≤t acceptance is a significant predictor in all transitions lanted patie days 90 , age a ( ) HR:1.032;95%CI , while the effect of year is only significant on the mortality intensity in pa ut transplant () HR:0.716;95%CI:0.571 0.899− and the effect of surgery only on the mortality intensity in transplanted patients :1.011 1.053− tients witho ( ) 0.968 . For the odel with the sulting surv ost survival unction HR:0.131;95%CI:0.018− second interval (time >90 days ), however, the only significant covariate was age at acceptance () HR:1.061;95%CI:1.008 1.116−. Finally, turn e-dependent Cox m four mu ) shows the re ival functions, St, obtained from f odels to Stanford data. The range of time has been restricted to 90 days to em ifferences between the four functions. In view of results of Table 3.7, the upperm model. All the other proposed models ulting survival f for the entire range of time. ing to the comparison of the tim lti-state models described here, Figure 3.5(a  () itting all these m phasize the d it is not surprising that curve corresponds to the homogeneous Markov produce similar survival curves. Figure 3.5(b) shows the res Chapter 3. Multi-state Markov models _____________________________________________________________________________________________________________________________________________ 77 Chapter 4 Multi-state non-Markov models Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 79 these models may carry severe limitations to use a semi-Markov assumption in which future of e process does not depend on the current time but rather on the duration in the current ate. Semi-Markov models (Andersen et al., 2000) are also called “clock reset” models, ecause each time the patient enters a new state time is reset to 0. The aim of the present odels and propose a new non-Markov or which the transition intensities are allowed to depend not only on the current le, whe conside sity e me of t 4.1. Introduction. Traditionally, statistical methods for analyzing multi-state models depend on the Markov assumption. Under the Markov assumption, the transition intensities depend on the current time and the state currently occupied; they do not depend on the patient history (length of stay in the current state; patient characteristics measured before, etc.). By ignoring the disease history behaviour, which can make the model inappropriate. It is a fact that the future health of recently diseased individuals may be different from those who have been diseased for a long time. One alternative approach is th st b chapter is to review Cox semi-Markov m approach, f time, but also on the time of transition to the current state. For examp n ring the illness-death model (see Figure 2.4), we obtain a non-Markovian model by allowing the mortality inten after the occurrence of the disease () transition 2 3→ to depend on the current time as well as th ti ransition to the diseased state. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 80 4.2. S Transition intensities can depend on the current state, states previously visited, the entry time in each state, the time since the entry into the current state, covariates, etc. As mentioned before, these models can rapidly become very complicated nd therefore simpler models are necessary. If the process is not Markov, then semi-Markov models can be a wise choice. Semi-Markov models make the assumption that the history of the process ariates. .2.1 Cox semi-Markov models. Assuming again the illness-death model, the transition intensities, emi-Markov models. , a depends only on the state currently occupied, the elapsed time since the entry into the current state and the cov 4 () 1jtZ α , he difference between Markov and semi-Markov models are expressed on transition, . According to semi-Markov models, future 2,3j=, are modelled as in [3.5]. T 23→ evolution not only depends on the current state, but also on 12 t ( ) 12 tt≤, the time that the individual remains in that same state. Therefore, these transitions intensities will be expressed as () () ( ) T 23 12 230 12 23 αα −=−exp ii ttZ tt Zβ i Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 81 The transition probabilities depend on the history of the process through the time 12 t, () 12 ,, hj p stZt , and can be estimated by,  ()  ( ) () 22 12 23 12 1 <≤ sut =− − ∏ ,| , p st Zt dA u t Z , )  ( )  () ( ) 12 23 12 01+−−−, |  ()  () (  12 11 12 00=−,| , | p tZ p t ZdA tZ p t Z dA t t Z , and  ()  () =−, 23 22 12 1 | , | , p st Z p st Zt . The rest of the transition probabilities are estimated as in the Cox Markov model of section 3.3.2. The survival probability ( ) StZ  can now be estimated by ()  ( )  () 11 12 0, 0,StZ p tZ p tZ=+ with changes on its values at observed failure times. .3. Non-Markov models. an-Meier estimator for the estimation of 4 There are few references to non-Markov multi-state models in the literature, some exceptions are the works of Strauss and Shavelle (1998), Aalen, Borgan and Fekjaer, (2001), Datta and Satten (2001), and Glidden (2002). Strauss and Shavelle (1998) developed an extension of the Kapl transition probabilities, avoiding the Markov assumption. The proposed estimator is constructed by partitioning the survival probability in proportion to the number of live and uncensored patients in each state. Aalen et al. (2001) and Datta and Satten (2001) Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 82 estimator occu t Mark ators of process under independent censoring. Later, Glidden (2002) developed robust confidence bands for those event curves. ving th to compare estimators d stimate ti ). For this purpose, w initially consider the illness-death model. Later, we consider the extension of this methodology for more complex multi-state models. 4.3.1 The illness-death model. Assume that we have an illness-death model. Let studied the performance of the Aalen-Johansen of stage pancy probabilities when the process is no ovian. These authors established the consistency of Aalen-Johansen estim the prevalence functions in a non-Markov Our research on non-Markov models has two main goals. The first goal is to develop an approach based on less restrictive assumptions than those based on the Markov property by completely remo e Markov assumption. The second goal is eveloped here for the transition probabilities with Aalen- Johansen e s (derived under the Markov assump on e ( )() { } 001,, Xt t X ≥= denote non-homogeneous stochastic process. We assume that a finite number of independent istories from the process are observed. We then represent the stochastic behaviour of e process by a random vector where is the potential sojourn time in h prior to transition to state j. An individual in the healthy state is exposed to two mutually exclusive events, e” and “death”, which do not occur simultaneously, and so only the first (small) e therefore divide patient’s history (“course”) into two a h () 12 13 23 ,, TTT hj T th state “diseas of these events is observed. W Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 83 groups according to whether the disease occurred ( ) 123→→ or not () 13→. The history belongs to the second group if 13 12 TT < and then a transition from state 1 into 12 3 TT a transition from state 1 to state 2 occurs at time 12 T , and later at time 12 23 TT a transition from state 2 into state 3 takes place. Of course, several issues influence the observation of the variables hj T. Rightcensoring may appear due to time limitation in the following-up, l llow-up cases, 13 12 TT state 3 is observed. Otherwise ≤ ost to fo etc. On the other hand, whenever () 1 + < , one gets a right-censored value of 12 T () 13 with censoring time T, and no information on is available. In the s e ts a righ d value of 23 Tam way, whenever 12 13 TT≤, one ge t-censore 13 T ( ) 12 with censoring time T. ter, we assume that comple not be observed due to right-censoring. Therefore, individual’s times are at risk of being right-censored by some censoring variable (denoting the potential censoring time). In this way, individual’s times can be censored before t Throughout this chap te information for each individual may C he illness () () 12 13 i.e., <min , CTT or after the illness () ( ) 13 in , TC 12 23 12 i.e., and <+ ≤m CT T T . We will assume that C is independent of ( ) 12 13 ,, TTT 23 . Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 84 1.1. N nparametric estimation of transition probabilities. In this section, our main target is the estimation of the transition probability, 4.3. o () () () () , hj p st X t jX s h===stP, < , which gives the probability of being in state j at time t enough to consider the estimation of the transition probabilities , conditionally on being in state h at time s. In multi-state modelling these are one of the most important targets to be estimated. For the illness-death model, it is ( ) () 11 , p st , 12 , p st and () 22 , p st . All the others can be obtained from these since () () () 13 11 12 ,1 , , p st p st p st=− − and ( ) ( ) 23 22 ,1 , p st p st=− . Therefore we have, () ( ) 11 12 13 12 13 ,,, p st T tT tT sT s=>>>>P, [4.1] () ( ) ,,, , 12 12 12 13 12 23 12 13 p st T T T tT sT s= +>>>P, .2] T tT≤≤ [4 () ( ) 22 12 12 13 12 23 12 12 13 12 23 ,,, ,, p st T tT T T T tT sT T T T s=≤≤+>≤≤+>P. [4.3] hese quantities [4.1]-[4.3] are determined by the joint distribution , T . T of 13 23 , TT () 12 Specifically, knowledge of the distribution of ( ) 13 min 12 , Z TT is enough for the , = recovery of () 11 p st : () ( ) () 11 , Z t pst Z s>P while expectations of type > =P, () ( ) ( ) 12 12 23 12 13 ,STTTTT φφ =+≤ ⎡ ⎤ ⎣ ⎦ EI () 12 , arise when handling p st () ()( ) () , ,, , st uusutt φφ ==<≤>vvvI, () () () , 12 ,st S pst Z s φ =>P, Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 85 and () 22 , p st ()  ()( ) () , ,, , st uuust φφ ==≤>Ivv v, ()  ( )  () , 22 , s t , s s φ S pst S φ =. The relevance of s from the fact that they do not rely on the Ma i-state models. Under the Markov assumption, the conditional probability of moving from one state to another only depends on the state currently occupied. These probabilities are not influenced by the states prev ited and he illnessidua t ry time, , a fact that may be unrealistic in me a assumption the proposed estimation methods come rkov assumption, typically used in mult iously vis the times of transition among them. In the case of t death model, this assumption implies that the survival prognosis of an indiv l being in state 2 does no depend on the ent 12 T so pplications. From a mathematical perspective, the Markov restricts the distribution of ( ) 12 13 23 ,, TTT to those distributions satisfying, for any st<, the equation 12 23 12 23 12 12 13 ⎡⎤ 12 23 12 23 12 0 12 13 0 0 +> +> ⎡⎤=+>+>=≤∀<≤ ⎣⎦ ,, TT tTT sT sTT ss P that is, the conditional independence between ≤ ≤ = ⎣⎦ ,, TT tTT sT sTT P [4.4] 12 23 TT + and given that s, TT≤ ition probab the current time and on the time of entry into state 2, . Suppose a sample of n individuals under study followed from an entry time to death or censoring. We denote the sample information as 12 T12 23 TT+> 12 13 , 12 . Proposed non-Markov estimators allow the trans ilities to depend on Ts≤ 12 T ( ) () 1 δ δδ ρ δ η −,,., , iiiiii i i UV , 1in ≤ ≤ Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 86 h are assumed t allwhic o be independent and identic y distributed copies of () () 1,,., , UV δ δδρ η −, where δ () 12 13 min , , UTT = C is the observed sojourn time in state 1; C is an indicator of whether a transition occurs; n , () () 12 13 δ =≤min , TT I 12→ () 23 12 mi V TCT − is the observed soj ≤− CT =ourn time in state 2; () 23 2 ρ = T I , so that δ ρ is an indicator o hether a tran ition 23 1f w s occurs; η =≤I → () ( ) 1 δ η − 13 TC , so that is an indicator of whether a transition 13 occurs. Note that, u → 0 δ =, no information on ( ) ,V ρ is available; while, if 1 δ = , nder then the observation of η is not possible. ble data for some individual will be, ( ) 13 0001,,,, T if in the Therefore, the availa ( ) process a direct transition from state 1 into state 3 occurs; 12 23 110,, ,, TT if a transition state 1 into state 2 occurs and 12 , T − if the individuals transfers from state 1 to state 2 at time , and afterwards have a ensored sojourn time in state 2; and finally from afterwards a transition to state 3; 12 100,, , TC () 12 T ( ) 0000,,,, C if the individual does not iding a censored sojourn time in state 1. a on of the transition probability [4 , we need to make inference c leave state 1, prov For the estim ti .1] on Z , whose distribution function we denote by ( ) ( ) = ≤ Ht Z t P. In such cases, since ZC= and we observe U () min , ( ) ( ) 1 Z C γδ δη =+−=≤I tor (Kaplan and Meier, 1958) based on the pairs , we consider the Kaplan- Meier product-limit estima () γ , ii U . Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 93 2, suggests that ( ) 13 , d p st can be e Lemma stimated by the sample average,  ()  () () ( ) ( ) () () 13 1 1 1 1 ,iiii iii n d i UsU Vt GU V nHs pst ii δ δρ δ − = >+≤ −+ − =∑ II . For the estimation of the expectation ( ) ( ) ( ) 13 13 12 13 ,DTTT φφ T>= ⎡ ⎤ ⎣ ⎦ EI we need the following Lemma, Lemma 3. For each function φ we have () ()() () ,1 1 UU DGU φ δη φ − ⎡ ⎤ − ⎢ ⎥ =− ⎢ ⎥ ⎣ ⎦ E. Proof. Write )( ) () () 13 () ( 13 DT φφ = ⎨ 13 12 13 13 12 13 ,, TC TTT TCTT ⎧⎫ ≤ ⎪⎪ >⎬ ⎡ ≤⎤ ⎪⎪ ⎣ ⎦ I I EI . ow, since C is independent of the process, ⎩⎭ E N () () ( ) ,1TCTT GT 13 12 13 13 − ⎡≤ ⎤=− ⎣⎦ EI On the other hand, . () ( )( ) 13 12 13 1TT TC δ η <≤=−II . Since, under , we have () 11 δη −= 13 UT = , we get () ()() () ()() () 13 13 13 ,1 ,1 11 TT UU DGT GU φ δη φ δη φ −− ⎡⎤⎡ −− ⎢⎥⎢ == −− ⎢⎥⎢ ⎣⎦⎣ EE . ⎤ ⎥ ⎥ ⎦  Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 94 Then, () ( ) ( ) 13 13 12 13 ,DTTTT φφ =>⎡ ⎣ EI⎤ ⎦ may be estimated by the sample average,  ()  ()( ) ()()  () 13 13 1 , iii sU DTTTT 12 13 1 1 1i n i t nGU δ η − = − − φφ <≤ => ⎡ ⎤ ⎣ ⎦ and ( , nd =∑I EI , ) 13 p st can be estimated by  ()  () () ( ) ( )  () 13 1 1 1 1 1 ,ii i n nd i sU t GU nHs pst i δ η − = <≤ − − − =∑I . [4.11] Then one possible estimator for ( ) 13 , p st is given by the following expression,  ()  () () ( ) () () 13 1 , 1 1 , iii n i i UsU Vt nHs pst δ 1 iiii GU V δ ν − −+ . [4.12] • = >+≤ − =∑I Note that, since () ( ) 13 , p st T tZ s≤> then, the above estimator [4.12] can be seen as an estimatio lving the =P, n problem invo joint distribution of ( ) , Z T. This problem is discussed below. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 95 the asymptotic results for the estimated transition probabilities 4.3.1.2. Some asymptotic results. In this section our main interest is  () , hj p st . Special attention will be paid to the following issues: (i) consistency; (ii) convergence in distribution to a normal; and (iii) (limit) variance ion. or tim estimat We first introduce some results f an es ator of the general expectation () ( ) ()()()( 12 12 23 12 13 13 13 12 13 , , , LZT TT T T T TT T T ϕϕ ϕϕ =⎡ ⎤ ⎣⎦ =+≤+ > ⎡⎤⎡ ⎣⎦⎣ E EIEI ) ⎤ ⎦ [4.13] ote that () 13 , p st is related to this expectation, through the choice () ()( N ) == or , ,, , st uuust ϕφ >≤Ivv v. Note that f st ≤ we get () () () ,,, s tstst LSD φ φφ = Also+. , for ( ) ( ) ( ) , , s st ϕφ , , tu u s t = =≤>Ivv and , we get st≤ () () ,, s ts LS t φ φ = (because ( ) ,,0 st uu φ = , the second expectation on the right-hand side of [4.13] is zero). ent shows that The same argum  ( )  () ,, s ts LS t φ φ =. Hence, estimators of () L ϕ are of much practical interest in the scope of the illness-death model. Now, note that, rather than the pair ( ) , Z T we observe where () ,UY YU V δ =+ . Besides, it is known whether the value of is an uncensored observation of T or not . So the problem of estimating Y () 1 ν = () 0 ν = ( ) L ϕ has to do with the tegrals in the presence of covariates (here, the role of the ill be played by the Z). In contrast with the existing related literature es (e.g. Stute, 1993), the Z may be censored too; however, U will offer an uncensored value of Z whenever evaluation of Kaplan-Meier in covariate w developed in the nineti 1 ν = , and so any Kaplan-Meier integral based Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 96 an be computed regardless of the referred issue. This is explained elow. In i= and YU V on the weights of  K c b troduce the Kaplan-Meier integral  () [] () () , K in i i LWUY ϕϕ =∑, n 1 where are the ordered values of the () ( ) ( ) 12 ... n YY Y≤≤≤ ´s i Yiiii δ =+ ; is the concomitant of the ith order statistic , and [] i U () i Y K in attached to () i Y under  K. Explicitly, W stands for the weight which is [] [] 1 1 1 11 i ij K in j Wni n j νν − = ⎡⎤ =− ⎢⎥ −+ − + ⎢⎥ ⎣⎦ ∏ , herew [] i ν is the censoring indicator of . Note that () i Y 0 K in W = for any censored . Since whenever () i Y [] [] ii UZ= [] 1 i ν =, we may write Y. A number of results are availabl or this  () [] () () 1 , n K in i i i LWZ ϕϕ = =∑  () L ϕ e f , in the case of complete observation on Z. When right-censored values of Z correspond to censored observations f T (as in our case), the results are easily adapted. The needed identifiably conditions are satisfied in this application (they are consequences of the independence between C and ,,TTT which of course implies the independence between C and o () 12 13 23 () , Z T ), since: (a) T and C are independent random variables; (b) ν is independent of Z onditionally on T. c Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 97 Define, ZTT () L τ ( )( ) , ϕ ϕτ ≤⎤ ⎦ I, =⎡ ⎣ E where τ stands for the upper bound of the support of . Specifically, we have (Stute, 1993): () min ,YU V TC δ =+ =  () () LL τ ϕ ϕ → w. p. 1 as , provided that the limit exists. In the previous section, we obtained the following estimators for the transition robabilities in the illness-death model: n→∞ p )  ()  ()  ( 11 , pst 1 1 − =− Ht Hs ,  ()  ()  () 12 1 φ =− , , st L pst Hs ,  ()   ()   () 22 φ =, , st L pst , φ , ss L  ()  ()   ( )   () , 23 22 , ,1 ,1 s t s s L pst pst L φ φ =− =− ,  ()  ()  ()  ()  ()  ( )  () , 13 11 12 1 ,1 , ,1 11 st L Ht pst pst pst Hs Hs φ − =− − =− − −− . As mentioned before, an alternative estimator for ( ) 13 , p st (see [4.11]) is now expressed as  ()  ()  () , 13 ,1 st L pst Hs φ •=−, herew () ( ) ,,, st uv u s t φ =>≤Iv. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 98 Estimators  () 1, j p st , as introduced in [4.6], [4.8] and [4.12], do not satisfy the natural restriction ,1 j pst  () 3 1 1 j= = ∑. This problem can be overcome by considering, when estimating the distribution H of Z, empirical integrals based on the weights associated to , the estimator of the survival time T. Certainly, introduce ≤ , note now that,  K  () [] () [] () :: 11 nn KK in in in in ii Ht W U t W Z t ∗ == =≤= ∑∑ II  ()  () t L Ht ϕ ∗ = where ( ) ( ) , tuut ϕ = ≤v. Further fol e general convergence I more,  () () t LHt ϕ →lows from th result above. Introduce now the following estimators:  ()  ()  () 11 1 1 ϕ ϕ ∗− =− , t s L pst L ,  ()  ()  () 12 1 φ ϕ ∗=− , st L L , s pst ,  ()  ()  () , 13 ,1 st s L pst L φ ϕ ∗=−; for whi ϕ = =≤ ∑I, i= , ch  () [] () 1 n K tin i i LWZt  () [] () () ,, n K st in i i LWsZtYt φ =<≤> ∑I 1 Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 99  () () 1 , n K i LWsZYt φ = =<≤ ∑I. Because [] () i i [] () , st in i i Z Y≤ whenever K in W  () [] () () , 1 , n K st in i i i LWsZtY φ = t = <≤ ≤ ∑I0, then ≠ . herefore, T  ()   ()  ( ) [] () [] () () [] () () ,,st st s LL L ϕφ φ ++= [] () [] () [] ()  () KK in in ii K in i t WZs WsZt WZt L ϕ =≤+<≤ ≤ = ∑∑ II I from which it f 11 1 11 1 , , nn n KK K in in ii ii i ii i nn ii n i W Z s Z tY t W s Z tY t == = == = = ≤ <≤≤+<≤> = ∑∑ ∑ ∑ II I ollows that in s W+  () 3 1 1 ,1 j j pst ∗ = = ∑. Note that in the case , we obtain 0s=  ()  ( )  () 0, 13 0, t p tL Kt φ ∗==, which means distribution function of the survival time of the process, K, may be estimated by that the  () 0,t L φ . The main goal in the remainder of this section will be establishing the asymptotics for  () , hj p st and  () , hj p st ∗. Asymptotics for  () 11 , p st First, consider the empirical weighted average U which can be regarded as an estimator for  ()  () () 1 n H in i i RW ϕϕ = =∑  ( )  () R Z ϕϕ ⎡ ⎤ = ⎣ ⎦ E. For simplicity we assume ontinuity in the following. Introduce ( ) ( ) ( ) 0 0 min ,Hz Z z τ τ =≤P, where 0 τ denotes c Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 100 urse, the upper bound of the support of U. Of co 0 H τ equals whenever the support of Z is contained in that of C (otherwise, we can not expect consistency in the right tail of the distribution). Thus if we define H  ()   ()( ) 00 0 RdHZZ ττ ϕ ϕϕ τ ⎡⎤ == ≤ ⎣⎦ ∫, 93): EI we have (Stute and Wang, 19   ()  ( ) 0 RR τ ϕ ϕ → w. p. 1, provided that the limit exists. Let N denote the distribution function of U, and ored ons, respectively; that is, , 0,1 j Nj=, for the subdistributions pertaining to the censored and uncens subpopulati ( ) ( ) ,, 0,1 j Nx x jj γ ==. We also have 1995): U=≤P (Stute,  ()  () ()  ()( )  ()  () 12 1 11 1 n ii ii in i i U RU nGU ϕϕ ϕγ Ur ϕ ηγη = ⎧⎫ ⎪⎪ =+−− ⎨⎬ − ⎪⎪ ⎩⎭ ∑ ϕ + where  () () ()  ( ) () ( ) ) ( ()  () () () ()  ()( ) 0 1 1 1 1 1 , 1 xw w xd Nx N xZ Z Z Nx ϕ ϕ η ϕγ ϕτ < =− ⎡⎤ < =− ⎣⎦ ⎡⎤ << ⎣⎦ − ∫I EI I 1 1 1 1 Nw Gw xU U x GU − ⎢⎥ − ⎢⎥ = I E  () ( )  () () () 10 2, 1 x xd N ϕ ϕ η η < =− ∫Ivv v v N and  () () 12 , n ron ϕ − =P Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 101 rovided that the foll ing condition hold: C1.0 p ow s  () ( ) ()  () () ()  () () () () 0 22 ZZ U z ϕτϕγϕ 2 1 22 111 dNz GZ GU G z ⎡⎤ ⎡⎤ ⎢⎥ <∞ −⎢⎥ −− ⎢⎥ ⎣⎦ ⎣⎦ , C2.0 ≤== ⎢⎥ ∫ I EE  () () ()  () () () 0 0 12 12 00 ZC Z Z zC z dH z τ ϕτϕ ⎡⎤ ≤= < ⎣⎦ ∞ ∫ EI where () 0 Cz is defined as () () () () () () 0 0 . 11 zdG x = Cz Nx Gx −− ∫ his result, together with the Central Limit Theorem (CLT), immediately leads to T   ()  ()  ( ) ( ) 0 12 0 0,nR R N τ ϕϕ σ ϕ ⎡⎤ −→ ⎣⎦ in law, where  ()  ( ) 2 00 var ϕ σ ϕξ = and  () ()   ()( )  () 012 1.UU GU ϕϕ ξηγη +−− − tions we need to consider transfor s such as 1 U ϕ ϕγ = Moreover, in our applica mation  ()  ()  () () 12 2 ,,gR 1R θ ϕϕ ϕ = where typically ϕ ( ) ,gxy xy=, ( )( )( ) ,11g=xy x y− − r )() ( ,1gxy x y=−. Then, the multivariate delta method (Serfling, age 124) ensure CLT and the o 1980, p  () ()() () () ( ) 00 12 12 1 2 0 12 ,, 0,ngRRN ττ , θϕϕ ϕ ϕ σ ϕϕ ⎡⎤ −→ ⎣⎦ in law, where () () () () () () () () () () () () () 00 00 00 12 12 12 2 2 22 2 012 0,1 0,2 ,, ,, 0,12 ,, , 2 xy R R xy R R xy R R gg xy gg xy ττ ττ ττ ϕϕ ϕϕ ϕϕ σϕϕ σ σ σ == = ⎡ ⎤ ∂∂ ⎡⎤ =+ ⎢ ⎥ ⎢⎥ ∂∂ ⎣⎦ ⎣ ⎦ ∂∂ +∂∂ and where Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 102 () 22 0,1 0 1 σ σ ϕ =, () 22 0,2 0 2 σ σ ϕ =, ( ) 12 0,12 0 0 cov , . ϕϕ σξξ = We can easily see that can be simplified as 0,12 σ () () () () () ( 12 12 0,12 1 1 2 1U GU ϕϕ )( ) () ( ) 00 1 2 1 UU U R R ττ ϕϕγ σ ηη ⎢⎥ ⎡ =− ⎣ ⎢⎥ − EE γ ϕϕ ⎡⎤ ⎤ −− ⎦ ⎣⎦ , hich can be used for computation of the limiting variances and covariances w 2 0,1 σ , 2 0,2 σ , 0,12 σ . Besides, the first summand is also written as () () () () () ()( ) () 012 12 2. 1 1 UU ZZZ GZ GU ϕϕγ ϕϕ τ ⎡⎤ ⎡⎤ ≤ ⎢⎥ =⎢⎥ − ⎢⎥ −⎣⎦ ⎣⎦ I EE  ()   ( ) ( )   ( ) ( ) 11 ,1 1 ts pst R R ϕϕ =− − where  () ( ) tuu ϕ =≤It and Now note that  ( ϕ ) ( ) suus=≤I. The existence of  ()( ) 0 tZZ ϕ τ ⎡ ⎤ ≤ ⎣ ⎦ EI is guaranteed since the nction  t ϕ is bounded (similarly for  s ϕ ). On the other hand,  ()( ) () fu 0 t τ ≤ . Hence, the estimator ,  () 0 t Z ZHt ϕτ ⎡⎤ ≤= ⎣⎦ EI for each 11 p st converges almost surely to ( ) 11 , p st provided that 0 t τ ≤ . We also obtain, for 0 t τ ≤ ,  () () () ( ) 11 11 0 ,, 0,,n p st p st N st σ ⎡⎤ 12 −→ where ⎣⎦ in law, () () () () ( ) () () () () () 2 22 2 00,1 0,2 0,12 24 11 1 ,2 11 1 Ht Ht st Hs Hs Hs σσ σ σ −− =+− −− − 3 , () () () ( ) () () () () () 2 2 2 0,1 2 1, 11 Ht HU Zt Ht GZ NU γ σ ⎡⎤ −− ⎡⎤ ≤⎢⎥ =− − ⎢⎥ −⎢⎥ − ⎣⎦ ⎣⎦ I EE () () () ( ) () () () () () 2 2 2 0,2 2 1, 11 Hs HU Zs Hs GZ NU γ σ ⎡⎤ −− ⎡⎤ ≤⎢⎥ =− − ⎢⎥ −⎢⎥ − ⎣⎦ ⎣⎦ I EE Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 109 The distribution could be estimated by [4.5], whereas the expectation involving H () , Z T is estimated by Y, where  () [] () () ,, 1 , n K st in st i i i LWZ ψψ = =∑ () ( ) ,,, st uus ψ =≤≤Ivvt. This would lead to  ()  ()   ()  () () , 22 , , s t s s pst Hs L ψ =− . However, here we investigate the asymptotics for Hs L ψ −  ()  ()  ()  ()  ()   ( )   () , , , , , 22 s t st s s H ∗ For this estimator, the theory reviewed in the preceding section will be used. First, we have ss L Hs L pst L s L φ ψ φ ψ ∗− == −.  ()  ()( )  ()( ) , 22 , , ,, st ss ZT T pst ZT T φ τ φ τ ⎡⎤ ≤ ⎣⎦ →⎡⎤ ≤ ⎣⎦ EI EI w. p. 1 where τ denotes the upper bound of the support of Y. Hence, if the support of T is contained in that of C, we get consistency. Assume ( ) 1T τ ≤ =P in what follows. We have  () () () ( ) 12 22 22 ,,0,n p st p st N st π ⎡⎤ −→ ⎣⎦ , in law, e wher ()  ()  ( )  ()  ( )  () 2 ,, 22 2 12 12 24 ,, 1 ,2 3 , s ts ss ss ss LL st LL L φφ ππ π π φφ =+− t φ , () () ()() () ()  () 22 1 2 12 ,,1 1 ZsTt CstY L GT ν π ,, 1st MY φ ⎡⎤ ⎡⎤ ≤> ∨ − ⎢⎥ − ⎣⎦ ⎣⎦ =− − ⎢⎥ −⎢⎥ I EE Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 110 ) () (()() () ()  () 22 1 , 22 ,,1 , 11ss sT CssY L GT MY νφ ⎡⎤ ⎤ ≤> ∨ − ⎢⎥ − − −⎢⎥ − ⎣⎦ 2Z s π ⎡ = ⎢⎥ I EE ⎣⎦ and () () ()()() () ()  ()  () 11 ,, 12 2 ,,,1 . 11 s sst ZsTt CstYCssY LL GT MY ν πφ ⎡⎤ ⎡⎤ ≤> ∨ ∨ − ⎢⎥ =− − ⎢⎥ −⎢⎥ − ⎣⎦ ⎣⎦ I EE φ For the uncensored case, we simply obtain ()  ( )  ()  ()  () , ,, 3 , st ss st ss LL L φφφ ⎡⎤ −. 2, st L π φ = ⎣⎦ Again, using analogous procedures it is possible to estimate the limiting variance of () 2,st π  () 22 , p st by plug-in methods:  ()    ()    ()   ()    ( )   () 2 ,, 22 2 12 12 24 ,, 1 ,2 3 , s ts ss ss ss LL st LL L φφ ππ π π φφ =+− t φ , where  ()()  () ()  () ( ) ()  () ()   () () 2 2 1 2 , 122 1 ,1 1, 11 nii iii st iii CstY UsYt L nGY MY ν ν πφ = ⎧⎫ ∨− ≤> ⎪⎪ =− ⎨⎬ ⎪⎪ − − ⎩⎭ ∑II −  ()()  () ()  () ( ) ()  () ()   () () 2 1 2 , 222 1 1 1, 11 ii iii ss iii UsYs L nGY MY ν ν πφ = ⎫ − ≤> ⎪⎪ =−− ⎨⎬ ⎪⎪ −− ⎩⎭ ∑II and 2 , nCssY ⎧ ∨ ()  ()  () ()  ()  ()()  () ()   ()   () 11 ,, 12 22 1 ,,1 1. 11 iii i ii n s sst ii Y t CstYCssY LL nGY MY νν πφφ = ⎫ > ∨ ∨− i U s ⎧ ≤ ⎪⎪ =− − ⎨⎬ ⎪⎪ −− ⎩⎭ ∑I Regarding technical conditions C1 and C2, we mention that I  , s t φ may or may not verify these assumptions. Some care will be needed in applications. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 111 Asymptotics for  () , 12 p st ∗ Our goal here is to investigate the asymptotic proper or the estimator pertaining ties f to () () () , 12 ,1 st S pst Hs φ =− where ) t. Recall that s ()( ,,, st usut φ =<≤>Ivv ince ( ) ,,0 st u φ =u we have , () , () s t Lst S φ φ =. In order to avoid problems with censoring in the right tail of the distribution, we can write () () ( ) ( ) ( ) ,,, st LsZtTtHtHssZtT φ =<≤>=−−<≤≤t ⎡ ⎤⎡⎤ ⎣ ⎦⎣ EI EI ⎦ . Here, the expectation involving () , Z T is estimated by Y where   ()  [] () () ,, 1 , n K st st in i i i LWZ φφ = =∑  ()( ,,, st usut φ =<≤≤Ivv ) t. This would lead to  ()  ()  ()   ( )  () , 12 ,1 s t Ht Hs L pst Hs φ −− =− . However, here we investigate the asymptotics for  ()  ()  ()   ( )  ()  ()  () ,, 12 ,1 1 st st s Ht Hs L L φφ pst L Hs ϕ ∗∗ −− ∗ ∗ == − −. First, we have  () ( ) () ()( ) , 12 ,1, s pst ZT T , st ZT T φτ ϕ τ ∗→−≤ ⎡⎤ ⎣⎦ EI w. p. 1 ⎡⎤ ≤ ⎣⎦ EI Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 112 where τ denotes the upper bound of the support of Y. Hence, if the support of T is contained in that of C, we get consistency. Assume ( ) 1T τ ≤ =P in what follows. Introduce ) . We also have ()( 2,, ,Csttx sZ tT tx∨= <≤ >∨P  () ()  () ( ) 12 12 12 ,,0,npstpst N st π ∗ ⎡⎤ −→ ⎢⎥ ⎣⎦ , in law, where  ()  () ()  ( ) () ()  ( ) () () 2 22 2 ,, 12 12 24 1 ,2 11 1 st st ss LL st LL L φφ ππ π π ϕϕ =+− −− − 3 s ϕ ,  () () ()() () () () 2 22 2 1, 2 ,,,1 , 11st sZtTt CsttY L GT MY ν πφ ⎡⎤ ⎡⎤ <≤ > ∨ − ⎢⎥ =− ⎢⎥ −⎢⎥ − ⎣⎦ ⎣⎦ I EE − () () ()( ) () () () 2 22 1  21GT π = ⎢ − ⎣ EE 2 ,1 , 1s Zs CsY L MY νϕ ⎡⎤ ⎡⎤ ≤− ⎢⎥ − − ⎥ ⎢⎥ − ⎦ ⎣⎦ I nd a  ()()() () () () () 21 12 , 2 ,, , 1 . 1 s st CsttYCsY LL MY ν πϕ ⎡⎤ ∨− ⎢⎥ =− − ⎢⎥ − ⎣⎦ E φ For the uncensored case, we simply obtain )  () () () ( () () () () 2, ,, 3 ,1 4 st sst sst s L st L L L L φ πϕφϕφ ϕ ⎡⎤ =−−+ ⎣⎦ . Again, it is possible to estimate th variance 1L −  () 2,st π of  () 12 , p st ∗ e limiting by plug-in methods:   ()    () ()    ( )  () ()    ( )  () () 2 22 2 ,, 12 12 24 1 ,2 11 1 st st ss LL st LL L φφ ππ π π ϕϕ =+− −− − 3 s ϕ , Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 113 where ()()  ( () )  () ( ) ()  () ()    () () 1, 22 , iii st i L Y πφ =−− ⎨⎬ ∑ 2 2 2 2 1 ,, 1 1 11 nii ii CsttY sU t Y t nGMY ν ν = ⎧⎫ ∨− <≤ > ⎪⎪ ⎪⎪ −− ⎩⎭ II   ()  () ()  () () ()  () ()  () () 2 2 1 2 222 1 ,1 1, 11 nii ii s iii CsY Us L nGY MY ν ν πϕ = ⎧⎫ − ≤ ⎪⎪ =− − ⎨⎬ ⎪⎪ −− ⎩⎭ ∑I and    ()  ()( )  () ()  ()  () 21 12 , 2 1 ,, , 1 1. 1 n iii s st ii CsttYCsY LL nMY ν πϕ = ⎧⎫ ∨− ⎪⎪ =− − ⎨⎬ ⎪⎪ − ⎩⎭ ∑ φ Note that the above functions , s t φ and t ϕ may or may not verify the assumptions C1 and C2. Asymptotics for  () 12 , p st  () 12 , p st We now come back to as defined above. We have then,  ()  ()  ()   ( )  ()  ()  () , 12 , ,1 s t st Ht Hs L pst Hs L s R 1 φ φ ∗∗ −− =− = ave ϕ − First, we h  () () ( )  ()( ) 0 , 12 , ,1 st s ZT T pst ZZ φ τ ϕ τ ⎡⎤ ≤ ⎣⎦ →⎡⎤ −≤ ⎣⎦ EI EI w. p. 1 Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 114 where 0 τ denotes the upper bound of the support of U whereas τ denotes de upper bound of the support of Y. Hence, if the support of T is contained in that of C, we get consistency. Assuming in what follows, we have () 0,1ZT ττ ≤≤=P  () ()  () ( ) 12 12 12 ,,0,n p st p st N st σ ⎡⎤ −→ ⎣⎦ , in law, e wher  ()  () ()  ( ) () ()  ( ) () () 2 22 2 ,, 12 12 24 1 ,2 11 1 st st LL st Hs Hs Hs φφ σσ σ σ =+− −− − 3 ,  () () ()() () () () 22 2 1, 2 ,,,1 st sZtTt CsttY L ν σφ 2 11 GT MY ⎡⎤ ⎡⎤ <≤ > ∨ − ⎢⎥ =− ⎢⎥ − ⎣⎦ ⎣⎦ − ⎢⎥ − I EE ,  () () () ( ) () () () () () 22 22 1 1 Hs HU Zs Hs GZ σ ≤⎢⎥ 2 1NU γ ⎡⎤ −− ⎡⎤ =− − ⎣⎦ ⎣⎦ and − ⎢⎥ −⎢⎥ I EE ,  () ( ) ( ) () ( ) ( ) () () () () () () 2,, 1 1CsttY Hs HU νγ σφ ⎡⎤ ∨−−− =− − ⎢⎥ 12 , 11 st LHs NU MY −− ⎢⎥ ⎣⎦ F E. or the uncensored case we simply obtain  () () () () () () () 2, ,, 314 st st st LHs L HsL φ σφφ 1Hs ⎡⎤ =−−+ ⎣⎦ −. Again, it is possible to estimate the limiting variance  () ,st σ of  () 12 , 2 p st by plug-in methods:   ()    () ()    ( )  () ()    ( )  () () 2 22 2 ,, 1 ,2 11 1 st st LL st Hs Hs Hs φφ σσ σ σ =+− −− − , where 12 12 24 3 Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 115   ()()  () ()  () ( ) ()  () ()  () () 2 2 2 1, 22 1 ,, 1 , 11 ii iii st ii CsttY L nGY MY ν σφ = ⎫ ∨− ⎪ =−− ⎨⎬ ⎪⎪ −− ∑ 21nsU t Y t ν ⎧ <≤ > ⎪ II i⎩⎭   ()  () ()  ()  () () ()  () ()  () 2 2 21 1nii ii i Hs HU Us n γ ν ⎧⎫ −− ≤ ⎪⎪ I () 222 1 , 11 ii Hs GY NU σ = =− − ⎨⎬ ⎪⎪ −− ⎩⎭ ∑ and  ()  () ()  ()()()  () ()  () ()    ()  () 2 1iiii s HU ⎪⎪ a 12 , 1 ,, 1 1 . 11 n st iii CsttY H LHs nNU MY νγ σφ = ⎧⎫ ∨−−− =− − ⎨⎬ −− ⎪⎪ ⎩⎭ ∑ In this section, we have introduced and investigated up to six different estim tors for the three different transition probabilities: ( ) , 11 p st , ( ) , 12 p st and () , 22 p st . As shown before, more different estimators for these parameters could have been introduced, and no general answer regarding their relative efficiency can be given in principal. Two main estimators were considered here: using estima rs  () , to 11 p st ,  ()  ()  () 22 , p st 12 , p st 11 , p st ∗ and (see, [4.6], [4.8] and [4.9]); or using estimators ,  () 12 , p st ∗ Kaplan-Meier in not on T . For obtaining the asymptotic results, we have used existing theory devoted to tegrals. When the expectation to be estimated depends only on Z (but ), general theory in Stute and Wang (1993) and Stute (1995) is to be used. However, the estimation of expectations involving ( ) , Z T requires using the extended theory in Stute (1993, 1996). In the following sections, a simulation study and an application are performed to illustrate differences between Aalen-Johansen estimators and non-Markovian estimators. To emphasize the differences between both approaches, we shall use in both cases the transition probability ( ) 22 , p st and the corresponding  () 22 , p st . non-Markov estimator Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 116 4.3.1.3. Simulation study. model. Three settings were considered here, differing on transition from state 2 to state 3. In the first setting, times In this section, simulation studies were undertaken to compare Markovian (Aalen-Johansen) and non-Markovian estimators for non-homogeneous illness-death ( ) 12 13 23 ,, TTT () 12 13 23 ,, TTT are exponential distributed and mutually independent, yielding a homogeneous Markov process (see Lemma in the appendix A). For the second setting, times are mutually independent, and are exponentially distributed, but for s were generated from a Weibull distribution with shape parameter 2 ter 0.05. In the third setting, independent exponentially distributed tim and , while . For all settings, the rate parameter used here for distributions and were 0.039 that 12 T 13 T 23 T time and scale parame es were considered for 12 13 23 12 13 and 0.026 so T T 12 17. TT =× T T ( ) 13 12 06.≥= TT P. In the first setting we used the rate parameter 0.065 for hat settin nd 3 do fulfill the Markov assumption. perform uncensored case and two different levels of censoring). For each subject, a censoring time was generated from an exponential distribution, with rate parameter 0.013 and 0.035 ( proportion of individuals with censored observations (either in state 1 or state 2) can be found through the following expression: 23 T. Note t gs 2 a not The ance of both Markovian (Aalen-Johansen) and non-Markovian estimators was examined under the presence of three levels of censoring (the determining the heaviness of the censoring). Subjects alive at this censoring time were then censored at that time. The Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 117 ( ) () () 12 13 12 13 12 23 12 13 Proportion of censoring =⎡< ⎤+ ⎣⎦ ⎡ ⎤ +⎡ ≤ ⎤× < + ≤ ⎣⎦ ⎣ ⎦ P PP min , min , min , . CTT TTCCTTTTC In this section we compare Markovian and non-Markovian approaches through the estimation of the transition probabilities ( ) 22 , p st for some fixed value of s, 10 s = serve that fixing s, () and 25 s =. We choose to compare both approaches through this transition probability to emphasize the differences between them. Ob ( ) ( ) 22 , 23 1, p st p st =− is a function of time t alone. For this simulation study,  () , 22 p st is the non-Markovian estimator examined. For the evaluation of the simulations, we shall borrow the work of Couper and Pepe (1997), which contain er things, a compreh lation study. As absolute bias, integrated variance and the integrated MSE. For each setting we derived the analytic expression of s among oth ensive simu Couper and Pepe (1997), we use the integrated () 22 , p st (given in the appendix A) so that the bias and the MSE of the estimator, could be For each configuration, examined. 100 K = data sets were generated, with two different sample sizes and . k 100 N =200 N = Let 22 ,  () p st denote the estimated transition probability at time t for some s fixed on the kth generated data set. For each time t considered, we may obtain the mean for all generated data sets,  ()  () 22 22 1 1 =∑ ,, Kk = k p st p st . K We then define, for some fixed value of s, the pointwise estimates of the bias, variance and Mean Square Error (MSE) as:  () ( )  () 22 22 bias =−,, tpstpst , Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 118   () ()  ()  () () 2 1 1 1= =− −∑ var , , , Kk k 22 22 22 p st p st p st K , and   () ()  () () () 2 22 22 22 1 MSE =− ∑ ,,, Kk 1 = k p st p st p st K . To summarize the results, we also calculate the integrated absolute bias, integrated variance and the integrated MSE, defined in Table 4.1. Tables 4.2, 4.3 and 4.4, present the summary statistics for the two estimators (Markovian and non-Markovian) for settings 1, 2 and 3, respectively. Table 4.1. Summary statistics measuring integrated bias, integrated varian integrated mean square error. Statistic Definition Estimator ce and Integrate () 1 bias t s tdt ∫  () 1bias t ts t = ∆ ∑ d absolute bias Integrated variance  () () 1 ∫var , t s 22 p st dt   () () 1∆ ∑var , t ts pst 22 = Integrated MSE 1 22 MSE ∫, t s  () () p st dt 1 22 MSE = ∑, t ts pst   () () ∆ 1 As shown in Tables 4.2-4.4, the non-Markovian estimator revealed itself to be is a fixed value; 50 with =0.2 and 1 250; ,..., . ststsl l =+ =+∆ ∆ = signific ile the Marko (and espe last antly more precise in settings 2 and 3 while the Markovian estimator obtains better results whenever the process is Markov (setting 1). The non-Markovian estimator obtains in all settings a small bias, wh vian estimator is grossly biased in the last two settings. It is clear that the bias clearly dominates the performance of the Markovian estimator in the last two settings cially in the one), leading to a larger MSE. In fact, in some cases (see Tables 4.3 and 4.4) the performance of the Markovian estimator is disappointing despite having obtained a smaller variance than the non-Markovian estimator. The reason for Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 125 trate the differences between the Aalenator and proposed non-Markov estimators of the transition probabilities. Here we present some figures to illus Johansen estim In Figure 4.3, we present estimated transition probabilities,  () , 22 p st , for a fixed value 200s= and (days). ansition 500s= Figure 4.3. Estimated tr probabilities for  () 22 , p st , 200s = (top) and 500s= (bottom): Aalen-Johansen estimator (solid line) and non-Markov model (bold solid line). 0 0,2 200 400 600 800 1000 1200 1400 P2 t) 0,4 0,6 0,8 1 2(200, 1 0 0,2 0,4 Tim ) P22(5 t) 0,6 0,8 500 700 900 1100 1300 1500 e (days 00, Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 126 4.3.1.5. Conclusions. When using multi-state models for the analysis of survival data, some care is necessary in choosing an appropriate model. Traditionally, statistical methods for nalyzing such models depend on the Markov assumption. We have shown through our simulation study that unless the process is, in fact, Markov, the non-Markovian the Aalendepartures from the true transition probabilities when vival data sets, our estimator always rovides good results. Another quantity of interest in multi-state modelling is the cause-specific umulative incidence function, as defined by Kalbfleisch and Prentice (1980). In fact, if ur interest is the estimation of failure probabilities, cumulative incidence functions are . Note that a estimator here proposed is a wise choice. Simulations suggest that Johansen estimator is highly susceptible to the process is not Markov. Moreover, for large sur p 4.3.1.6. Cumulative incidence functions. c o the appropriate quantity to use. We denote the cumulative incidence for disease as () () 12 12 12 13 =≤≤P, CI t T t T T ( ) 12 CI t should not be confused with the ansition probability, () 12 0, p t one individual’s bein tr . The cumulative incidence function represents the probability of g or having been in the diseased state at time t, while () 12 0, p t is the conditional probability of a healthy individual’s being at entry time in the diseased state at later time t. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 127 at, Note th )( ( ) ( ) ( ) 12 12 13 12 12 23 12 13 ϕ =≤≤=⎡ + ≤⎤ 12 ⎣ ⎦ ,, t t T t T T T T T T T PEI, where () ) tv t ϕ =≤I. In the same spirit as above, this quantity can be e by  () CI ,u u stimated ( ()  () () 12 11 δρ δ 1 − = ≤ −+ iii iii i t GU V n , and it is consistent under the assumption of independence between =∑I n U CI t ( ) 12 13 23 ,,TTT and C. The formal details are omitted. For the uncensored case, it is simplified by  () () 12 1 δ 1 = i n =≤ ∑I n ii CI t U t . lative incidence function for death is the probability of one individual’s being dead within time t. Because death is an absorbing state, the cumulative incidence for death is the transition probability from state 1 into state 3, The cumu ( ) ( ) 13 13 0, CI t p t =. 1.7. Inference on the joint distribution function of ) TT . ) purpose, we will assume (besides the independence between the censoring variable C and the process) that and re independent random variables. Let be the distribution function of . For the 4.3. ( 12 T 13 23 ,, Finally, we consider the problem of estimating the distribution functions of and in the illness-death model. For this 13 T ( 12 23 ,TT 13 T () 12 23 ,TT 13 F13 T a Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 128 simply compute the Kaplan-M it ,1U estimation of 13 F, we eier product lim , based on )  13 F () ( δ η −. The consistency of this estimator follows from the independence of the variables . On the other hand, the joint distribution function of can be written as an expectation: ) for 12 T and 13 T () 12 23 ,TT ()( 12 23 12 23 ,,TsTt TT ψ ≤≤= ⎡⎤ ⎣⎦ PE ( ) ( ) ,,uust ψ = ≤≤Ivv. Of course, for  () ( ) ,,uuu ψ ψ =−vv we have ()  ()  ()()  ()( 12 23 12 12 23 12 12 23 12 13 12 12 23 12 13 ,, , , TsTt TTT TT T T T TT T T T ψ ψψ ⎡⎤ ≤≤= + ⎣⎦ ⎡⎤⎡ =+≤++> ⎣⎦⎣ PE EIEI ) ⎤ ⎦ and hence it is seen that   () S ψ as previously defined would underestimate the target. Note that ()  ()  () ()( ) ()( ) 12 23 12 12 23 12 13 12 23 12 12 23 12 13 12 23 12 23 ,, , , TsTt TTT TT TTC TT T TT TTCTT ψ ψ ⎡⎤ ≤≤= + ⎣⎦ ⎡⎤ ≤+≤ =+ ⎢⎥ ⎡ ≤+≤ ⎢⎥ ⎤ ⎣ ⎦ ⎣⎦ PE II E EI I Now, in order to introduce an estimator for the above expectation, besides the sumpas tion that the censoring variable is independent of the process we will need the following assumptions: H1: 12 13 independent of TT; H2: () 12 13 23 12 independent of conditionally on TT T T≤I. Of course, H1 and H2 are also valid if we assume independence between T and () ,TT . Using these assumptions we have, 13 12 23 ()( ) ( ) ( ) () () () 12 13 12 23 12 23 ⎣⎦ 12 13 12 23 12 23 12 23 13 12 12 23 ,, , 1 1 TT TTCTT TTTT TTCTT FT GT T− − ⎡≤ +≤ ⎤=⎡≤ ⎤×⎡+≤ ⎤ ⎣⎦⎣⎦ ⎡⎤ ⎡⎤ =− ×− + ⎣⎦ ⎣⎦ EI I EI EI Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 129 On the other hand, () ( ) 23 TC 12 13 12 TT T δ ρ ≤=. Since, under ≤+II 1 δ ρ =, we have T= and 12 U23 VT δ =. Then we obtain () () () ()  ()  () 12 23 12 12 23 ,, , TsTt TTT UU V ψ ψδδρ ⎡⎤ ≤≤= + ⎣⎦ + ⎢⎥ =⎢⎥ 13 11 FU GU V δ − − ⎡⎤ ⎡ ⎤ ⎡⎤ −×−+ ⎢⎥ PE E ⎣⎦ ⎣ ⎦ ⎣⎦ which can be estim he sample average:  () ated through t  ()  () () () 13 111 . i iiii iFU GU V n ρ δ − − =−−+ ⎡⎤ 12,23 , 1 , ii iii nUU V Fst δδ ψ + = ⎡ ⎤ ⎣ ⎦ ⎣⎦ ∑  () 12,23 ,Fst is of The formal derivation of the consistency (and further asymptotics) of much interest, and it will be considered in future research. 4.3.2. Extension to mor ls. Methods presented in the previous sections regarding the scope of the illnessdeath model can be extended to more complex multi-state models. Here we will consider, merely as an example, the bivariate model and the progressive four-state model. Note that the progressive three-state model is a particular case of the dea y mod aplan-Meier estimates. s a la meth bec e complex multi-state mode illnessth model when 13 T=∞. For the mortalit el, for the survival analysis, and under our approach we obtain K Note however that if the state structure ha rge number of states, then these ods rapidly ame complicated. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 130 4.3.2 Model. The bivariate model, depicted in Figure 2.8 (page 34), is the multi-state model 1 dead’, ‘individual 2 dead’ a is model, the sto astic b aviour o proces that C is independent of we may now denote the sample information as .1 The bivariate for bivariate parallel data, with states ‘both alive’, ‘individual nd ‘both dead’. This model is described in detail in Hougaard (2000). For th ch eh f the s is represented by a random variable () 12 13 24 34 ,,,TTTT . Assuming , () 12 13 24 34 ,,,TTTT () δ ηρ ,,,, iiii i UV , 1in ≤ ≤ which are assumed to be independent and identically distributed copies of ( ) δ ηρ ,,, , UV , where: C is the observed sojourn time in state 1; C is an indicator of whether a transition occurs; , so that () 12 13 min , , UTT = () () 12 13 δ =≤min , TT I 12→ () 13 η =≤ TC I ( ) 1 δ η − is an indicator of whether a transition occurs; 13→ () ( ) () ( ) () ( ) 24 12 12 13 34 13 13 12 =−≤ +−<IImin , min , min , min , V TCT T TC TCT T TC is the observed sojourn time in the intermediate state (state 2 if 1 δ = , or state 3 if 0 δ = and 1 η =); () () () ( ) () ( ) 24 12 12 13 34 13 13 12 ρ =≤− ≤ +≤− <min , min , TCT T TC TCT T TC II II , so that δ ρ is an indicator of whether a transition 24→ occurs; and ( ) 1 δ ηρ − is an in of whether a transition 34→ occurs; dicator 0 δ =, information on ( ) ,V ρ Note that, under is still available provided that 1 η =; while if 1 δ =, then observation of η is not possible. Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 131 Again, we consider ( ) 12 13 min , Z TT= , and that . Furthermore, we assume that H is the distribution function of Z  H is the Kaplan-Meier product-limit estimator of H. As for the illness-death model,  H is based on the ( ) ( ) ,1U δ δη +− . For the bivariate m del, it is enough to consider the estimation of the transition probabilities o () 11 , p st , () 12 , p st , ( ) 22 , p st , ( ) 13 , p st and ( ) 33 , p st . All the others can be obtained from these since ( ) ( )()() 14 11 12 13 ,1 , , , p st p st p st p st=− − − , () () 24 22 ,1 , p st p st=− and () ( ) 34 33 ,1 , p st p st=− . Therefore, we have () ( ) 11 12 13 12 13 ,,, p st T tT tT sT s=>>>>P, () ( ) 12 12 12 13 12 24 12 13 ,,, , p st T tT T T T tT sT s=≤≤ +>>>P, () ( ) 22 12 12 13 12 24 12 12 13 12 24 ,,, ,, p st T tT T T T tT sT T T T s=≤≤+>≤≤+>P, () ( ) 13 13 12 13 13 12 13 ,,, , 34 p st T tT T T tT sT s=≤> +>>>P, and T () ( ) 33 13 12 13 13 34 13 12 13 13 34 ,,, ,, p st T t T T T tT sT T T T s=≤>+>≤>+>P. These quantities are determined by the joint distribution of . Specifically, knowledge of the distribution of T () 12 13 24 34 ,,,TTTT ( ) 12 13 min , Z TT= is enough for the ery of () 11 , p st recov , () ( ) () 11 , Z t pst Z s > => P P, while expectations of type () ( ) ( ) 12 12 24 12 13 ,STTTT φφ =+≤T ⎡ ⎤ ⎣ ⎦ EI arise when handling () 12 , p st () ()( ) () , ,, , st uv u s u t t φφ ==<≤>vvI, Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 132 () () () , ,st S pst 12 Z s >P and φ =, () 22 , p st ()  ()( ) () ,st I,, ,uuust φφ ==≤>vv v, ()  ( )  () , 22 , s t , s s S S pst φ φ =. In the same spirit, expectations of type ( ) ( )( ) 13 13 34 12 13 ,DTTTT φφ =+>T ⎡ ⎤ ⎣ ⎦ EI are needed for () 13 , p st () ( ) ( ) () , ,, , st uusut φφ ==<≤>vv vI t () () () , 13 ,st D pst Z s >P φ =; nd for () 33 , p st a,  () ()  () , 33 , , s t s s pst D φ =. Using Lemma 2, expectations D φ ( ) S φ , can be estimated using sample averages [4.7], yielding the following estimates for the transition probabilities,  ()  ()  () 11 1 1 − =− , Ht pst Hs ,  ()  () () ( ) ( )  () () 12 1 11 1 δ ρ − <≤ +> =∑II , iiiii nsU t U V t pst , () () =−+ − ii iGU V nHs and  ()( ) δ ρ () ()( ) () () 1 1 22 1 1 δ ρ − − = = −+ ≤+> −+ =∑ ∑   II , ii iiiii ii i n i GU V UsUVs GU V st . ≤+>II iiiii n UsUVt p Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 133 imate of C, denoted by , is now ased on Note that the Kaplan-Meier product-limit est  G () ,1UV ν +− with ( ) 1 ν δρ δ ηρ =+− b. For the expectations () D φ , we need the following Lemma: Lemma 4. For each function φ , we have () ()() () () ,1 1 UU V D GUV φ δηρ φ − ⎡ ⎤ +− ⎢ ⎥ = ⎢ ⎥ −+ ⎢ ⎥ ⎣ ⎦ E. Proof. Analogous to Lemma 2. Then, estimators for the rest of the transition probabilities may be introduced:   ()  () () ( ) ( ) ( )  () () 13 1 1 1 1 1 δ ηρ = − =∑II , iii iii i t n pst , ) () () <≤ +> nsU t U V − −+ − ii GU V Hs and ()( )( 1 δ ηρ () () ()  () ( )() 1 1 33 1 1 1 δ ηρ ≤+ − UsU s . − − = = −+ > −+ =∑ ∑   II , iii ii iii iii ii n i n i GU V V GU V pst Estimation of other parameters under the bivariate model can be introduced following ideas similar to those discussed for the illness-death model. The formal derivation of the asymptotics pertaining to the introduced estimators is left for future research. ≤+>−II iii UsUVt Chapter 4. Multi-state non-Markov models _____________________________________________________________________________________________________________________________________________ 134 The four-state progressive model is fully characterized by three transition intensities. We therefore represent the stochastic behaviour of the process by a random variable , assuming that C is independent of the process. We denote the mple information as 4.3.2.2 The progressive four-state model. () 12 23 34 ,,TTT sa ( ) 12 12 12 23 12 23 δ δδδδδ ρ ,,., , iiiiii iii UVW , 1in ≤ ≤ hich are assumed to be independent and identically distributed copies of ) w ( 12 12 12 23 12 23 δ δδδδδρ , U ,., , VW , where C is the observed sojourn time in state 1; () 12 =min , UT () 12 12 δ =≤ TC I is an indicator of whether a transition 12→ () 23 12 =−min , occurs; V TCT is the observed sojourn time in state 2; () 23 23 12 δ =≤− TCT I, so that 12 23 δ δ is an indicator of whether a transition occurs; m 23→ () 34 23 −in , W = TCT is the sojourn time in state 3; () δ =≤−− TCTT I, so that 34 34 12 23 12 23 34 δ δδ is an indicator of whether a transition occurs. Figure 4.4. The four-state progressive model. 34→ State 1 State 2 State 3 State 4State 1 State 2 State 3 State 4State 1 State 2 State 3 State 4 Bibliography _____________________________________________________________________________________________________________________________________________ 237 Aalen, O. (1976). Nonparametric inference in connection with multiple decrement models. Scandinavian Journal of Statistics 3, 15-27. Aalen, O., Borgan, O., and Fekjaer, H. (2001). Covariate adjustment of event histories estimated from Markov chains: the additive approach. Biometrics 57, 993-1001. Aalen, O. and Johansen, S. (1978). An empirical transition matrix for nonhomogeneous Markov chains based on censored observations. Scandinavian Journal of Statistics 5, 141-150. Aalen, O., Vernon, T.F., De Angelis, D., Day, N.E. and Gill, O.N. (1997). A Markov model for HIV disease progression including the effect of HIV diagnosis and treatment: Application to AIDS prediction in England and Wales. Statistics in Medicine 16, 2191-2210. Alioum, A. and Commenges, D. (2001). MKVPCI: A computer program for Markov models with piecewise constant intensities and covariates. Computer Methods and Programs in Biomedicine 64, 109-119. Andersen, P.K. (1988). Multistate models in survival analysis: a study of nephropathy and mortality in diabetes. Statistics in Medicine 7, 661-670. Andersen, P.K., Abildstrom, S.Z. and Rosthøj, S. (2002). Competing Risks as a Multi- State Model. Statistical Methods in Medical Research 11, 203-215. Andersen, P.K., Borgan, O., Gill, R.D. and Keiding, N. (1993). Statistical Models Based on Counting Processes. Springer, New York. Andersen, P.K., Esbjerj, S. and Sorensen T.I.A. (2000). Multistate models for bleeding episodes and mortality in liver cirrhosis. Statistics in Medicine 19, 587-599. Andersen, P.K., Hansen, L.S. and Keiding, N. (1991). Assessing the influence of reversible disease indicators on survival. Statistics in Medicine 10, 1061-1067. Bibliography _____________________________________________________________________________________________________________________________________________ 238 Andersen, P.K., Hansen, L.S. and Keiding, N. (1991b). Non- and semi-parametric estimation of transition probabilities from censored observations of a nonhomogeneous Markov process. Scandinavian Journal of Statistics 18, 153-167. Andersen, P.K., Horowitz, M.M., Klein, J.P., Socie, G., Stone, J.V. and Zhang, M.J. (1999). Modelling covariate adjusted mortality relative to a standard population. Statistics in Medicine 18, 1529-1540. Borgan O. (1998). Aalen-Johansen estimator. Encyclopedia of Biostatistics (eds. P. Armitage and T. Colton), vol 1, pp 5-10, Wiley: Chichester. Borgan O. (1998b). Kaplan-Meier estimator. Encyclopedia of Biostatistics (eds. P. Armitage and T. Colton), vol 1, pp 2154-2160, Wiley: Chichester. Breslow, N. (1974). Covariance analysis of censored survival data. Biometrics 30, 89–99. Chen, H.H., Duffy, S.W. and Tabar, L. (1996). A Markov chain method to estimate the tumour progression rate from preclinical to clinical phase, sensitivity and positive predictive value for mammography in breast cancer screening. Statistician 45, 307–317. Chevret, S., Leporrier, M. and Chastang, C. (2000). Measures of treatment effectiveness on tumour response and survival: a multi-state model approach. Statistics in Medicine 19, 837-848. Commenges, D. (2002). Inference for multi-state models from interval-censored data. Statistical Methods in Medical Research 11, 167-182. Commenges, D., Joly, P., Letyenneur, L. and Dartigues J.F. (2004). Incidence and mortality of Alzheimer’s disease or dementia using an illness-death model. Statistics in Medicine 23, 199-210. Bibliography _____________________________________________________________________________________________________________________________________________ 239 Couper D, Pepe MS. (1997). Modelling prevalence of a condition: chronic graft-versus- host disease after bone marrow transplantation. Statistics in Medicine 16, 1551- 1571. Cox D.R. (1972). Regression models and life tables (with discussion). Journal of the Royal Statistical Society Series B 34, 187-220. Cox, D.R. and Miller, H.D. (1965). The theory of stochastic processes. Chapman and Hall, London. Cox, D.R. and Oakes, D. (1984). Analysis of survival data. Chapman and Hall, London. Crowley, J. and Hu, M. (1977). Covariance Analysis of Heart Transplant Survival Data. Journal of the American Statistical Association 72, 27-36. Datta, S. and Satten, G.A. (2001). Validity of the Aalen-Johansen estimators of stage occupation probabilities and Nelson Aalen integrated transition hazards for non- Markov models. Statistics and Probability Letters 55, 403-411. De Boor, C. (1978). A practical guide to splines, Berlin: Springer. Duffy, S.W. and Chen, H.H. (1995). Estimation of mean sojourn time in breast cancer screening using a Markov chain model of entry to and exit from preclinical detectable phase. Statistics in Medicine 14, 1531–1543. Durrleman, S. and Simon, R. (1989). Flexible regression models with cubic splines. Statistics in Medicine 8, 551–561. Efron B. (1977). The efficiency of Cox’s likelihood function for censored data. Journal of the American Statistical Association 72, 557-565. Eilers, P.H.C. and Marx, B.D. (1996). Flexible smoothing with B-splines and penalities. Statistical Science 11, 89-121. Esbjerg, S., Keiding, N. and Koch-Henriksen, N. (1999). Reporting delay and corrected incidence of multiple sclerosis. Statistics in Medicine 18, 1691-1706. Bibliography _____________________________________________________________________________________________________________________________________________ 240 Escolano, S., Golmard, J.L., Korinek, A.M. and Mallet, A. (2000). A multi-state model for evolution of intensive care unit patients: prediction of nosocomial infections and deaths. Statistics in Medicine 19, 3465-3482. Fleming, T. R. and Harrington, D. P. (1984). Nonparametric estimation of the survival distribution in censored data. Communications in Statistics: Theory and Methods, 13, 2469–2486. Fleming, T. R. and Harrington, D. P. (1991). Counting processes and survival analysis. Wiley, New York. Frydman, H. (1995). Nonparametric estimation of a Markov illness-death process from interval-censored observations, with application to diabetes survival data. Biometrika 82, 773-789. Gentleman, R. and Crowley, J. (1991). Local full likelihood estimation for the proportional hazards model, Biometrics 47, 1283-1296. Gentleman, R.C., Lawless, J.F., Lindsey, J.C. and Yan, P. (1994). Multi-state Markov models for analysing incomplete disease history data with illustrations for HIV disease. Statistics in Medicine 13, 805-821. Glidden D. (2002). Robust inference for event probabilities with non-Markov event data. Biometrics 58, 361-368. Gottschau, A. and Hogh, B. (1995). Interval censored survival data and multistate compartmental models in the analysis of first appearance of Plasmodium falciparum parasites in infants. Statistics in Medicine 14, 2727-2736. Gray R.J. (1994). Spline-based tests in survival analysis. Biometrics 50:640– 652. Greenwood, M. (1926). The natural duration of cancer. Reports on public Health and Medical Subjects. London: Her Majesty’s Stationery Office 33: 1-26. Bibliography _____________________________________________________________________________________________________________________________________________ 241 Hansen, B.E., Thorogood, J., Hermans, J., Ploeg, R.J., Van Bockel, J.H. and Van Houwelingen, J.C. (1994). Multistate modelling of liver transplantation data. Statistics in Medicine 13, 2517-2529. Harrel, F.E. and Lee, K.L. (1986). Verifying assumptions of the Cox proportional hazards model. Proceedings of the eleventh international conference of the SAS user’s group. Atlanta Georgia, February 9-12, 823-828. Hastie, T.J. and Tibshirani, R.J. (1986). Generalized additive models. Statistical Science 1, 297–318. Hastie, T.J. and Tibshirani, R.J. (1987). Generalized additive models: some applications. Journal of the American Statistical Association 82, 371-386. Hastie, T.J. and Tibshirani, R.J. (1990). Exploring the nature of covariate effects in the proportional hazards model. Biometrics 46, 1005-1016. Helms F., Czado C. and Gschlobl S. (2004). Calculation of LTC Premiums based on direct estimates of transition probabilities. 10 June 2005 <http://wwwm4.ma.tum.de/Papers/Czado/paper393.ps>. Hougaard P. Analysis of Multivariate Survival Data. Springer, New York, 2000. Hui-Min, W., Ming-Fang, Y. and Hsiu-Hsi, C. (2004). SAS macro program for nonhomogeneous Markov process in modelling multi-state disease progression. Computer Methods and Programs in Biomedicine 75, 95-105. Jackson, C.H., Sharples, L.D. (2002). Hidden Markov models for the onset and progression of bronchiolitis obliterans syndrome in lung transplant recipients. Statistics in Medicine 21, 113-128. Joly, P., Commenges, D., Helmer, C. and Letenneur, L. (2002). A penalized likelihood approach for an illness-death model with interval censored data: application to age-specific incidence of dementia. Biostatistics 3, 433-443. Bibliography _____________________________________________________________________________________________________________________________________________ 242 Juang, B.H. and Rabiner, L.R. (1991). Hidden Markov models for speech recognition. Tecnometrics 33, 251-272. Kalbfleisch, J.D. and Lawless, J. (1985). The analysis of panel data under a Markov assumption. Journal of the American Statistical Association 80, 863-871. Kalbfleisch, J.D. and Prentice, R.L. (1980). The Statistical Analysis of Failure Time Data. Wiley: New York. Kaplan, E.L. and Meier, P. (1958). Nonparametric estimation from incomplete observations. Journal of the American Statistical Association 53, 457-481. Kay, R. (1986). A Markov model for analysing cancer markers and disease states in survival studies. Biometrics 42, 855-865. Keiding N. (1991). Age-specific incidence and prevalence: A statistical perspective. Journal of the Royal Statistical Society Series A 154, 371–396. Keiding, N., Klein, J.P. and Horowitz, M. (2001). Multi-state models and outcome prediction in bone marrow transplantation. Statistics in Medicine 20, 1871-1885. Klein, J.P., Klotz, J.H. and Grever M.R. (1984). A biological marker model for predicting disease transitions. Biometrics 40, 927-936. Klein, J.P., Keiding, N. and Copelan, E.A. (1994). Plotting summary predictions in multistate survival models: probabilities of relapse and death in remission for bone marrow transplantation patients. Statistics in Medicine 13, 2315-2332. Klotz, J.H. and Sharples, L.D. (1994). Estimation for a Markov heart transplant model. Statistician 43, 431– 436. Lagakos, S.W., Barraj, L.M. and De Gruttola, V. (1988). Nonparametric analysis of truncated survival data, with applications to AIDS. Biometrika 75, 515-523. Bibliography _____________________________________________________________________________________________________________________________________________ 243 Longini, I.R., Clark, W.S., Byers, R.H., Ward, J.W., Darrow, W.W., Lemp, G.F. and Hethcote, H.W. (1989). Statistical analysis of the stages of HIV infection using Markov model. Statistics in Medicine 8, 831– 843. Longini, I.M., Clark, W.S., Gardner, L.I. and Brundage, J.F. (1991). The dynamics of CD4+ T-lymphocyte decline in HIV infected individuals: a Markov modelling approach. Journal of Acquired Immune Deficiency Syndromes 4, 1141–1147. Mantel, N. and Byar, D.P. (1974). Evaluation of Response-Time Data Involving Transient States: An Illustration Using Heart-Transplant Data. Journal of the American Statistical Association 69, 81-86. Marshall, G. and Jones, R.H. (1995). MARKOV: a computer program for multi-state Markov models with covariables. Computer Methods and Programs in Biomedicine 47, 147-156. Marshall, G. and Jones, R.H. (1995) Multi-state models and diabetic retinopathy. Statistics in Medicine 14, 1975-1983. Nelson,W. (1969). Hazard plotting for incomplete failure data. Journal of Quality Technology 1, 27–52. O´Sullivan, F. (1988). Nonparametric estimation of relative risk using splines and crossvalidation, SIAM journal of Scientific and Statistical Computation 9, 531-542. Paes, A. and Lima, A. (2004). A SAS macro for estimating transition probabilities in semiparametric models for recurrent events. Computer Methods and Programs in Biomedicine 75, 59-65. Pérez-Ocón, R., Ruiz-Castro, J.E. and Gámiz-Pérez, M.L. (2001). A piecewise Markov process for analysing survival from breast cancer in different risk groups. Statistics in Medicine 20, 109-122. Bibliography _____________________________________________________________________________________________________________________________________________ 244 Peto, R. (1972). Contribution to the discussion of the paper by D. R. Cox. Journal of the Royal Statistical Society Series B, 34, 205–207. aint-Pierre, P., Combescure, C., Daurès, J.P. and Godard, P. (2003). The analysis of asthma control under a Markov assumption with use of covariates. Statistics in Medicine 22, 3775-3770. atten, G.A. and Longini, I.M. (1996). Markov chains with measurement error: Estimating the ‘true’ course of a marker of the progression of human immunodeficiency virus disease. Journal of the Royal Statistical Society Series C 45, 275–295. erfling, R.J. (1980). Approximation theorems of mathematical statistics. New York, Wiley. trauss, D. and Shavell, R. (1998). An extended Kaplan-Meier estimator and its applications. Statistics in Medicine 17, 971-982. tute, W. (1993). Consistent estimation under random censorship when covariables are present. Journal of Multivariate Analysis 45, 89-103. tute, W. (1995). The central limit theorem under random censorship. Annals of Statistics 23, 422-439. tute, W. (1996). Distributional convergence under random censorship when covariables are present. Scandinavian Journal Statistics 23, 461-471. Stute, W. and Wang, J.L. (1993). The strong law under random censorship. Annals of Statistics 21, 1591-1607. Su, C.T., Wu, S.C. and Chang, C.C. (2000). Multiaction maintenance subject to actiondependent risk and stochastic failure. European Journal of Operational Research 125, 133-148. S S S S S S S Bibliography _____________________________________________________________________________________________________________________________________________ 245 Turnbull, B.W., Brown, B.W. and Hu, M. (1974). Survivorship Analysis of Heart Transplant Data. Journal of the American Statistical Association 69, 74-80. Bibliography _____________________________________________________________________________________________________________________________________________ 246