Evaluating the Interplay between Transmissibility and Virulence of SARS-CoV-2 by Mathematical Modeling

ABSTRACT

During the first months of 2020 SARS-CoV-2 spread to all continents, virtually reaching all countries. In the subsequent months, new variants emerged in different regions of the world. A mathematical model based on the Covid-19 natural history encompassing the age-dependent fatality was applied to evaluate SARS-CoV-2 transmissibility and virulence. Transmissibility was assessed by calculating the basic reproduction number R ₀ and virulence by counting the proportion of severe Covid-19 cases and deaths. The model parameters were adjusted against the data observed in the state of São Paulo, Brazil, considering two different levels of virulence. The severe Covid-19 cases and deaths were three times higher and R ₀ was 25% lower when the more virulent SARS-CoV-2 variant was compared to the less virulent variant. However, under the high-virulence scenario the number of transmitting individuals is 25% lower, mainly due to the isolation of symptomatic individuals. The corollary that transmission increases in the low virulence scenario is also true. The estimated parameters, using data from São Paulo up to May 13, 2020, showed that the Covid-19 epidemic predicted with low virulence SARS-CoV-2 transmission matched the observed data just before the beginning of the relaxation, which occurred by mid-June 2020. The assessment of the interplay between transmissibility can be applied to explain in somehow the appearance of gamma and omicron variants of concern in São Paulo.

Keywords:
mathematical model; SARS-CoV-2; basic reproduction number; asymptomatic individuals; Covid-19 data

1 INTRODUCTION

The severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), an RNA virus, spreads mostly through air by means of aerosols ¹, but can be also transmitted by droplets produced by coughing, sneezing or talking, or by fomites ^{2) (3}. The coronavirus disease 2019 (Covid-19) caused by SARS-CoV-2 infection was declared a pandemic by the World Health Organization on March 11, 2020. In general, the fatality rate is considerably higher in senior patients (60 years or more) ⁴.

SARS-CoV-2 is an RNA-virus presenting relatively high mutation rate (See ^{Appendix A} for details). Briefly, the original SARS-CoV-2 variant was first detected in Wuhan, China by December 2019, and the omicron variant, that appeared in South Africa in November 2021, is the most phylogentically divergent variant of concern. The omicron variant showed high transmissibility but less virulence than the original SARS-CoV-2.

We developed a mathematical model to evaluate the interplay between transmissibility and virulence during the Covid-19 epidemic in the state of São Paulo, Brazil. We hypothesized that transmissibility is measured by the transmission rate (ultimately, by the basic reproduction number R ₀), and that virulence is assessed by the ratio between asymptomatic and symptomatic individuals. It is worth stressing that the additional mortality due to Covid-19 results in an asymptotic value, which is not an equilibrium value at steady sate. Hence, a special method is required to assess the steady state to obtain R ₀. Additionally, as individuals with severe Covid-19 are hospitalized (isolated), and they must not be considered as transmitting agents of SARS-CoV-2, rather only the asymptomatic or pre-symptomatic or mild Covid-19 individuals are spreading the disease ¹⁹. In section 2, we present and analyze a model that describes SARS-CoV-2 transmission, in section 3, using the data for São Paulo State, we estimate the basic reproduction number R ₀ considering two different virulence levels. The discussion and conclusion are given in sections 4 and 5.

2 MATERIAL AND METHODS

One of the main aspects of Covid-19 is the high fatality rate of older people. For this reason, we divided the population in two groups, one composed of young (60 years old or less, denoted by subscript y) and the other containing individuals of 60 years old or more, denoted by subscript o. The community’s vital dynamic is described by the percapita rates of birth (ϕ) and mortality (µ), and φ is the aging rate, that is, the flow from the young subpopulation y to the senior subpopulation o. Another aspect is SARS-CoV-2 transmission by asymptomatic and pre-symptomatic individuals, that is, individuals that never develop symptoms and those before the onset of the disease ²⁰. Hence, for each subpopulation j( j = y, o), individuals are divided into seven compartments: susceptibles S _j , isolated (quarantine) Q _j , exposed and incubating E _j , asymptomatic A _j , pre-symptomatic (or pre-diseased) individuals P _j , symptomatic individuals with mild Covid19 M _j , and individuals with severe Covid-19 D _j . All young and senior individuals in classes A _j , M _j , and D _j enter to the same recovered class R. Hence, the SQEAPMDR model has 15 compartments. The compartments irrespective to age are E = E _y + E _o (exposed), A = A _y + A _o (asymptomatic), P = P _y + P _o (preymptomatic), M = M _y + M _o (mild Covid-19), and D = D _y + D _o (severe Covid-19). Table 1 summarizes the model compartments (variables).

Table 1:
Summary of the model variables, for j = y (young) and j = o (senior).

Symbol	Meaning
S j	Susceptible individuals
Q j	Isolated individuals
E j	Exposed and incubating SARS-CoV-2 individuals
A j	Asymptomatic individuals
P j	Pre-diseased (pre-symptomatic) individuals
M j	Mild (non-hospitalized) Covid-19 individuals
D j	Severe (hospitalized) Covid-19 individuals
R	Recovered (immune) individuals

The susceptible individuals are either in isolation (quarantine) or circulating. To express isolation (the flow from S to Q), we consider a unique pulse in isolation at time t = τ, described by u _j S _j δ(t − τ), with j = y, o - The fraction of isolated people is u _j and δ(x) is the Dirac delta function, that is, δ(x) = ∞, if x = 0, otherwise, δ (x) = 0, with ${\int }_0^{\infty }\delta \left(x\right)dx=1$. The equations that describe isolation are as follows

