Assessing the impact of non-pharmaceutical interventions (NPIs) and BCG vaccine cross-protection in the transmission dynamics of SARS-CoV-2 in eastern Africa.

OBJECTIVE: The outbreak of the novel coronavirus disease 2019 (COVID-19) is still affecting African countries. The pandemic presents challenges on how to measure governmental, and community responses to the crisis. Beyond health risks, the socio-economic implications of the pandemic motivated us to examine the transmission dynamics of COVID-19 and the impact of non-pharmaceutical interventions (NPIs). The main objective of this study was to assess the impact of BCG vaccination and NPIs enforced on COVID-19 case-death-recovery counts weighted by age-structured population in Ethiopia, Kenya, and Rwanda. We applied a semi-mechanistic Bayesian hierarchical model (BHM) combined with Markov Chain Monte Carlo (MCMC) simulation to the age-structured pandemic data obtained from the target countries.
RESULTS: The estimated mean effective reproductive number (Rt) for COVID-19 was 2.50 (C1: 1.99-5.95), 3.51 (CI: 2.28-7.28) and 3.53 (CI: 2.97-5.60) in Ethiopia, Kenya and Rwanda respectively. Our results indicate that NPIs such as lockdowns, and curfews had a large effect on reducing Rt. Current interventions have been effective in reducing Rt and thereby achieve control of the epidemic. Beyond age-structure and NPIs, we found no significant association between COVID-19 and BCG vaccine-induced protection. Continued interventions should be strengthened to control transmission of SARS-CoV-2.

Introduction

The emergence of COVID-19 pandemic was expected to have devastating consequences in Africa due to the weak healthcare systems [1-3]. However, fatalities have remained low and most cases are asymptomatic [4]. This is attributed to previous exposure to epidemics such as Ebola, demographic factors, host genetics and environmental factors [2, 5]. Apart from Africa’s young population, Bacillus Calmette − Guérin (BCG) vaccine against tuberculosis was proposed to reduce the severity of COVID-19 [6-8].

Most countries implemented NPIs to limit human-to-human transmission of SARS-CoV-2 and therefore lower the reproduction number (R0)—the number of secondary infections acquired from a primary case [9-11]. It was imperative to quantify enforced NPIs in terms of their efficacy and appropriate use to influence and improve public health policy. Indeed, several models have been used to unravel COVID-19 [11-13].

The aim of this study was to examine the association between age-structure and BCG vaccine-induced protection from COVID-19 and to assess the impact of NPIs implemented in Ethiopia, Kenya, and Rwanda.

Main text

Methods

Data sources

COVID-19 data (2020–2021) were obtained from the World Health Organization (WHO) and the Johns Hopkins University (JHU) repositories [14, 15]. Population data were sourced from United Nations (UN) records [16] while BCG vaccination records from 1980–2019 (both sexes combined) were obtained from the WHO [17, 18].

BCG vaccination data were segmented into 10-year age-groups, and the mean percentage vaccination coverage (pvc) was calculated, assigning zero coverage to age-groups above 40 years [19], (Additional file 1: Table S1). pvc was used to infer the number of BCG-vaccinated individuals (N) in country m (Eq. 11).

COVID-19 data were split into two age-groups, 0–39 and 40 and above. This was based on the fact that BCG vaccination was introduced in EACs in the early 1980s, and therefore, only individuals aged below 39 years were assumed to be vaccinated by 2019 [19]. Finally, implementation dates of NPIs were obtained from the respective government websites and media houses (Additional file 1: Table S2).

Model formulation

At the onset of the pandemic, the Imperial College London (ICL) proposed a BHM that uses observed deaths to infer the true number of infections [11]. Deaths were expressed as a function of infection-fatality-ratio, infection-to-onset, and onset-to-death distributions [11]. Infections were expressed as a product of the time-varying reproduction number (Rt) with a discrete convolution of previous infections weighted by an infection-to-onset distribution specific to SARS-CoV-2 [11]. Rt was inferred from the initial R0 before interventions and the effect sizes from the interventions [11]. The ICL model has been applied in several studies [11, 20–22] and its’ adapted structure used in this study is shown in Additional file 1: Figure S1.

