A vaccination model for COVID-19 in Gauteng, South Africa.

The COVID-19 pandemic provides an opportunity to explore the impact of government mandates on movement restrictions and non-pharmaceutical interventions on a novel infection, and we investigate these strategies in early-stage outbreak dynamics. The rate of disease spread in South Africa varied over time as individuals changed behavior in response to the ongoing pandemic and to changing government policies. Using a system of ordinary differential equations, we model the outbreak in the province of Gauteng, assuming that several parameters vary over time. Analyzing data from the time period before vaccination gives the approximate dates of parameter changes, and those dates are linked to government policies. Unknown parameters are then estimated from available case data and used to assess the impact of each policy. Looking forward in time, possible scenarios give projections involving the implementation of two different vaccines at varying times. Our results quantify the impact of different government policies and demonstrate how vaccinations can alter infection spread.

Introduction

With the declaration of the COVID-19 pandemic by the World Health Organization on March 11, 2020, many nations began to consider the implications of the spread of this disease within their borders and across the world (Cheng et al., 2020; Garba et al., 2020; COVID-19 Outbreak in Countries from, 2021). Whilst waiting for the development of an effective and safe vaccine, governments focused attention on other interventions in an attempt to curb the spread of infection. Broadly speaking, these can be split into preventative actions, such as social distancing and stay-at-home mandates, and surveillance actions, such as contact tracing and testing. Nonpharmaceutical interventions include masking, hand-cleaning, sanitizing, and air filtering. The COVID-19 pandemic provides a rare opportunity to explore the impact of these interventions on a novel infection and thereby link them to early-stage outbreak dynamics.

By the end of 2020, safe and effective vaccines were starting to become available in some countries. The vaccines have been shown to improve outcomes for infected individuals and to reduce risk of infection and onward transmission from vaccinated to unvaccinated individuals (Katella, 2021). Under this scenario of vaccine availability, we are able to scrutinize the data to determine the impacts of vaccination and the effects of delaying the vaccine roll-out beyond the start of 2021 in certain geographical areas.

Our interest lies not simply in understanding the data, but in using mathematical modeling to explore a range of alternative scenarios in which the preventative and surveillance actions are initiated at different time points. In doing this, we can gauge the impacts of the actions and provide a systematic framework for this approach.

South Africa had the largest COVID-19 burden in Africa, with its cumulative cases contributing about 36% of the total in Africa and at a rate three times more than the second-place country (COVID-19 Outbreak in Countries from, 2021; Africa Centers for Disease Control, 2021). Data on COVID-19 cases, policy changes on management of the disease, and vaccination rates and protocols are readily available in the public domain (COVID-19 Resources for Republic of, 2021). South Africa was also one of the first African countries to roll out COVID-19 vaccination (Africa Centers for Disease Control, 2021). Whilst it would be possible to consider infection dynamics across the country, given the restrictions on movement (COVID-19 Resources for Republic of, 2021) and the need for individuals to interact in order for infection to spread, we focus on Gauteng, the single most-populous region of South Africa. Gauteng contains several large population centers including the cities of Johannesburg, Ekurhuleni, and Tshwane. Although it is the smallest province in the country, covering less than 2% of the total land area, Gauteng accounts for more than 25% of South Africa's population. This densely populated region has registered 1,095,360 confirmed cases of COVID-19 since April 2020 (Statista and Confirmed Coronavirus, 2021). We investigate how the various changes in guidelines over time related to closures and social distancing affect the transmission and spread of COVID-19, with our model parameterized for Gauteng. Changes in government guidelines and human behavior can lead to a time-varying transmission rate, due to changing contact rates and possibly due to new variants. We formulate a system of ordinary differential equations (ODEs) to represent the transmission of COVID-19, including contingencies for vaccination.

