Spatial variability of Middle East respiratory syndrome coronavirus survival rates and mortality hazard in Saudi Arabia, 2012–2019

MERS-CoV data

All publicly available laboratory-confirmed cases of MERS-CoV reported between 13 June 2012 and 30 April 2019 were compiled by the Saudi Ministry of Health via their Command and Control Center (CCC) website (Ministry of Health of Saudi Arabia, 2019) and WHO (WHO, 2019). With the emergence of MERS-CoV in Saudi Arabia, the CCC was launched in 2014 to take immediate proactive action against future public health challenges by way of continuous coordination and monitoring. For the current study, the MERS-CoV data was thoroughly reviewed, and a range of checks was performed to ensure the data were consistent, complete, and fit for the study. A data dictionary was developed for the variables for each case of MERS-CoV, including gender, age, underlying comorbidity, exposure to camels, healthcare status, outcome, date of symptom onset, date of death, city of residence, province and region. The variables were categorized according to gender (male/female), age (<20, 20–40, 41–60, or ≥61), comorbidity (yes/no), exposure to camels (yes/no), HCW (yes/no), regions (five categories as shown below), provinces (13 categories as shown below) and outcome (alive/dead). The data contained no reported comorbidity variable for 58 patients and no exposure to camel variable for 939 patients (713 were unknown or missing a source of infection and 226 were primary or secondary cases without camel exposure status). The MERS-CoV data had been reported by Saudi city, so the data were aggregated for the following 13 provinces, comprising five regions: Riyadh and Qassim (Central region); Tabuk, Jouf and Hail (North region); Assir, Baha, Jazan and Najran (South region); the Eastern province (East region); and Makkah and Madinah (West region). A spatial database of MERS-CoV incidence and mortality in Saudi Arabia was created using geographical information systems (GIS).

Statistical analysis

A comparative epidemiological analysis using a Chi-square test was undertaken to determine the significant differences between the demographic and clinical characteristics of the survival and mortality rates of MERS-CoV cases at the regional level in Saudi Arabia. The variables used in the analysis included age and gender and whether the patients were HCWs, demonstrated underlying comorbidities, and had been exposed to camels. The statistical analyses were performed using MedCalc software, Version 19.0.6 and a P-value < 0.05 was considered statistically significant.

Survival analysis and hazard modeling

Survival analysis aims to model and analyze time-to-event data, that is, a dataset of the time that an event occurred as its endpoint. In the current study, two MERS-CoV parameters, the time of symptom onset and the date of death, were evaluated for 14 days and 45 days after the appearance of symptoms. The incubation period of MERS-CoV, the elapsed time between contracting the virus and displaying symptoms, ranged from 2 to 14 days (Ahmed, 2018). The 14-day period is medically important because it indicates the point at which most MERS-CoV patients typically experience symptoms. The 45-day period was chosen to exhibit the progression of MERS-CoV for an advanced period. The 14-day and 45-day MERS-CoV survival rates were estimated using the Kaplan–Meier estimator at both the national and regional levels in Saudi Arabia. The log-rank test was used to compare survival curves according to the demographic and clinical factors affecting MERS-CoV patients. However, the log-rank test delivers a statistical but not clinical appraisal of the effect of the factors (Bradburn et al., 2003). Thus, the Cox proportional-hazards (CPH) model was applied to assess the survival of MERS-CoV patients regarding several concurrent factors and estimate the extent of the influence of each constituent factor. The CPH model is broadly applicable and is the most widely used multivariate approach for survival analysis (Bradburn et al., 2003). The CPH model examines the relationships of covariates to the time-to-event as expressed by the hazard function. The CPH model is interpreted using HRs, which are expressed as the ratio of the estimated hazard function under two different values of a covariate (George, Seals & Aban, 2014). An HR above one indicates that the event (in this case, death resulting from MERS-CoV) has a higher likelihood of occurring and that the covariate is positively related to the event probability and negatively related to the duration of survival. An HR less than one indicates that an event has a lower probability of occurring, and an HR of one indicates that the covariate does not influence the hazard of an event.

The probability of the endpoint (i.e., time to death) is termed the hazard, and the hazard is modeled as a group of predictor variables. The predictors (risk factor) are dichotomous and coded as 1 if present (died) and 0 if absent (alive). The value of the hazard ratio can then be inferred as the instant relative risk of an event at any time for a patient with the risk factor present compared to a patient with the risk factor absent, given both patients are equivalent for all other predictors. A variable is removed from the model if its related significance level is larger than the P-value (0.05). Survival time is the variable of the time to reach the event of interest (the outcome, in our case, death). The endpoint variable represents the occurrence of the event (i.e., death) and comprises patient cases that have reached the endpoint (coded as 1) or cases that have not reached the endpoint (coded as 0) as they are removed from the study (in our case, at 14 and 45 days). The predictor variables are the variables we expect to predict survival time (gender, age, underlying comorbidity, exposure to camels, healthcare status, region and province) (MedCalc, 2020). When building our model, we applied a filter for particular regions and provinces.

