Quantification of the natural history of visceral leishmaniasis and consequences for control

Background Visceral leishmaniasis has been targeted for elimination as a public health problem (less than 1 case per 10,000 people per year) in the Indian sub-continent by 2017. However, there is still a high degree of uncertainty about the natural history of the disease, in particular about the duration of asymptomatic infection and the proportion of asymptomatically infected individuals that develop clinical visceral leishmaniasis. Quantifying these aspects of the disease is key for guiding efforts to eliminate visceral leishmaniasis and maintaining elimination once it is reached. Methods Data from a detailed epidemiological study in Bangladesh in 2002–2004 was analysed to estimate key epidemiological parameters. The role of diagnostics in determining the probability and rate of progression to clinical disease was estimated by fitting Cox proportional hazards models. A multi-state Markov model of the natural history of visceral leishmaniasis was fitted to the data to estimate the asymptomatic infection period and the proportion of asymptomatic individuals going on to develop clinical symptoms. Results At the time of the study, individuals were taking several months to be diagnosed with visceral leishmaniasis, leading to many opportunities for ongoing transmission. The probability of progression to clinical disease was strongly associated with initial seropositivity and even more strongly with seroconversion, with most clinical symptoms developing within a year. The estimated average durations of asymptomatic infection and symptomatic infection for our model of the natural history are 147 days (95 % CI 130–166) and 140 days (95 % CI 123–160), respectively, and are significantly longer than previously reported estimates. We estimate from the data that 14.7 % (95 % CI 12.6-20.0 %) of asymptomatic individuals develop clinical symptoms—a greater proportion than previously estimated. Conclusions Extended periods of asymptomatic infection could be important for visceral leishmaniasis transmission, but this depends critically on the relative infectivity of asymptomatic and symptomatic individuals to sandflies. These estimates could be informed by similar analysis of other datasets. Our results highlight the importance of reducing times from onset of symptoms to diagnosis and treatment to reduce opportunities for transmission. Electronic supplementary material The online version of this article (doi:10.1186/s13071-015-1136-3) contains supplementary material, which is available to authorized users.


Multi-state Markov model of natural history of visceral leishmaniasis
A general continuous-time multi-state Markov model of disease progression consists of R disease states and M individuals, each of whom is in one of the R states at any particular time. The state occupied by the i th individual at time t is denoted by S i (t) and the movement of individuals between the states is governed by a set of transition intensities q rs (t, z) (r, s = 1, . . . , R), which may depend on time and a set of (potentially individual-specific) explanatory variables z. The transition intensity q rs represents the instantaneous risk of moving from state r to state s for r = s, q rs (t, z) = lim δt→0 P r(S i (t + δt) = s|S i (t) = r)/δt, and q rr := − s =r q rs . The intensities form an R × R matrix, Q, whose rows sum to zero. Fitting the multi-state model to observations of individuals' disease states enables Q to be estimated.
The Markov assumption is that the future state of the system depends only on its current state and not its history, i.e. q rs (t, z, F t ) = q rs (t, z), where F t is the history of the process before time t. This is equivalent to assuming that all individuals in a state have the same expected outcome regardless of their previous states; for example, an individual with asymptomatic infection will progress to KA with the same probability regardless of whether this is their first or second infection.

Likelihood for multi-state model
The likelihood for the multi-state model given the data is calculated from the transition probability matrix P (u, t), whose (r, s) th entry p rs (u, t) is the probability of being in state s at time t > u, given that the state at time u was r (note that the process may have passed through other states between times u and t). P (u, t) is calculated from the forward Kolmogorov equations [1]: For time-homogeneous Markov processes such as the models we consider (in which the transition intensities q rs are independent of t), P (u, t) = P (t − u), and P (t − u) may be calculated by the matrix exponential which is the solution, in matrix form, of (A1). If the data for the i th individual consists of a series of n i observation times (t i1 , t i2 , . . . , t ini ) and corresponding disease states (S i (t i1 ), S i (t i2 ), . . . , S i (t ini )), the contribution to the likelihood from each pair of successive observed states is which is the (S i (t ij ), S i (t i,j+1 )) th entry of the transition probability matrix P (t) evaluated at t = t i,j+1 − t ij . The full likelihood L(Q) is given by the product of all such L i,j over all M individuals and all transitions The transition intensity matrix Q is estimated by maximising the likelihood in (A4) as a function of Q.

