key: cord-103291-nqn1qzcu
authors: Chapman, Lloyd A. C.; Spencer, Simon E. F.; Pollington, Timothy M.; Jewell, Chris P.; Mondal, Dinesh; Alvar, Jorge; Hollingsworth, T. Deirdre; Cameron, Mary M.; Bern, Caryn; Medley, Graham F.
title: Inferring transmission trees to guide targeting of interventions against visceral leishmaniasis and post-kala-azar dermal leishmaniasis
date: 2020-02-25
journal: nan
DOI: 10.1101/2020.02.24.20023325
sha: 
doc_id: 103291
cord_uid: nqn1qzcu

Understanding of spatiotemporal transmission of infectious diseases has improved significantly in recent years. Advances in Bayesian inference methods for individual-level geo-located epidemiological data have enabled reconstruction of transmission trees and quantification of disease spread in space and time, while accounting for uncertainty in missing data. However, these methods have rarely been applied to endemic diseases or ones in which asymptomatic infection plays a role, for which novel estimation methods are required. Here, we develop such methods to analyse longitudinal incidence data on visceral leishmaniasis (VL), and its sequela, post-kala-azar dermal leishmaniasis (PKDL), in a highly endemic community in Bangladesh. Incorporating recent data on infectiousness of VL and PKDL, we show that while VL cases drive transmission when incidence is high, the contribution of PKDL increases significantly as VL incidence declines (reaching 55% in this setting). Transmission is highly focal: >85% of mean distances from inferred infectors to their secondary VL cases were <300m, and estimated average times from infector onset to secondary case infection were <4 months for 90% of VL infectors, but up to 2.75yrs for PKDL infectors. Estimated numbers of secondary VL cases per VL and PKDL case varied from 0-6 and were strongly correlated with the infector's duration of symptoms. Counterfactual simulations suggest that prevention of PKDL could have reduced VL incidence by up to a quarter. These results highlight the need for prompt detection and treatment of PKDL to achieve VL elimination in the Indian subcontinent and provide quantitative estimates to guide spatiotemporally-targeted interventions against VL.

. PKDL has therefore been recognised as a major 49 potential threat to the VL elimination programme in the ISC 50 (10), which has led to increased active PKDL case detection. 51 Nevertheless, the contribution of PKDL to transmission in 52 field settings still urgently needs to be quantified. 53 Although the incidence of asymptomatic infection is 4 to 54 17 times higher than that of symptomatic infection in the 55 ISC (21), the extent to which asymptomatic individuals con- (Fig. 1A) . The data from this study are fully de-105 scribed elsewhere (8, 31) . Briefly, month of onset of symptoms, 106 treatment, relapse, and relapse treatment were recorded for 107 VL cases and PKDL cases with onset between 2002 and 2010 108 (retrospectively for cases with onset before 2007), and year of 109 onset was recorded for VL cases with onset before 2002. There 110 were 1018 VL cases and 190 PKDL cases with onset between 111 January 2002 and December 2010 in the study area, and 413 112 VL cases with onset prior to January 2002. 113 Over the whole study area, VL incidence followed an epi-114 demic wave, increasing from approximately 40 cases/10,000/yr 115 in 2002 to ≥90 cases/10,000/yr in 2005 before declining to 116 <5 cases/10,000/yr in 2010 (Fig. 1B) . PKDL incidence fol-117 lowed a similar pattern but lagging VL incidence by roughly 118 2yrs, peaking at 30 cases/10,000/yr in 2007. However, VL 119 and PKDL incidence varied considerably across paras (aver-120 age para-level incidences: VL 18-124 cases/10,000/yr, PKDL 121 0-31 cases/10,000/yr, Table S5 ) and time (range of annual 122 para-level incidences: VL 0-414 cases/10,000/yr, PKDL 0-120 123 cases/10,000/yr, Fig. S15 ). ú CI = credible interval, calculated as the 95% highest posterior density interval † risk of subsequent VL/asymptomatic infection if susceptible ‡ based on assumed infectiousness § in the absence of background transmission and relative to living directly outside the case household.