A number of authors have focused their attention on building and analyzing mathematical models to explore the impact of COVID-19 in Africa. These models provide important context for our work. For instance, they highlight the pandemic disruption caused to HIV treatment programs (Jewell et al., 2020) and the impact of response strategies (Taboe et al., 2020; Van Zandvoort et al., 2020). More generally, mathematical models and statistical analysis have been used to predict spread and scale of outbreaks at the start of the pandemic (Atangana & Araz, 2021; Musa et al., 2020), to understand the continental heterogeneity in COVID-19 impact (Musa et al., 2022), and to explore the comparatively low transmission and mortality rates that have been observed across the continent (Bouba et al. Kong).

In the specific context of South Africa, two studies have been published that investigate the impact of interventions using systems of ODEs. In the first, the impact of social distancing on the number of COVID-19 cases was explored for the period March–May 2020 (Nyabadza et al., 2020); our modeling study extends the period of consideration and increases the model complexity to include additional population compartments to reflect more-recent knowledge about infection transmission. The second model (Mukandavire et al., 2020) focuses again on the early stages of the outbreak and considers the utility of a vaccination in containing disease. Again, our work builds on this by considering specific vaccines and implementation no earlier than January 2021 when they became available in some countries (but not South Africa).

Based on the multitude of policy shifts during the COVID-19 outbreak, and our knowledge of how such shifts affect parameter values from (Burton et al., 2021), we planned for specific parameters to vary based on certain time periods over the course of the pandemic. Specific time points where values changed correspond to inflection points in the data and recorded policy changes. Vaccine distribution was at a very low level during the time period under consideration; thus the parameters for the model without vaccination were estimated using incident data, with only two parameters changing as time-varying step functions. After estimating our parameters from the data, including appropriate time-varying rates, we present the fit of our simulated model to our data. We then consider various scenarios involving varying transmission and vaccination rates to illustrate the course of the infection if certain behaviors continued.

Understanding the data

Data for this project was obtained from a repository of South Africa COVID-19 data, maintained by the Data Science for Social Impact research group at the University of Pretoria (Katella, 2021; Data Science for Social Impact Research Group @ University of Pretoria, 2021). It consists of weekly cumulative confirmed cases for the province of Gauteng, South Africa and is shown in Fig. 1 (Coronavirus Disease, 2021).

Fig. 1

Depiction of weekly cumulative confirmed case data for the province of Gauteng, South Africa (Katella, 2021; Data Science for Social Impact Research Group @ University of Pretoria, 2021; Coronavirus Disease, 2021). Approximate behavior-change dates are highlighted in red. We have labelled the time intervals between changes as T, i = 1, 2, …, 5, as described in the text.

The epidemic curve demonstrates a number of key features: initial exponential growth in case numbers due to the novel form of the virus; reduction in the rate of increase in July 2020 associated with state policy intervention; and subsequent repeated changes in the rate of increase in cases due to periodic changes in policy and intervention which lead to increased rates of infection (when interventions were relaxed) or reduced rates of infection (when interventions were tightened).

We translate the changes in government policy and individual behaviors into our model by assuming that two key model parameters—namely, the infection transmission coefficient and the rate of identifying positive COVID-19 cases through testing—vary over time.

To correctly identify the time-dependent nature of these two parameters, we combined information from two sources—the epidemic curve and the government guidelines—in the following way. From the data on cumulative cases, we were able to estimate points of inflection. These points separate regions in which the epidemic is worsening (rate of cumulative cases increasing over time) from the periods in which it is improving (rate of cumulative cases decreasing over time). Our hypothesis was that these points of inflection should broadly correspond to times where government guidelines changed. Clearly there are likely to be time-delays with any implementation, but our choices of time-dependent parameters align well with changes taking place rapidly (for example, a lockdown will reduce the transmission parameter β).

Interrogation of government guidelines provided the following timeline of policy changes:

March 30 to July 12, 2020: This was a learning and adjustment period. Having declared a national disaster, the government imposed a three-week lockdown where citizens were told to stay at home and only essential services and businesses were able to continue to work at large. No alcohol or cigarette sales were permitted and citizens were not allowed to travel or attend any form of gatherings. Phased relaxation of this lockdown was introduced using a five-level alert system. By early June 2020, the country was at alert level 3, in which alcohol sales were permitted, schools and universities were reopened, and limited air travel was allowed. South Africa recorded its highest number of infections in one day (13,674) within this first wave. This prompted some renewed strict restrictions such as bans on alcohol and family visits, and night curfew being restored.

