Estimates of the risk of large or long-lasting outbreaks of Middle East respiratory syndrome after importations outside the Arabian Peninsula.

We quantify outbreak risk after importations of Middle East respiratory syndrome outside the Arabian Peninsula. Data from 31 importation events show strong statistical support for lower transmissibility after early transmission generations. Our model projects the risk of ≥10, 100, and 500 transmissions as 11%, 2%, and 0.02%, and ≥1, 2, 3, and 4 generations as 23%, 14%, 0.9%, and 0.05%, respectively. Our results suggest tempered risk of large, long-lasting outbreaks with appropriate control measures.

Introduction

Clusters of patients infected with Middle East respiratory syndrome (MERS) coronavirus continue to occur in countries throughout the Middle East, where the virus is thought to be endemic in camels (Kayali and Peiris, 2015). While rare, countries elsewhere in the world experience importations from infected individuals traveling from the endemic region (Carias et al., 2016). Most identified importations of MERS from travelers have not resulted in documented transmissions in the destination country (Nishiura et al., 2015); however, the recent large cluster of 186 infected patients stemming from a single introduction in the Republic of Korea (ROK) (Korea Centers for Disease Control and Prevention, 2015) demonstrated that explosive outbreaks are possible.

The ROK outbreak, combined with a non-negligible likelihood of further exportations of MERS from Middle Eastern countries (Carias et al., 2016), is cause for continued concern for importation of MERS to other countries. For public health officials requiring quantitative assessment of the risk posed by incoming infected travelers, it is important to have a nuanced understanding of the full spectrum of possible outcomes, especially when they are highly variable (Fisman et al., 2014); modeling studies can play an important role in this regard.

Recent studies (Nishiura et al., 2015, Kucharski and Althaus, 2015, Chowell et al., 2015) have quantified the variability implied by different data sets of MERS cluster sizes resulting from importation of cases. These analyses found that the data are potentially consistent with high transmission variability associated with the occurrence of superspreading events, similar to what was observed during severe acute respiratory syndrome (SARS) outbreaks in 2003 (Lloyd-Smith et al., 2005). These studies quantified transmission probabilities using a negative binomial offspring distribution within a branching process outbreak model, assuming that every infected individual transmits with an average of R 0 transmissions and dispersion parameter k, where k  < 1 implies high over-dispersion (Lloyd-Smith et al., 2005).

In this paper, we extend the results of the above work to allow the reproductive number R to vary across subsequent generations of transmissions during an outbreak. The ROK outbreak consisted of a large number of transmissions from the initial traveler and from a few patients in the next transmission generation. Then, once local officials determined that a MERS outbreak was occurring and implemented control measures in response, there was an extremely rapid decrease in transmissions such that the entire outbreak was extinguished after three total generations of transmission following the introduction (Korea Centers for Disease Control and Prevention, 2015). This type of differential transmissibility before vs. after implementation of control measures has also been observed during localized outbreaks of SARS (Lloyd-Smith et al., 2005, Wallinga and Teunis, 2004) and Ebola (Toth et al., 2015, Shuaib et al., 2014).

A simple way to model a post-control change in average transmissibility is to use one parameter for the reproduction number in early generations (R 0) and another for later generations (R c, or post-control reproductive number), as assumed in several previous modeling studies of observed outbreaks and public health response for different diseases (Lloyd-Smith et al., 2005, Wallinga and Teunis, 2004, Toth et al., 2015, Chowell et al., 2004). We hypothesized that a model allowing this type of switch would produce a substantially better fit to the data from outbreak clusters caused by MERS importations. Given results from our previous work assessing Ebola importation risk (Toth et al., 2015), we also hypothesized that this model might produce substantially different results for the risk of a very large outbreak compared to a model assuming a single reproductive number across all transmission generations.

Data