Based on the relative infectiousness of VL and the di erent 151 types of PKDL from the xenodiagnostic data, in the absence 152 of any other sources of transmission, the estimated probability 153 of being infected and developing VL if living in the same 154 household as a single symptomatic individual for 1 month 155 following their onset was 0.018 (95% CI: 0.013, 0.024) for VL 156 and ranged from 0.009 to 0.023 (95% CIs: (0.007,0.013)-(0.018, 157 0.031)) for macular/papular PKDL to nodular PKDL. Living 158 in the same household as a single asymptomatic individual, 159 the monthly risk of VL was only 0.00037 (95% CI: 0.00027, 160 0.00049), if asymptomatic individuals are 2% as infectious as 161 VL cases. 162 The risk of infection if living in the same household as an 163 infectious individual was estimated to be more than 10 times 164 higher than that if living directly outside the household of an 165 infectious individual (hazard ratio = 12.0), with a 95% CI 166 well above 1 (8.3, 16.7). The estimated spatial kernel (Fig. 167  S16 ) around each infectious individual shows a relatively rapid 168 decay in risk with distance outside their household, the risk of 169 infection halving over a distance of 84m (95% CI: 71, 99m).

mission. We assess the contribution of di erent infectious 172 groups to transmission in terms of their relative contribu-173 tion to the transmission experienced by susceptible individuals 174 ( Fig. 2A and Fig. S17 ). The contribution of VL cases was 175 fairly stable at around 75% from 2002 to the end of 2004 176 before decreasing steadily to 0 at the end of the epidemic, 177 while the contribution of PKDL cases increased from 0 in 178 2002 to ≥73% in 2010 (95% CI: 63, 80%) (Fig. S17) . Only a 179 small proportion of the total infection pressure on susceptible 180 individuals, varying between 9% and 14% over the course of 181 the epidemic, was estimated to have come from asymptomatic 182 and pre-symptomatic individuals. Reconstructing the Transmission Tree. By sampling 1,000 193 transmission trees from the joint posterior distribution of 194 the transmission parameters and the unobserved data (as de-195 scribed in Materials and Methods), we can build a picture of 196 the most likely source of infection for each case and how infec-197 tion spread in space and time. Fig. 3 shows the transmission 198 tree at di erent points in time in part of the south-east cluster 199 of villages. Early in the epidemic and at its peak (Figures 3A 200 and 3B), most new infections were due to VL cases. Towards 201 the end of the epidemic, some infections were most likely due 202 to PKDL cases and there was some saturation of infection 203 around VL cases (Fig. 3C) . The inferred patterns of trans-204 mission suggest that disease did not spread radially outward 205 from index cases over time, but instead made a combination 206 of short and long jumps around cases with long durations of 207 symptoms and households with multiple cases. . Arrows show the most likely source of infection for each case infected up to that point in time over 1,000 sampled transmission trees, and are coloured by the type of infection source and shaded according to the proportion of trees in which that individual was the most likely infector (darker shading indicating a higher proportion). Asymptomatic infections are not shown for clarity. S/A = susceptible or asymptomatic, E = pre-symptomatic, I = VL, R = recovered, D = dormantly infected, P = PKDL (see SI Text). GPS locations of individuals are jittered slightly so that individuals from the same household are more visible. An animated version showing all months is provided in SI movie 1. There is considerable heterogeneity in the estimated contri- by each VL/PKDL case is typically less than 1 (Fig. S19A ).

The times after onset of symptoms in the infector at which 336 secondary VL cases become infected are typically longer for 337 PKDL infectors than for VL infectors (Fig. 4B) detected after a longer delay than subsequent cases and there 356 will be some delay in mounting a reactive intervention, such 357 as active case detection and/or targeted IRS around the index 358 case(s), interventions will need to be applied in a large radius 359 (up to 500m) around index cases to be confident of capturing 360 all secondary cases and limiting transmission.