Exactly observed transition times
Equation (A4) gives the likelihood for the model when all observations of individuals' disease states are at arbitrary times. However, some observations in the dataset can be regarded as being at exact transition times with no transitions having occurred since the last observation, such as the date of the end of KA treatment, which is preceded by the observation of the onset of symptoms. For such observations the contribution to the likelihood is ] and then enters state

Censored observations
The contribution to the likelihood from each pair of successive states given in (A3) applies when individual i is known to be in state S i (t i,j+1 ) at time t i,j+1 . However, for some observations in the dataset the individual is known only to be in a certain set of states, since either the ELISA or LST test is missing. For such an observation S i (t i,j+1 ), known to be in the set of states C, the contribution to the likelihood is

Covariates
The dependence of infection and disease progression on characteristics of individuals such as sex and age is potentially very important in designing control interventions for VL. Individual-specific covariates can be incorporated into the multi-state model in the form of a proportional hazards model, by replacing the transition intensity matrix elements by where q (0) rs are the baseline transition intensities, z i is the value of the covariate z for the i th individual, and β rs represents the effect of the covariate (exp(β rs ) is the hazard ratio for a unit increase in z). The modified transition intensities are then used to determine the likelihood, and the likelihood is maximised over the baseline intensities q rs (0) and log-hazard ratios β rs .

Classification of states
The full classification of observations in the data into the different states in the 5-state Markov model, including censored states for missing rK39 ELISA and LST readings is given in Table A1.
Key: + positive for this marker, − negative for this marker, ? missing test Results tables for 5-state Markov model

Model fit
To assess the goodness of fit of the multi-state model without covariates, we calculated the observed and expected numbers and percentages of individuals (prevalences) in each state for the study period from 2002-2004 at 3-month intervals (see Table A7). The expected number of individuals in each state at time t was calculated by multiplying the number of individuals under observation at time t by the initial proportion of individuals in state r and the transition probability matrix for the time interval t, P (t). Table A7 and Figure A3 show that the observed and expected numbers and prevalences match closely for susceptible individuals, KA patients and recovered/dormant individuals (states 1, 3 and 4), but that the predicted number of asymptomatically infected individuals (state 2) is underestimated by the model and the number of deaths (state 5) is over-estimated. The discrepancy in the number of asymptomatic individuals is likely due to variation in the transition rates with time (e.g. the infection rate q 12 changing as asymptomatic infection and KA incidence changed during the study period), which is not accounted for in DOI 10.1186/s13071-015-1136-3 the model; while the discrepancy in the number of deaths is likely due to the large number of individuals that died who are not included in the model fitting (37 out of 64) as their only observation is their date of death.  Comparison with model with separate states for asymptomatic and pre-symptomatic infection To assess whether there is a difference in the duration of asymptomatic infection for individuals who progress to KA (referred to as pre-symptomatics) and those who recover without developing symptoms (referred to as asymptomatics from here on), and whether this affects the estimate of the proportion of infected individuals that develop KA, we fitted the 6-state model shown in Figure A4 with separate states for pre-symptomatics (state 2) and asymptomatics (state 6) to the data. The state space of the model was otherwise the same as in the 5-state model. The transition intensity matrix for this model is   The fitted transition intensity matrix for this model is  Table A8, and are also very similar to those for the 5-state model (Table A2). Pre-symptomatic individuals appear to progress to KA more quickly on average than asymptomatic individuals recover from infection (in 135 days, 95% CI 109-167 days, compared to 159 days, 95% CI 138-183 days), but there is considerable overlap in the 95% confidence intervals for these estimates, so it is not clear that the difference is significant. The mean duration of KA is shorter for the 6-state model than the 5-state model (127 days, 95% CI 113-143 days, as opposed to 147 days, 95% CI 123-160 days), but the 95% confidence intervals are again overlapping. The 6-state model has a negative log-likelihood of 1717.1 and an Akaike information criterion (AIC) value of 3456.1, which is much lower than for the 5-state model (AIC=3541.1), implying that separating pre-symptomatics and asymptomatics into two states does yield a significant improvement in the fit of the model to the data. However, the similarity in the model estimates suggests that grouping pre-symptomatic and asymptomatic individuals together is a reasonable modelling assumption, and does not significantly bias the results of the model fitting.