July 13 to October 10, 2020: Temporary school closures were announced around July 23 for four weeks. The risk-adjusted alert system level was reduced to 2 in mid-August before further relaxation of restrictions went down to level 1. Almost all economic activities opened with increased mobility and crowding within and outside households. Towards the end of this period, a new variant, beta, emerged in the coastal regions of South Africa with a transmissibility rate twice that of the original variant, alpha.

October 11 to December 27, 2020: Differentiated regional alert levels were used as a refined version of interventions that were less necessary in some areas than others. A surge in cases up to 145% was observed in some coastal regions; this was attributed to interactions in educational settings, the festive season, and other big events. Spread of infection was seen to occur mainly in the younger age group of 15–19 year-olds. Some restrictions were put in place during the festive season to restrict spread.

December 28, 2020 to March 30, 2021: The national alert level 3 prohibited indoor and outdoor gatherings for 14 days from December 28. Additional restrictions in mid-January until mid-February involved closing land ports of entry. Easing of restrictions to level 3 opened public spaces and economic activities. By the end of February, the national alert level was at 1. Vaccination of health care workers began.

March 31 to June 6, 2021: Vaccination of health workers continued and began with people of age 60 years or more. The delta variant emerged in Gauteng, which caused the alert response to move from level 1 to level 2.

For details of the alert levels, refer to (Lockdown has Pushed COVID, 2020; President Cyril Ramaphosa, 2021a; Level 1 Lockdown in Numbers, 2020; President Cyril Ramaphosa, 2021b; COVID-19 Special Public Health Surveillance, 2021; Disaster Management Act, 2021a; Coronavirus COVID, 2021; President Cyril Ramaphosa, 2020; Disaster Management Act, 2021b; Disaster Management Act, 2021c; Disaster Management Act, 2020a; Disaster Management Act, 2020b).

Combining this information with the points of inflection found in the cumulative cases data, we estimated the following four dates for behavioral change: July 13, 2020, October 12, 2020, December 28, 2020, and March 30, 2021. These dates correspond to the times at which certain parameters change values in our model, resulting in five time periods which we label , i = 1, 2, …, 5, as shown in Fig. 1. More details on our process are given in Section 4.

A number of pharmaceutical and biotech companies have developed vaccines for COVID-19, each of which differs in the biotechnology used, efficacy, and geographic availability. In South Africa, two commonly available vaccines are those developed by Pfizer and Johnson & Johnson. Based on the information presented in (Katella, 2021), the Pfizer vaccine can be assumed to be 95% effective at preventing individuals from contracting the disease and 95% effective at preventing symptomatic disease (vaccinated individuals becoming ill and developing symptoms). The vaccine developed by Johnson & Johnson claims to be 64% effective at preventing disease specifically in South Africa and 81% effective at preventing symptomatic disease (vaccinated individuals becoming ill and developing symptoms) (Katella, 2021). We use these values to choose vaccine-related parameters when we consider the implementation of vaccination programs in South Africa.

Model structure

Our full model system comprises nine compartments in a coupled system of nonlinear ODEs. The compartments are taken to represent: S susceptible individuals, V fully vaccinated individuals, E exposed individuals, E exposed after vaccination individuals, A asymptomatic individuals, I pre-symptomatic and symptomatic infected individuals, Q individuals with confirmed infections, H individuals in the hospital with confirmed infections, and R recovered individuals. We do not consider partial vaccination—consequently, we assume that individuals remain in the unvaccinated class until the vaccine is effective (typically two weeks after the final vaccine dose administered).

