In regions where cholera is endemic, the considerable year-to-year variation in the size of its seasonal outbreaks has motivated the development of models to retrospectively explain multiannual cycles in disease incidence and prospectively forecast the occurrence of large cholera events (1–3). Previous studies of endemic cholera have relied largely on data aggregated over relatively large spatial scales (2, 4, 5) and statistical models for recurrent outbreaks in urban centers have so far achieved prediction at lead times of ≤1 mo (6, 7), too short to be of use in early warning systems. Studies have also shown that the El Niño southern oscillation (ENSO) influences multiannual cycles in disease incidence (4, 2, 8, 9) but because such global climate drivers have been expected to act over large spatial scales and to synchronize epidemiological responses across space (the Moran effect) (10, 11) and because the physical and ecological linkages between global drivers and local transmission remain obscure (12–14, 16), the fine-scale spatiotemporal dynamics of endemic cholera have been largely ignored.

However, the rapid growth and pronounced demographic and socioeconomic heterogeneity of the emerging megacities of the modern era (17, 18) lead to significant fine-scale spatial structure with unexplored consequences for disease dynamics. In the case of fecal–oral diseases such as cholera, for which sanitary conditions and access to clean water are determinants of disease risk, heterogeneities on very local scales can be important, a fact that has been appreciated since Snow's mapping of the association between cholera risk and water contamination (19). On the other hand, a body of recent work has established a link between global climate drivers—in particular the ENSO—and cholera dynamics. The existence of key predictors on these two vastly different scales begs the following questions: (i) What linkages mediate the action of global climate drivers at the local scale? (ii) Can we combine data at these two scales to obtain a useful cholera early-warning system?

We formulated a probabilistic model for cholera dynamics within the city of Dhaka and fit it to O1 El Tor cholera hospitalization data broken down into administrative subdivisions called thanas (Fig. 1). The monthly data span the period 1995–2008 and include all cases since the last disappearance of the O1 classical and O139 strains. The cholera case data used here were collected by International Centre for Diarrhoeal Disease Research, Bangladesh (ICDDR,B) hospital located in Dhaka city. This hospital provides free treatment to ~110,000 diarrhea cases annually. A systematic 4% sample (every 25th patient) from 1984 to 1995 and 2% sample (every 50th patient) from 1996 onward were taken from all patients that reported to Dhaka hospital irrespective of disease severity. Fecal specimens were collected from these patients and were cultured for diagnosis of cholera following standard procedures. Inspection of the data suggested the existence of two distinct regions within Dhaka: one comprising the older core districts and the other, the newer urban periphery (Fig. 1). Cholera attack rates in the core districts are higher than in the periphery (Fig. 2 B and C and Movie S1).

Open in a separate window

Fig. 1.

The thanas (administrative subdivisions) of Dhaka. We divided thanas into two groups: older, core thanas (orange) and newer, peripheral thanas (blue).

Open in a separate window

Fig. 2.

Cases of cholera per 10,000 in Dhaka (in black) superimposed to an index for the El Niño southern oscillation (ENSO) (in blue) and annualized flood extent at the level of the country as a whole (red). (Top) The whole city of Dhaka; (Middle) the core, older region of the city; (Bottom) the peripheral, new regions. The interannual variability in the cholera rate is coherent with that of the climate covariates, especially for the older thanas and for the very large El Niño of 1998.

Our model assigns a probability to transitions in cholera rate for each thana each month. We categorize the cholera attack rate in each thana into “no cholera,” “low cholera,” and “high cholera,” using a criterion described in section S1. The model is thus a finite-state Markov chain (20, 21). We allow the transition probabilities to depend on temporal and spatial covariates as well as the state of neighboring thanas. Because all thanas’ transitions must be considered simultaneously, the model is a multidimensional inhomogeneous Markov chain (MDIMC), schematized in Fig. 3. Each month, each thana is in one of three states. In the next month, each thana probabilistically changes its state according to transition probabilities that depend on temporal and spatial covariates as well as the state of neighboring thanas. We enforced the total-probability constraint through choice of parameterization and use of the “barrier method”; we verified that the conclusions summarized below are all robust to alternative parameterization (see section S2 and section S3 for a description of the barrier method and alternative parameterization analysis, respectively).

Open in a separate window

Fig. 3.