Our results demonstrate the importance of accounting for 362 spatial clustering of infection and disease when modelling 363 VL transmission. Previous VL transmission dynamic models 364 (23, 41-43) have significantly overestimated the relative con-365 tribution of asymptomatic infection to transmission (as up 366 to 80%), despite assuming asymptomatic individuals are only 367 1-3% as infectious as VL cases, by treating the population 368 as homogeneously mixing, such that all asymptomatic indi-369 viduals can infect all susceptible individuals via sandflies. In 370 reality, asymptomatic individuals do not mix homogeneously 371 with susceptible individuals as they are generally clustered 372 together around or near to VL cases (25, 28), who are much 373 more infectious and therefore more likely to infect suscepti-374 ble individuals around them, even if they are outnumbered 375 by asymptomatic individuals. Asymptomatic infection also 376 leads to immunity, and therefore local depletion of suscep-377 tible individuals around infectious individuals. Hence, for 378 the same relative infectiousness, the contribution of asymp-379 tomatic individuals to transmission is much lower when spatial 380 heterogeneity is taken into account.

Nonetheless, our results suggest that asymptomatic indi-382 viduals do contribute a small amount to transmission and 383 that they can "bridge" gaps between VL cases in transmission 384 chains, as the best-fitting model has non-zero asymptomatic 385 relative infectiousness. Superficially, this appears to conflict 386 with preliminary results of xenodiagnosis studies in which 387 asymptomatic individuals have failed to infect sandflies ac-388 cording to microscopy (44). However, historical (12, 45) and 389 experimental (46) data show that provision of a second blood 390 meal and optimal timing of sand fly examination are criti-391 cal to maximizing sensitivity of xenodiagnosis. These data 392 suggest that recent xenodiagnosis studies (11, 44), in which 393 dissection occurred within 5 days of a single blood meal, may 394 underestimate the potential infectiousness of symptomatic 395 and asymptomatic infected individuals. Occurrence of VL in 396 isolated regions where there are asymptomatically infected 397 individuals, but virtually no reported VL cases (27, 47), also 398 seems to suggest that asymptomatic individuals can generate 399 VL cases. However, it is possible that some individuals who de-400 veloped VL during the study went undiagnosed and untreated, 401 and that we have inferred transmissions from asymptomatic 402 individuals in locations where cases were missed. We will in-403 vestigate the potential role of under-reporting in future work. 404

The analysis presented here is not without limitations. As 405 can be seen from the model simulations (Fig. S20) , the model is 406 not able to capture the full spatiotemporal heterogeneity in the 407 observed VL incidence when fitted to the data from the whole 408 study area, as it underestimates the number of cases in higher-409 incidence paras (e.g. paras 1, 4 and 12). There are various 410 possible reasons why the incidence in these paras might have 411 been higher, including higher sandfly density, lower initial lev-412 els of immunity, variation in infectiousness between cases and 413 within individuals over time, dose-dependence in transmission 414 (whereby flies infected by VL cases are more likely to create 415 VL cases than flies infected by asymptomatic individuals (22)), 416 

where K(d) = e ≠d/is the spatial kernel function that determines 7 . CC-BY 4.0 International license It is made available under a author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

is the (which was not peer-reviewed) The copyright holder for this preprint .