Asymptomatic individuals remain asymptomatic throughout the disease process, while pre-symptomatic individuals are not presently symptomatic but will later become symptomatic. Note that pre-symptomatic individuals transmit the infection at a rate similar to symptomatic individuals, and thus the pre-symptomatic and symptomatic individuals are combined in the I class. Asymptomatic individuals transmit the infection at a much lower rate than pre-symptomatic and symptomatic individuals, and this difference is indicated in the force of infection term which has a larger coefficient c on the I term than the A term. The transitions between compartments are illustrated in Fig. 2. The parameters (with units) are described in Table 1.

Fig. 2

Flow diagram of compartmental model with vaccination.

Table 1

Parameters in the model, including their definitions and units.

Symbol	Interpretation	Units
β	transmission rate	per day
1/α	length of exposure period	days
ε	proportion of asymptomatic out of infectious individuals	unitless
μ	COVID-19 death rate in I and Q	per day
μH	COVID-19 death rate in H	per day
κ	testing rate resulting in isolation	per day
ρ	hospitalization rate	per day
c	scaling factor of infected compartment	unitless
γ1	recovery rate of asymptomatic individuals	per day
γ2	recovery rate of infected individuals	per day
γ3	recovery rate of individuals with confirmed cases	per day
γ4	recovery rate of hospitalized compartment	per day
θ	vaccination rate	per day
π1	scaling factor on the force of infection for V	unitless
ε1	proportion of asymptomatic out of infectious vaccinated	unitless

We assume that the behavior and policy changes discussed in Section 2 altered the transmission rate β at which the disease spreads as well as the testing infrastructure and capacity rate κ. For this reason we allow these two parameters to vary over time using the inflection points and the policy mandates identified in Section 2. Equations for our compartmental model are given by:where all model parameters are non-negative and are defined in Table 1. The force of infection for susceptibles S is given by

and the vaccinated individuals V have a modified force of infection,

due to having a reduced chance of becoming exposed (with 0 < π1 < 1) (Katella, 2021). Note that the denominator of the force of infection only includes the compartments that are available to be in contact with the susceptibles S i.e., we assume that the two compartments Q and H are isolated and are not in contact with the susceptibles.

Having split the time frame for our modeling into five distinct intervals T, i = 1, …, 5, the time-dependent parameters β(t) and κ(t) are assumed to take the form:

Since an infectious pre-symptomatic or symptomatic individual is more likely to transmit COVID-19 than an asymptomatic individual, we scale the contact rate with I by a factor of c within the force of infection (Centers for Disease Control and Prevention, 2021). Due to the lack of widespread testing of asymptomatic individuals in this province, we do not include a route for asymptomatic individuals to move into the Q class. We assume that individuals in I, Q, and H compartments could die from a COVID-19 infection. Individuals in the A, I, Q, and H compartments can recover from the disease. The proportion of E that become asymptomatic A is larger than without vaccination, (Katella, 2021) which corresponds to

As with many modeling problems involving infectious diseases, COVID-19 is an infection where the value of the basic reproduction number has practical significance. It corresponds to the initial outbreak and can be used to see how the infection was growing at that time. Using the Next Generation Matrix Method (Diekmann et al., 2010; van den Driessche & Watmough, 2002) with β and κ held constant at the values estimated in interval T1, we calculate the basic reproduction number for this model to be

We can interpret the two terms as representing the two routes of transmission through the A (first term) and I (second term) compartments. Our derivation can be found in the Appendix. Using parameter estimates from our time period T1 (see Section 4 below), and Equation (2) we estimate that, for Gauteng Province, at the start of the outbreak

Field studies for the same period in Gauteng Province estimate to fall between 1 and 2.19 (The Initial Daily COVID, 2020), which gives us confidence in our parameter estimates.

Parameter estimation

We obtained estimates for many of the model parameters from the literature; these are given in Table 3 together with their references. The remaining parameters were estimated using our data from Gauteng (Coronavirus Disease, 2021). To account for inconsistent reporting rates over the course of a given week, we transformed the data that was reported as cumulative confirmed daily cases into weekly data. We estimate unknown parameters by formulating a least-squares optimization problem with the goal of minimizing the difference between the recorded number of weekly cumulative confirmed cases in Gauteng and our model's output. The objective function to be minimized iswhere the vector WQ contains the weekly cumulative number of confirmed infections from the model and the vector WQ∗ contains the corresponding values from the data.

