Background: A sequence of annual measles epidemics has been observed from January 2013 to April 2017 in the South West Shoa Zone of the Oromia Region, Ethiopia. We aimed at estimating the burden of disease in the affected area, taking into account inequalities in access to health care due to travel distances from the nearest hospital. Methods: We developed a dynamic transmission model calibrated on the time series of hospitalized measles cases. The model provided estimates of disease transmissibility and incidence at a population level. Model estimates were combined with a spatial analysis to quantify the hidden burden of disease and to identify spatial heterogeneities characterizing the effectiveness of the public health system in detecting severe measles infections and preventing deaths. Results: A total of 1819 case patients and 36 deaths were recorded at the hospital. The mean age was 6.0 years (range, 0-65). The estimated reproduction number was 16.5 (95% credible interval (CI) 14.5-18.3) with a cumulative disease incidence of 2.34% (95% CI 2.06-2.66). Three thousand eight hundred twenty-one (95% CI 1969-5671) severe cases, including 2337 (95% CI 716-4009) measles-related deaths, were estimated in the Woliso hospital’s catchment area (521,771 inhabitants). The case fatality rate was found to remarkably increase with travel distance from the nearest hospital: ranging from 0.6% to more than 19% at 20 km. Accordingly, hospital treatment prevented 1049 (95% CI 757-1342) deaths in the area. Conclusions: Spatial heterogeneity in the access to health care can dramatically affect the burden of measles disease in low-income settings. In sub-Saharan Africa, passive surveillance based on hospital admitted cases might miss up to 60% of severe cases and 98% of related deaths.
This study was conducted in the South West Shoa Zone of the Oromia Region in Ethiopia (Fig. 1a), with an estimated population of 1,341,702 inhabitants in 2014, of whom 50.3% were men and 49.7% were women. The main hospital is located in Woliso town, 114 km southwest of the capital Addis Ababa, representing the nearest hospital for 521,771 individuals living within an area of 30 km radius from Woliso town (53,065 inhabitants). The hospital has 200 beds with an annual average bed-occupation rate of 84%; single-patient air-borne infection isolation rooms are not available in the hospital. Epidemiological evidences: a Study area and spatial distribution of woredas. b Age distribution of measles patients hospitalized at the Woliso hospital between January 2013 and April 2017. The inset shows the estimated measles seroprevalence by age, as obtained on the basis of model estimates. c Time series of case patients recorded during the study period, overall, and in most affected woredas. The inset shows the cross correlation in the timing of epidemics in Woliso and most rural areas. d Cumulative incidence of hospitalizations per 10,000 individuals (h) by woreda/kebele and distance from Woliso hospital (d). The solid line represents estimates obtained by the negative binomial regression model; the shaded area represents 95% CI Data on age, sex, residence at woreda (i.e., district) and kebele (i.e., neighborhood) level, date of hospital admission, and death/discharge of measles case patients from 2013 to 2017 were obtained from the registers of Woliso hospital. Incidence of hospitalizations by woreda and kebele were calculated by assuming population projections for the 2014, based on the 2007 census conducted by the Central Statistical Agency of Ethiopia (Table 1) [15]. Travel distances to the Woliso hospital for different kebeles and woredas were obtained from administrative hospital records on distances of all health posts and largest villages distributed in the main hospital’s catchment area (see Table 1). The case fatality rate (CFR) for hospital admitted cases was calculated as the percentage of fatal cases among measles patients recorded. Routine vaccination coverage for this area was derived from administrative records: on average, 88% of children are routinely vaccinated against measles at 9 months of age. Two immunization campaigns were conducted in the area from May 29 to June 5, 2013, and from March 13 to March 20, 2017, targeting children 9–59 months of age [16]; the achieved vaccination coverage is unknown. In 2016, the vaccination status of case patients was assessed for 295 children in the age group 9 months to 5 years. Measles cases patients. Epidemiological characteristics of measles cases admitted to Woliso hospital (South West Shewa Zone, Oromia Region, Ethiopia) from January 1, 2013, to April 9, 2017 Patients’ records related to different illness conditions recorded at the Woliso hospital between 2014 and 2016 were considered to estimate hospitalization incidence over time and to assess differences in the access to health care and related outcomes with respect to travel distances from the hospital. Collected data consisted of routine health data and medical records, were encrypted and anonymous, and did not contain any information that might be used to identify individual patients; therefore, the study did not require informed consent. Synchrony in the timing of epidemics across different woredas was assessed by calculating the cross-correlation of time series at different time lags. The aim of this analysis is twofold: (i) to evaluate whether the observed seasonal pattern is an artifact of averaging asynchronous local epidemics and (ii) to support the hypothesis that observed measles cases were the result of a unique synchronous epidemic with similar epidemiological characteristics across different woredas. The baseline analysis combines results of a dynamic transmission model, calibrated on the time series of hospitalized measles cases occurring between 2013 and 2017, with a spatial regression analysis, providing estimates of the measles hospitalization rate at different distances from the Woliso hospital. We restricted the analysis to measles cases from Woliso, Wonchi, Ameya, and Goro woredas, which represent the main hospital catchment area, consisting of 521,771 inhabitants and accounting for 83.1% of recorded case patients. Under the assumption of homogeneous mixing transmission, the baseline model provided estimates of the basic reproductive number (R0), the age-specific immunity profile, and the average measles incidence in the considered area. The estimated total number of infection cases in the population was disaggregated into smaller spatial units (woredas and kebeles), by assuming the same incidence rate across all spatial units and proportionally to the population size of each spatial unit. A regression model was applied to counts of observed hospitalized cases in each spatial unit to estimate the corresponding hospitalization rate; distance from the hospital was used as the independent variable and the estimated total number of cases in each spatial unit as offset. Obtained results were used to quantify the hidden burden of measles disease. In the rest of this section, we detail the dynamic transmission model, the performed spatial analysis, how we calculated the hidden burden of disease, and the performed sensitivity analyses. Measles transmission dynamics between 2013 and 2017 is simulated through a deterministic, non-stationary, age-structured transmission model. In the model, the population is stratified in 86 1-year age classes, according to available data on the age distribution of the Ethiopian population in 2013 [17]. The crude birth rate of the population is 0.0325 years−1; individuals die according to age-specific mortality rates as reported between 2013 and 2015 and reflecting a crude mortality rate of 0.0083 days−1 [17]. The population of any age a is divided into five epidemiological classes: individuals protected by maternal antibodies (Ma), susceptible individuals (Sa), exposed individuals (Ea), infectious individuals (Ia), and individuals who acquired immunity against measles through either vaccination or natural infection (Ra). We assume that newborn individuals are protected against measles infection for 6 months on average by the passive transfer of maternal immunity [1], after which they become susceptible to the infection. Susceptible individuals can acquire infection after contact with an infectious individual under the assumption of homogeneous mixing and become exposed without symptoms; at the end of the latent period, lasting 7.5 days on average, infectious individuals can transmit the infection for 6.5 days on average; the resulting generation time is 14 days [18]. After recovery, individuals are assumed to gain lifelong immunity. Newly infected individuals are hospitalized with a certain, age-independent, probability ph, representing the average hospitalization rate in the main hospital catchment area. Seasonal variations in the transmission rate are considered: during school holidays, overlapping with the rainy season [14], the transmission rate is decreased by a factor r. Routine vaccination of children is simulated at 9 months of age [3] with homogenous coverage across woredas at 88%. The latter estimate was obtained by administrative records on infant immunization occurring between 2013 and 2016 in the main hospital catchment area. Vaccine efficacy at the first dose of routine administration is assumed at 85% [19]. The follow-up campaigns conducted in 2013 (from May 29 to June 5) and in 2017 (from March 13 to March 20), targeting children 9–59 months of age [16], are also considered. The coverage of the 2013 supplementary immunization activities (SIAs), cS, was estimated among free model parameters. Vaccine efficacy during SIAs is assumed to be 95% [19]. Epidemiological transitions are described by the following system of ordinary differential equations: where t represents time and a the individuals’ chronological age; b(t) and d(t,a) are the crude birth and the age-specific mortality rates at time t; 1/μ is the average duration of protection provided by maternal antibodies; 1/ ω and 1/γ are the average duration of the latent and the infectivity periods; cR(t, a) and cS(t, a) are the coverage associated with the first-dose routine vaccination and SIAs for individuals of age a, at time t; εR and εS represent the vaccine efficacy associated with routine vaccination of infants and SIAs. Specifically, cS denotes the vaccinated fraction of individuals who were not yet immunized by natural infection or routine programs. N(t) and H(t) represent the total population of the hospital main catchment area and the cumulative number of hospitalized measles cases at time t; ph is the fraction of measles infections that are hospitalized, and β(t) is the measles transmission rate defined as follows: At the end of the year, the chronological age of individuals is incremented by 1. The number of hospitalized measles cases in a time interval [t1,t2] is computed as H(t2) − H(t1). Model estimates were obtained by simulating measles transmission between January 1, 2013, and March 20, 2017. Simulations are initialized on January 1, 2013. As the result of past natural infection and immunization campaigns, only a fraction s0 of the population is assumed to be susceptible to the infection. The age distribution of susceptibles at the beginning of 2013 was assumed to mirror the age distribution of hospitalized cases between January 2013 and March 2017. Specifically, the initial fraction of susceptible and immune individuals in each age group are Sa(0) = Nas0Za/∑a=085Za and Ra(0) = Na − Sa(0), respectively, where Na is the number of individuals of age a at the beginning of 2013 in Woliso, Ameya, Goro, and Wonchi [17] and Za is the observed total number of hospitalized measles cases of age a. Free model parameters (s0, β, rβ, ph, cS) were calibrated using a Markov Chain Monte Carlo (MCMC) approach based on the negative binomial likelihood of observing the weekly number of hospitalized case patients reported between January 1, 2013, and the beginning of the 2017 SIA. The scale parameter defining the negative binomial distribution was jointly estimated with other free parameters within the MCMC procedure. Details are provided in the Additional file 1. The fundamental quantity regulating disease dynamics is the basic reproduction number (defined as R0 = 〈β〉/γ, where 〈β〉 is the average of β(t) over the year), which represents the average number of secondary infections in a fully susceptible population generated by a typical index case during the entire period of infectiousness. The larger the R0, the higher the disease transmissibility. If R0 > 1, the infection will be able to spread in a population. Otherwise, the infection will die out. For endemic diseases like measles, R0 provides insights into the proportion p of population to be successfully vaccinated to achieve disease elimination; the equation p = 1–1/R0 is widely accepted (e.g., [5, 18, 20]). For instance, if R0 = 10, at least 90% of children have to be routinely immunized to eliminate the disease. A negative binomial regression was used to study the relationship between incidence of hospitalization by kebeles/woredas and distance from Woliso hospital. Specifically, the observed number of hospitalized cases from each spatial unit is the response variable, the distance from the hospital is the independent variable, and the estimated total number of measles cases in each spatial unit (as estimated by the transmission model) is used as the offset. Detailed origin of patients at the kebele level was used to better identify the travel distances for patients living within the Woliso woreda, where the hospital is located (Table 1). In the negative binomial regression, we assume that counts of hospitalized cases hi (the response variable) associated with a given location i are distributed as a negative binomial of mean μi determined by the number of infection in the location ci (the offset) and the distance of location from the hospital di (the regressor) as follows: where b1, b2 are unknown parameters that are estimated from the observed hospitalized cases hi. In order to take into account the uncertainty on incidence estimates obtained with the dynamic model, 10,000 draws from the posterior distribution of incidence estimates associated with 10,000 samples of the posterior distribution of free model parameters were considered to generate a distribution of regression model fits. Obtained results therefore account for the combined uncertainty due to the regression model and the dynamic transmission model. We investigate the spatial variation in the incidence of hospitalized patients in the population as a consequence of different illness conditions. The aim is to characterize the relationship between hospitalization and distance from the hospital. The relative risk of being hospitalized at different distances from the hospital was computed by considering the incidence of hospitalization in each kebele/woreda divided by the incidence of hospitalized cases from Woliso town. The relative risk was fitted by an exponential function using distance as the independent variable (i.e., by fitting a linear model to the logarithm of the relative risk without intercept). Finally, a proportional test was used to assess possible statistical differences in the case fatality rate at hospital between cases coming from different sites. Persons living in Woliso town do not have distance barriers to access to the Woliso hospital. The probability of severe disease after measles infection was therefore computed as the fraction of measles patients from Woliso town that have been hospitalized for two nights or more among all measles infections estimated by the transmission model for this spatial unit. For severe cases, we indicate here those cases that from a clinical point of view are physiologically unstable and require supportive care (fluid resuscitation, oxygen, etc.) that can be provided only inside a well-resourced hospital. The resulting probability of developing severe measles illness ps was used in combination with the estimated number of measles infections at different kebeles and woredas ci to estimate the potential number of severe cases occurring at different distances from the hospital as psci. For each considered spatial unit i, missed severe cases were computed as the difference between the estimated number of severe cases and the number of patients recorded at the hospital, namely mis=psci−hi. Missed severe cases were considered untreated and counted as additional deaths. The overall number of deaths caused by measles was estimated as the sum of missed deaths and measles deaths observed among hospital admitted patients. Averted deaths due to hospital treatment were estimated by considering all severe cases psci as counterfactual deaths that would have occurred in the absence of adequate treatment. A variety of sensitivity analyses were conducted to evaluate to what extent some crucial assumptions made in the above described analysis may affect the obtained results. We evaluated whether the assumption of decreased transmissibility during school holidays (or rainy season) is necessary to explain the observed pattern, by fitting a model with constant transmission rate against the time series of measles hospitalized cases. Since the fraction of immunized individuals during the SIA in 2013 is unknown, we also considered two alternative models with cS = 0 (SIA not conducted in 2013 in the considered area) and cS = 0.92 (the highest coverage reported for past campaigns, namely 92% [3]). We explored whether the assumption of homogeneous mixing, consisting in applying the same transmission rate to all age groups, can affect the model ability in reproducing the observed epidemiological patterns. To do this, we fitted the time series of cases with a transmission model encoding age-specific contacts as recently estimated for Ethiopia by Prem et al. [21]. In this case, increased mixing in schools corresponds to higher transmission rate among school-age children. Models’ performances were assessed through the Deviance Information Criterion (DIC). A sensitivity analysis was also conducted by fitting a transmission model to the time series of measles cases observed in Woliso, Wonchi, Ameya, and Goro separately. Specifically, a single epidemic was simulated in the four woredas simultaneously, by assuming the same initial conditions and by assuming that populations from different locations mix homogeneously. All epidemiological parameters were assumed to be equal across different woredas, but a different hospitalization rate was considered for each woreda. An additional sensitivity analysis was performed to test whether estimates on the spatial variation of the hospitalization rates change when patients recorded from all woredas of the South West Shoa Zone are considered or when patients’ sex is considered. Finally, estimates on the overall number of measles deaths and on the overall case fatality rate were estimated by relaxing the assumption that all missed/untreated severe measles cases die. Details are provided in Additional file 1.