Study of Disease Progression and Relevant Risk Factors in Diabetic Foot Patients Using a Multistate Continuous-Time Markov Chain Model

The diabetic foot is a lifelong disease. The longer patients with diabetes and foot ulcers are observed, the higher the likelihood that they will develop comorbidities that adversely influence ulcer recurrence, amputation and survival (for example peripheral arterial disease, renal failure or ischaemic heart disease). The purpose of our study was to quantify person and limb-related disease progression and the time-dependent influence of any associated factors (present at baseline or appearing during observation) based on which effective prevention and/or treatment strategies could be developed. Using a nine-state continuous-time Markov chain model with time-dependent risk factors, all living patients were divided into eight groups based on their ulceration (previous or current) and previous amputation (none, minor or major) status. State nine is an absorbing state (death). If all transitions are fully observable, this model can be decomposed into eight submodels, which can be analyzed using the methods of survival analysis for competing risks. The dependencies of the risk factors (covariates) were included in the submodels using Cox-like regression. The transition intensities and relative risks for covariates were calculated from long-term data of patients with diabetic foot ulcers collected in a single specialized center in North-Rhine Westphalia (Germany). The detected estimates were in line with previously published, but scarce, data. Together with the interesting new results obtained, this indicates that the proposed model may be useful for studying disease progression in larger samples of patients with diabetic foot ulcers.


Introduction
Macro-and microvascular complications (peripheral arterial disease (PAD), ischaemic heart disease (IHD) and chronic renal failure (CRF)) increase both the morbidity and mortality of diabetic patients. One of the most complex and economically challenging complications of these patients is foot ulceration. About 15% of all patients with diabetes experience at least one episode of this complication during their lifespan [1]. Without provision of preventative measures such as adapted shoes and insoles, and regular podiatric care, following the initial "healing", a diabetic patient with a history of a foot lesion will have more than one relapse per year [2][3]. Thus, diabetic foot syndrome should be considered a lifelong disease with patients being "in remission" rather than "cured" after the initial healing of an index diabetic foot wound [4]. More than 80% of all diabetes-related amputations are preceded by an ulcer [5]. Hence every new episode potentially increases the risk of lower extremity amputation [6], increases the associated costs [7] and the risk of premature death [8,9]. Recent publications suggest that the prognosis of patients with diabetic foot ulcers is similar to many cancers [10,11]. Amputationor ulcer-free survival may therefore reflect the effectiveness of diabetic foot ulcer management strategies more effectively [12].
However, there are few quantitative data available about the dynamics of disease progression and the effects on this of demographic variables such as age or gender, or of additional risk conditions, prevalent either at baseline or appearing during the course of the disease. (e.g. PAD, IHD, CRF; definitions: see Table 1). Reliable knowledge about the effect of those risk factors could result in more effective strategies for the prevention and management of patients with diabetic foot ulcers.
Multistate Markov chain models assume that the time course of chronic diseases can be described by assigning every individual patient to a distinct state at any defined time and that the next state depends on the current but not the preceding states. So transition characteristics and the influence of risk factors (predictors) can be estimated using the maximum likelihood or the Bayesian estimator.
In the present study we therefore developed a nine-state continuous-time Markov chain model for quantifying disease progression and the time-dependent influence of relevant risk factors for diabetic foot patients and then tested this model on real world data. Results obtained in this paper allow us to better understand complex interplay of active and inactive episodes of diabetic foot disease.

Materials and Methods
The study data The long-term data of 260 consecutive patients presented between June 1998 and December 1999 with a new diabetic foot ulcer to a single specialized diabetic foot center in North-Rhine Westphalia (Germany) were used for analysis. Patients were followed either until December 31, 2014 or death. The date of the last contact with each patient (right censoring) was documented. For all patients, data which were relevant to the course of the diabetic foot disease (outcomes as well as predictors) were documented in the treating center with exact dates. All patients gave informed consent for participation in the study. They agreed to be contacted personally or allowed the investigators to obtain information on their outcomes from their relatives or their family physicians. A part of this cohort (247 patients without previous unilateral major amputation at study inclusion) has been described in detail elsewhere [13]. All living patients were divided into eight groups on the basis of their ulcer history (previous or current) and amputation status (none, minor or major), adding death as an absorbing ninth state (see Fig 1). Patients in states 1, 3, 5, and 7 had by definition an active ulcer, while those in stages 2, 4, 6, and 8 were considered "in remission" [4]. Patients could switch between the groups several times during the follow-up period. Potential risk factors used for the parameter estimation included time-independent (gender) and time-dependent (age, PAD, CRF, and IHD) factors.