Table 3

Parameter values with a citation or estimated value from data. Note that T1 is for the first time period, T2 is the second time period, and so on, which correspond to Fig. 1.

Parameters from Literature
Symbol	Value	Symbol	Value	Symbol	Value
α	0.25 (Johansson et al., 2021; Sanche et al., 2020)	ε	0.3 (Centers for Disease Control and Prevention, 2021)	c	1.25 (Centers for Disease Control and Prevention, 2021)
γ1	0.14 (Cheng et al., 2020; Johansson et al., 2021)	γ2	0.14 (Cheng et al., 2020; Johansson et al., 2021)	γ3	0.1 (Johansson et al., 2021; Renardy et al., 2020)
π1 Pfizer)	0.95 (Katella, 2021)			ε1 (Pfizer)	0.95 (Katella, 2021)
π1 (J & J)	0.64 (Katella, 2021)			ε1 (J & J)	0.81 (Katella, 2021)
			Constant Estimated Parameters
ρ	0.142	μ	0.014	μH	0.094
θ	0.003	γ4	0.067
			Estimated Parameters that Change Over Time
Symbol	T1	T2	T3	T4	T5
β	0.539	0.169	0.277	0.199	0.419
κ	0.739	0.339	0.098	0.506	0.499

To aid the optimization process we bounded each parameter using information from (Parameter Estimates for 2019 Novel, 2019) and from the sources listed in Table 2.

Table 2

Bounds imposed on our parameters in the optimization problem.

Parameter Bounds
Symbol	[Lower Bound, Upper Bound]	Sources
ρ	[1/12, 1/3]	Zhang et al., 2020, 2021
μ	[0, 1/50]	Mizumoto & Chowell, 2020; Zhang et al., 2021
μH	[1/14.9, 1/6.4]	Linton et al., 2020; Sanche et al., 2020
γ4	[1/17.3, 1/8]	Sanche et al., 2020
βi	[0.00001, 1]	Parameter Estimates for 2019 Novel, 2019
κi	[1/20, 1]	Zhang et al., 2020

Our parameter estimation simulations begin on March 30, 2020 and take daily time steps until the date our data ends, which is June 6, 2021. Initial values for our simulations include the population of Gauteng for S0 and the exact value from the data for Q0. The initial values of E0, A0, and I0 are chosen using the information about the average lengths of time spent in those classes and the asymptomatic infected being about one-third of the total infected (Centers for Disease Control and Prevention, 2021). The values of H and R are taken as zero because of the paucity of information linked to the onset of the outbreak. Initial values are summarized below:

We solve our least-squares optimization problem in MATLAB version R2021a using the MultiStart and fmincon functions (Burton et al., 2021; Edholm et al., 2019; Levy et al., 2017). Given a starting point for our objective function J, the fmincon algorithm outputs a local minimum on the surface of J. To help find the global minimum, MultiStart allows us to exhaustively test different starting values throughout our bounded range. We used 10,000 different starting points, each of which converged to a unique local minimum on the surface of J. The smallest objective function value obtained was J = 0.01. To align with what took place in Gauteng, one final condition we have on our parameter values is that the resulting simulation must produce two infection peaks, the second of which should be greater than or approximately equal to the first.

Table 3 depicts the parameter estimates that we obtained either from fitting the model to data or from the literature. The only exception is our estimate for θ which we chose to correspond to a vaccination capacity of around 40,000 individuals per day (Department: Health, 2021).