Bangladesh (protocol #2007-003) and the Centers for Disease Con- 

• recovered (i.e. treated for primary VL, VL relapse or PKDL, or self-resolved from PKDL, or recovered from asymptomatic

Upon infection, individuals either develop pre-symptomatic infection with probability pI or asymptomatic infection with probability 1 ≠ pI (see Table S2 for values of fixed parameters used in the model 

2 of 37 . CC-BY 4.0 International license It is made available under a author/funder, who has granted medRxiv a license to display the preprint in perpetuity. whereis the rate constant for spatial transmission between infected and susceptible individuals; K(dij) is the spatial kernel function that scales the transmission rate by the distance dij between individuals i and j; " (Ø 0) is a rate constant for additional within-household transmission; 1ij is an indicator function for individuals living in the same household, i.e. between an exponentially decaying spatial kernel and a Cauchy-type kernel in our previous study (1) (the exponential kernel 73 gave a marginally better fit), we use the exponential kernel here: and asymptomatic individuals, we take the relative infectiousness of pre-symptomatic individuals, h0, to be the same as that of 93 asymptomatic individuals (i.e. h0 = h4).

One-hundred and thirty-eight of the 190 PKDL cases underwent one or more examinations by a trained physician to 95 determine the type and extent of their lesions (Table S1 ). Data from a recent xenodiagnosis study in Bangladesh (20) 

Pre-symptomatic -- Thus, individual j's infectiousness at time t is given by

[3]

113 Incubation period. Following previous work (1), we model the incubation period as negative binomially distributed NB(r, p) 114 with fixed shape parameter r = 3 and 'success' probability parameter p, and support starting at 1 (such that the minimum 115 incubation period is 1 month):

We estimate p in the MCMC algorithm for inferring the model parameters and missing data (see below).

VL onset-to-treatment time distribution. Several VL cases with onset before 2002 have missing symptom onset and/or treatment times (only their onset year is recorded), and may therefore have been infectious at the start of the study period. In order to be able to infer the onset-to-treatment times of these cases, OT Õ j = min(R Õ j , D Õ j ) ≠ I Õ j (j = 1, . . . , nI 0 ), in the MCMC algorithm 121 (see below) we model the onset-to-treatment time distribution as a negative binomial distribution NB(r1, p1) and fit to the 122 onset-to-treatment times of all VL cases for whom both onset and treatment times were recorded ( Figure S2A ):

124 to obtain r1 = 1.34 and p1 = 0.38 (corresponding to a mean onset-to-treatment time of 3.2 months). assume that all cases not recorded as having immediate recurrence of symptoms su ered treatment relapse and that the time to relapse follows a geometric distribution Geom(p4) with PMF:

where fitting to the recorded gaps gives p4 = 0.13 (corresponding to a mean time to relapse of 7.9 months). Relapse cases are 156 assumed to be uninfectious from their treatment month to their relapse time and their duration of symptoms upon relapse 157 is assumed to follow the same distribution as the onset-to-treatment time for a first VL episode (Eq. (5)). We assume all 158 relapse cases were treated for relapse before the end of the study, since the latest treatment time for primary VL in a case that 168

169 while the probabilities of pre-symptomatic or asymptomatic infection in month t given susceptibility up to month t ≠ 1 are, respectively:

[12]

Model for initial status of non-symptomatic individuals. As there was transmission and VL in the population before the start The probabilities of each non-symptomatic individual initially present (i.e. with Vj = 0) being susceptible, asymptomatically infected, or recovered from asymptomatic infection at time t = 0 can then be found by calculating the probability of avoiding infection in every month from their birth to the start of the study, summing over the probabilities of being infected in one of the months between their birth and the start of the study and recovering after the start of the study, and summing over the 6 of 37 probabilities of being infected in a month before the start of the study and recovering before the start of the study, respectively: pS 0 (aj) := P(Aj > 0, Rj > 0) = e ≠⁄ 0 a j [13] pA 0 (aj) := P(Aj AE 0,

where aj is the age of individual j in months at t = 0. Since we assume that non-symptomatic individuals who are born, or 183 who immigrate into the study area, after the start of the study (with Vj > 0) are susceptible, for notational convenience we 184 define the probabilities for these individuals as pS 0 (aj) = 1, pA 0 (aj) = pR 0 (aj) = 0. We estimate the historical asymptomatic infection rate, ⁄0, by fitting the model to age-prevalence data on leishmanin skin 186 test (LST) positivity amongst non-symptomatic individuals from a cross-sectional survey of three of the study paras conducted 187 in 2002 (28) (see Figure S4 ). We assume that entering state R corresponds to becoming LST-positive, as LST positivity is 188 a marker for durable, protective cell-mediated immunity against VL (28, 29), and estimate ⁄0 by maximising the binomial With these definitions, the complete data likelihood for the augmented data Z = (Y, X) given the model parameters ◊ = (-, -, ', ", p) is composed of the products of the probabilities of all the di erent individual-level events over all months: the joint posterior distribution of the model parameters ◊ = (-, -, ', ", p) and the missing data X given the observed data Y

[17] 223 We do this using a Metropolis-within-Gibbs MCMC data augmentation algorithm in which we iterate between sampling from 224 the conditional posterior distribution of the parameters given the observed data and the current value of the missing data, 

t=V j qj(t) + qj(T + 1) is a normalising constant to account for the fact that we know that j was not pre- and ", which are non-negative, since there is little information available with which to construct informative priors (Table S3) .

The mean of the prior distribution foris chosen as 100m based on our previous findings (1) . A beta distribution, Beta(a, b) , 266 is chosen as a conjugate prior for the incubation period parameter p, since it is a probability (p oe where -= (-, -, ', "), so p can be updated e ciently in the MCMC by drawing from this full conditional distribution rather than using a random walk Metropolis-Hastings update. 1 ≠ Aj,0) for Rj,0 = T + 1. by repeating the following steps. Note that throughout the following we suppress notation of conditional dependencies in the 320 likelihood terms where they are obvious to maintain legibility. The algorithm also accounts for the fact that some individuals 321 were born or migrated or died during the study when updating the unknown pre-symptomatic infection times and asymptomatic 322 infection and recovery times (using the birth/migration/death times as bounds on the proposed unobserved times), but we 323 omit these details from the following description for simplicity. 

(b) Accept the infection time move with probability

where 

(c) i. If Aj = 0:

if Rj = 0,

if Rj > 0. ii. If Aj = T + 1:

A. If A Õ j = T + 1, then R Õ j = T + 1, so accept immediately as the likelihood does not change.

Step 4(c)iC, except with

. 365 iii. If Aj oe [1, T ]:

A. If A Õ j = 0, follow Step 4(c)iA, but with Q replaced by

Step 4(c)iB except with

Step 4(c)iC but with 6. Update missing treatment times of VL cases during the study: 382 Update the treatment time of the VL case whose treatment time is missing but whose onset time is known, conditional 383 on the treatment time being before their PKDL onset: 384 (a) Propose a new treatment time as

Update the onset and treatment times of all cases who potentially had active VL at the start of the study (t = 1) who N(0, 4) ) " = 0 p1) . 

Update the treatment times of cases who potentially had active VL at the start of the study whose treatment times were 390 not recorded but whose onset times are known, one by one. For each case j: p1) .

9. Update whole relapse period of cases missing both relapse and relapse treatment times: 392 Update the relapse and relapse treatment times of all VL cases who su ered relapse during the study who are missing Step 4 in the above algorithm may appear complicated, but essentially consists of proposing a new asymptomatic infection 

is the empirical covariance of the last k ≠ f (k) + 1 samples offrom the chain, with

the mean of the last k ≠ f (k) + 1 samples; 0 is the initial guess for the covariance matrix, and k0 determines the rate at 442 which the influence of 0 on k+1 decreases (the weight of 0 halves after the first 2k0 iterations). We use k0 = 1000 here.

If f (k) = f (k ≠ 1) (i.e. if k is odd with f (k) chosen as above), an additional observation is added to the estimate of the 444 covariance matrix

if k is even), the new observation replaces the oldest

It has been shown that N (k , 2.38 2 /n -), where is the covariance matrix of the posterior distribution, is the optimal 453 proposal distribution for rapid convergence and e cient mixing of the MCMC chain for symmetric product-form posterior 454 distributions as nae OE, and leads to an acceptance rate of 23.4% (38, 39) . This corresponds to a scaling of c k = 1 in Eq. (23).

However, we are in a context with a large amount of missing data, which is strongly correlated with some of the transmission 456 parameters (see Parameter estimates below), so the posterior distribution is not symmetric, and this scaling is not optimal. 457 We therefore follow (36) and scale c k adaptively as the algorithm progresses to target an acceptance rate of approximately 458 23.4% for updates to -. We do this by rescaling c k by a factor of x k > 1 every time an acceptance occurs and by a factor of 459 x ‹/(‹≠1) k < 1 every time a rejection occurs such that the acceptance rate ‹ approaches 23.4% in the long run,

if proposal is rejected.

In order to satisfy the 'Diminishing Adaptation' condition (40), which is necessary to ensure the Markov chain is ergodic and converges to the correct posterior distribution, it is required that c k tends to a constant as k ae OE. So that the adaptation 463 diminishes as k increases, we use the sequence 465 where m0 is the number of iterations over which the scaling factor x k decreases from 2 to 1.5. Here, we use m0 = 100.

Model comparison 467 We compare the goodness of fit of models with di erent asymptomatic and pre-symptomatic relative infectiousness (between 0% and 2% of that of VL cases), with and without additional within-household transmission, to test di erent assumptions about how infectious asymptomatic and pre-symptomatic individuals are, using DIC (41). DIC measures the trade-o between model fit and complexity and lower values indicate better fit. Since some variables were not observed, we use a version of DIC appropriate for missing data from (42), which is based on the complete data likelihood L(◊; Z) = P(Y, X|◊). This is equivalent to the standard version of DIC for fully observed data except that it is averaged over the missing data:

where D(◊) is the deviance (the measure of model fit), given (up to an additive constant dependent only on the data) by .

[26]

The relative contribution of state X to the infection pressure on the ith VL case at their infection time, i.e. the probability 483 that i's infection source is X (Fig. 2B in the main text), is:

, X oe {A, E, I, P}.

[27]

The probability that the ith VL case is infected from the background transmission is ' ⁄i(Ei ≠ 1)

.

Reconstructing the epidemic Reconstructing the transmission tree. We reconstruct the transmission tree following the 'sequential approach' described 488 in (44). We draw N samples (◊ k , X k ) (k = 1, . . . , N) from the joint posterior distribution from the MCMC, calculate the 489 probability that infectee i was infected by individual j conditional on their infection time Ei and uncertainty in the parameter values and missing data (over the posterior distribution). We use N = 1000 here.

Calculating transmission distances and times. The mean infector-to-VL-infectee distance and mean infector-onset-to-VLinfectee-infection time for each VL and PKDL infector ( Figures 4A and 4B in the main text) are calculated from the sample of N transmission trees by averaging the distances and times from each infector to their VL infectees within each tree, and then averaging these quantities over all the trees in which that VL/PKDL case is an infector:

where · 

[33]

The absolute contribution of each infectious state to the e ective reproduction number at time t is:

where X oe {A, I} denotes the infectious state, and, as described above, in the main text we split the numbers of secondary 517 infections (Rj) arising from VL and PKDL for cases that had both.

To assess the fit of the model and simulate hypothetical interventions against PKDL, we create a stochastic simulation version of 520 the individual-level spatiotemporal transmission model described above. We follow standard stochastic simulation methodology 521 for discrete-time individual-level transmission models (47), converting infection event rates into probabilities in order to 522 determine who gets infected in each month. We assume that an individual's progression through di erent infection states 523 following infection occurs independently of the rest of the epidemic (i.e. is either governed by internal biological processes or 524 random external processes of detection), which enables the simulation of an individual's full infection history from the point of 525 infection.

So that we can simulate durations of PKDL infectiousness, we fit a negative binomial distribution NB(r5, p5) to the observed 527 PKDL onset-to-treatment times and onset-to-resolution times for self-resolving PKDL cases in the data: 

Given these pieces of information, the simulation algorithm proceeds as follows: negative binomially distributed (48). The PMF of a size-biased negative binomial random variable X ú corresponding to 542 X ≥ NB is: and assign a PKDL infectious by drawing from Cat({h1, h2, h3, hu}, p) . 562 ii. Else the individual recovers without developing PKDL, so draw a recovery time: Figure   575 S5 and Table S4 respectively. Based on the deviance distributions and DIC values, the best-fitting model is the model with 576 additional within-household transmission and the highest level of relative pre-symptomatic and asymptomatic infectiousness 577 (both 2% as infectious as VL). Hence, we focus on the output of this model in the main text and below. and ") and incubation period distribution parameter p for the di erent models are shown in Table S4 . The parameter estimates 580 are very similar across the di erent models and vary in the way expected -the spatial transmission rate constantand 581 background transmission rate ' are lower for models with additional within-household transmission (" > 0) and decrease with 582 increasing relative asymptomatic infectiousness h4, and the mode foris slightly larger for models with " > 0 (since a flatter 583 kernel shape compensates for the extra within-household transmission). The posterior distributions for the incubation period 584 distribution parameter p correspond to a mean incubation period of 5.7-6.9 months (95% HPDIs (4.8,6.6)-(6.0,7.8) months). The log-likelihood trace and posterior distributions for the parameters for the best-fitting model are shown in Figure S6 .

The parameters are clearly well defined by the data, as the posterior distributions di er significantly from the weak prior 587 distributions.

The corresponding autocorrelation plots are shown in Figure S7 . The high degree of autocorrelation evident for all the 589 parameters is due to strong correlation between the transmission parameters and the missing data, in particular between the 590 spatial transmission rate constantand the asymptomatic infection times. Figure S8 shows thatis strongly negatively 591 correlated with the mean asymptomatic infection timeĀ. This is expected since a higher overall transmission rate leads to show that there is some negative correlation betweenand ', -and ', and " and p. These correlations are not surprising: the 603 more transmission that is explained by proximity to infectious individuals (the higher -), the less needs to be explained by the 604 background transmission (the lower '); the flatter the spatial kernel (the larger -), the fewer infections need to be explained by 605 the background transmission; and the more infections are accounted for by transmission within the same household (the higher 606 "), the longer the incubation period (the lower p) needs to be (due to long times between onsets of cases in the same household).

The acceptance rate for the transmission parameter updates (Step 1 in the MCMC algorithm) was 23. demonstrate that the data augmentation algorithm works as expected. Figure S10 shows the incidence curve of VL and PKDL 615 cases for the whole study area and the inferred incidence curve of asymptomatic infections (averaged over the MCMC chain).

The number of asymptomatic infections increases and decreases with the number of VL cases as expected given the assumption 617 that the incidence ratio of asymptomatic to symptomatic infection is fixed.

The posterior probabilities that individuals were asymptomatically infected during the study (shown in Figure S11 , with is higher). This is as expected given the structure of the model (the decrease in the risk of infection with distance from an 622 infectious individual encoded in the spatial kernel) and the estimates of the transmission parameters. The examples shown in Figure S12 demonstrate that non-symptomatic individuals' asymptomatic "infection" time posterior distributions (red). Note that asymptomatic "infection" in months 0 and T + 1 = 109, represent asymptomatic infection before the study and no asymptomatic infection before the end of the study, respectively. (A) Individual who migrated into a house with an active VL case from outside the study area in month 53 and therefore had a high initial probability of asymptomatic infection, followed by further peaks in asymptomatic infection risk in months 64 and 69 with the PKDL and VL onsets of two other household members in months 63 and 68 respectively. (B) Individual born in month 2 with a high probability of having avoided asymptomatic infection for the duration of the study. (C) Individual who was 23-years-old at the start of the study with a moderately high risk of having been asymptomatically infected before the study and a small peak in asymptomatic infection risk when a fellow household member had VL onset in month 45. 

The role of case proximity in transmission of visceral leishmaniasis in a highly endemic village in 649

Visceral Leishmaniasis in the Indian Subcontinent: Modelling Epidemiology and Control

Feasibility of eliminating visceral leishmaniasis from the Indian subcontinent: explorations with a set 653 of deterministic age-structured transmission models

Elimination of visceral leishmaniasis in the Indian subcontinent: a comparison of predictions from 655 three transmission models

Risk factors for kala-azar in Bangladesh

Increasing Failure of Miltefosine in the Treatment of Kala-azar in Nepal and the Potential Role of Parasite 658

Drug Resistance, Reinfection, or Noncompliance

Single locus genotyping to track Leishmania donovani in the Indian subcontinent: Application in Nepal

Seasonal and nocturnal landing/biting behaviour of Phlebotomus 662 argentipes (Diptera: Psychodidae)

Studies on seasonal man sandfly (phlebotomus argentipes) contact at night

Biting rhythm & biting activity of phlebotomid sandflies

Seasonal distribution of phlebotomine sand flies-vector of visceral leishmaniasis

Use of RK39 for diagnosis of post kala-azar dermal leishmaniasis in Nepal

Bait preference of phlebotomus argentipes (ann. & 672 brunn.)

The importance of systems ecology: implications of vector-based field studies and quantitative analysis of 674 vector control programs

Spatial Distribution of Phlebotomus argentipes (Diptera: Psychodidae) in Eastern 677

Case Study Evaluating Multispatial Resolution Remotely Sensed Environmental Evidence and Microclimatic 678

Relative abundance of Phlebotominae sandflies with emphasis on vectors of kala-azar

National Vector Borne Disease Control Programme, Guidelines on Vector Control in Kala-Azar Elimination (year?)

Policy recommendations from transmission modelling for the elimination of visceral leishmaniasis in 683 the Indian subcontinent

Infectiousness of visceral leishmaniasis, post kala-azar dermal leishmaniasis and asymptomatic subjects 685 with their capacity to infect sand fly vectors by xenodiagnosis

Quantifying the infectiousness of post-kala-azar dermal leishmaniasis towards sandflies

Clinical and immunological aspects of post-kala-azar dermal leishmaniasis in Bangladesh

Quantification of parasite load in clinical samples of leishmaniasis patients: Il-10 level correlates with 692 parasite load in visceral leishmaniasis

Study of parasite kinetics with antileishmanial drugs using real-time 694 quantitative PCR in Indian visceral leishmaniasis

Relapse after treatment with Miltefosine for visceral Leishmaniasis is associated with increased infectivity of 696 the infecting Leishmania donovani strain

Monitoring of Parasite Kinetics in Indian Post-Kala-azar Dermal Leishmaniasis

Quantification of the natural history of visceral leishmaniasis and consequences for control. Parasites 700 & Vectors

Characterization and identification of suspected counterfeit miltefosine 702 capsules

Loss of leishmanin skin test antigen sensitivity and potency in a longitudinal study of visceral leishmaniasis 704 in Bangladesh

Age trends in asymptomatic and symptomatic Leishmania donovani infection in the Indian 706 subcontinent: a review of data from diagnostic and epidemiological studies

Estimating parameters in stochastic compartmental models using markov chain methods

Bayesian inference for partially observed stochastic epidemics

Bayesian Analysis for Emerging Infectious Diseases

Linking Time-Varying Symptomatology and Intensity of 715 Infectiousness to Patterns of Norovirus Transmission

Coupling and ergodicity of adaptive Markov chain Monte Carlo algorithms

Accelerating adaptation in the adaptive Metropolis Hastings random walk algorithm

Weak convergence and optimal scaling of random walk Metropolis algorithms. The 722 Annals Appl

Optimal scaling for various Metropolis-Hastings algorithms

Examples of Adaptive MCMC

Deviance information criteria for missing data models

The deviance information criterion: 12 years on (with discussion)

Methods to infer transmission risk factors in complex outbreak data

Di erent epidemic curves for severe acute respiratory syndrome reveal similar impacts of control 733 measures

A Bayesian Approach to Quantifying the E ects of Mass Poultry Vaccination upon the Spatial and 735

Temporal Dynamics of H5N1 in Northern Vietnam

Dynamics of the 2001 UK Foot and Mouth Epidemic: Stochastic Dispersal in a Heterogeneous 737 Landscape

The MathWorks, Inc

Julia: A fresh approach to numerical computing

Julia v1.0.5. The Julia Project

 3  37  2  3  NW  998  830  18  3  24  4  4  NW  2185  1768  197  35  124  22  5  NW  934  774  15  2  22  3  6  NW  1640  1363  94  22  77  18  7  SE  604  493  42  7  95  16  8  SE  585  496  8  0  18  0  9  SE  969  809  33  6  45  8  10  SE  701  595  17  3  32  6  11  SE  1300  1080  28  9  29  9  12  SE  2807  2388  102  25  47  12  13  SE  933  816  36  7  49  10  14  SE  446  391  23  3  65  9  15  SE  660  574  15  3  29  6  16  NW  905  762  26  9  38  13  17  NW  2080  1764  75  13  47  8  18  NW  3212  2653  61  11  26  5  19  NW  774  647  62  18  107  31  Total  24781  20798  1018  190  54