We developed a data set of cluster sizes from MERS importations to countries entirely outside of the Arabian Peninsula (Table 1 ); we excluded data from Jordan, the Kingdom of Saudi Arabia, Kuwait, Oman, Qatar, the United Arab Emirates, and Yemen, countries where it was not always clear whether the initial or subsequent cases within clusters acquired infection from exposure to MERS cases or animals (camels). The data were extracted from World Health Organization reports (World Health Organization, 2015) as well as published accounts of individual clusters (Yavarian et al., 2015, Puzelli et al., 2013, Abroug et al., 2014, The Health Protection Agency U. K. Novel Coronavirus Investigation team, 2013). Our data set consists of 31 importation events, of which 23 resulted in no confirmed or suspected transmissions (clusters of size 1) and the other 8 resulted in clusters of size 2–186. Following Nishiura et al. (2015), we also recorded the total number of generations of transmission that occurred after the introduction.

Table 1

Cluster data from reported Middle East respiratory syndrome importations outside the Arabian Peninsula.

Country	Cluster sizea	Transmission generationsb
Algeria	1	0
Algeria	1	0
Austria	1	0
China	1	0
Egypt	1	0
France	2	1
Germany	1	0
Germany	1	0
Germany	1	0
Greece	1	0
Iran	7	3
Iran	2	1
Italy	3	1
Lebanon	1	0
Malaysia	1	0
Netherlands	1	0
Netherlands	1	0
Philippines	1	0
Philippines	1	0
Philippines	1	0
Republic of Korea	186	3
Spain	1	0
Thailand	1	0
Tunisia	2	1
Tunisia	1	0
Turkey	1	0
United Kingdom	3	1
United Kingdom	1	0
United States	2	1
United States	1	0
United States	1	0

Each row represents a unique individual infected traveler to the indicated country.

Cluster size includes the initial infected traveler and any subsequent infected persons epidemiologically linked to that traveler; a cluster of size 1 indicates no known transmission from the traveler in the destination country.

Transmission generations are the maximum number of transmission links from an infected person in the cluster back to the initial traveler.

Methods

For each generation of transmission, we assumed a negative binomial offspring distribution with parameter set , where i is the generation of transmission (i  = 0 from the initial traveler). This assumption results in the following equations.

First, the probability that x independent cases in generation i produce a total of y cases in generation i  + 1 is

Next, given n independent introductions (generation 0), the joint probability of a cluster of total size j consisting of exactly g generations of transmission, under parameter set , is

We used the above equations to evaluate ten different models. In Model 0, we assumed constant parameter values across all generations of transmission, i.e.,. In Models 1a, 1b, and 1c, we assumed the initial patient transmitted with reproductive number R 0 and dispersion parameter k 0, and all subsequent patients transmitted with a post-control reproductive number R c and dispersion parameter k c i.e., . Because we found that allowing k c to range freely in the optimization scheme resulted in wide uncertainty (due to few multi-generation clusters in the data), we chose to test three different assumptions for this parameter. In Model 1a, we assumed that k c  =  k 0; in Model 1b, we assumed k c  = 1, a special case in which the negative binomial distribution reduces to the geometric distribution; and in Model 1c, we assumed infinite k c, another special case in which the negative binomial distribution reduces to the Poisson distribution.

In Models 2a, 2b, and 2c, we assumed that the reproductive number and dispersion parameter switched from R 0 to R c and k 0 to k c after two generations of transmission, i.e., , and made the same three assumptions regarding k c as described above. In summary,

For each of these three parameterizations , we quantified the likelihood of observing the 31 clusters of size j extinguished after g generations using the formula

We compared the maximum likelihood fits under the three models using the Akaike information criterion (AIC), which evaluates model parsimony in determining statistical support for the hypothesized difference in transmission across outbreak generations (Blumberg et al., 2014).

We also developed a model extension to test the robustness of our results against the possibility that there were additional MERS exportations outside the Arabian Peninsula causing clusters that were not detected. If undetected clusters exist, the data set in Table 1 might be biased toward larger cluster sizes, as smaller clusters presumably would be more likely to go undetected.

To quantify the implications of undetected clusters, we made the following assumptions for this part of the analysis. Let N be the number of undetected clusters, and u be the probability that an individual infected patient outside the Arabian Peninsula goes undetected. We assumed that if any one patient in a cluster was detected with MERS, then the entire cluster was detected, due to the intensive contact tracing that would be initiated after the first detection. Under those assumptions, the probability that a cluster of size j would go undetected is . We also assumed that transmission among patients in an undetected cluster was governed by the R 0, k 0 parameters only, under every model, because the presumed mechanism for shifting to R c, k c (implementation of transmission control measures) would only be relevant if detection occurred.