The estimated mean values of the generation interval, serial interval, incubation period, and latent period have varied geographically and temporally since the beginning of the COVID-19 outbreak in Hubei, China (Griffin et al., 2020). These values affect the estimate of and decisions about intervention actions (Ferretti et al., 2020; Griffin et al., 2020; Lipsitch et al., 2003; Wallinga & Lipsitch, 2007); as we indicated earlier, our prediction for falls within the bounds estimated from Gauteng daily reported cases which suggests that we have reasonable parameter estimates for the onset of COVID-19 in the region. More recent work on the mean generation intervals (time lag between infections of a primary case and its secondary case) indicates that the sum of the mean latent period (Xin et al., 2021) and mean infectious period may be shorter than predicted earlier in the epidemic (Tang et al., 1262). This would affect our parameters, α and γ2, but for the purpose of this work, we fix the value of these parameters using estimates from the initial period of outbreak.

Results

In this section, we present results from our numerical exploration of different responses to public health policy changes. We consider dynamics both with and without vaccination implementation since vaccines were not readily available in South Africa until June 2021.

Numerical solution of our model system (1) in the absence of any vaccination intervention using the estimated parameters given in Table 3 is shown in Fig. 3. The figure indicates that we have established a good model fit to the cumulative cases of COVID-19 in Gauteng, and so we label this as the baseline scenario from which we explore hypothetical, alternative scenarios. Also see the two peaks in confirmed cases from this baseline scenario in Fig. 4.

Fig. 3

Simulation output using our estimated parameters is plotted alongside the confirmed cumulative case data and exhibits a good fit to the shape of the infection curve. Our values for β and κ change four times at the locations indicated by ∗ (for further details on these change times, see Section 2).

Fig. 4

Simulation output using our estimated parameters is plotted alongside the weekly confirmed case data and exhibits a good fit to the shape of the infection curve.

In Section 5.1, we explore how different responses to the public health interventions would have impacted cumulative COVID-19 cases in the absence of vaccination. In Section 5.2 we compare key metrics output from the model for the two vaccines commonly available in South Africa, assuming different dates when they became generally available.

Infection dynamics without vaccination

Fig. 5 shows model predictions against the actual cumulative case data assuming no additional policy or behavior changes at the end of each of the time intervals T, i = 1, 2, 3, 4. This means that the parameters in place at the end of a specific time interval were used for the simulation for the remainder of the entire simulation run, and these are denoted as Scenarios 1–4 in Table 4. When this happened at the end of a time period where case numbers were increasing, our model predictions show large, unsurprising, increases in COVID-19 infections compared with the actual data. Similarly, at the end of a time interval where restrictions were imposed and had caused lower infection rates, if those behaviors continued into the next time period, we see a significant reduction in the number of new cases. We explore these impacts further in Table 4, to see how the interventions and relaxation of restrictions affected the global maximum number of infected individuals in the population (calculated as A + I + Q + H), the number of cumulative confirmed cases, and how many infected people were hospitalized (calculated as H) when that number was at a peak. These results are shown in Table 4 as Scenarios 1–4, respectively; our choice of metrics was motivated by the data that has been regularly reported throughout the pandemic as critical for public health decision making.

Fig. 5

This plot simulates what could have taken place if infection dynamics that result from specific values of β and κ had persisted rather than changed to new values.

Table 4

Summary of key information from simulating scenarios. The table column (Global) Max Infected also records the date at which the largest number of infected individuals occurs. Row 1 summarizes our baseline simulation as seen in Fig. 3. Scenarios 1–4 summarize the simulations shown in Fig. 5, continuing the parameters to the end from the first, second, third, and fourth changes. The information found in Scenarios 5–8 was obtained by holding parameters from each given group constant for the duration of the simulation.

Scenario	Max Infected (date)	Cumulative Confirmed Cases	Total in Hospital	1st Hospital Peak (date)	2nd Hospital Peak (date)
Baseline simulation using estimated parameters
Baseline	70,322 (7/16/20)	474,648	257,670	16,123 (7/20/20)	14,922 (1/8/21)
Summarizing simulations inFig. 5
1	1,402,803 (9/18/20)	6,906,308	3,822,977	385,847 (9/22/20)	N/A
2	70,323 (7/16/20)	226,514	125,276	16,123 (7/19/20)	N/A
3	836,901 (3/18/21)	2,823,107	1,557,400	16,123 (7/20/20)	1,118,056 (3/25/21)
4	70,323 (7/16/20)	416,728	230,551	16,123 (7/20/20)	14,922 (1/8/21)
Fixing parameters from a given section for duration
5	456 (3/30/20)	724	287	N/A	N/A
6	937,647 (10/23/20)	2,804,247	1,552,217	125,081 (10/30/20)	N/A
7	456 (3/30/20)	750	300	N/A	N/A
8	837,389 (11/14/20)	5,302,811	2,935,303	214,753 (11/19/20)	N/A