Model description
On the supposition of full observability, the nine-state model can be decomposed into eight competing risks survival submodels (the absorbing state death cannot be an origin state). An additional assumption of the proportional hazards allows us to estimate the influence of risk factors on transition intensities using Cox-like regression. All unknown parameters were estimated using the maximum likelihood method.
Since all transitions are fully observable (the dates of departure from the origin state and arrival into the destination state or censoring were known), we separately analyzed the groups of transitions having common origin states. For any origin state i we defined the set of possible destination states by a subset D i & D 0 = {1,2,3,4,5,6,7,8,9}. Under the settings of the competing risks model with proportional hazards [14], [15] the intensity of transition (hazard function) from state i to state j, j 2 D i , at age t was defined as where β i,j T is transposed vector of the Cox-regression coefficients, u t is the vector of timedependent explanatory factors (covariates), and λ 0,i,j (t) stands for the baseline hazard functions. If the destination state is a unique element of D i then the full hazard function can be calculated as follows Here |D i | denotes the number of possible transitions from state i. Note that is the conditional probability of the transition from state i in destination state j at time t.
For an individual that arrives in the origin state i at time t 0,i we can calculate the conditional (left truncated) survival function using the formula The conditional probability density function corresponding to the transition from state i to state j at time t is equal to Parameter estimation To estimate the vector of unknown parameters ω (this vector can include Cox-regression coefficients β and parameters defining the baseline hazard functions λ i,j (t)) we maximized the loglikelihood function. The Cox-regression parameters that significantly differed from zero and corresponded to predictive variables were chosen using step-wise regression with backwards elimination. Parameter estimation of the continuous-time Markov chain models with observed covariates in the case of partially observable data have been discussed elsewhere [16], [17]. Since mean times to transition from states 1-8 were relatively small, for parameter estimation we used the model with constant (age-independent) baseline hazard functions and included an additional covariate-the age of patient at arrival in the origin state. Possible time-dependency on the binary (switching) covariates PAD, IHD and CRF have been also taken into account in formulas (4)- (5). Statistical analyses were performed with R version 3.0.1. The estimates of the constant baseline hazard functions refer to the baseline group-70-year-old-women at arrival in the origin state not suffered from PAD, CHD, and IHD.