The new likelihood L for a given test value of N undetected clusters is comprised of the product of the joint probabilities that each of the 31 clusters was of the given size and number of generations and was detected, times the probability that N clusters were unobserved; this latter factor includes the probabilities for undetected outbreaks of any size.

We estimated the infinite sum using a partial sum that had converged to six decimal places. The likelihood was maximized for N  = 31 and N  = 93, representing scenarios where 50% and 75% of importation clusters were undetected, respectively, over the parameters u, R 0, k 0, and R c if applicable, for Models 0, 1b, and 2b.

Results

Model 0 produced an MLE of R  = 0.87 (95% CI: 0.46–1.90) and k  = 0.035 (0.016–0.069). Model 1a, assuming a change in the reproductive number after the first generation of transmission, produced an MLE of R 0  = 5.2 (1.7–29.9), R c  = 0.19 (0.05–0.53), k 0   =   k c  = 0.068 (0.031–0.14). Model 1b produced R 0  = 5.5 (1.7–35.4), R c  = 0.14 (0.05–0.29), and k 0  = 0.061 (0.025–0.13). Model 1c produced R 0  = 5.5 (1.7–35.7), R c  = 0.14 (0.05–0.27), and k 0  = 0.061 (0.025–0.13).

Model 2a, assuming the change occurred after the second generation of transmission, resulted in estimates of R 0  = 2.0 (1.0–6.7), R c  = 0.064 (0.007–0.27), and k 0   =   k c  = 0.078 (0.034–0.16). Model 2b produced R 0  = 2.2 (1.0–6.8), R c  = 0.060 (0.008–0.18), and k 0  = 0.076 (0.032–0.16). Model 2c produced R 0  = 2.2 (1.0–6.8), R c  = 0.060 (0.008–0.17), and k 0  = 0.076 (0.032–0.16).

Each version of Models 1 and 2 produced an MLE with substantially higher likelihood and lower AIC than Model 0, the two-parameter (R, k) model previously implemented (Nishiura et al., 2015, Kucharski and Althaus, 2015, Chowell et al., 2015). Of these, models 2b and 2c produced the lowest AIC value (Table 2 ); we chose Model 2b to represent an optimal model under this criterion.

Table 2

Results of fitting models to the cluster data.

	Control generation	Parametersa	log likelihoodb	AIC valuec
Model 0	None	(R, k) = (0.87, 0.035)	−50.6	105.2
Model 1a	1	(R0, k0, Rc, kc) = (5.2, 0.068, 0.19, 0.068)	−45.1	96.3
Model 1b	1	(R0, k0, Rc, kc) = (5.5, 0.061, 0.14, 1)	−44.7	95.5
Model 1c	1	(R0, k0, Rc, kc) = (5.5, 0.061, 0.14, ∞)	−44.8	95.5
Model 2a	2	(R0, k0, Rc, kc) = (2.0, 0.078, 0.064, 0.078)	−44.2	94.5
Model 2b	2	(R0, k0, Rc, kc) = (2.2, 0.076, 0.060, 1)	−44.0	94.0
Model 2c	2	(R0, k0, Rc, kc) = (2.2, 0.076, 0.060, ∞)	−44.0	94.0

For Model 0, the reproductive number R is the average number of transmissions from each individual regardless of the transmission generation; for Models 1a, 1b, and 1c, the initial reproductive number R0 and dispersion parameter k0, apply to the initial traveler only (generation 0), and the post-control reproductive number Rc and dispersion parameter kc apply to any infected persons in generations ≥1; for Models 2a, 2b, and 2c, R0 and k0 apply for both generations 0 and 1, and Rc and kc apply for generations ≥2.

Parameters were optimized according to the shown maximal log likelihood.

AIC = Akaike information criterion, used to determine the optimal model (Model 2b represents an optimal model, with lowest AIC value).