Infection model specification

The infected population (c) on day (t) was modelled as a discrete renewal process. The model was initialized by a serial interval distribution (g) with density g(τ), whereby g is gamma distributed with a mean of 6.5 and standard deviation of 0.62 (Eq. 1) [23]. g is shared across all the countries [11].

The number of infections (c) on day t, in country, m, were approximated by a discrete convolution function (Eq. 2).

Daily g was discretized by the serial interval () (Eq. 3).

The current number of infections were determined by infections in the previous days, weighted by .

is a function of interventions (I) imposed at time t in country m (Eq. 4) [11].

The implemented intervention (I) is denoted by k which is 1 if k is enforced in country m at time t, and 0 otherwise. Exponentiation of Eq. 4 constrains R0, to be positive. Further, α1, …7 determines the impact of each intervention on R. Prior distributions of α are Gamma distributed, . R0 assumes a prior distribution specified below (Eq. 5).

Death model specification

Daily deaths (dt, m) for days t ∈ {1, …, n} and countries m ∈ {1, …, p} were projected using a function dt, m = E[dt, m] whereby dt, m represents daily deaths, and it follows a negative binomial distribution with mean = dt, m and variance = , where ψ1 follows a positive half-normal distribution (Eq. 6a and 6b) [11].

Observed deaths are associated with cases using the infection-fatality-ratio (ifr, probability of death given infection) of 0.1% and the infection-to-death (π) distribution [20]. The model applies an adjusted ifr (ifr) that incorporates the attack rate and the population size [20]. Therefore, where is the age-group-specific attack rate. where, is the number of infections in age-group α, and the population size. The infection-to-death (π) distribution consists of infection-to-onset (π′) and onset-to-death distributions. π was initialized using values from Verity et al.[11, 20]. Infection-to-onset is Gamma distributed with a mean of 5.1 days and coefficient of variation of 0.86 while onset-to-death is also Gamma distributed with a mean of 18.8 days and a coefficient of variation of 0.45 (Eq. 7) [11]. The expected number of deaths dt, m, on day t, in country m was estimated by Eq. 8. where, is the number of new infections on day τ in country m. π is discretized via Eq. 9

BCG vaccine-induced protection

To assess vaccine-induced protection from COVID-19, the number of BCG-vaccinated individuals (N) in country m was assumed to have anti-SARS-CoV-2 antibodies. N was applied as a scaling factor to estimate susceptible individuals (S) on day t, in country m (Eq. 10).

The number of infections (i) on day t, in country m, was estimated by a discrete convolution function (Eq. 11).

The daily g was discretized by distribution (Eq. 3). Similarly, we computed weighted by the interventions (I) at time t in country m (Eq. 4). The expected number of deaths (dt, m) on day t, for country m was estimated by Eq. 12. where, is the number of new infections on day τ in country m. πt-τ is discretized via Eq. 9.

The SEIR and eSIR compartmental models

We explored an extension of the susceptible-exposed-infectious-recovered (SEIR) model to compute Rt [24]. Ordinary differential equations (ODE) that describe the dynamics of this model, are extensively described in [24]. Moreover, we used the extended SIR (eSIR) model to predict Rt values [25, 26].

Model implementation

Modelling was implemented in R (version 4.0.4) using ICL covid19model version 10 [11]. The ICL model was run in Stan R package using 500 iterations [11, 27]. Computation of the Root Mean Square Error (RMSE) and Mean Absolute Error (MAE) was executed in the ehaGoF R package [28] while the SEIR model was executed using the SEIR-fansy R package [24].

Model validation and comparison