Schematic representation of MDIMC model. Thanas of the same hue are governed by the same rules. Each thana can be in one of three states (0, 1, or 2 for no cholera, low cholera, or high cholera, respectively) as indicated by the different shades within each of the two groups. Probability of transitions between states from one month to the next can vary according to season, district group, state of neighboring districts, and climate covariates.

The model captures the inter- and intraannual dynamics of cholera in each thana. Likelihood-ratio tests can be used to test the statistical significance of covariates. Although binning the original numerical data into categories discards some information, it focuses attention on levels of variation that are relevant to public health. Comparison of the binned and original data shows that the percentiles of these two data sets match well and that loss of differentiation is present only for large outbreaks and rare extreme outbreaks (see Figs. S1 and S2 and Movie S2 for an animation of the binned data). A detailed description of the model and specific analyses performed can be found in sections S2–S4, Figs. S3–S5, and Tables S1 and S2.

We identified significant spatial effects at two different scales. The first effect corresponds to the core–periphery distinction shown in Fig. 1. Cholera transmission dynamics (as indicated by model transition matrices) differ significantly between core and peripheral thanas (P < 0.0001). The second effect is more local: We found a significant spatial effect between neighboring thanas (P < 0.05), implying propagation of outbreak risk from thana to thana.

The evidence for two distinct regions within the city raises the question of the interaction between climate forcing and the spatial dimension of disease dynamics. Consistent with previous studies for an earlier period (4, 8), we find that ENSO per se is a significant driver of cholera outbreaks in Dhaka (P < 0.0001). Warming in the Pacific precedes an effect on disease dynamics by 11 mo as indicated by the lag with highest likelihood. Importantly, we find a significantly greater effect of ENSO in the core relative to the periphery (P < 0.05). Thus, even though the city covers only 160 km², there is important heterogeneity in cholera's sensitivity to this global driver.

Previous studies of temporal cholera patterns have suggested that flooding may mediate the effect of ENSO locally (2, 22–24) via a dynamical linkage between sea surface temperature (SST) in the Pacific and regional precipitation in Bangladesh (16). Flooding events occur each summer in this region, vary in magnitude from year to year, and can be extensive in Dhaka. We identify a significant correlation between flooding and cholera cases both during and after the monsoon season (P < 0.01). We further investigated the interaction between the effect of flooding and region of the city and found the correlation varied significantly spatially (P < 0.01). For example, flooding correlates strongly with the midmonsoon season cholera level within the core regions of the city (r = 0.764 for July), but is essentially uncorrelated within the peripheral regions (r = −0.009 for July). Moreover, in the 2 y with most extensive flooding—1998 and 2007—there was earlier onset of the fall cholera season, which suggests that flooding may not only increase cholera attack rates but also hasten the postmonsoon epidemic peak. Statistical evidence for a seasonal shift is limited, however, by the low number of extreme events in these data.

Importantly, the remote forcing of Dhaka cholera by SST in the Pacific at a relatively long lead time (11 mo) affords an opportunity for anticipating outbreaks using model-based forecasts. To assess our model's prediction performance, we took a cross-validation approach, sequentially removing 1 y of data, refitting the model to the remaining data, and then simulating 11 mo into the future to assess forecasting accuracy. Fig. 4 shows that the predicted rate of cases is highly coherent with the observed rate. We checked that the model fulfills the basic requirement for forecasting skill, namely that it outperforms the simpler model given by the seasonal means (25); our model's 11-mo predictions' sum of squared error (SSE) is 79% that of the seasonal mean model. When we remove months of the year with consistently low cholera (December through March), our model's SSE is 75% that of the seasonal model. We also considered a more stringent benchmark based on ENSO. Specifically, we adjusted the seasonal prediction according to the linear relationship between ENSO and the cholera anomaly. Not surprisingly, this benchmark improves on the strictly seasonal model: Its SSE is 84% of that of the seasonal mean model. Our model outperforms the stricter benchmark, with an SSE 12% smaller during months when cholera is present.

Open in a separate window

Fig. 4.

Eleven-month predictions (blue dots) with empirical estimate of 95% confidence interval (yellow shaded region) vs. the observed data (black line). (Left Inset) Distribution of predictions for October 1998. (Right Inset) Distribution of predictions for May 2003.