We compared the risk assessment implications of the optimal model against those of other models. The optimal model produces a higher probability of smaller outbreaks across one or two generations of transmission, but a much lower probability of very large outbreaks or of outbreaks exceeding several transmission generations (Table 3 ).

Table 3

Risk assessment implications of each model.

	Control generation	Probability of >(10, 100, 500, 1000) total transmissions	Probability of >(1, 2, 3, 4, 5) generations of transmission
Model 0	None	(3.9%, 1.0%, 0.3%, 0.1%)	(11%, 4.5%, 2.6%, 1.7%, 1.2%)
Model 1a	1	(12%, 1.5%, 0.007%, 0.00002%)	(26%, 3.6%, 0.6%, 0.12%, 0.02%)
Model 1b	1	(11%, 1.6%, 0.011%, 0.00006%)	(24%, 3.3%, 0.5%, 0.07%, 0.01%)
Model 1c	1	(11%, 1.6%, 0.011%, 0.00006%)	(24%, 3.3%, 0.4%, 0.06%, 0.009%)
Model 2a	2	(11%, 1.6%, 0.008%, 0.00002%)	(23%, 14%, 0.8%, 0.05%, 0.003%)
Model 2b	2	(11%, 2.0%, 0.018%, 0.00011%)	(23%, 14%, 0.9%, 0.05%, 0.003%)
Model 2c	2	(11%, 2.0%, 0.018%, 0.00011%)	(23%, 14%, 0.8%, 0.05%, 0.003%)

Probabilities of exceeding selected numbers of total transmissions/generations of transmission after a single importation of Middle East respiratory syndrome, under three different models. Model 3 was the optimal model given the data in Table 1, according to criterion summarized in Table 2.

The results under the assumption of undetected clusters (Table 4 ) show that Model 2b is still optimal according to AIC, although the change in AIC compared to Model 0 becomes smaller as the number of assumed undetected clusters increases. Also, as the number of assumed undetected clusters increases, the optimal model’s estimate of “worst-case” outbreak sizes at the 0.1% or 0.01% probability level move closer to those of the simpler Model 0 (Fig. 1 panels A, C, E). However, the optimal model still produces much lower estimates of the probability of outbreaks lasting several generations across all assumptions for undetected clusters (Fig. 1 panels B, D, F).

Table 4

Sensitivity analysis – results of fitting models to the cluster data given that portion of importation clusters were undetected.

Undetected Fraction	Model	Control generation	Parameters*	log likelihood	AIC value
50%	Model 0	None	(R, k, u) = (0.76, 0.028, 0.46)	−91.9	189.8
	Model 1a	1	(R0, k0, Rc, kc, u) = (2.5, 0.038, 0.23, 0.038, 0.46)	−89.5	187.0
	Model 1b	1	(R0, k0, Rc, kc, u) = (2.7, 0.032, 0.14, 1, 0.46)	−88.7	185.3
	Model 2a	2	(R0, k0, Rc, kc, u) = (0.96, 0.041, 0.078, 0.041, 0.44)	−88.6	185.2
	Model 2b	2	(R0, k0, Rc, kc, u) = (1.5, 0.042, 0.063, 1, 0.46)	−87.7	183.5
75%	Model 0	None	(R, k, p) = (0.62, 0.022, 0.22)	−119.2	244.3
	Model 1a	1	(R0, k0, Rc, kc, p) = (1.1, 0.024, 0.30, 0.024, 0.22)	−118.6	245.2
	Model 1b	1	(R0, k0, Rc, kc, p) = (1.4, 0.019, 0.15, 1, 0.23)	−117.4	242.7
	Model 2a	2	(R0, k0, Rc, kc, p) = (1.4, 0.028, 0.075, 0.028, 0.22)	−117.1	243.0
	Model 2b	2	(R0, k0, Rc, kc, p) = (1.4, 0.026, 0.065, 1, 0.22)	−116.9	241.8

Fig. 1