Using MedCalc software, the CPH model was applied to estimate the HRs and 95% confidence intervals (CIs) for the predictors of the two primary endpoints: 14-and 45-day mortality. Both univariate and multivariate analyses were applied using the CPH model. The former describes the MERS-CoV survival rate regarding the factor under examination but disregards effects from any other factors; whereas, the latter simultaneously considers several covariates that may affect patients’ prognoses to enable better clinical assessment. The variables used in the CPH model included a patient’s age and gender and whether the patient was an HCW, had a pre-existing illness, and had been exposed to camels.

The MERS-CoV mortality risk factors for the 14-and 45-day post-symptom periods were estimated for national, regional and provincial levels in Saudi Arabia, and the spatial variations were assessed according to region, province and mortality hazard factors using both univariate and multivariate CPH model analyses. For each post-symptom period, ten CPH models were built: national, each region and province separately, incorporating regions as a covariate, and incorporating provinces as a covariate for both the univariate and multivariate model analyses. The regional and provincial factors characterized the geographical representation, following a previously used approach (Ayele, Zewotir & Mwambi, 2017) for a survival analysis of mortality using the CPH model.

Spatial analysis

The spatial variations of the risk factors associated with MERS-CoV mortality were assessed using univariate local spatial autocorrelation and multivariate spatial clustering analyses. Spatial autocorrelation measures the association within values across locations (Getis, 2008). It is the dependance of a particular variable’s values on the values of the same variable recorded at nearby locations (Fortin & Dale, 2005). We applied the Getis–Ord Gi* statistic, a local spatial autocorrelation that detects the locations of statistically significant hot spots (the clustering of high values) and cold spots (the clustering of low values) within the setting of nearby features (Ord & Getis, 1995) using Esri ArcGIS 10.7.1. Given a set of provinces and associated risk factor values, the Gi* statistic assesses whether the expressed pattern is high-value cluster or low-value cluster (Mitchell, 2005). To consider the spatial interaction between the outcome and covariates examined in this study, we applied the Bivariate Local Moran Index via GeoDa package. In principle, this index identifies the association between the value for one variable at a location and the average of the neighboring values for another variable, that is, its spatial lag (Anselin, Syabri & Kho, 2006).

Clustering is an unsupervised machine learning approach that involves determining natural groupings within data (Duque, Ramos & Suriñach, 2007), and the K-means algorithm is one of the most popular methods for clustering data (Spielman & Thill, 2008). We used multivariate spatial clustering analysis to identify clusters of Saudi provinces according to patients’ characteristics variables for MERS-CoV mortality using the K-means algorithm via Esri ArcGIS. The resulting maps created distinct clusters of variables (demographic, occupational and clinical patient characteristics for MERS-CoV mortality); whereby, the provinces that were part of a cluster were as similar as possible and the provinces between clusters were as dissimilar as possible. Each variable’s contribution, evaluated using r² values, was used to differentiate the provinces, which revealed how much of the variation in the original data was kept after the clustering process: the higher the r² value for a specific variable, the better that variable was for discerning the provinces (Pimpler, 2017).

Characteristics of MERS-CoV patients

The analysis of the relationship of MERS-CoV mortality to patients’ characteristics is shown in Table 1. A total of 2,037 laboratory-confirmed cases of human infection with MERS-CoV were reported in Saudi Arabia between 13 June 2012 and 30 April 2019. Among these, 760 cases died (37.30%), with an overall crude mortality rate of 2.37 per 100,000 population. A higher number of deaths resulted from MERS-CoV cases involving males (n = 558; 73.42%) compared to females (n = 202; 26.58%); a statistically significant difference (p = 0.0007). The overall median age of the deceased cases was 62 years (interquartile range = 50–73), ranging from 2 to 114 years, while the mean age of the cases was 60.36, with a standard deviation of ±17.99.

The greatest number of deaths occurred among elderly patients (age ≥ 61), representing 57.13% (n = 434) of the total deceased cases, which was statistically significant (p < 0.0001). MERS-CoV mortality was higher among non-HCWs (n = 729; 95.93%) than among HCWs (n = 31; 4.07%), and this difference was statistically significant (p < 0.0001). MERS-CoV patients with underlying comorbidities had a higher risk for mortality compared to those without (n = 696 (94.18%) vs. n = 43 (5.82%); p < 0.0001). A lower mortality rate was found among those exposed to camels compared to those not exposed, although this difference was not statistically significant (p = 0.3651). The Central region reported 366 deaths resulting from MERS-CoV, accounting for 48.16% of the total number of MERS-CoV deaths, followed by the West region (n = 220, with 28.95% mortality) and the East region (n = 98, with 12.89% mortality), but the variation across regions was not statistically significant (p = 0.7706).