Results
A brief description of the study of patients at first presentation is given in Table 1. Mean patient age at first presentation was 69.0±10.9 years, 59.2% of the patients were male, 56.2% of the patients had a history of previous ulceration, 21.9% of the patients had a minor and 5.0% (n = 13) a major amputation before study inclusion. PAD, CRF and IHD were present at first presentation in 58.1%, 23.5% and 20.0% of the patients, respectively.
There was a higher prevalence of male patients in the studied cohort (p = 0.002, exact binomial test) and more than 50% of the male subjects presented with a recurrent diabetic foot ulcer (p = 0.03). The mean time of observation was 6.3±5.3 years (range 0.0-16.7). During the observation period, 16.5%, 29.7% and 35.1% of the patients developed PAD, CRF, and IHD, respectively, and 211 patients died (87%).
The number of transitions (or right-censored events on the diagonal) from origin to destination states are given in Table 2.
The maximum likelihood estimates of baseline transition intensities and their confidence intervals for models with/without selection of the covariates are presented in Table 3.
As expected, the highest baseline mortality λ 0,8,9 was observed in diabetic foot patients with bilateral major amputation (0.451, CI = (0.262-0.777) in the model with backward elimination and 0.292, CI = (0.001-165.986) in the model without selection of risk factors. Estimates of transition intensities 1!6, 3!6, 3!8, 5!8, 6!7 and 6!9 obtained using an approach without selection of risk factors have abnormally wide confidence intervals. This can indicate that sample size is not large enough or the number of transitions is too small. The mean time to transition from the origin state varied between 0.335 and 2.871 years in the model with backward elimination. The highest probabilities of transition for patients with diabetic foot ulcer corresponded to transitions where there was no change in the level of amputation (none, minor, or major; transitions 1!2, 3!2, 5!4, 7!6).
To understand disease progression in diabetes foot patients better, we calculated the marginal (not conditional on any covariates) and adjusted (based on given covariate profile) transition probabilities P ij (0,t) from state i to state j in age interval (0,t] subject to independent right censoring and possibly left truncation using the R package mstate [18]. The resulting plot is given in Fig 2. In each panel of Fig 2 are shown three profiles of transition probabilities-marginal (solid line), for patients from the baseline group (all covariates are equal to zero, dashed line), and for patients with only one risk factor (all but one covariates are equal to zero, dotted line). The first profile refers to the non-parametric model. The second and the third profiles refer to the semi-   . Presence of peripheral arterial disease (PAD) reduces probability of ulcer healing without amputation in patients without any amputation (transition 3!2). PAD and ischaemic heart disease (IHD) increase risk of death in patients in state 2, 3, and 4. These results are in accordance with estimates of the Cox-regression coefficients obtained by using the parametric approach (see Table 4 below).
In Table 4 the relative risks for potential risk factors and their confidence intervals for models with/without selection of the covariates are shown.

Discussion
Although a large proportion of all lower extremity amputations is preceded by diabetic foot ulcers, little is known about the details of progression of the diabetic foot syndrome and the influence of relevant risk factors on its course. Studies on the long-term prognosis of patients with diabetic foot disease show that more than a third of the patients incur a new ulceration within the first 12 months after their index lesion has healed, and 70% experience at least one recurring episode within the next five years [19]. A previous analysis of the cohort described in this paper generated comparable results (1-, 3-, and 5-year cumulative incidences of at least one recurring ulcer episode: 35, 63, and 77%) [20]. Using a multistate Markov chain model to study the course of this chronic condition, we were able to show the high frequency of transitions between active (1, 3, 5 and 7) and inactive stages of disease (2, 4, 6 and 8). We have quantified mean times (in years) to transition for all origin states and conditional probabilities of transition for all possible transitions. We have found that patients with first and with reulceration but without any previous amputation have similar conditional probability of healing without amputation. The later have lower conditional probability of healing with amputation (ν 3,4 + ν 3,6 + ν 3,8 < ν 1,4 + ν 1,6 ). Conditional probability of death increases with stage in patients with active ulcer (ν 1,9 < ν 3,9 < ν 5,9 < ν 7,9 ). Conditional probability of death is higher in patients with previous major amputation. Mean time to transition varies between 1.967 and 2.871 years in healed patients and between 0.335 and 0.427 years in patients with active ulcer. The preliminary results obtained in this study are supported by other work previously carried out in this area and are therefore biologically plausible. For example, patients with active ulcers and a previous major amputation had the highest death rates and the presence of other comorbidities had an influence on the probability of healing and survival. Previous studies have not taken into account possibility of developing comorbidities during the observation period. These studies investigated the impact of comorbidities that were present at baseline. With our model we were also able to investigate the impact of comorbidities that were not present at baseline but that developed during the observational period. This enabled us to study the differing impact at various stages of the disease, or more precisely, on the transitions between those stages. For example, ischaemic heart disease (IHD) was a predictor for transition to death from all stages except "diabetic foot ulcer without any previous ulcer" (stage 1), "diabetic foot in remission without prior amputation (stage 2)", and "bilateral major amputation" (state 8). This observation might explain why some previous studies find an impact of IHD on patients with diabetic foot ulcers while others do not. In our own analysis of this previous cohort [3] and a report from Faglia et al. [20], the presence of IHD was only significant in univariate analysis (HR = 1.59 (1.11-2.28; p<0.05) and HR = 1.892 (1.077-3.324; p = 0.026, respectively), but lost significance in multivariate analysis: On contrast, Izumi et al. reported ischaemic heart disease (IHD) to be the only significant predictor for death in their cohort [21]. The major difference between the studies was, that the latter analyzed first-time amputees with diabetes while the former two followed diabetic foot ulcer cases from all stages (with and without amputation). Other examples are the diverse impacts of age, gender and the presence of peripheral arterial disease (PAD) on specific transitions.
One of the major strengths of our study is the almost complete follow-up over a 17 year long observation period. Among censored patients only 18 patients (6.9%) were unavailable for follow-up before death or the end of the observation period. Furthermore, all data were available with exact dates, enabling us to estimate the time-dependent impact of the potential risk conditions. However, since some transitions occurred rarely (for example transitions 58, 38), we could not obtain the reliable estimates for all Cox-regression coefficients. A further limitation of this study is the fact that all patients in the underlying cohort had experienced at least one active diabetic foot ulcer episode at study inclusion. Therefore, we cannot provide information on transition intensities and relative risk factors for patients with non-ulcerated high risk feet without a previous ulcer episode ("state 0"), which leaves the picture of disease progression in the diabetic foot syndrome somehow incomplete.
To better understand the role of relevant risk factors in disease progression and to increase the power of statistical testing, we need further studies based on a wider range of risk factors and possibly larger data sets. In addition, we plan to validate the model in a population outside that from which it was drawn.

Conclusions
The diabetic foot is one of the most complex and economically challenging complications of diabetic patients. However, little is known about the dynamics of disease progression and the influence of associated factors. We proposed a nine-state continuous-time Markov chain model for quantifying disease progression and the influence of time-dependent factors. We applied our model to the long-term data of 260 patients with diabetic foot ulcer collected in a single specialized center in Germany. Our study combines two important properties that have not been observed in previous investigations of the diabetic foot patients-the long observation period (17-year follow up) and a detailed subdivision of the patients in groups with respect to their amputation (none, minor or major) and health (ulceration or remission) status. It has allowed us to carry out a more detailed analysis of the disease progression and the influence of relevant risk factors in these patients. We have investigated the differing impact of comorbidities (PAD, IHD, CRF) at various stages of diabetic foot. The model is therefore helpful to understand the complex interplay of active and inactive episodes of diabetic foot disease and will now be used to comprehensively analyze larger prospective datasets from patients with diabetes and lower extremity complications.