Projected outbreak risk from a single infected traveler outside the Arabian Peninsula. Model-derived probabilities of an outbreak exceeding a given total number of transmissions (A, C, E) or total number of transmission generations (B, D, F). Dashed = Model 0, assuming no change in transmission parameters (reproductive number and dispersion parameter) across generations; dotted = Model 1b, assuming the transmission parameters change after one generation of transmission; solid = Model 2b (the optimal model), assuming the transmission parameters change after two generations of transmission. Panels A and B use MLE parameters derived assuming there were no undetected importation clusters beyond those listed in Table 1. Panels C and D assuming 50% of importation clusters were undetected. Panels E and F assuming 75% of importation clusters were undetected.

Discussion

We have considered a simple method to assess the statistical support for differential transmission in earlier versus later generations after a new introduction of MERS, based only on outbreak data for the sizes of transmission clusters and total number of transmission generations that produced them. This method demonstrated strong statistical support for assuming a higher reproductive number in earlier generations after a MERS introduction in a non-endemic area.

Projections from the optimal model have important implications for assessing the risk posed by new introductions of MERS. Compared to previous assessments (Nishiura et al., 2015, Kucharski and Althaus, 2015, Chowell et al., 2015) that were similar to those from our Model 0, our optimal Model 2b suggests a higher probability of moderately sized outbreaks (e.g., on the order of 10 total transmissions) across one or two total generations of transmission, but a much lower probability of outbreaks significantly larger than the one in ROK or of outbreaks of any size lasting several generations. These conclusions are robust to assuming that only 50% of MERS importations outside the Arabian Peninsula have been detected. If the non-detection rate is much higher than 50%, then our optimal model would produce closer estimates to previous models for the probability of very large outbreaks, but the conclusion that outbreaks are less likely to last several generations than previous predictions is robust to high rates of non-detection.

The results from all the models we fit to the data suggest very high transmission variability from the index patient (and perhaps also from subsequent patients, depending on the model), as the MLE for the parameter k 0 was less than 0.1 for each model, which indicates even higher over-dispersion than what was estimated for SARS (Lloyd-Smith et al., 2005). The MLE value of k 0 was even lower in the analyses assuming there were undetected clusters, as undetected clusters were likely small, making the ROK outbreak even more extreme compared to the average. The implications of very high initial variability are 1) a high probability of no transmissions from the index patient, even if R 0  > 1; and 2) a relatively high probability of a superspreading event, i.e., an unusually large number of transmissions, if any do occur. For example, using the MLE (R 0, k 0) from our optimal Model 2b (Table 2) there would be 77% chance of no transmissions from the initial traveler, but a 5% chance of more than 12 transmissions and a 1% chance of more than 40 transmissions from the initial traveler.

For public health officials in countries anticipating further introductions of MERS-CoV from travelers, it is important to anticipate the non-negligible possibility of an explosive outbreak in early generations of transmission driven by superspreading. There are several reasons that superspreading might occur from an infected individual, including unusually high levels of viral shedding, long length of infectious period, or high numbers of person-to-person contacts, particularly when numerous contacts coincide with peak timing of infectiousness and/or if contacts have unusual susceptibility, such as hospital patients. Investigations of the MERS superspreading events in ROK suggest that patient symptoms (frequent and vigorous coughing) during close proximity with many others in crowded hospital areas contributed to unusually high numbers of transmissions from certain individuals (Oh et al., 2015).

While the potential for superspreading exists, our results also suggest that a prompt public health response in the early stages of a new outbreak, with efforts to prevent further transmission similar to what has been implemented previously, would most likely reduce the risk of a very large or long-lasting outbreak to negligible levels. Compared to projections from our optimal model, previously published models extrapolate higher probability of MERS outbreaks that are larger or longer-lasting than what occurred in ROK, but those models did not fully incorporate the rapid decline in transmission rate that was achieved in later generations of the ROK outbreak once it had been identified. Nonetheless, any model-based extrapolation beyond the data is subject to potentially wide uncertainty and should be interpreted with caution.

Regardless of the true risk posed by infected travelers, the key elements of a coordinated strategy to mitigate new outbreaks of MERS, as with any emerging infection, are continued awareness, targeted surveillance strategies based on importation risk from travelers, appropriately detailed travel histories of ill patients, pre-positioned availability of laboratory diagnostics, and a strong public health response once a potential case is suspected or recognized.