Spatial pattern of MERS-CoV survival

Figure 1 displays the results of the spatial pattern of MERS-CoV survival in Saudi regions. Approximately 59.30–65.25% of MERS-CoV cases across the Saudi regions were not fatal, but 34.75–40.70% of the cases resulted in death, as shown in Fig. 1A. The distribution of demographic, occupational, and clinical characteristics of MERS-CoV survival is shown in Figs. 1B–1F. MERS-CoV mortality was more frequent in males compared to females across all regions (ranging from 1.42% to 12.47% for males and 0.20–5.50% for females) out of the total number of MERS-CoV cases, as shown in Fig. 1B. MERS-CoV mortality was significantly higher among elderly patients aged ≥61 for all regions (ranging from 0.79% in the North region to 9.97% in the Central region), as shown in Fig. 1C. MERS-CoV mortality was higher among non-HCWs (1.57–17.33%) than HCWs (0.10–0.64%) across all regions (see Fig. 1D). A greater number of patients with underlying comorbidities (1.52–17.53%) died compared to those without (0.10–0.96%) across all regions (see Fig. 1E). A higher mortality rate was found among MERS-CoV patients who had not been exposed to camels (2.46–23.32%) than those who had been exposed (0.36–6.56%) for all regions, as shown in Fig. 1F.

The spatial distribution of the total MERS-CoV mortality cases varied across the provinces, with most cases reported from Riyadh and Makkah followed by the Eastern and Qassim provinces (see Fig. 2A). The spatial variation patterns of the risk factors, based on the demographic, occupational, and clinical characteristics of the MERS-CoV patients, also indicated spatially varied effects and were mainly elevated in Riyadh, Makkah, Eastern and Qassim provinces but were degraded in the northern and southern provinces (see Figs. 2B–2H and Appendix Figs. 1A–1E). The univariate local spatial autocorrelation Getis–Ord Gi* analysis suggested that there were statistically significant hot spots in Riyadh province for all risk factors except for age <20, for which there were no hot spots across the entire country. Madinah province was a hot spot for the age 20–40 factor, and Riyadh and Madinah provinces were hot spots for the exposed to camel factor. Statistically significant cold spots were observed in Jouf province for all risk factors and detected for age <20, age 20–40, and patients without comorbidities in Jazan province. The bivariate spatial analysis with Bivariate Local Moran Index showed no statistically significant spatial autocorrelation between the examined covariates and MERS-CoV mortality attributable to <20 years and non-comorbidity in all provinces. However, the identification of such an autocorrelation for all the other covariates was observed in Riyadh province. The multivariate spatial clustering analysis applied here optimized within-cluster similarities and between-cluster differences for Saudi provinces based on the risk factors with subcategories (12 variables) for demographic, occupational and clinical patient characteristics of MERS-CoV mortality. The results identified three spatial clusters of provinces, with the r² for variables ranging between 0.85 and 0.97. Clusters were found in Riyadh province (Cluster 1); Eastern, Makkah and Qassim provinces (Cluster 2); and the remaining provinces in the north and south of the country (Cluster 3), see Appendix Fig. 1F.

Survival analysis of MERS-CoV

The mortality rate was 25.52% at 14 days, 32.35% at 45 days and 37.30% overall. The survival rates for MERS-CoV at 14 and 45 days were estimated using the Kaplan–Meier method. These estimated survival rates were 70.20% (at 14 days), 57.70% (at 45 days) and 56.50% (overall). Figs. 3 and 4 present the Kaplan–Meier curves for regional, demographic, occupational and clinical factors for 14-and 45-day mortalities, respectively. The Central region (69.30%) had a lower 14-day survival rate than the East (70.10%), north (70.20%), west (71.40%) and south (72.70%) regions, as shown in Fig. 3A. For 45 days, the East region (56%) had a lower survival rate than the Central (57.30%), west (58.70%), south (60.30%) and north (70.20%) regions (see Fig. 4A). However, the differences in survival rates across all Saudi regions were not statistically significant (p = 0.9045 and 0.9279 for 14-day and 45-day survival, respectively). Lower survival probabilities were found among males compared to females (see Figs. 3B and 4B) (14 days = 67.80% vs. 75%, respectively, and 45 days = 54.10% vs. 65.20%, respectively). Similarly, lower survival probabilities were found among the elderly compared to younger people, see Figs. 3C and 4C. Survival rates were also lower among non-HCWs than HCWs (see Figs. 3D and 4D): 14 days = 63.60% vs. 96.20%, respectively, and 45 days = 50.30% vs. 91.30%, respectively. Survival rates were lower among patients with comorbidities than among those without (see Figs. 3E and 4E): 14 days = 60.60% vs. 93%, respectively, and 45 days = 47.50% vs. 89.60%, respectively, and for patients who were not exposed to camels compared to those who were exposed (see Figs. 3F and 4F): 14 days = 69.40% vs. 73.10%, respectively, and 45 days = 54% vs. 58.20%, respectively.