More relevant to public health is the model's skill at predicting the chance of a large outbreak. We evaluated the performance of our model at this task. Specifically, we used the model to predict the chance of monthly cholera incidence exceeding its 75th percentile in 11 mo. Thus, for example, Fig. 4, Left Inset shows the distribution of model forecasts for October 1998 based on November 1997 data. Although the mean prediction differs from the observation, almost all (>99%) model simulations resulted in large events. Thus, our model would have, with high confidence, predicted a large outbreak in October 1998 due to the large preceding El Niño event. Fig. 4, Right Inset repeats this hindsight analysis for May 2003 (using data from June 2002). Here, we find a reduced, but still large (~87%) probability of large outbreak. By construction, our model tends to underestimate the size of outbreaks: We find that an outbreak probability >26% (and not 50%) should be interpreted as an indication of risk (see section S5 and Figs. S6 and S7).

Using the 75th percentile definition of a large event, we find false positive and false negative rates of ~25%, but that the chance of observing an outbreak when none is predicted is only 10% (10/99 mo). Section S6 contains further discussion of our probabilistic forecasts and their interpretation.

An 11-mo lead time makes a cholera early-warning system based on the MDIMC model and ENSO practically operational. A short lead-time model with an accurate flooding index would allow updating of 11-mo predictions immediately before and during the postmonsoon season.

Because it is annualized, the flooding index considered here is of limited practical application to prediction. However, because most of the flooding occurs during the monsoon season, this index is well correlated with cholera throughout the city in the monsoon season (r = 0.91; Fig. S8). Our results suggest that a complementary short lead-time model based on an accurate monthly flooding index would allow updating of long-range predictions 1–2 mo before the onset of the postmonsoon cholera season. The evidence that the effect of flooding on cholera is fully felt only 4–7 wk after initial flooding (23) further supports this suggestion. Our results suggest too that efforts to forecast flooding itself are likely to prove valuable for cholera early warning.

Because our modeling framework permits assessment of the effects of spatial and temporal variables on disease in a straightforward and effective way, it should prove useful for investigations of other diseases in other large cities. It complements previous approaches to the larger-scale spatiotemporal dynamics of infectious diseases inspired by metapopulation models in ecology (26, 27). Indeed, the approach described here can be viewed as an extension of presence/absence models (26) to more than two states. As with the earlier approaches, the MDIMC framework is simple, flexible, and straightforward to fit by maximum likelihood. The likelihood-ratio test affords a convenient means of discounting model complexity. Our modeling framework also complements more recent efforts to model the detailed spatial spread of cholera in Haiti at the shorter timescales of an individual outbreak following the introduction of the disease into a fully susceptible population (28–30). Models incorporating spatial effects of “gravity” type have proved useful in this context (31, 32); their utility in the context of urban and endemic cholera remains to be seen. More generally, spatial structures more complex than the nearest-neighbor effect we have considered here can be incorporated and tested within the MDIMC framework. On another front, we are currently exploring more systematic methods of identifying spatial heterogeneities of the type described here. Specifically, we are developing statistical approaches for the objective identification of regions sharing similar dynamics. Although more sophisticated models yielding more skillful predictions are without doubt possible for cholera, our results—the most skillful forecasts to date at lead times of practical utility—make unequivocally clear that such forecasts must account appropriately for heterogeneities at local scales.

Interestingly we found that the action of the global driver, ENSO, on cholera in Dhaka is locally differentiated: The sensitive city core acts to couple climate variability to interannual variation of disease and propagates the risk of disease to the rest of the city. The partition of the city into two dynamically distinct regions maps well onto the distribution of demographic parameters closely linked to socioeconomic conditions: The core districts are largely those with the highest population density, the highest number of the poorest kind of housing, and the greatest reliance on municipal tap water (as opposed to wells; Fig. S9). The existence of a sensitive core is consistent, too, with the well-established roles of sanitary conditions and access to clean water as key risk factors for a fecal–oral disease such as cholera. For these reasons, we hypothesize that the flooding-induced breakdown of sanitary conditions is likely the principal mediator of the effect of climate on the postmonsoon peak in the core districts of Dhaka. More generally, our results show that spatial structure within a city is central not only in the fine-scale structure of transmission classically illustrated by Snow's famous map of London's 1854 cholera epidemic (19), but also to the disease's dynamical response to climate forcing at the megacity scale.