Table of Contents
Model structure
To model the transmission dynamics of COVID-19 in the first COVID-19 epidemic wave in South Africa a stochastic compartmental transmission SEIR model was used hereafter called the “ARI COVID-19 SEIR” model. Figure 1 shows the structure of the ARI COVID-19 SEIR model. The ARI COVID-19 SEIR Model was constructed in a Macro-Enabled Microsoft (Ms) Excel File for user-friendliness (visual interaction with parameters) and Database Query Support. The Model had a Visual Basic Application (VBA) code for Sensitivity and Variable Analysis.
The ARI COVID-19 SEIR Model structure
Figure 1 shows that the ARI COVID-19 SEIR model had the following population classes based on the assumed clinical diagnosis of COVID-19 within the population:
-
Susceptible (S)-Individuals within the population of the model who can incur the disease however have not been infected yet.
-
Exposed (E)-Individuals within the population of the model in an incubation period who are not yet infectious.
-
Asymptomatic (Ia)-Individuals within the population of the model who are infected and are infectious that are transmitting the disease to others however are not showing any symptoms throughout their infectiousness.
-
Pre-symptomatic Infectious (Ip)-Individuals within the population of the model that are transmitting the disease during their incubation period.
-
Infected with Mild and Moderate Symptoms (Ib)-Individuals within the population of the model with mild and moderate symptoms who are infectious.
-
Infected with Severe and Critical Symptoms (Ic)-Individuals within the population of the model with severe and critical symptoms who are infectious however have not yet been hospitalised.
-
Hospitalised COVID-19 Cases (H)-Individuals within the population of the model with severe and critical symptoms who have been hospitalised.
-
Death due to COVID-19 (D)-Individuals who have died due to COVID-19 or indirect consequences of the COVID-19 epidemic.
-
Recovered (R)-Individuals within the population of the model who have recovered from the disease with immunity or partial immunity.
Model equations
The transmission and severity of COVID-19 within-population classes in the model were simulated using ordinary differential equations (ODEs). The differentiation in the ODEs was conducted using Euler’s method with 1-day estimation steps. The total population (N) in the model is represented by Eq. 1 based on the conservation of mass. Vitals (new births and non-COVID-19 deaths) were not considered in the model due to the relatively small annual growth rate of 1.40% in the South African population [39] and the low incidence of COVID-19 in neonatal [40]. The ARI COVID-19 SEIR model instead used the 2020 South African population estimates from the United Nations (UN) World Population Prospects [41]. The total infections are given by Eq. 2.
$$N=S+E+Ia+Ip+Ib+Ic+H+R+D$$
(1)
$$Iabc=Ia+Ip+Ib+Ic+H$$
(2)
The change in the susceptible population class is given by Eq. 3 where β is the Effective daily contact rate, this is the average number of adequate contacts per infective per day. The product of S and I in Eq. 3 is referred to as the mass incident term.
$$\frac{\partial S}{\partial t}=-A\beta S\frac{Iabc}{N}$$
(3)
$$A\;=\;-\frac{\displaystyle\frac{\displaystyle\frac{\sqrt{Country\;Area}}\pi}N}{Effective\;Social\;distance}\;+\;1$$
(4)
A is the Population density factor given by Eq. 4. Where the Effective Social distance is the minimum distance between the infector and infectee which prevents infection. For COVID-19 a distance of 2 m was assigned [42]. As the average distance between individuals tends towards the effective social distance, the Population Density Factor (A) tends towards 0. The Population Density Factor assumes a uniform distribution of the population within a confined area. The change in the exposed class is given by Eq. 5.
$$\frac{\delta E}{\delta t}= A\beta S\frac{Iabc}{N}-{k}_{1}{\xi }_{1}E-{k}_{2}{\xi }_{2}E$$
(5)
Where K1 and K2 are the rates at which exposed individual moves to the infected class. K1 is inversely proportional to the average incubation period of COVID-19 Asymptomatic Cases (Tinc,1) and K2 is inversely proportional to the average incubation period of COVID-19 Symptomatic Cases (Tinc,2) in the population. ξ1, ξ2, ξ3, ξ4 are the proportions of the exposed and pre-symptomatic who will be Asymptomatic (ξ1), Symptomatic (ξ2), Develop Mild and Moderate Symptoms (ξ3) and Severe and Critical Symptoms (ξ4), respectively.
$$\xi_1\;+\;\xi_2\;=\;1$$
(6)
$${\xi }_{3}+{\xi }_{4} =1$$
(7)
The change in the infectious class is given by Eqs. 8, 9, 10, 11 and 12.
$$\frac{\delta Ip}{\delta t}={k}_{2}{\xi }_{2}E-{\theta {\xi }_{3}I}_{p}-{\theta {\xi }_{4}I}_{p}$$
(8)
Equation 8 reduces to Eq. 9 by substituting Eq. 7.
$$\frac{\delta Ip}{\delta t}={k}_{2}{\xi }_{2}E-{\theta I}_{p}$$
(9)
$$\frac{\delta Ia}{\delta t}={k}_{1}{\xi }_{1}E-{\Upsilon }_{1}Ia$$
(10)
$$\frac{\delta Ib}{\delta t}={\theta {\xi }_{3}I}_{p}-{\Upsilon }_{2}Ib$$
(11)
$$\frac{\delta Ic}{\delta t}={\theta {\xi }_{4}I}_{p}-hIc-{\mu }_{o}Ic$$
(12)
The change in the Hospitalised, Death and Recovered class is given by Eqs. 13, 14, and 15, respectively.
$$\frac{\delta H}{\delta t}=hIc-{\Upsilon }_{3}H-\mu H$$
(13)
$$\frac{\delta D}{\delta t}=\mu H+{\mu }_{o}Ic$$
(14)
$$\frac{\delta R}{\delta t}={\Upsilon }_{1}Ia+{\Upsilon }_{2}Ib+{\Upsilon }_{3}H$$
(15)
Where, ϒ1, ϒ2 and ϒ3 are the daily recovery rates of individuals with Asymptomatic, Mild and Moderate Symptoms and Severe and Critical Symptoms, respectively. θ is the rate at which pre-symptomatic individuals develop symptoms. h is the rate at which individuals who have developed severe and critical cases are hospitalised. μo is the daily death rate due to direct and indirect effects of COVID-19 in individuals with severe and critical symptoms who have not been hospitalised. \(\mu\) is the daily death rate due to COVID-19 in hospitalised individuals.
$$\frac{\delta Iabc}{\delta t}=\frac{\delta Ia}{\delta t}+\frac{\delta Ip}{\delta t}+\frac{\delta Ib}{\delta t}+\frac{\delta Ic}{\delta t}+\frac{\delta H}{\delta t}$$
(16)
$$\frac{\delta Iabc}{\delta t}={k}_{1}{\xi }_{1}E-{\Upsilon }_{1}Ia+{k}_{2}{\xi }_{2}E-{\theta I}_{p}+{\theta {\xi }_{3}I}_{p}-{\Upsilon }_{2}Ib+{\theta {\xi }_{4}I}_{p}-hIc-{\mu }_{o}Ic+hIc-{\Upsilon }_{3}H-\mu H$$
(17)
$$\frac{\delta R}{\delta t}{=\Upsilon }_{1}Ia+{\Upsilon }_{2}Ib+{\Upsilon }_{3}Ic$$
(18)
$$\frac{\delta R}{\delta t}{=\Upsilon }_{1}\left(Ia+Ib+Ic\right)={\Upsilon }_{1}(Iabc)$$
(19)
Model basic reproductive number and herd immunity
For the feasible region in the octant of the mathematic model, there exists an equilibrium in which there is no disease called the Disease-Free Equilibrium (DFE). This condition is satisfied by the stability of the rate of change in the population. Thus:
$$\frac{\delta N}{\delta t}=\frac{\delta S}{\delta t}+\frac{\delta E}{\delta t}+\frac{\delta Ia}{\delta t}+\frac{\delta Ip}{\delta t}+\frac{\delta Ib}{\delta t}+\frac{\delta Ic}{\delta t}+\frac{\delta H}{\delta t}+\frac{\delta R}{\delta t}+\frac{\delta D}{\delta t}$$
(20)
$$\frac{\delta N}{\delta t}=-A\beta S\frac{Iabc}{N}+A\beta S\frac{Iabc}{N}-{k}_{1}{\xi }_{1}E-{k}_{2}{\xi }_{2}E+{k}_{1}{\xi }_{1}E-{\Upsilon }_{1}Ia+{k}_{2}{\xi }_{2}E-{\theta I}_{p}+{\theta {\xi }_{3}I}_{p}-{\Upsilon }_{2}Ib+{\theta {\xi }_{4}I}_{p}-hIc-{\mu }_{o}Ic+hIc-{\Upsilon }_{3}H-\mu H+\mu H+{\mu }_{o}Ic+ {\Upsilon }_{1}Ia+{\Upsilon }_{2}Ib+{\Upsilon }_{3}H$$
(21)
$$\frac{\delta N}{\delta t}=0$$
(22)
$$\frac{\delta N}{\delta t}\ge 0\ge S+E+Iabc+R+D$$
(23)
At Disease Free Equilibrium (DFE), E = 0, I = 0, R = 0, D = 0. Therefore, substituting DFE values into Eq. 23 gives Eq. 24.
$$\frac{\delta N}{\delta t}\ge 0\ge S$$
(24)
At each point in time in the model, there also exists an equilibrium in which there is a maximum/minimum for each class. This equilibrium is called the Endemic Equilibrium (EE). Thus, taking Eq. 3 and Eqs. 9, 10, 11, 12, 13, 14 and 15 where the rate of change is 0 gives Eqs. 25, 26, 27, 28, 29, 30, 31, 32 and 33
$$\frac{\partial S}{\partial t}=-A\beta S\frac{Iabc}{N}=0$$
(25)
$$\frac{\delta E}{\delta t}= A\beta S\frac{Iabc}{N}-{k}_{1}{\xi }_{1}E-{k}_{2}{\xi }_{2}E=0$$
(26)
$$\frac{\delta Ip}{\delta t}={k}_{2}{\xi }_{2}E-{\theta I}_{p}=0$$
(27)
$$\frac{\delta Ia}{\delta t}={k}_{1}{\xi }_{1}E-{\Upsilon }_{1}Ia=0$$
(28)
$$\frac{\delta Ib}{\delta t}={\theta {\xi }_{3}I}_{p}-{\Upsilon }_{2}Ib=0$$
(29)
$$\frac{\delta Ic}{\delta t}={\theta {\xi }_{4}I}_{p}-hIc-{\mu }_{o}Ic=0$$
(30)
$$\frac{\delta H}{\delta t}=hIc-{\Upsilon }_{3}H-\mu H=0$$
(31)
$$\frac{\delta D}{\delta t}=\mu H+{\mu }_{o}Ic=0$$
(32)
$$\frac{\delta R}{\delta t}={\Upsilon }_{1}Ia+{\Upsilon }_{2}Ib+{\Upsilon }_{3}H=0$$
(33)
Substituting Eq. 27–31 into Eq. 2 gives Eq. 34.
$$Iabc=Ia+Ip+Ib+Ic+H= \frac{E{k}_{2}{\xi }_{2}}{\theta }+\frac{E{k}_{1}{\xi }_{1}}{{\Upsilon }_{1}}+\frac{{\theta {\xi }_{3}I}_{p}}{{\Upsilon }_{1}}+\frac{{\theta {\xi }_{4}I}_{p}}{(h+{\mu }_{o})}+\frac{hIc}{({\Upsilon }_{3}+\mu )}$$
(34)
Substituting Eqs. 27 and 30 into Eq. 34 gives Eq. 35.
$$Iabc=\frac{E{k}_{2}{\xi }_{2}}{\theta }+\frac{E{k}_{1}{\xi }_{1}}{{\Upsilon }_{1}}+\frac{{E{{\xi }_{2}\xi }_{3}k}_{2}}{{\Upsilon }_{2}}+\frac{{E{{\xi }_{2}\xi }_{4}k}_{2}}{(h+{\mu }_{o})}+\frac{h{E{{\xi }_{2}\xi }_{4}k}_{2}}{(h+{\mu }_{o})({\Upsilon }_{3}+\mu )}$$
(35)
From Eq. 26, making S the subject of the equation gives Eq. 36:
$$S=\frac{NE {(k}_{1}{\xi }_{1}+{k}_{2}{\xi }_{2})}{A\beta Iabc}$$
(36)
Substituting Eq. 35 into Eq. 36 gives Eq. 37.
$$S=\frac{N{(k}_{1}{\xi }_{1}+{k}_{2}{\xi }_{2}) \theta {\Upsilon }_{1}{\Upsilon }_{2}(h+{\mu }_{o})({\Upsilon }_{3}+\mu )}{A\beta ({k}_{2}{\xi }_{2}{\Upsilon }_{1}{\Upsilon }_{2} \left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{1}{\xi }_{1}\theta {\Upsilon }_{2 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{2}{{\xi }_{2}\xi }_{3}\theta {\Upsilon }_{1 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+ {k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 }\left({\Upsilon }_{3}+\mu \right)+h{k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2}}$$
(37)
The relative critical point for the model is when DFE = EE. Using Eq. 24 and Eq. 36 we derive Eqs. 38, 39 and 40.
$$S=\frac{N{(k}_{1}{\xi }_{1}+{k}_{2}{\xi }_{2}) \theta {\Upsilon }_{1}{\Upsilon }_{2}(h+{\mu }_{o})({\Upsilon }_{3}+\mu )}{A\beta ({k}_{2}{\xi }_{2}{\Upsilon }_{1}{\Upsilon }_{2} \left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{1}{\xi }_{1}\theta {\Upsilon }_{2 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{2}{{\xi }_{2}\xi }_{3}\theta {\Upsilon }_{1 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+ {k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 }\left({\Upsilon }_{3}+\mu \right)+h{k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 }}\le N$$
(38)
$$\frac{N{(k}_{1}{\xi }_{1}+{k}_{2}{\xi }_{2}) \theta {\Upsilon }_{1}{\Upsilon }_{2}(h+{\mu }_{o})({\Upsilon }_{3}+\mu )}{A\beta ({k}_{2}{\xi }_{2}{\Upsilon }_{1}{\Upsilon }_{2} \left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{1}{\xi }_{1}\theta {\Upsilon }_{2 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{2}{{\xi }_{2}\xi }_{3}\theta {\Upsilon }_{1 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+ {k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 }\left({\Upsilon }_{3}+\mu \right)+h{k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 }}\le 1$$
(39)
$$1 \le {R}_{0}=\frac{A\beta ({k}_{2}{\xi }_{2}{\Upsilon }_{1}{\Upsilon }_{2} \left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{1}{\xi }_{1}\theta {\Upsilon }_{2 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+{k}_{2}{{\xi }_{2}\xi }_{3}\theta {\Upsilon }_{1 }\left(h+{\mu }_{o}\right)\left({\Upsilon }_{3}+\mu \right)+ {k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 }\left({\Upsilon }_{3}+\mu \right)+h{k}_{2}{{\xi }_{2}\xi }_{4}\theta {\Upsilon }_{1 }{\Upsilon }_{2 })}{{(k}_{1}{\xi }_{1}+{k}_{2}{\xi }_{2}) \theta {\Upsilon }_{1}{\Upsilon }_{2}(h+{\mu }_{o})({\Upsilon }_{3}+\mu )}$$
(40)
Equation 40 is what is defined as the basic reproductive number (R0) for the ARI COVID-19 SEIR model. The basic reproductive number is the number of secondary infections that one infected person would produce in a fully susceptible population through the entire duration of the infectious period. R0 provides a threshold condition for the stability of the disease-free equilibrium point [1]. If R0 is greater than 1 then there is an endemic equilibrium thus there will be an epidemic. If R0 is less than 1 then the disease will die out and remain at a relatively low level to the population size. As can be seen from Eq. 40, R0 can be summarised as a ratio of the daily contact number over the daily recovery rate. It can also be defined by Eq. 41:
$$Ro=(Number\,of\,contacts\,per\,time)\times\,(Probability\,of\,transmission\,per\,contact)\times\,(Duration\,of\,Infection)$$
(41)
The basic reproductive number assumes a completely susceptible population however the infection productiveness of the population changes as infections increase. The effective reproductive number, Re, is the number of secondary infections that one infected person would produce through the entire duration of the infectious period. It can be estimated based on the susceptible class given by Eq. 42:
$$Re={R}_{0}\times S$$
(42)
Herd immunity is an important concept in epidemiology. Herd immunity is when the population has enough people immune such that the disease will not spread if it was suddenly introduced randomly into the population. Consider if α is the fraction immune due to vaccination/acquired immunity. Then Heard Immunity is given by Eq. 43:
$$Herd\,Immunity\,(\mathrm{\alpha })=(1-\frac{1}{Ro})$$
(43)
Model parameters
Determining the hospital discharge rate in South Africa
The average Daily Hospital Discharge Rate (ϒ3) in South Africa was calculated based on clinical information from admitted patients with laboratory-confirmed COVID-19 in selected hospitals in South Africa under the National Institute for Communicable Diseases (NICD) DATCOV surveillance system. The NICD sentinel hospital surveillance system was designed to monitor and describe trends of COVID-19 hospitalizations and the epidemiology of hospitalized patients in South Africa [43]. The number of hospitals reporting in the NICD DATCOV surveillance system increased in the reporting period. Initially, 204 Facilities were reported, and this increased to 434 Facilities by 4 September 2020 [44]. Therefore, caution was taken when taking averages between reporting case dates. The average Daily Hospital Discharge Rate (ϒ3) for COVID-19 patients was calculated using the Number of Discharged Alive and Admitted patients Data in the NICD DATCOV surveillance system from 24 May to 01 October 2020. Equation 44 was used to calculate the Daily Discharge:
$${Daily\,Discharge \left(n\right)}_{i+1}={Discharged\,Alive\,(n)}_{i+1}-{Discharged\,Alive\,(n)}_{i}$$
(44)
Where n is the number of patients and i is the reported case date.
The average Daily Hospital Discharge Rate (ϒ3) was then calculated using Eq. 45.
$$\gamma_3\;=\;\frac{\overline{Daily\;Discharge{\;\left(n\right)}_{i+1}}}{Admitted_{i+1}}$$
(45)
Determining the death rate in hospitalised cases in South Africa
The average Daily Death Rate or Daily Case Fatality Rate (CFR) for COVID-19 patients (μ1) was calculated using the Number of Daily Deaths and Admitted patients Data in the NICD DATCOV surveillance system from 24 May-01 October 2020. The WHO guideline in estimating the CFR in [45] was followed. Equation 46 was used to calculate the Daily Deaths:
$${Daily\,Deaths\,\left(n\right)}_{i+1}={Died\,(n)}_{i+1}-{Died\,(n)}_{i}$$
(46)
where n is the number of patients and i is the reported case date.
The Daily Case Fatality Rate (CFR) for COVID-19 patients (μ) was then calculated using Eq. 47.
$$\mu\;=\;\frac{\overline{Daily\;Deaths\;{\left(n\right)}_{i+1}}}{{Admitted}_{i+1}}$$
(47)
Determining the death rate of unreported severe and critical cases in South Africa
Excess mortality is a count of deaths from all causes relative to what would normally have been expected. Excess mortality/deaths allow for accounting for miscounted or underreported COVID-19 Deaths and indirect Deaths related to the COVID-19 pandemic. National statistical agencies publish weekly deaths and averages of past ‘normal’ deaths [46]. In South Africa, the South African Medical Research Council (SAMRC) published the Excess deaths from 29 December 2019 to 01 October 2020 using information obtained from the National Population Register [47, 48]. The Unreported Excess Deaths (Natural) to COVID-19 Death Ratio was calculated for the period, 25 March to 01 October 2020 with data from [47] using Eq. 48:
$$Unreported\,Excess\,Deaths\,(Natural)\,to\,COVID-19\,Death\,Ratio=\frac{Excess\,Deaths\,(Natural)i-{Weekly\,Reported\,COVID-19\,Deaths }_{i}}{{Weekly\,Reported\,COVID-19\,Deaths }_{i}}$$
(48)
where i is the Weekly Reported Date. The daily death rate due to direct and indirect effects of COVID-19 in individuals with severe and critical symptoms who have not been hospitalised (μo) was then calculated using Eq. 49:
$${\upmu }_{o}=\stackrel{-}{Unreported\,Excess\,Deaths\,\left(Natural\right)\,to\,COVID-19\,Death\,Ratio}\times\,{\upmu }_{1}$$
(49)
Determining the COVID-19 patient admission status in South Africa
The average admission status for COVID-19 patients was calculated using the Number of patients Currently in Hospital (n), General Ward (n), High Care (n), Intensive Care Unit (n), Isolation Ward (n), On Oxygen (n) and On Ventilator (n) Data in the NICD DATCOV surveillance system from 24 May-1 October 2020. The average admission Status was calculated using Eqs. 50, 51, 52, 53, 54 and 55:
$$\mathrm{General\,Ward }(\mathrm{\%})=\overline{\frac{{\mathrm{General\,Ward\,}\left(\mathrm{n}\right)}_{i}}{{\mathrm{Currently\,in\,Hospital\,}(\mathrm{n})}_{i}}}\times 100$$
(50)
$$\mathrm{High\,Care }(\mathrm{\%}) =\overline{\frac{{\mathrm{High\,Care\,}(\mathrm{n}) }_{i}}{{\mathrm{Currently\,in\,Hospital\,}(\mathrm{n})}_{i}}}\times\,100$$
(51)
$$\mathrm{Intensive\,Care\,Unit\,}(\mathrm{\%})=\overline{\frac{{\mathrm{Intensive\,Care\,Unit\,}(\mathrm{n}) }_{i}}{{\mathrm{Currently\,in\,Hospital\,}(\mathrm{n})}_{i}}}\times 100$$
(52)
$$\mathrm{Isolation\,Ward\,}(\mathrm{\%})=\overline{\frac{{\mathrm{Isolation\,Ward\,}(\mathrm{n})}_{i}}{{\mathrm{Currently\,in\,Hospital\,}(\mathrm{n})}_{i}}}\times 100$$
(53)
$$\mathrm{On\,Oxygen }(\mathrm{\%}) =\overline{\frac{{\mathrm{On\,Oxygen\,}(\mathrm{n})}_{i}}{{\mathrm{Currently\,in\,Hospital\,}(\mathrm{n})}_{i}}}\times 100$$
(54)
$$\mathrm{On\,Ventilator\,}(\mathrm{\%}) =\overline{\frac{{\mathrm{On\,Ventilator\,}(\mathrm{n})}_{i} }{{\mathrm{Currently\,in\,Hospital\,}(\mathrm{n})}_{i}}}\times 100$$
(55)
The admission status was then calculated using the Hospitalised Cases (H) in the model with Eq. 56:
$$\mathrm{Admitted }\left(\mathrm{n}\right) =Admitted\,(\%)\times H$$
(56)
where the Admitted (%) is the General Ward (%), High Care (%), Intensive Care Unit (%), Isolation Ward (%), On Oxygen (%) and Ventilator (%) respectively. A summary of the model parameters used in the ARI-COVID-19 Model is given in Table 1:
Model NPI scenarios and seeding
Modelling periods and reported case data
To model the impact of NP1s in South Africa, the National Lockdown Alert Level 5, 4, and 3 policies implemented by the South African government were modelled as scenarios in the ARI COVID-19 SEIR model, respectively. A “No lockdown” scenario where there was no policy response from the South African government was also modelled in the ARI COVID-19 SEIR model.
In each scenario, the ARI COVID-19 SEIR model was seeded using COVID-19 deaths with the Unreported Excess Deaths (Natural) to COVID-19 Death Ratio used to account for excess deaths in the National Lockdown Alert Level 3 scenario. For the “No lockdown” Scenario, reported active cases were used to seed the model due to no reported COVID-19 and Excess (Natural) Deaths in this period. Table 2 shows the model classification of the South African COVID-19 policy and the seeding period used in the model.
Cumulative Daily Reported COVID-19 Case, Recovery and Death Data for South Africa were obtained from the Johns Hopkins University (JHU) Center for Systems Science and Engineering (CSSE) COVID-19 Database [56] for the period of 22 January 2020 to 1 October 2020. Since reported case data come after a period of clinical diagnosis, reported case dates thus are lagged from the “real-time” date of infection. Therefore, an Average Date of infection date was estimated using Eq. 57:
$$Average\,Date\,of\,Infection=Reported\,Case\,Date-Average\,Time\,for\,Clinical\,Diagnosis$$
(57)
where the Average time for Clinical diagnosis (Ttesting) is the average time taken for an infected person to be diagnosed and the diagnosis outcome to be classified and reported as a COVID-19 case. The COVID-19 Active Cases were determined using Eq. 58:
$$Active\,Cases=Sum\,of\,Cases-Sum\,of\,Recovered\,Cases-Sum\,of\,Deaths$$
(58)
Regression and statistical analysis
To seed the model, a non-linear regression analysis was conducted between the Reported COVID-19 Deaths and Death due to COVID-19 (D) from the model for the model seeding periods stated in Table 2. Seeding the Models with points that are oversensitive results in large deviations between modelled data versus reported case data. This deviation or noise introduced by “over-sensitive” data points creates a significant error in the model results. Thus, it was important to decide which Data points can be used in the regression analysis. For the modelling of NPIs, this was particularly important at the start of the pandemic when data values are relatively small. To decide on seeding data points, a Data Point Sensitivity term described by Eq. 59 was used:
$$\mathrm{Data\,Point\,Sensitivity}=\frac{\mathrm{Active\,Case\,or\,Reported\,COVID}-19\mathrm{\,Deaths }}{Lowest\,Possible\, Data\,Unit}$$
(59)
where the Lowest Possible Data unit is the lowest possible unit of measurement for that data. In this case, it is 1 COVID-19 Reported Case. Data points with a Data Point Sensitivity greater than 5% were ignored in the regression analysis. To conduct the regression analysis, the Residual and Normalised Errors were determined using Eq. 60 and Eq. 61:
$$Residual=Modelled\,Data- \mathrm{Active\,Case\,Data}/\mathrm{Reported\,COVID}-19\mathrm\,{ Deaths}$$
(60)
$$Normalised\,Error=\frac{Residual}{Active\,Cases\,or\,Reported\,COVID-19\,Deaths}$$
(61)
To allow the goodness of fit of Modelled data to Active Case or Reported COVID-19 Deaths Data, the Average Normalised Error of all Data points used in the Regression Analysis was reduced to 0 by changing the Effective Daily Contact Number (β) using the What-If Analysis Function in MS Excel. The pooled sample variance (s2) was then calculated using Eq. 62:
$${s}^{2}=\frac{\sum_{i=1}^{n1}{(xi-\overline{x1 })}^{2}+\sum_{i=1}^{n1}{(xj-\overline{x2 })}^{2}}{n1+n2-2}$$
(62)
Where xi is the Reported Case Data points, \(\overline{x1 }\) is the Reported Case Data Points Mean, xj is the Model Data points, \(\overline{x2 }\) is the Model Data Points Mean, n1 and n2 are the sample sizes. The t-value: Two-Sample Assuming Equal Variance was used to calculate the t-value using Eq. 63
$$t=\frac{\overline{x1 }-\overline{x2}}{\sqrt{{s }^{2} (\frac{1}{n1}+\frac{1}{n2})}}$$
(63)
For the No lockdown scenario, since reported active cases were used to seed this model scenario, a sensitivity analysis was run to determine the fraction of reported cases that are pre-symptomatic, mild & moderate, asymptomatic, and critical & severe. This was done using a VBA code following the computational steps outlined in Fig. 2.
ARI COVID-19 SEIR Model Scenario Statistical Analysis Computational Steps
The combination of the Φ1 Φ2 Φ3 Φ4 that resulted in the lowest T-value (most significance) between model data and reported case data for the No Lockdown Scenario are shown in Table 3 and were chosen to seed the No Lockdown Scenario.
For other Model Scenarios (South Africa National Lockdown Alert Level 5, 4 and 3), Reported COVID-19 Deaths were used to seed the Models. However, for the National Lockdown Alert Level 3 the reported COVID-19 deaths were adjusted to include unreported COVID-19 Deaths using Eq. 64:
$$Reported\ COVID-19\,Deaths\ \left(Seed\right)=Reported\ COVID-19\ Deaths+Reported\ COVID-19\ Deaths\times \stackrel{-}{Unreported\ Excess\ Deaths\ \left(Natural\right)\ to\ COVID-19\ Death\ Ratio}$$
(64)
Validation and limitations of model
The functionality and data produced by the ARI COVID-19 SEIR model were validated using the following:
-
For model functionality, the model ODEs were validated by setting the model scenario reproductive number (Ro) to 1 by changing the Effective Daily Contact Number (β) using the What-If Analysis Function in MS Excel. At Ro = 1, there should be no transmission of the population between model classes and model class values should be either at initial values or 0.
-
For model data, a comparison between COVID-19 Admitted data in the NICD DATCOV surveillance system and Hospitalised COVID-19 Cases (H) in the model. Comparison between Excess (Natural) deaths data and Death due to COVID-19 (D) in the model. Comparison between the reported date of peak Active cases and Total infections (Iabc) in the model. Comparison between the ARI COVID-19 SEIR Model results and results from the South Africa National COVID-19 Epi Model (NCEM), the Center for Disease Dynamics, Economics & Policy (CDDEP) COVID-19 Model and the Imperial College London COVID-19 Model.
The following limitations were identified during the modelling:
-
A non-linear regression analysis could not be conducted for the National Lockdown Alert Levels 2 and 1 due to these policies being implemented in the negative exponential phase of the first COVID-19 epidemic wave in South Africa. The limitation of non-linear regression analysis in the negative exponential phase of an epidemic is that the decrease in cases is not only due to a reduction in contact but also due to the reduction in the susceptibles and increase in the recovered cases in the population (decrease in the mass incident term) at disease-free equilibrium.
-
There was a limitation in determining the CFR in South Africa during the no-lowdown period before the first COVID-19 epidemic wave via regression analysis due to the low number of reported COVID-19 and excess deaths. The CFR could have been adjusted using functions that relate to active COVID-19 cases however such functions are limited. There was also a limitation in using the CFR to seed the model for the no lock down scenario hence active COVID-19 cases were used instead.
-
In this study, bootstrapping was not performed to estimate the variance of the ARI COVID-19 SEIR Model Effective Daily Contact Number (β) due to the small sample size used in the model seeding period.
-
Pharmaceutical interventions were not accounted for in the model.
-
The COVID-19 NPIs investigated in this study were limited to the National Lockdown Alert Levels.
-
The NPIs investigated in this model were a collection of measures of the National Lockdown Alert Levels. Individual NPI measures in the National Lockdown Alert Levels could not be modelled individually due to limited data.