The reliability of the ICL model was assessed by comparing model predictions against the observed data between 06/16/20 and 04/11/2021 using RMSE and MAE metrics [29, 30]. Additional validation was performed using an approach suggested by Flaxman et al. through an importance sampling leave-one-out cross validation scheme [11, 21, 28]. Moreover, we compared the ICL model with compartmental models using the predicted Rt and case-death-recovery counts [24-26].

Results

Scenario analysis of COVID-19 trends

We evaluated the effectiveness of NPIs under two scenarios: presence or absence of an age-structured population and BCG vaccination. Assuming that the population is homogenous and not structured by age, the ICL model estimated Rt values were 2.50 (C1: 1.99–5.95), 3.51 (CI: 2.28–7.28) and 3.53 (CI: 2.97–5.60) in Ethiopia, Kenya and Rwanda respectively (Table 1).

Table 1

Comparison of predicted case-death counts and time-varying reproduction number (Rt)

Country	Cases	Deaths	ICL (Rt)1	eSIR (Rt)2	SEIR (Rt)3	RMSE4	MAE5
Ethiopia (ET)	363,714	6412	2.50 (1.9–5.95)	2.75	2.98	4.257	0.027
ET + BCG	363,714	6412	1.67 (1.5–3.19)			4.138	0.026
ET(40 +) + BCG	134,716	5408	5.25 (3.3–8.16)			3.558	0.026
Kenya (KE)	252,938	5266	3.51 (2.8–7.28)	2.70	2.51	4.652	0.033
KE + BCG	252,938	5266	5.34 (3.5–7.99)			31.014	0.219
KE(0–39) + BCG	159,259	821	5.18 (3.8–7.87)			33.642	1.449
KE(40 +) + BCG	93,679	4445	5.15 (3.2–7.67)			3.701	0.031
Rwanda (RW)	99,559	1322	3.53 (2.7–5.60)	3.10	2.03	0.932	0.051
RW + BCG	99,559	1322	6.32 (4.53–13.34)			0.962	0.053
RW(0–39) + BCG	62,693	180	5.21 (3.5–8.62)			0.151	0.057
RW(40 +) + BCG	36,866	1142	5.93 (4.4–9.97)			0.834	0.052

 1ICL (Rt) - Imperial College London (ICL) model estimates of the time-varying reproduction number (Rt). 2eSIR (Rt) - the extended susceptible-infected-removed (eSIR) compartmental model estimates of the time-varying reproduction number. 3SEIR(Rt) - susceptible-exposed-infectious-recovered (SEIR) compartmental model estimates of the time-varying reproduction number. 4RMSE measures the model (ICL) prediction accuracy against the observed data in a regression analysis. It is the Root of the Mean of the Square of Errors between the predicted and the observed COVID-19 cases and deaths. 5MAE measures the accuracy of the model fit in terms of performance in its predictions - the Mean of Absolute value of Errors between the predicted and the observed COVID-19 cases and deaths. The mean Rt values projected by the ICL model overlapped with the SEIR and eSIR models. However, the ICL model tends to overestimate Rt values while the SEIR model had less variability (Table 1) [31]

We observed a good model fit between the predicted and the reported cases across East Africa. Larger RMSE values, particularly in Kenya, indicate a wider divergence between the predicted and observed values. Similarly, the computed MAE values ranged between 0.026–1.449 (Table 1). In general, lower RMSE and MAE values provide better support for the model fit.

Additionally, our results indicate that lockdowns and curfews profoundly reduced the Rt in Kenya (Fig. 1A) while in Ethiopia, the declaration of emergency and regional lockdowns reduced human-to-human transmissions (Fig. 1B). The dusk-to-dawn curfews in Rwanda had the most effect in lowering Rt. Beyond age-structure, we found no significant association between COVID-19 and BCG-vaccine induced protection (Fig. 2 and Additional file 1: Figure S2, S3).

Fig. 1