Scenario 1 confirms that if no intervention had taken place at the end of the time period T1, then the public health outcome would have been significantly worse. The healthcare system would have been overloaded with several hundred thousand patients. In some sense, Scenario 2 is the best-case scenario if the public health interventions implemented in period T2 had been continued at the same level through June 6, 2021. Of course, this does not take into consideration any economic or social welfare issues which would certainly have been negatively impacted by the more-extreme interventions. Having said that, this scenario gives an indication of the levels of infection and hospitalization that could have occurred. Scenarios 3 and 4 mimic 1 and 2 but add insight into the size of hospital peaks, in particular, the second one.

Scenarios 5–8 present a more-theoretical exploration of the model predictions since they take the parameter estimates from each of the time intervals T, i = 2, 3, 4, 5, and assume these hold for the entire period of simulation starting in March 2020. This idealization would clearly not have been possible (since it requires retrospective implementation), but it highlights the importance of considering the largest possible time interval for parameter estimation to avoid focusing on a single interval which may produce different parameter values and infection outcomes. As an example, in Scenario 6 we see that if the parameters estimated for the time interval T2 were used for the model simulation from the outset, the number of COVID cases and hospitalizations would have been predicted to be significantly lower. Meanwhile, for Scenario 8—which has high parameter values from T5—the outputs are unrealistically large reflecting the possibility of unchecked spread. Whilst there has always been potential to allow model parameters to vary throughout numerical experiments, the explicit public health interventions seen globally for COVID have brought this feature to the fore.

Impact of vaccination

In this analysis we ignore the differences in biotechnology used to deliver the dose, cost of implementation, and the number of doses required to become fully vaccinated, and instead focus on the effectiveness of each vaccine. Using the associated parameter estimates for T5 given in Table 3, we undertook a series of numerical experiments to explore how delays in rolling out the vaccination program and vaccine efficacy impact the key public health metrics.

We considered two start dates for vaccination: January 4, 2021 represents an aspirational date corresponding to the roll out of programs in several countries including the UK and USA; and June 7, 2021 which reflects a realistic start date for the program of vaccination in South Africa. In both cases, we simulated the model system from March 30, 2020 to an end date April 4, 2022.

We present the results of our exploration in Table 5 and note that they appear consistent with an intuitive understanding of what vaccination would achieve. The program reduces COVID numbers and hospitalizations; the more effective the vaccine, the better it does; and delays to implementing a vaccination strategy can significantly increase the number of infections and subsequently people needing hospital care.

Table 5

Key metrics from vaccination simulations. The top portion considers vaccination beginning on 1/4/21 while the bottom portion begins vaccination on 6/7/21. The end date in both cases is 4/4/22. The Cumulative Confirmed Cases is the number of cases from the start of the simulation for vaccination until April 4, 2022.

V start 1/4/21	Max Infected	Cumulative Confirmed Cases	Total in Hospital	3rd H Peak	Total V	Total EV
Pfizer	60,296	97,896	65,447	330	10,503,215	908
J & J	60,296	774,180	436,284	11,055	10,155,676	618,777
None	626,811	4,671,953	2,597,767	160,981	0	0
V start 6/7/21	Max Infected	Cumulative Confirmed Cases	Total in Hospital	3rd H Peak	Total V	Total EV
Pfizer	201,280	1,653,919	920,419	51,170	6,835,309	46,521
J & J	375,935	2,960,281	1,643,341	89,092	5,951,945	735,231
None	626,811	4,513,906	2,503,479	160,981	0	0