$$\begin{gathered} \left\{\begin{array}{l}\frac{d}{dt}{S}_y=\varphi N-\left(\phi +\mu \right)S_y-\lambda S_y-u_y{S}_y\delta \left(t-\tau \right) \\ \frac{d}{dt}{S}_o=\phi S_y-\mu S_o-\lambda \psi S_o-u_o{S}_o\delta \left(t-\tau \right)\\ \frac{d}{dt}{Q}_j=u_j{S}_j\delta \left(t-\tau \right)-\mu Q_j.\end{array}\right. \end{gathered}$$ (2.1)

On the other hand, the natural history of Covid-19 is the same for the young (j = y) and senior (j = o) subpopulations in circulation. We assume that individuals in the asymptomatic (A _j ), pre-symptomatic (P _j ), and a fraction z _j of mild Covid-19 (M _j ) class are transmitting the virus. Other infected classes ((1 − z _j ) M _j and D _j ) are under voluntary or forced isolation. Susceptible individuals are infected at a rate λ _j S _j (known as the mass action law ²¹), where λ _j is the per-capita incidence rate (or force of infection) defined by λ _j = λ (δ _jy + ψδ _jo ), with λ being

$\lambda =\frac{\epsilon }{N}\left({\beta }_{1y}{A}_y+{\beta }_{2y}{P}_y+{\beta }_{3y}{z}_y{M}_y+{\beta }_{1o}{A}_o+{\beta }_{2o}{P}_o+{\beta }_{3o}{z}_o{M}_o\right),$ (2.2)

where δ _ij is the Kronecker delta, with δ _ij = 1 if i = j, and 0, if i ≠ j; and β _1j , β _2j and β _3j are the transmission rates - the rates at which a virus encounters a susceptible people and infects him/her. The parameter ψ ≥ 1 measures more susceptibility among seniors than youngs. The parameter ε ≤ 1 diminishes the transmission rates - the protection factor ε decreases viral transmission by individual (face mask, hygiene, etc.) and collective (social distancing) protective measures.

Susceptible individuals are infected at rate λ _j and enter into class E _j . After an average period 1/σ _j in class E _j , where σ _j is the incubation rate, exposed individuals develop into asymptomatic A _j (with probability l _j ) or pre-symptomatic class P _j (with probability 1 − l _j ). After an average period 1/γ _j in class A _j , where γ _j is the recovery rate of asymptomatic individuals, asymptomatic individuals acquire immunity (recovered) and enter the recovered class R. Asymptomatic individuals may manifest symptoms at the end of this period, and a fraction 1 − χ _j becomes mildly ill (class M _j ). As for the symptomatic individuals, after an average period 1/γ _1j in class P _j , where γ _1j is the infection rate of pre-symptomatic individuals, these individuals become either severely ill (class D _j ) (with probability 1 − k _j ) or mildly ill (class M _j ) (with probability k _j ). Individuals in class D _j acquire immunity after a period 1/γ _2j , where γ _2j is the recovery rate of severe Covid-19, and enter the recovered class R or die under the disease-induced (additional) mortality rate α _j . Individuals in the mild Covid-19 class M _j acquire immunity after a period 1/γ _3j , where γ _3j is the recovery rate of mild Covid-19, and enter the recovered class R. Table 2 summarizes the model parameters. The description of the assigned values to the Covid-19 natural history parameters can be found in ²² (the parameters marked with * are estimated considering two SARS-CoV-2 virulence levels).

Table 2:
Summary of the model parameters (j = y, o) and values (rates in days ⁻¹ and proportions are dimensionless) from ²². (*) Values fitted considering two levels of the SARS-CoV-2 virulence.

Symbol	Meaning	Value
µ	Natural mortality rate	1/(78.4 × 365)
ϕ	Birth rate	1/(78.4 × 365)
φ	Aging rate	6.7 × 10−6
σy (σ o )	Incubation rate	1/5.8 (1/5.8)
γy (γ o )	Recovery rate of asymptomatic individuals	1/12 (1/14)
γ1y (γ 1o )	Infection rate of pre-symptomatic individuals	1/4 (1/4)
γ2y (γ 2o )	Recovery rate of severe Covid-19 individuals	1/12 (1/21)
γ3y (γ 3o )	Infection rate of mild Covid-19 individuals	1/13 (1/16)
τ	Time of the introduction of isolation	March 24, 2020
zy (z o )	Proportion of transmission by mild Covid-19 individuals	0.5 (0.2)
ψ	Increased susceptibility among seniors	1.17
χy (χ o )	Proportion of remaining as asymptomatic individuals	0.98 (0.95)
ly (l o )	Proportion of asymptomatic individuals	∗
ky (k o )	Proportion of mild (non-hospitalized) Covid-19	∗
ε	Protection factor	∗
uy (u o )	Proportion in isolation	∗
αy (α o )	Additional mortality rate	∗
β1y (β 1o )	Transmission rate due to asymptomatic individuals	∗
β2y (β 2o )	Transmission rate due to pre-symptomatic individuals	∗
β3y (β 3o )	Transmission rate due to mild Covid-19 individuals	∗

The SARS-CoV-2 transmission model is described by the system of ordinary differential equations, named the SQEAPMDR model considering young (y) and senior (o) subpopulations. The dynamic equations for S, Q, E, A, P, M, D, and R are obtained through the balance of inflow and outflow in each compartment. We drop out the pulses in equation (2.1) describing quarantine individuals and transfer them to the pulse conditions at time of quarantine. Hence, the equations for the susceptible individuals are

$\left\{\begin{array}{l}\frac{d}{dt}{S}_y=\varphi N-\left(\phi +\mu \right)S_y-\lambda S_y\\ \frac{d}{dt}{S}_o=\phi S_y-\mu S_o-\lambda \psi S_o,\end{array}\right.$ (2.3)

for isolated individuals Q _j and infected individuals, with j = y, o,

$\left\{\begin{array}{l}\frac{d}{dt}{Q}_j=-\mu Q_j\\ \frac{d}{dt}{E}_j=\lambda \left({\delta }_{jy}+\psi {\delta }_{jo}\right)S_j-\left({\sigma }_j+\mu \right)E_j\\ \frac{d}{dt}{A}_j=l_j{\sigma }_j{E}_j-\left({\gamma }_j+\mu \right)A_j\\ \frac{d}{dt}{P}_j=\left(1-l_j\right){\sigma }_j{E}_j-\left({\gamma }_{1j}+\mu \right)P_j\\ \frac{d}{dt}{M}_j=\left(1-{\chi }_j\right){\gamma }_j{A}_j+k_j{\gamma }_{1j}{P}_j-\left({\gamma }_{3j}+\mu \right)M_j\\ \frac{d}{dt}{D}_j=\left(1-k_j\right){\gamma }_{1j}{P}_j-\left({\gamma }_{2j}+\mu +{\alpha }_j\right)D_j,\end{array}\right.$ (2.4)

and for recovered individuals,

$\frac{d}{dt}R={\chi }_y{\gamma }_y{A}_y+{\gamma }_{3y}{M}_y+{\gamma }_{2y}{D}_y+{\chi }_o{\gamma }_o{A}_o+{\gamma }_{3o}{M}_o+{\gamma }_{2o}{D}_o-\mu R,$ (2.5)

where N _j = S _j + Q _j + E _j + A _j + P _j + M _j + D _j , and N = N _y + N _o + R obeys

$\frac{d}{dt}N=\left(\varphi -\mu \right)N-{\alpha }_y{D}_y-{\alpha }_o{D}_o,$ (2.6)

with the initial number of population being N(0) = N ₀ = N _0y + N _0o , where N _0y and N _0o are the size of the young and senior subpopulations at t = 0.

The initial conditions (simulation time t = 0 corresponding to the calendar time when the first case was confirmed: February 26 for São Paulo State) supplied to the system of equations (2.3), (2.4), and (2.5) are, for j = y, o,

$S_j\left(0\right)=N_{0j},X_j\left(0\right)=n_{X_j},\text{where}{X}_j=Q_j,E_j,A_j,P_j,M_j,D_j,R,$ (2.7)

where $n_{X_j}$ is a non-negative number. For instance, $n_{E_y}=n_{E_o}=0$ implies that there was not any exposed person (young or senior) at the beginning of the epidemic. The isolation implemented at τ = 27 (corresponding to calendar time March 24, 2020) is described by the pulse conditions at time τ

$\left\{\begin{array}{l}{S}_j\left({\tau }^{+}\right)=S_j\left({\tau }^{-}\right)\left(1-u_j\right)\text{and}{Q}_j\left({\tau }^{+}\right)=Q_j\left({\tau }^{-}\right)+S_j\left({\tau }^{-}\right)u_j,\\ X_j\left({\tau }^{+}\right)=X_j\left({\tau }^{-}\right),\text{where}{X}_j=E_j,A_j,P_j,M_j,D_j,R,\end{array}\right.$ (2.8)

with τ ^- = lim_t→τ t (for t < τ), and τ ⁺ = lim_τ←t t (for t > τ). For São Paulo State, at t = 0 we use N _0y = 37.8 million and N _0o = 6.8 million ³⁶.

In this paper, we use, at t = 0,

$\left\{\begin{array}{l}\left(S_y\left(0\right)=N_{0y},E_y\left(0\right)=30,A_y\left(0\right)=24,P_y\left(0\right)=6,M_y\left(0\right)=6,D_y\left(0\right)=0\right)\\ \text{and}\\ \left(S_o\left(0\right)=N_{0o},E_o\left(0\right)=20,A_o\left(0\right)=16,P_o\left(0\right)=4,M_o\left(0\right)=3,D_o\left(0\right)=1,\right)\end{array}\right.$ (2.9)

plus Q _y (0) = Q _o (0) = R(0) = 0. The initial conditions’ setup are found in ²². Briefly, at t = 0, one senior individual was diagnosed with severe Covid-19. For other variables, using the ratio 4 : 1 for asymptomatic:pre-symptomatic and 3 : 1 for mild:severe (non-hospitalized:hospitalized) Covid-19, if we assume that there is one person in D _o (the first confirmed case), then there are 3 persons in M ₀. The sum 4 is the number of persons in class P _o , implying that there are 16 in class A _o , and the sum 20 is the number of persons in class E _o . For young population, we used D _y (0) = 0 and arbitrarily M _y (0) = 6. Finally, at the beginning of the epidemic, there are not isolated and immunized persons.

The system of equations (2.3), (2.4), and (2.5) result in an asymptote, which is not an equilibrium value (or steady state of the system) due to the additional equation (2.6). In ^{Appendix B}, the steady-state of the system of equations in terms of fractions corresponding to equations (2.3), (2.4) and (2.5) was analyzed to obtain the basic reproduction number R ₀, which is given by equation (B.5) in ^{Appendix B}. We consider this equation as an approximate value to the timevarying population model, and write the fractions as $s_y^0=N_{0y}/N$ and $s_o^0=N_{0o}/N$, resulting in

$R_0=R_{0y}+R_{0o},\text{with}{R}_{0y}=\left(R_{1y}+R_{2y}\right)\frac{N_{0y}}{N_0}\text{and}{R}_{0o}=\left(R_{1o}+R_{2o}\right)\frac{N_{0o}}{N_0},$ (2.10)

where R _0y and R _0o are the partial reproduction numbers for young and senior subpopulations, and N _0y and N _0o are the initial numbers for young and senior subpopulations, with N ₀ = N _0y + N _0o , and

$\left\{\begin{array}{l}{R}_{1y}=l_y\frac{{\sigma }_y}{{\sigma }_y+\varphi }\frac{{\beta }_{1y}}{{\gamma }_y+\varphi }+\left(1-l_y\right)\frac{{\sigma }_y}{{\sigma }_y+\varphi }\frac{{\beta }_{2y}}{{\gamma }_{1y}+\varphi }\\ R_{1o}=l_o\frac{{\sigma }_o}{{\sigma }_o+\varphi }\frac{{\beta }_{1o}\psi }{{\gamma }_o+\varphi }+\left(1-l_o\right)\frac{{\sigma }_o}{{\sigma }_o+\varphi }\frac{{\beta }_{2o}\psi }{{\gamma }_{1o}+\varphi }\\ R_{2y}=\left[l_y\frac{{\sigma }_y}{{\sigma }_y+\varphi }\left(1-{\chi }_y\right)\frac{{\gamma }_y}{{\gamma }_y+\varphi }+\left(1-l_y\right)\frac{{\sigma }_y}{{\sigma }_y+\varphi }{k}_y\frac{{\gamma }_{1y}}{{\gamma }_{1y}+\varphi }\right]\frac{z_y{\beta }_{3y}}{{\gamma }_{3y}+\varphi }\\ R_{2o}=\left[l_o\frac{{\sigma }_o}{{\sigma }_o+\varphi }\left(1-{\chi }_o\right)\frac{{\gamma }_o}{{\gamma }_o+\varphi }+\left(1-l_o\right)\frac{{\sigma }_o}{{\sigma }_o+\varphi }{k}_o\frac{{\gamma }_{1o}}{{\gamma }_{1o}+\varphi }\right]\frac{z_o{\beta }_{3o}\psi }{{\gamma }_{3o}+\varphi }.\end{array}\right.$ (2.11)

Letting z _y = z _o = 0 (R _2y = R _2o = 0), we retrieve the basic reproduction number obtained in ²³. We chose two key parameters to evaluate the interplay between transmissibility and virulence.

• (1) Transmissibility - The force of infection λ given by equation (2.2) is the per-capita incidence rate, and λS _y (λS _o ) is the total number of the new cases in the young (or senior) subpopulation per unit of time. Notwithstanding, the intensity of infection is proportional to the amount of virus released by infectious individuals and the capacity to infect susceptible individuals. The transmission rates depend on other characteristics not considered here (social network, movement, demography, genetic, nutritional and health conditions, etc.). For this reason, we assume that the basic reproduction number R ₀ calculated by the estimated transmission rates is an approximate measure for the transmissibility of SARSCoV-2. To estimate the transmission rates β _1j , β _2j , and β _3j , j = y, o, the proportion of individuals in isolation u, and the protection factor ε, we used the accumulated number of severe Covid-19 cases Ω = Ω_y + Ω_o given by

$\frac{d}{dt}{\Omega }_y=\left(1-k_y\right){\gamma }_{1y}{P}_y\text{and}\frac{d}{dt}{\Omega }_o=\left(1-k_o\right){\gamma }_{1o}{P}_o,\text{with}{\Omega }_y\left(0\right)=0\text{and}{\Omega }_o\left(0\right)=0,$ (2.12)

• which represents the moving from class P to class D.

• (2) Virulence - By enhancing the virus’ capacity to infect target cells, an increased number of cells are infected. Consequently, the possibility of manifesting mild or severe disease increases, increasing the number of symptomatic Covid-19 cases. For this reason, we equate virulence with the ratio between asymptomatic and symptomatic individuals. It is worth noticing that another effect of increased virulence is the high amount of virus released by infectious individuals in the environment, which is likely to increase transmission ⁸. We consider two broadly separate values to express the proportion of asymptomatic individuals and estimate the transmission rates to calculate the basic reproduction number R ₀.

Besides the severe Covid-19 cases, the number of Covid-19 fatalities are registered. Based on these data, we can estimate the additional mortality rates α _j , j = y, o, using the number of deaths due to severe Covid-19 cases given by Π = Π_y + Π_o , where

$\frac{d}{dt}{\Pi }_y={\alpha }_y{D}_y\text{and}\frac{d}{dt}{\Pi }_o={\alpha }_o{D}_o,\text{with}{\Pi }_y\left(0\right)=0\text{and}{\Pi }_o\left(0\right)=0.$ (2.13)

In the estimation of additional mortality rates, we must bear in mind that the time at which new cases and deaths were registered does not have direct correspondence; instead, they are delayed by ∆ days, that is, Π(t + ∆) = αD(t). From those individuals leaving compartment D(γ _2j , µ, and α _j ), we are counted only the deaths caused by severe Covid-19.

3 RESULTS

We studied the transmissibility and virulence of SARS-CoV-2 in the state of São Paulo, Brazil ²⁴. At the beginning of the epidemic, we had two data sets: severe Covid-19 cases (those in hospitals who were tested and confirmed) and deaths. During the period of the lack of mass testing (PCR and serology), the epidemic curve of SARS-CoV-2 can be characterized by the severe Covid-19 curve D.

To evaluate the interplay between transmissibility and virulence, we assumed that the basic reproduction number R ₀ measures transmissibility, and the proportions of asymptomatic individuals l _y and l _o represent virulence. Letting l = l _y = l _o , we consider two values: l ₁ = 0.8 (lower virulence, with the ratio between asymptomatic and symptomatic individuals being 4 : 1) and l ₂ = 0.2 (higher virulence, with the ratio 1 : 4).

To estimate the transmission rates, the proportion of individuals in isolation, and the protective factor (x stands for one of these parameters), we calculate

$\mathrm{Sum}=\sum _{i=1}^n{\left[\Omega \left(x,t_i\right)-{\Omega }^{ob}\left(t_i\right)\right]}^2,$

where Ω(x, t _i ) is the accumulated severe Covid-19 cases calculated from equation (2.12) and Ω^o (t _i ) is the accumulated severe Covid-19 registered cases at day t _i . We searched for the value of x by minimizing the Sum. To evaluate Ω, the numerical solutions of equations (2.3), (2.4) and (2.5) are obtained using the initial conditions given by equation (2.9). To estimate the additional mortality rates, we substitute Ω(x, t _i ) by Π(α, t _i ) given by equation (2.13) and Ω^o (t _i ) by Π^o (t _i ).

The observed data from February 26 to November 30 ²⁴ is partitioned into two sets - The first set is used to estimate the model parameters (input data set), and the estimated model is then confronted with the second set to assess its prediction ability (test data set) ²⁵. We estimate the model parameters using severe Covid-19 (Ω^ob ) and deaths (Π^ob ) data from February 26 to May 13, 2020 (input data set). Then, we use the estimated parameters to predict the epidemic under interventions (partial quarantine and protective measures) from May 14 to November 30 (test data set) and compared the outcomes with the observed Ω^ob and Π^ob data.

3.1 Input data set - February 26 to May 13, 2020

We estimated the model parameters using the severe Covid-19 data Ω^ob from February 26 to May 13 following the estimation procedure presented in ²². Briefly, the input data set was split into three periods aiming the estimation of the basic reproduction number (without any interventions), the proportion isolated by quarantine, and the effect of protective measures. For low virulence (l ₁ = 0.8), we borrowed the estimated parameters from ²². For high virulence (l ₂ = 0.2), we obtain:

• (1) Transmission rates (data from February 26 to April 3) - Letting u = 0 and ε = 1, the estimated values are β _y = 0.69 and β _o = 0.79 (both in days ⁻¹), resulting in the basic reproduction number R ₀ = 7.02 using equation (2.10). The partial reproduction numbers are R _0y = 6.14 and R _0o = 0.88. This period, at the beginning of the Covid-19 epidemic without any kind interventions, characterizes the natural epidemic, allowing the estimation of the basic reproduction number.

• (2) Proportion of individuals in partial quarantine (data from March 24 to April 13) - Using the previously estimated β _y and β _o , and letting ε = 1, the estimated value is u = 0.48.

• (3) Protective measures (data from April 4 to May 13) - Using the previously estimated β _y , β _o , and u, the estimated value is ε = 0.53.

• Considering the previously estimated parameters β _y , β _o , u, and ε based on Ω^ob , the Covid19 fatality data Π^ob allow us to estimate:

• (4) Additional mortality rates (data from March 16 to May 28) - The estimated values for high virulence are α _y = 0.0018 and α _o = 0.0071 (both in days ⁻¹) letting ∆ = 15 days. The first death due to Covid-19 was on March 16, and we considered 15 more observed data in the input data set due to ∆. The estimates for low virulence are found in ²².

The estimated parameters for l ₁ = 0.8 and l ₂ = 0.2 are summarized in Table 3. The probability of seniors becoming mildly ill k _j is lower than youngs ³⁷. The ratio R _0o /R _oy is 0.195 for low, and 0.143 for high virulent variants.

Table 3:
The estimated parameters ( j = y, o) considering two SARS-CoV-2 virulence levels (rates in days ⁻¹ and proportions are dimensionless).

Parameter	Low (l 1 = 0.8) 22	High (l 2 = 0.2)
ky (k o )	0.92 (0.75)	0.94 (0.78)
ε	0.5	0.53
uy (u o )	0.53 (0.53)	0.48 (0.48)
αy (α o )	0.00185 (0.0071)	0.00185 (0.0071)
β1y (β 1o )	0.78 (0.90)	0.69 (0.79)
β2y (β 2o )	0.78 (0.90)	0.69 (0.79)
β3y (β 3o )	0.78 (0.90)	0.69 (0.79)
R0y (R 0o )[R 0]	7.73 (1.51) [9.24]	6.14 (0.88) [7.02]

Estimations of the epidemic parameters were done considering São Paulo state’s data from February 26 to May 13, 2020, for Ω; and from March 16 to May 28, 2020, for Π. Notice that since April 4, the observed data correspond to the epidemic with interventions (quarantine u and protective measures ε).

Figure 1 shows the estimated curve of accumulated Covid-19 cases Ω and the observed data Ω^ob (a), and the accumulated deaths Π and the observed Π^ob (b) for the high (dashed line) and low (continuous line) virulence. Both low and high virulence levels provided quite similar estimated curves for the input data set.

Figure 1:
The estimated accumulated curve Ω and the observed data Ω^ob (a), and the accumulated deaths Π and the observed data Π^ob (b) for the higher (dashed line) and lower (continuous line) virulence. The natural epidemic corresponding to the low virulence is shown. Estimations were done considering data from São Paulo State.

To assess the contributions of young and senior subpopulations in the overall epidemic, Figure 2 shows the estimated curve of accumulated Covid-19 cases Ω_y , Ω_o , and Ω (a), and the accumulated deaths Π_y , Π_o , Π (b) for the high (dashed line) and low (continuous line) virulence. Notice that Ω = Ω_y + Ω_o and Π = Π_y + Π_o are those shown in Figure 1. For low virulence (high virulence), the numbers of severe Covid-19 cases on May 13 are 35408 (32251) and 22719 (24475) for young and senior subpopulations, with the ratio being 1.56 (1.32). For low virulence (high virulence), the numbers of Covid-19 fatality cases on May 28 are 988 (910) and 3123 (3406) for young and senior subpopulations, with the ratio being 0.32 (0.27).

Figure 2:
The estimated accumulated curves Ω_y , Ω_o , and Ω (a), and the accumulated deaths Π_y , Π_o , and Π (b) for the higher (dashed line) and lower (continuous line) virulence.

As shown in Figures 1 and 2, both low and high virulence levels of SARS-CoV-2 provided a very good estimation of the epidemic under interventions. We must, however, choose one of them to characterize the Covid-19 epidemic in São Paulo, which task is addressed by using the test data set.

3.2 Test data set - May 14 to November 30, 2020

We aim the evaluation of the estimated parameters’ prediction ability for two different virulence levels. The model is solved numerically during the test data sets using the parameters given in Table 2 and those fitted using the input data set given in Table 3. The results are compared with the observed data Ω^ob and Π^ob from São Paulo.

Figure 3 shows the extension of Figure 1 until November 30, 2020: the accumulated Covid-19 cases Ω and the observed data Ω^ob (a) and deaths Π and the observed data Π^ob (b) for high (dashed line) and low (continuous line) virulence sets. The estimated curve of Ω separates from the observed trend Ω^ob around on July 1 for low and on May 7 for high virulent variants; and the estimated curve of Π separates from the observed trend Π^ob around on July 21 for low and on May 16 for high virulent variants. The curves of Ω and Π for low virulent variant accompanied the trend of observed data in more 55 and 66 days than for higher one. It is worth stressing that, after splitting, the estimated curves of Ω and Π situate below the observed Ω^ob and Π^ob data for the lower virulent variant, but above the observed data for the higher virulent one.

Figure 3:
The extended curves until November 30, 2020: the accumulated Covid-19 cases Ω the observed data Ω^ob (a) and deaths Π the observed data Π^ob (b) for the higher (dashed line) and lower (continuous line) virulence. The curves correspond to the prediction of epidemic under quarantine and protective measures.

To assess the contributions of young and senior subpopulations in the overall epidemic, Figure 4 shows the extension of Figure 2 until November 30, 2020: accumulated Covid-19 cases Ω_y , Ω_o , and Ω (a), and the accumulated deaths Π_y , Π_o , Π (b) for high (dashed line) and low (continuous line) virulence sets. For low virulence (high virulence), the numbers of severe Covid-19 cases on May 13 are 242544 (742418) and 142104 (516293) for young and senior subpopulations, with the ratio being 1.56 (1.32). For low virulence (high virulence), the numbers of Covid-19 fatality cases on May 28 are 5254 (16080) and 18418 (66767) for young and senior subpopulations, with the ratio being 0.29 (0.24).

Figure 4:
The extended curves until November 30, 2020: the accumulated Covid-19 cases Ω_y , Ω_o , and Ω (a), and the accumulated deaths Π_y , Π_o , Π (b) for the higher (dashed line) and lower (continuous line) virulence. The curves correspond to the prediction of epidemic under quarantine and protective measures.

Figures 1 to 4 were obtained using R ₀ = 9.24 for low and R ₀ = 7.02 for high virulent variants, which illustrates the low transmissibility of high virulent one. Also, the number of deaths in high virulent variant epidemic is very huge compared to the low virulent one. Let us evaluate the ratio between asymptomatic and symptomatic individuals disregarding the age.

Figure 5 shows the curves of asymptomatic (A) and pre-symptomatic (P) individuals (a), and the estimated curves of mild (M) and severe (D) Covid-19 cases (b) for the high (dashed line) and low (continuous line) virulence. Figure 5 also shows the curve of carriers capable of transmitting the virus (A + P) and the curve of apparent Covid-19 (M + D).

Figure 5:
Estimated curves of asymptomatic (A) and pre-symptomatic (P) individuals (a), and the estimated curves of mild (M) and severe (D) Covid-19 cases (b) for the higher (dashed line) and lower (continuous line) virulence. The sum A + P and M + D are shown.

For the low virulence virus, the peaks of asymptomatic (A), pre-symptomatic (P), mild (M), and severe (D) Covid-19 are 2230000, 195000, 586000, and 66000, and for the high virulence virus, the peaks are 485000, 664000, 1910000, and 189000. (In the low virulence scenario, the peak of A increased by 460%, and the peaks of P, M, and D were reduced by 340%, 326%, and 286% compared to the high virulence values.) Asymptomatic and pre-symptomatic individuals are carriers, which makes the control of infection more difficult. Comparing the lower and higher virulence, the ratio between the peaks of carriers (A + P) is 212%, and the ratio between peaks of apparent Covid-19 (M + D) is 31%.

Additionally, from the numerical simulations (Figure 5), the peaks of SARS-CoV-2 transmitters (A + P + z _y M _y + z _o M _o ) for the low and high virulence scenarios are 2680000 and 1990000, occurring on June 19 and July 21. The ratio of low and high virulence peaks for all infectious individuals is 135%. The peaks of exposed (E), releasing virus (A + P + M + D), and harboring virus (E + A + P + M + D) individuals in the low virulence set were 1430000, 3060000, and 4410000, respectively, and for the high virulent variant, the peaks were 1210000 (85%), 3210000 (105%), and 4360000 (99%), respectively. The percentage shown in parentheses is the ratio high to low virulence. The number of individuals releasing the virus are slightly higher in the high virulence set, but the number of individuals carrying the virus is nearly equal under low or high virulence transmission. Notice that exposed individuals (new infection cases) are higher in low virulence due to a higher R ₀.

4 DISCUSSION

During the first months of the Covid-19 outbreak (input data set), the absence of mass testing of SARS-CoV-2 infection and clinical follow-up of infected individuals precluded the detection of asymptomatic and pre-symptomatic individuals. Considering this initial period, we explained São Paulo’s Covid-19 data using two broadly different values representing the ratio between asymptomatic and symptomatic individuals (4 : 1 and 1 : 4). For each proportion l, we estimated the transmission rates β _y and β _o , and calculated the basic reproduction number R ₀. It is important to remember that the low and high virulent SARS-CoV-2 variants provided quite similar fitted Ω and Π curves (see Figure 1).

Figure 1 showed the parameters estimation. For high virulent variant, the reproduction number was R ₀ = 7.02, with the partial numbers R _0y = 6.14 and R _0o = 0.88 (R _0o /R _0y = 0.143); while for low virulent variant, the reproduction number was R ₀ = 9.24, with the partial numbers R _0y = 7.73 and R _0o = 1.51 (R _0o /R _0y = 0.195). The ratio between senior and young subpopulations is N _0o /N _0y = 0.180. The quotients between low and high virulent variants for R ₀, R _0y , and R _0o are 132%, 126%, and 171%, respectively. In an high virulent variant epidemic, the senior subpopulation is more protected when in interaction with young subpopulation than in a low variant epidemic. The effect of age on the epidemic was shown in Figures 2 and 4.

Remember that Figure 1 was estimated considering quarantine, and Figure 3 was extension of the estimated curve until November 30, 2020. In a direct observation, the predicted curve of severe Covid-19 cases for low variant epidemic matched the observed data 55 more days than the high virulent one (for fatality cases, more 66 days).

Let us discuss the model predictive ability by comparing the extended curves in Figure 3 with the observed data for the low and high virulence scenarios. However, it can be observed in Figure 3 that the low virulence epidemic fits very well into the subset of test data from May 13 to June 30 (Ω) or July 20 (Π), which is the last day of the quarantine. The first observed point detaching from the low virulence curve Π occurs on July 21, while the severe Covid-19 curve Ω splits on July 1, 20 days earlier. Hence, this finding corroborates our assumption ∆ = 15 days in Π(t + ∆) = αD(t) to estimate the fatality rates α _y and α _o . The effects of the relaxation appeared on July 1 and since then the observed number of deaths was systematically above the estimated curve of Π considering the quarantine only (relaxation began by mid-June 2020, and its effect on the epidemic curve D appeared 9 days later). Hence, the subset of test data from July 1 must be described by a model incorporating intermittent pulses of relaxation ²³, which is not dealt with here. In ³⁷, the relaxation incorporated in the model fitted very well the test data set considering the low virulent variant. We concluded that this epidemiological scenario can be explained by the spread of less virulent SARS-CoV-2 variant during the quarantine and relaxation (test data set).

Figure 5 showed enhanced virulence increasing the deaths due to Covid-19 by around 300%; however, the transmissibility measured by R ₀ was reduced to ∼ 76%. The higher virulence increased by around 300% the number of mild and severe Covid-19 cases, potentially more transmissible but in isolation; however, the number of asymptomatic, pre-symptomatic, and a fraction of mild Covid-19 infected individuals transmitting the virus is reduced by 74%. The rate of SARS-CoV-2 transmission in the population is reduced, and the number of deaths goes up as virulence increases. Initially, the model variables (compartments) have similar values, but the peaks of the more virulent SARS-CoV-2 occur around one month later. In particular, the number of deaths in the more virulent epidemic is three times higher, although both variants presented similar fatality rates in the first three months.

However, during the test data set, besides the relaxation, the emerging of variants of SARS-CoV2 may occur. The persistent transmission of SARS-CoV-2 should result in a plethora of new mutations and variants. Some of these variants acquired enhanced ability to infect cells and/or higher levels of replication, which may in turn increase virulence and, consequently, the risk of death. Moreover, as more viruses are released in the environment, the level of transmission also increases. However, individuals with severe symptoms are isolated in hospitals, and mild cases are tested and detected and advised to self-isolate. Conversely, asymptomatic and presymptomatic individuals (carriers) circulate freely and may get in close contact with susceptible individuals, unlike individuals manifesting any suspicious symptoms (not necessarily Covid-19). This fact was exposed by our model, that showed that the presence of high virulence viruses results in more deaths but is less transmissible, i.e., display a lower R ₀. Therefore, while the epidemic curve increases slowly, the curve of deaths increases much more quickly than in the low virulence scenario.

Indeed, the gamma variant was first detected in Manaus by November 2020. At the end of February 2021, São Paulo still maintained part of the population isolated, but the accumulated number of deaths of 68904 on March 24, 2021, surpassed the 57300 predicted by the epidemic with a low virulence SARS-CoV-2 variant. This high number of deaths is likely a consequence of the appearance of more virulent variants (but not so high as l = 0.2 as shown in Figures 1 to 5) during this period ¹⁰. Based on the model’s analysis, we may assume that more virulent viruses were present in the pool of SARS-CoV-2 variants, resulting in more severe Covid-19 cases and deaths. It should be stressed out, however, that an increased number of severe cases and deaths caused by the highly virulent variant does not imply increased transmissibility, according to our model. The omicron variant appeared in South Africa in November 2021. Nevertheless, the omicron variant was showed to be more transmissible, but having less virulence ³⁵.

5 CONCLUSION

In conclusion, throughout the course of an epidemic fast-mutating RNA viruses acquire new mutations resulting in a pool of variants. The new variants may display lower or higher levels of virulence than the original SARS-CoV-2 strain. In a short-run epidemic, the spread of highly virulent variants increase the number of severe cases and deaths, but, paradoxically, reduce viral transmission. However, in the long-run highly-transmissible low virulence variants will prevail. The competition between two strains is left to further work (see ³¹ for the drug-sensitive and resistant M. tuberculosis transmission modeling). Therefore, for the short-run epidemic control, strict isolation of symptomatic (and perhaps also of suspected cases) should be adopted to avoid the spread of highly virulent SARS-CoV-2 variants which will otherwise result in more severe cases and deaths.

REFERENCES

• ¹
  Greenhalgh T, et al. Ten scientific reasons in support of airborne transmission of SARSCoV-2. Lancet 2021; 397 (10285); 1603-1605.
• ²
  Yuan L, et al. Aerodynamic analysis of SARS-CoV-2 in two Wuhan hospitals. Published April 27, 2020. Available at: Available at: https://doi.org/10.1038/s41586-020-22713 (accessed May 10, 2020).
  » https://doi.org/10.1038/s41586-020-22713
• ³
  Zhu N, et al. A novel Coronavirus from patients with pneumonia in China, 2019. N. Eng. J. Med. 2020; 381: 1-7. doi: 10.1056/NEJMoa2001017.
  » https://doi.org/10.1056/NEJMoa2001017
• ⁴
  WHO, Report of the WHO-China Joint Mission on Coronavirus Disease 2019 (Covid-19), 16-24 February 2020 (2020).
• ⁵
  Sender R, et al. The total number and mass of SARS-CoV-2 virions. PNAS 2021; 118 (25) e2024815118. doi: 10.1073/pnas.2024815118.
  » https://doi.org/10.1073/pnas.2024815118
• ⁶
  WHO, Coronavirus Disease (Covid-19): Virus Evolution. 30 December 2020. Available at: Available at: https://www.who.int/news-room/q-a-detail/sars-cov-2-evolution (accessed February 23, 2021).
  » https://www.who.int/news-room/q-a-detail/sars-cov-2-evolution
• ⁷
  CDC, SARS-CoV-2 Variants. Updated Jan. 31, 2021. Available at: Available at: https://www.cdc.gov/coronavirus/2019-ncov/cases-updates/variant-surveillance/variant-info.html (accessed February 23, 2021).
  » https://www.cdc.gov/coronavirus/2019-ncov/cases-updates/variant-surveillance/variant-info.html
• ⁸
  Korber B, et al. Tracking changes in SARS-CoV-2 spike: Evidence that D614G increases infectivity of the Covid-19 virus. Cell 2021; 182: 812-817. doi: https://doi.org/10.1016/j.cell.2020.06.043
  » https://doi.org/10.1016/j.cell.2020.06.043
• ⁹
  Davies NG, et al. Estimated transmissibility and severity of novel SARS-CoV-2 Variant of Concern 202012/01 in England. MedRxiv 2020. doi: 10.1101/2020.12.24.20248822.
  » https://doi.org/10.1101/2020.12.24.20248822
• ¹⁰
  Faria NR, et al. Genomics and epidemiology of the P.1 SARS-CoV-2 lineage in Manaus, Brazil. Science 2021; 372 (6544); 815-821.
• ¹¹
  Martin DP, et al. The emergence and ongoing convergent evolution of the N501Y lineages coincides with a major global shift in the SARS-CoV-2 selective landscape. MedRxiv 2021. doi: 10.1101/2021.02.23.21252268.
  » https://doi.org/10.1101/2021.02.23.21252268
• ¹²
  Naveca FG, et al. COVID-19 in Amazonas, Brazil, was driven by the persistence of endemic lineages and P.1 emergence. Nature Medicine 2021; 27; 1230-1238.
• ¹³
  Zhan Y, Yin H, Yin JY. B.1.617.2 (Delta) Variant of SARS-CoV-2: features, transmission and potential strategies. Int J Biol Sci 2022;18(5); 1844-1851.
• ¹⁴
  Planas D, et al. Reduced sensitivity of SARS-CoV-2 variant Delta to antibody neutralization. Nature 2021; 596; 276-280.
• ¹⁵
  Jansen L, et al. Investigation of a SARS-CoV-2 B.1.1.529 (Omicron) Variant Cluster - Nebraska, November-December 2021. Morbidity and Mortality Weekly Report 2021; 70 (5152);1782-1784.
• ¹⁶
  Kannan SR, et al. Omicron SARS-CoV-2 variant: Unique features and their impact on pre-existing antibodies. Journal of Autoimmunity 2022; 126; 102779. doi: 10.1016/j.jaut.2021.102779.
  » https://doi.org/10.1016/j.jaut.2021.102779
• ¹⁷
  Saxena SK, et al. Characterization of the novel SARS-CoV-2 Omicron (B.1.1.529) variant of concern and its global perspective. J Medical Virology 2021. http://doi.org/10.1002/jmv.27524
  » http://doi.org/10.1002/jmv.27524
• ¹⁸
  Kumar S, Thambiraja TS, Karuppanan K, Subramaniam G. Omicron and Delta variant of SARS-CoV-2: A comparative computational study of spike protein. J Medical Virology 2021. http://doi.org/10.1002/jmv.27526
  » http://doi.org/10.1002/jmv.27526
• ¹⁹
  Yang HM, Lombardi Junior LP, Campos AC. Are the SIR and SEIR models suitable to estimate the basic reproduction number for the Covid-19 epidemic? MedRxiv 2020. doi: https://doi.org/10.1101/2020.10.11.20210831
  » https://doi.org/10.1101/2020.10.11.20210831
• ²⁰
  Arons MM, et al. Presymptomatic SARS-CoV-2 infections and transmission in a skilled nursing facility. The New Engl. Jour. Medicine 2020; April 24, 2020. doi: 10.1056/NEJMoa2008457.
  » https://doi.org/10.1056/NEJMoa2008457
• ²¹
  Anderson RM, May RM. Infectious Diseases of Human. Dynamics and Control. Oxford, New York, Tokyo: Oxford University Press; 1991: 757 p.
• ²²
  Yang HM, Lombardi Junior LP, Castro FFM, Yang AC. Mathematical modeling of the transmission of SARS-CoV-2 - Evaluating the impact of isolation in São Paulo State (Brazil) and lockdown in Spain associated with protective measures on the epidemic of Covid-19. PLoS ONE 2021; 16(6): e0252271. https://doi.org/10.1371/journal.pone.0252271
  » https://doi.org/10.1371/journal.pone.0252271
• ²³
  Yang HM, Lombardi Junior LP, Castro FFM, Yang AC. Mathematical model describing Covid-19 in São Paulo State, Brazil - Evaluating isolation as control mechanism and forecasting epidemiological scenarios of release. Epidemiology and Infection 2020 148: e155. doi: 10.1017/S0950268820001600.
  » https://doi.org/10.1017/S0950268820001600
• ²⁴
  SEADE, SP contra o novo coronavírus - Boletim completo. Available from: Available from: https://www.seade.gov.br/coronavirus/?utm_source=portal&utm_medium= banner&utm_campaign=boletim-completo (accessed August 14, 2020).
  » https://www.seade.gov.br/coronavirus/?utm_source=portal&utm_medium= banner&utm_campaign=boletim-completo
• ²⁵
  Alvarenga MY, Sameshima K, Baccalá LA, Yang HM. Non-linear analysis of the rhythmic activity in rodent brain. Mathematical Biosciences 1999; 157 (1-2): 287-302.
• ²⁶
  Informações Covid-19. http://www.transparencia.am.gov.br/Covid19/monitoramento-Covid-19/ (accessed February 27, 2021).
  » http://www.transparencia.am.gov.br/Covid19/monitoramento-Covid-19/
• ²⁷
  Freitas ARR, et al. The emergence of novel SARS-CoV-2 variant P.1 in Amazonas (Brazil) was temporally associated with a change in the age and sex profile of COVID-19 mortality: A population based ecological study. The Lancet Regional Health Americas 2021; 2; 100064. https://doi.org/10.1016/j.lana.2021.100021
  » https://doi.org/10.1016/j.lana.2021.100021
• ²⁸
  Barbosa GR, et al. Rapid spread and high impact of the Variant of Concern P.1 in the largest city of Brazil. MedRxiv 2021. doi: 10.1101/2021.04.10.21255111.
  » https://doi.org/10.1101/2021.04.10.21255111
• ²⁹
  Francisco Junior RS, et al. Turnover of SARS-CoV-2 Lineages Shaped the Pandemic and Enabled the Emergence of New Variants in the State of Rio de Janeiro, Brazil. Viruses 2021; 13 (10); 2013. doi: 10.3390/v13102013
  » https://doi.org/10.3390/v13102013
• ³⁰
  Martins AF, et al. Detection of SARS-CoV-2 lineage P.1 in patients from a region with exponentially increasing hospitalisation rate, February 2021, Rio Grande do Sul, Southern Brazil. Eurosurveillance 2021; 26 (12). doi: 10.2807/1560-7917.es.2021.26.12.2100276
  » https://doi.org/10.2807/1560-7917.es.2021.26.12.2100276
• ³¹
  Raimundo SM, Yang HM, Venturino E. Theoretical assessment of the relative incidence of sensitive and resistant tuberculosis epidemic in presence of drug treatment. Math. Biosc. Eng. 2014; 11(4): 971-993.
• ³²
  Diekmann O, Heesterbeek JAP, Roberts MG. The construction of next-generation matrices for compartmental epidemic models. J. R. Soc. Interface 2010; 7: 873-885.
• ³³
  Yang HM. The basic reproduction number obtained from Jacobian and next generation matrices - A case study of dengue transmission modelling. BioSystems 2014; 126: 52-75.
• ³⁴
  Yang HM, Greenhalgh D. Proof of conjecture in: The basic reproduction number obtained from Jacobian and next generation matrices - A case study of dengue transmission modelling. Appl. Math. Comput. 2015; 265: 103-107.
• ³⁵
  Bálint G, Vörös-Horváth B and Széchenyi A. (2022). Omicron: increased transmissibility and decreased pathogenicity. Signal Transduction and Targeted Therapy 2022; 7(1): 1-3.
• ³⁶
  Fundação Sistema Estadual (SEADE) database (36 Fundação Sistema Estadual (SEADE) database (https://www.seade.gov.br ) (accessed 9 April 2020).
  » https://www.seade.gov.br
• ³⁷
  Yang HM, Lombardi Junior LP, Castro FFM, Campos AC. Modeling the transmission of the new coronavirus in São Paulo State, Brazil - assessing the epidemiological impacts of isolating young and elder persons. Mathematical Medicine and Biology: A Journal of the IMA 2021; 38 (2): 137-177. doi: https://doi.org/10.1093/imammb/dqaa015
  » https://doi.org/10.1093/imammb/dqaa015
• ³⁸
  Yang HM, Lombardi Junior LP, Castro FFM, Yang AC. Evaluating the impacts of relaxation and mutation in the SARS-CoV-2 on the COVID-19 epidemic based on a mathematical model: a case study of São Paulo State (Brazil). Comput. Appl. Math. 2021; 40(8): 272 doi: 10.1007/s40314-021-01661-w.
  » https://doi.org/10.1007/s40314-021-01661-w

APPENDIX A SARS-COV-2 AND VARIANTS OF CONCERN

Like other RNA-based viruses, SARS-CoV-2 mutation rate is relatively high (about 3 × 10⁻⁶ mutations per nucleotide per replication cycle ⁵) but not as high as in other RNA viruses such as HIV or influenza ⁶. “Viruses of a particular lineage that acquire one or more non-synonymous new mutations are referred to as variant of the original virus” ⁶. Some of these new variants are of public health concern (known as ‘variants of concern’ or VOC) either because they have an increased capacity of causing severe illness or because they are more efficiently transmitted than other variants. For instance, the SARS-CoV-2 VOCs alpha (B.1.1.7), gamma (P.1 or B.1.1.28), delta (B.1.617) and omicron, that emerged in the United Kingdom, Brazil, India and South Africa, respectively (see ⁷ and references therein). Some of these VOCs display high lethality, high levels of transmission and are able to escape the host immune system. Indeed, in many cases more virulent variants of SARS-CoV-2 are also more transmissible ⁸.

Given the relatively high SARS-CoV-2 mutation rate the emergence of new variants is expected in the course of a pandemic. However, only a few of them are expected to confer significant selective advantage to the virus. Some of these new variants are of public health concern either because they have an increased capacity of causing severe illness or because they are more efficiently transmitted than other variants.

The original SARS-CoV-2 variant was first detected in Wuhan, China by December 2019. The first VOC alpha, appeared in the United Kingdom less than one year later, and quickly spread throughout the country. The alpha variant carries 23 mutations, including 14 non-synonymous substitutions, 3 deletions, and 6 synonymous mutations, three of them targeting the spike protein. These mutations are related to a 50% increase in transmissibility (according to the CDC), making the virus highly contagious and difficult to contain ⁹. The gamma VOC was first detected in Manaus by November 2020 ¹⁰. From there, this variant disseminated throughout the country creating havoc and a pattern of high lethality. The gamma variant harbours a total of 25 mutations (synonymous and non-synonymous) compared with its most closely related ancestral, including 10 amino acid substitutions in the spike protein, some of them associated with a higher degree of transmissibility ^{10) (11}. In addition, the viral load in patients infected with the gamma variant was ∼10-fold higher than in those with other SARS-CoV-2 variants ¹². In Brazil alone, the gamma variant was responsible for the death of about 400,000 people, in a 6-month period. The delta variant carries 4 mutations of interest in the spike protein. It was first detected in India in late 2020. This VOC also displayed increased virulence and transmissibility ¹³, but its mutations have also enabled the virus to evade immune response, and for that reason most vaccines were less effective against this variant ¹⁴. Finally, the omicron variant, that appeared in South Africa in November 2021, is the most phylogenetically divergent VOC to date. It is associated with enhanced transmissibility, reduced vaccine-induced immunity and increased risk of re-infection. Omicron harbors more than 50 mutations, 34 of them in the spike protein: 30 amino acid substitutions, 3 deletions, and 1 small insertion ¹⁵. Half of these amino acid substitutions occur in the receptor binding domain, which might explain the high transmissibility of the this variant ^{16), (17), (18}.

APPENDIX B THE STEADY-STATE ANALYSIS OF THE MODIFIED SQEAPMDR MODEL

The system of equations (2.3), (2.4) and (2.5), time-varying population, approach an asymptote, which is not an equilibrium value. However, in term of fractions this system attains steady-state. Defining the fraction x _j = X _j /N, for j = y, o, with X _j = {S _j , Q _j , E _j , A _j , P _j , M _j , D _j , R} , we have, using equation (2.6),

$\frac{d}{dt}\frac{X_j}{N}=\frac{1}{N}\frac{d}{dt}{X}_j-\frac{X_j}{N}\frac{1}{N}\frac{d}{dt}N=\frac{1}{N}\frac{d}{dt}{X}_j-x_j\left(\varphi -\mu -{\alpha }_y{d}_y-{\alpha }_o{d}_o\right),$

and the system of equations (2.3), (2.4) and (2.5) in terms of fractions become, for susceptible individuals,

$\left\{\begin{array}{l}\frac{d}{dt}{s}_y=\varphi -\left(\phi +\varphi \right)s_y-\lambda s_y+s_y\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \frac{d}{dt}{s}_o=\phi s_y-\varphi s_o-\lambda \psi s_o+s_o\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right),\end{array}\right.$ (B.1)

for isolated and infected individuals,

$\left\{\begin{array}{l}\frac{d}{dt}{q}_j=-\varphi q_j+q_j\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \frac{d}{dt}{e}_j=\lambda \left({\delta }_{jy}+\psi {\delta }_{jo}\right)s_j-\left({\sigma }_j+\varphi \right)e_j+e_j\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \frac{d}{dt}{a}_j=l_j{\sigma }_j{e}_j-\left({\gamma }_j+\varphi \right)a_j+a_j\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \frac{d}{dt}{p}_j=\left(1-l_j\right){\sigma }_j{e}_j-\left({\gamma }_{1j}+\varphi \right)p_j+p_j\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \frac{d}{dt}{m}_j=\left(1-{\chi }_j\right){\gamma }_j{a}_j+k_j{\gamma }_{1j}{p}_j-\left({\gamma }_{3j}+\varphi \right)m_j+m_j\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \frac{d}{dt}{d}_j=\left(1-k_j\right){\gamma }_{1j}{p}_j-\left({\gamma }_{2j}+\varphi +{\alpha }_j\right)d_j+d_j\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right),\end{array}\right.$ (B.2)

and for recovered individuals,

$\frac{d}{dt}r={\chi }_y{\gamma }_y{a}_y+{\gamma }_{3y}{m}_y+{\gamma }_{2y}{d}_y+{\chi }_o{\gamma }_o{a}_o+{\gamma }_{3o}{m}_o+{\gamma }_{2o}{d}_o-\varphi r+r\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right),$ (B.3)

where λ is the force of infection given by equation (2.2) re-written as

$\lambda ={\beta }_{1y}{a}_y+{\beta }_{2y}{p}_y+{\beta }_{3y}{z}_y{m}_y+{\beta }_{1o}{a}_o+{\beta }_{2o}{p}_o+{\beta }_{3o}{z}_o{m}_o,$

letting ε = 1 (absence of protective measures), with

$\sum _{j=y,o}\left(s_j+q_j+e_j+a_j+p_j+m_j+d_j\right)+r=1.$

This new system of equations has steady-states (the sum of derivatives of all classes is zero).

The trivial (disease-free) equilibrium point P ⁰ of the new system of equations (B.1), (B.2) and (B.3) is given by

$P^0=\left(s_j^0,q_j^0=0,e_j^0=0,a_j^0=0,p_j^0=0,m_j^0=0,d_j^0=0,r^0=0\right),$

for j = y and o, where

$\left\{\begin{array}{l}{s}_y^0=\frac{\varphi }{\varphi +\phi }\\ s_o^0=\frac{\phi }{\varphi +\phi },\end{array}\right.$ (B.4)

with $s_y^0+s_o^0=1$.

Let us assess the stability of P ⁰ by applying the next generation matrix theory considering the vector of variables x = (e _y , a _y , p _y , m _y , e _o , a _o , p _o , m _o ) ³². We apply method proposed in ³³ and proved in ³⁴. To obtain the basic reproduction number, diagonal matrix V is considered. Hence, the vectors f and v are

$f^T=\left(\begin{array}{c}\lambda s_y+e_y\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ l_y{\sigma }_y{e}_y+a_y\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \left(1-l_y\right){\sigma }_y{e}_y+p_y\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \left(1-{\chi }_y\right){\gamma }_y{a}_y+k_y{\gamma }_{1y}{p}_y+m_y\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \lambda \psi s_o+e_o\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ p_o{\sigma }_o{e}_o+a_o\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \left(1-p_o\right){\sigma }_o{e}_o+p_o\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\\ \left(1-{\chi }_o\right){\gamma }_o{a}_o+k_o{\gamma }_{1o}{p}_o+m_o\left({\alpha }_y{d}_y+{\alpha }_o{d}_o\right)\end{array}\right)and{v}^T=\left(\begin{array}{c}\left({\sigma }_y+\varphi \right)e_y\\ \left({\gamma }_y+\varphi \right)a_y\\ \left({\gamma }_{1y}+\varphi \right)p_y\\ \left({\gamma }_{3y}+\varphi \right)m_y\\ \left({\sigma }_o+\varphi \right)e_o\\ \left({\gamma }_o+\varphi \right)a_o\\ \left({\gamma }_{1o}+\varphi \right)p_o\\ \left({\gamma }_{3o}+\varphi \right)m_o\end{array}\right)$

where the superscript T stands for the transposition of a matrix, from which we obtain the matrices F and V (see ³²) evaluated at the trivial equilibrium P ⁰, which were omitted. The next generation matrix FV ⁻¹ is

$F{V}^{-1}=\left[\begin{array}{cccccccc}0& \frac{{\beta }_{1y}{s}_y^0}{{\gamma }_y+\varphi }& \frac{{\beta }_{2y}{s}_y^0}{{\gamma }_{1y}+\varphi }& \frac{{\beta }_{3y}{z}_y{s}_y^0}{{\gamma }_{3y}+\varphi }& 0& \frac{{\beta }_{1o}{s}_y^0}{{\gamma }_o+\varphi }& \frac{{\beta }_{2o}{s}_y^0}{{\gamma }_{1o}+\varphi }& \frac{{\beta }_{3o}{z}_o{s}_y^0}{{\gamma }_{3o}+\varphi }\\ \frac{l_y{\sigma }_y}{{\sigma }_y+\varphi }& 0& 0& 0& 0& 0& 0& 0\\ \frac{\left(1-l_y\right){\sigma }_y}{{\sigma }_y+\varphi }& 0& 0& 0& 0& 0& 0& 0\\ 0& \frac{\left(1-{\chi }_y\right){\gamma }_y}{{\gamma }_y+\varphi }& \frac{k_y{\gamma }_{1y}}{{\gamma }_{1y}+\varphi }& 0& 0& 0& 0& 0\\ 0& \frac{{\beta }_{1y}\psi s_o^0}{{\gamma }_y+\varphi }& \frac{{\beta }_{2y}\psi s_o^0}{{\gamma }_{1y}+\varphi }& \frac{{\beta }_{3y}{z}_y\psi s_o^0}{{\gamma }_{3y}+\varphi }& 0& \frac{{\beta }_{1o}\psi s_o^0}{{\gamma }_o+\varphi }& \frac{{\beta }_{2o}\psi s_o^0}{{\gamma }_{1o}+\varphi }& \frac{{\beta }_{3o}{z}_o\psi s_o^0}{{\gamma }_{3o}+\varphi }\\ 0& 0& 0& 0& \frac{p_o{\sigma }_o}{{\sigma }_o+\varphi }& 0& 0& 0\\ 0& 0& 0& 0& \frac{\left(1-p_o\right){\sigma }_o}{{\sigma }_o+\varphi }& 0& 0& 0\\ 0& 0& 0& 0& 0& \frac{\left(1-{\chi }_o\right){\gamma }_o}{{\gamma }_o+\varphi }& \frac{k_o{\gamma }_{1o}}{{\gamma }_{1o}+\varphi }& 0\end{array}\right]$

and the characteristic equation corresponding to FV ⁻¹ is

${\varkappa }^5\left[{\varkappa }^3-\left(R_{1y}{s}_y^0+R_{1o}{s}_o^0\right)\varkappa -\left(R_{2y}{s}_y^0+R_{2o}{s}_o^0\right)\right]=0,$

with the basic reproduction number R ₀ being given by

$R_0=\left(R_{1y}+R_{2y}\right)s_y^0+\left(R_{1o}+R_{2o}\right)s_o^0,$ (B.5)

where the initial fractions $s_y^0$ and $s_o^0$ are given by equation (B.4), and the partial basic reproduction numbers R _1y , R _2y , R _1o , and R _2o are given by equation (2.11) in the main text. The spectral radius ρ FV ⁻¹ is the biggest solution of a third-degree polynomial, not easy to evaluate. The procedure proposed in ³³ allows us to obtain the threshold R ₀ as the sum of coefficients of the characteristic equation. Hence, the trivial equilibrium point P ⁰ is locally asymptotically stable if R ₀ < 1.

Publication Dates

• Publication in this collection
  15 Dec 2023
• Date of issue
  Oct-Dec 2023

History

• Received
  04 May 2022
• Accepted
  22 Mar 2023