Country-level estimates of infections, deaths and Rt. Top: daily number of infections, brown bars are reported infections, blue bands are predicted infections, dark blue 50% credible interval (CI), light blue 95% CI. Bottom-left: daily number of deaths, brown bars are reported deaths, blue bands are predicted deaths. Bottom-right: time-varying reproduction number (R), dark-green 50% CI, light-green 95% CI

Fig. 2

The association between BCG vaccination, R0 and the case-death counts. A Association between R0 and case-death counts; B Association between age-structure, case-death counts and R0

Discussion

Parameters such as Rt and serial interval are estimators of the disease extent in a given country and they inform policy-makers about the most effective interventions [32]. Our findings show that, Ethiopia, Kenya and Rwanda were at a critical point in 2020, whereby Rt values remained above 2, however, while infections were high, fatalities remained low. Indeed, this is an African paradox where COVID-19 infections have consistently remained high while fatalities are low [7].

Beyond under-estimation of the disease extent, multiple factors have been associated with the low fatalities, namely, herd immunity due to anti-SARS-CoV-2 antibodies, climate, comorbidities and demographic structure [5]. These factors have not been studied conclusively to establish their association with COVID-19 [5].

While BCG vaccine offers cross-protection against other diseases, it has also been proposed to reduce the severity of COVID-19 [6, 33]. However, our findings show that there is no linkage between BCG vaccination and COVID-19 prevalence. In fact, the WHO did not find evidence of BCG vaccine-induced protection from COVID-19 [34].

Africa’s predominantly young population, with fewer comorbidities has been associated with the low prevalence of COVID-19 relative to other continents [35, 36]. Indeed, we observed a negative correlation between Rt and age-structured case-death counts. It is noteworthy that, in this study, the susceptible population was segmented according to age and BCG vaccination status prior to estimation of posterior parameters. Consequently, age had a confounding effect on Rt and case-death counts wherein the effect of BCG vaccination could not be separated from the effect of age-structure. Generally, majority of COVID-19 cases were identified to be in the age group 30–39 while most deaths comprised of those aged above 40 [14].

Limitations

The ICL model uses observed deaths to infer the true number of infections [11]. While this approach overcomes uncertainties associated with asymptomatic cases and low testing in the African context, inferring infections from deaths to estimate the burden of the disease is challenging given the low mortality recorded in Africa. Further, some parameters were set by assumption or used values from literature, which significantly affect the parameter estimation.

Additional file 1: Table S1. Country-level percentage coverage of Bacillus Calmette–Guérin (BCG) vaccination in 2019, by age group. Table S2. Country-level dates of implementation of non-pharmaceutical interventions. Figure S1. Schematic overview of the Imperial College London (ICL) model adopted in this study [11]. Figure S2. Country-level estimates of infections, deaths and Rt in Ethiopia. Scenarios: A) The population is not BCG vaccinates, homogenous and not structured by age; B) BCG vaccinated population aged 39 years and below; C) BCG vaccinated population aged 40 years and above. Top: daily number of infections, brown bars are reported infections, blue bands are predicted infections, dark blue 50% credible interval (CI), light blue 95% CI. Bottom-left: daily number of deaths, brown bars are reported deaths, blue bands are predicted deaths. Bottom-right: time-varying reproduction number (R), dark-green 50% CI, light-green 95% CI. Figure S3. Country-level estimates of infections, deaths and Rt in Kenya. Scenarios: A) The population is not BCG vaccinates, homogenous and not structured by age; B) BCG vaccinated population aged 39 years and below; C) BCG vaccinated population aged 40 years and above. Top: daily number of infections, brown bars are reported infections, blue bands are predicted infections, dark blue 50% credible interval (CI), light blue 95% CI. Bottom-left: daily number of deaths, brown bars are reported deaths, blue bands are predicted deaths. Bottom-right: time-varying reproduction number (R), dark-green 50% CI, light-green 95% CI.