Since we are undertaking this exploratory work using parameters from time period T1, the levels of infection which we obtain are larger than observed. Moreover, since we are extrapolating our simulation until April 2022, we are not surprised to see bigger infected levels than for the baseline case.

Discussion and conclusion

COVID-19 has provided mathematical modellers with access to detailed data on key infection metrics such as incidence of infection and infection-related deaths, together with data on public health interventions such as hospitalization and vaccination. There are certainly inaccuracies within the data sets, but they represent a rich data source to interrogate. The approach that we have taken here—identifying points of inflection in the data—does exactly that, and it led us to identify key dates at which public health interventions caused a change in the cases of COVID-19 in Gauteng Province, South Africa. This was undoubtedly helped by the daily reporting of COVID-19 cases and illustrates a methodology which could be used with other infectious disease data sets.

The existing published models for COVID-19 in South Africa focussed on the early stages of the epidemic. Building upon this, we have demonstrated the impact of public health interventions and have predicted what would have happened at each stage if those interventions and corresponding rates had continued. In part, we undertook these explorations to support any future public health decisions that may be needed, but also to highlight the impact of the interventions. Our focus on Gauteng Province, rather than all of South Africa, was intentional in this respect. Gauteng is the most-densely populated region in South Africa, and since COVID-19 infection is spread through close contact, we consider density as a proxy for likelihood of interaction between individuals. Moreover, travel restrictions that were widespread around the globe resulted in unprecedented levels of isolation between large urban areas. This meant that infection transmission between different cities in South Africa was likely to be less significant than transmission within cities.

Once we identified the distinct time periods T, i = 1, 2, …, 5, we were able to demonstrate the relative impact of the changing parameters by undertaking a series of numerical experiments in which the two time-dependent parameters β and κ were fixed to values estimated in one of the time periods. Results from these simulations showed that changes in our estimated parameter values lead to notably different infection dynamics and total incidents (see Fig. 5 and Table 4). This analysis indicates that when modeling COVID-19 one should carefully explore the data to extract all relevant information when estimating parameters. Additionally, since it is unlikely that parameter estimates from a limited time interval should be used to explore an extended period of the epidemic, future researchers should consider whether parameter estimates are realistic for the period of exploration or whether additional data should be sought.

Our simulations related to the impact of vaccination was informed by details of the efficacy of the Pfizer and Johnson & Johnson vaccines, since these are the vaccines most widely available in South Africa. We also considered the consequences of unequal availability of vaccines around the world by assuming two different start dates for the vaccination program. The start date of January 4, 2021 was chosen to be consistent with general roll-outs of programs in countries such as the UK, the USA, and some European Union countries and was compared to a start date of June 7, 2021, which aligns with the general availability of vaccines in South Africa. Our projections continue to April 2022 and so exact numbers must be treated with caution. However, the message is clear—vaccination provides an important element in the fight against COVID-19, and delay in delivering vaccines has a significant negative impact.

At the time of this paper's analysis, the Omicron variant had just emerged, and thus our model involved data and epidemiological properties of earlier variants. In the future, the model structure proposed here could either be modified to include multiple variants or adapted through parameter values to explore the dynamics of the dominant variant. Despite the dimensionality of our model (which has nine state variables), it is a simple representation of the dynamics of COVID-19. That was intentional, to allow us to address the questions we chose to explore; but it also allows the model to be modified in the future to address other important processes. For example, one could consider reinfection or wanning immunity from vaccination if the model covers a longer time period. Also one could include mobility data to better represent interactions and movement in models with spatial components, using ideas from (Potgieter et al., 2021). Besides mobility issues, the size distribution of households in a community may give more detail for possible types of interactions and for further information about individual actions (Wu et al., 2020). When continuing to explore vaccination scenarios, optimal control techniques may be used to choose allocation strategies.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

S0	E0	A0	I0	Q0	H0	R0
15,488,137	30	15	30	411	0	0