Overall, the 14-day Kaplan–Meier curves for the factors of gender, age, HCW, presence of comorbidities, and having been exposed to camels for each region showed relatively similar patterns to those at the national level but with variations in the survival probabilities (see Figs. 5A–5Y).

Risk assessment of MERS-CoV mortality

The MERS-CoV mortality risk factors for 14-and 45-day periods were estimated for the national level, each of the regions, and each of the provinces using both CPH modeled univariate and multivariate analyses. We also built models that incorporated regions and provinces as covariates.

The univariate risk factors for 14-and 45-day mortality are shown in Appendix Table 1. At the national level, during the 14-day mortality period, MERS-CoV survival was impaired according to the following HRs with varying statistical significance: males (1.34), age 41–60 (2.20), age ≥ 61 (5.83), non-HCWs (17.61) and underlying comorbidities (7.21). During the 45-day mortality period, similar risk factors were found compared to the 14-day period, but with relatively higher HRs for males and patients aged 41–60 and lower HRs for the other risk factors. At the regional level, during the 14-day period and in all regions except the North and South regions, patients aged ≥61 (HR = 2.56–11.39) and those with underlying comorbidities (HR = 4.30–10.95) were associated with MERS-CoV mortality. In the Central (HR = 9.47) and West regions (HR = 55.77), non-HCWs were linked to higher 14-day mortality rates, whereas being male (HR = 1.46) was only significant in the Central region. Overall, the univariate analysis for 45-day mortality identified similar risk factors as for the 14-day mortality, with relatively lower HRs. However, non-HCWs (HR = 5.93) in the East region, males (HR = 1.40) in the West region, and patients aged ≥61 (HR = 2.43) in the North region were at higher risk for the 45-day period but not for the 14-day period. Overall, the univariate analysis at the provincial level identified similar risk factors as for the regional level (Appendix Table 2). When regions and provinces were used as covariates in the CPH univariate models, only coming from Qassim province was significantly associated with MERS-CoV mortality at 14 days (HR = 2.12) (Appendix Tables 3 and 4).

The results of the multivariate analysis for the national level and each of the regions are shown in Table 2. The national model revealed that people aged ≥61 (HR = 2.10), non-HCWs (HR = 10.12), and those with underlying comorbidities (HR = 4.11) were independently related to higher 14-day mortality (p < 0.0001). The model for 45-day mortality identified similar risk factors to those of the 14-day model, but with an additional risk factor: patients aged 41–60 (HR = 1.44). The regions showed variations in the mortality hazard models. Those aged ≥61 (HR = 1.86–5.39) and patients with underlying comorbidities (HR = 2.69–6.69) were independently related to higher 14-day mortality in all regions except for the North and South regions. Non-HCWs were significantly associated with mortality in the West (HR = 29.34) and Central (HR = 4.55) regions, while those aged 41–60 were at risk only in the West region (HR = 2.44). During the 45-day mortality period, by region, similar risk factors were identified as for those in the 14-day models, but with relatively lower HRs. The spatial distribution of the HRs for the MERS-CoV risk factors for 14-and 45-day mortality based on multivariate analysis for each of the regions and provinces are shown in Fig. 6 and Appendix Fig. 2, respectively. The estimated HRs for 14- and 45-day mortality varied spatially by province. For the former, the highest HRs were found in Makkah (HR = 2.1 for age 41–60 and HR = 27.03 for non-HCWs), in Riyadh (HR = 8.24 for having comorbidity) and in Madinah (HR = 8.59 for age ≥ 61), while for the latter, the highest HRs were found in Makkah (HR = 3.52 for age ≥ 61 and HR = 6.15 for non-HCWs), in Riyadh (HR = 7.17 for having comorbidity) and in Madinah (HR = 11.14 for age 41–60) (p < 0.0001). Overall, the multivariate analysis at the provincial level identified similar risk factors as for the regional levels, as shown in Table 3, with the assistance to discover the specific risk factor for each province. The results of the multivariate CPH modeling that incorporated regions and provinces as covariates in addition to the other risk factors are shown in Tables 4 and 5. Regionality was not significantly associated with 14-and 45-day mortality, while for the provinces, coming from Makkah (HRs = 1.30 and 1.27) or Qassim (HRs = 1.77 and 1.70) provinces were independently related to higher 14-and 45-day mortality, respectively (p < 0.0001).