# Second quantization approach to COVID-19 epidemic

Leonardo Mondaini1, Bernhard Meirose2 and Felipe Mondaini3 1 Department of
Physics, Federal University of the State of Rio de Janeiro, Rio de Janeiro,
Brazil  2 Fysiska institutionen, Lunds universitet, Lund, Sweden  3 Centro
Federal de Educação Tecnológica Celso Suckow da Fonseca, CEFET-RJ/ Campus
Petrópolis, 25.620-003, Petrópolis, Brazil<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

In this article, a stochastic SIR-type model for COVID-19 epidemic is built
using the standard field theoretical language based on creation and
annihilation operators. From the model, we derive the time evolution of the
mean number of infectious (active cases) and deceased individuals. In order to
capture the effects of lockdown and social distancing, we use a time-dependent
infection rate. The results are in good agreement with the data for three
different waves of epidemic activity in South Korea.

*All authors contributed equally to this work.

,,

  * December 2021

Keywords: second quantization, infectious diseases, stochastic models

##  1 Introduction

The year 2020 will be remembered in history as the year in which humanity
suffered from the rapid spread of coronavirus disease 2019 (COVID-19), caused
by the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) and
responsible for the worst pandemic in the last hundred years. Governments
around the world were strongly impacted by the pandemic and, in addition to
efforts of a social and economic nature, it became evident the need to bring
together members of the international scientific community to analyze data on
the evolution of the disease and propose concrete actions for the control and
mitigation of its serious effects [law2020tracking]. In this sense, even
researchers from areas far from epidemiology were encouraged to collaborate
with their experience and specific knowledge in the study of possible models
for describing the evolution of the pandemic in different countries of the
world [thomas2020primer, zhuang2020preliminary, zhao2020preliminary,
neves2020predicting, medo2020contact], often serving as guides for the actions
of local governments. Theoretical particle physicists are used to dealing with
problems and proposing models to describe complex phenomena involving the
fundamental constituents of matter. In this work, we show how the standard
field theoretical language based on creation and annihilation operators
(building blocks of the second quantization method [Greiner]), usually used in
the construction of these elementary particle models, may also be adapted to
the construction of an epidemiological model for COVID-19 which can be
compared to the data available for several countries111M. Doi[doi1976] was the
first to point that, despite what one could infer from the name, the second
_quantization_ approach is not limited to _quantum_ systems. Indeed, this
approach may be introduced to the description of certain _classical_ many-
particle systems, such as the one we are dealing with in this work.. We
believe this can serve as an stimulus so that other researchers working on
these subjects can also contribute to this great scientific network. Indeed,
in a recent letter [MondainiBMST2015], we have already shown how this field
theoretical language may be used to derive closed master equations [Gardiner]
describing the population dynamics of multivariate stochastic epidemic models.
Since then we applied a simplified version of the theoretical framework
established in Ref. [MondainiBMST2015] to the description of tumour growth
[Mondaini_2017]. Nevertheless, no successful confrontation with real data was
possible on both works.

In what follows we adapt the methodology established in Ref.
[MondainiBMST2015] so as to derive an SIR-type stochastic model for COVID-19
epidemic, from which we obtain the time evolution of the mean number of
infectious (active cases) and deceased individuals in South Korea. Our choice
for South Korea was motivated by its mass testing and careful monitoring of
infected patients since appearance of their first symptoms
[dongarwar2021implementation]. This leads us to the first successful
description of real data and allows us to finally confirm the robustness and
versatility of the approach. Indeed, our main motivation comes from the fact
that, as remarked in Ref. [dodd2009many], for the types of models studied in
population biology and epidemiology, this field theoretical description is
notationally neater and more manageable than standard methods, often replacing
sets of equations by single equations with the same content. As an example of
this, we may stress that a single Hamiltonian function sums up the system
dynamics compactly and may be easily written down from a verbal description of
the transitions presented in our model. Besides the simplicity and notational
appeal, we may stress as another advantage of our SIR-type model, the use of a
time-dependent infection rate in order to capture the effects of a lockdown,
since the standard SIR model does not take these into account. This approach
is similar to the one adopted in Ref. [palladino2020modelling], although the
infection rate function we use differs from Ref. [palladino2020modelling] at
the initial stage of epidemic, i.e., before lockdown and social distance
measures are applied.

The rest of this work is organized as follows. In Section 2, we introduce the
basic aspects of our model which allow us to obtain differential equations
describing the time evolution of the mean number of infectious and deceased
individuals, while the time-dependent infection rate used for data comparison
is presented in Section 3. The differential equations are used to fit
available data (South Korea) for the time evolution of infectious (active
cases) individuals in Section 4 and deceased individuals in Section 5. Further
data comparisons for two other distinct waves of epidemic activity in South
Korea are presented in Section 6, whereas a discussion of the results is
presented in Section 7. Our conclusions are presented in Section 8.

##  2 Building the model

We start by considering interacting populations, whose total sizes are allowed
to change, composed of four types of individuals: susceptible, infectious,
recovered and deceased. Let us introduce $\mathcal{S}(t)$, $\mathcal{I}(t)$,
$\mathcal{R}(t)$ and $\mathcal{D}(t)$ as random variables which represent,
respectively, the number of susceptible, infectious, recovered and deceased
individuals at a given time instant $t$.

We then consider a multivariate process
$\\{(\mathcal{S}(t);\,\mathcal{I}(t);\,\mathcal{R}(t);\,\mathcal{D}(t))\\}_{t=0}^{\infty}$
with a joint probability function given by

$p_{(n_{S},n_{I},n_{R},n_{D})}(t)={\rm{Prob}}\\{\mathcal{S}(t)=n_{S};\,\mathcal{I}(t)=n_{I};\,\mathcal{R}(t)=n_{R};\,\mathcal{D}(t)=n_{D}\\}~{}.$
(1)

Our aim is to compute time-dependent expectation values of the observables
$\mathcal{S}(t)$, $\mathcal{I}(t)$, $\mathcal{R}(t)$ and $\mathcal{D}(t)$,
which may be defined in terms of the configuration probability according to

$\displaystyle\left<\mathcal{S}(t)\right>=\sum_{n_{S},n_{I},n_{R},n_{D}}n_{S}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)~{};$
$\displaystyle\left<\mathcal{I}(t)\right>=\sum_{n_{S},n_{I},n_{R},n_{D}}n_{I}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)~{};$
$\displaystyle\left<\mathcal{R}(t)\right>=\sum_{n_{S},n_{I},n_{R},n_{D}}n_{R}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)~{};$
$\displaystyle\left<\mathcal{D}(t)\right>=\sum_{n_{S},n_{I},n_{R},n_{D}}n_{D}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)~{}.$
(2)

The probabilistic state of the system may be represented by the vector

$\left|\nu_{S},\nu_{I},\nu_{R},\nu_{D}\right>_{t}=\sum_{n_{S},n_{I},n_{R},n_{D}}p_{(n_{S},n_{I},n_{R},n_{D})}(t)\left|n_{S},n_{I},n_{R},n_{D}\right>~{},$
(3)

with the normalization condition
$\sum_{n_{S},n_{I},n_{R},n_{D}}p_{(n_{S},n_{I},n_{R},n_{D})}(t)=1$.

The configurations are fully given in terms of occupation numbers
($n_{S},n_{I},n_{R},n_{D}$), thus a representation in terms of second-
quantized bosonic operators [Cardy_1998] seems natural. We introduce,
therefore, creation and annihilation operators for the susceptible,
infectious, recovered and deceased individuals, respectively, $a_{S}$ and
$a^{\dagger}_{S}$, $a_{I}$ and $a^{\dagger}_{I}$, $a_{R}$ and
$a^{\dagger}_{R}$, $a_{D}$ and $a^{\dagger}_{D}$.

These operators must satisfy the following commutation relations

$\displaystyle\left[a_{i},\,a^{\dagger}_{j}\right]=\delta_{ij}\,~{};$
$\displaystyle\left[a_{i},\,a_{j}\right]=\left[a^{\dagger}_{i},\,a^{\dagger}_{j}\right]=0\,~{},$
(4)

where $i,j=\\{S,I,R,D\\}$ and $\delta_{ij}$ is the Kronecker delta
($\delta_{ij}=1$ if $i=j$ and $\delta_{ij}=0$ if $i\neq j$). As usual in the
second quantization framework, we say that $a^{\dagger}_{S}$,
$a^{\dagger}_{I}$, $a^{\dagger}_{R}$ and $a^{\dagger}_{D}$ “create” ,
respectively, susceptible, infectious, recovered and deceased individuals when
applied over the reference (vacuum) state $\left|0,0,0,0\right>$. This allows
us to build our space from basis vectors of the form
$\left|n_{S},n_{I},n_{R},n_{D}\right>=\left(a^{\dagger}_{S}\right)^{n_{S}}\left(a^{\dagger}_{I}\right)^{n_{I}}\left(a^{\dagger}_{R}\right)^{n_{R}}\left(a^{\dagger}_{D}\right)^{n_{D}}\left|0,0,0,0\right>$.

This vacuum state has the following properties: $a_{S}\left|0,0,0,0\right>$
$=a_{I}\left|0,0,0,0\right>$ $=a_{R}\left|0,0,0,0\right>$
$=a_{D}\left|0,0,0,0\right>=0$ (from which “annihilation” operators) and
$\left<0,0,0,0|0,0,0,0\right>=1$ (inner product). Following the above
definitions, we also have

$\displaystyle a^{\dagger}_{S}\left|n_{S},n_{I},n_{R},n_{D}\right>=$
$\displaystyle\left|n_{S}+1,n_{I},n_{R},n_{D}\right>~{};$ (5) $\displaystyle
a^{\dagger}_{I}\,\left|n_{S},n_{I},n_{R},n_{D}\right>=$
$\displaystyle\left|n_{S},n_{I}+1,n_{R},n_{D}\right>~{};$ $\displaystyle
a^{\dagger}_{R}\,\left|n_{S},n_{I},n_{R},n_{D}\right>=$
$\displaystyle\left|n_{S},n_{I},n_{R}+1,n_{D}\right>~{};$ $\displaystyle
a^{\dagger}_{D}\,\left|n_{S},n_{I},n_{R},n_{D}\right>=$
$\displaystyle\left|n_{S},n_{I},n_{R},n_{D}+1\right>~{};$ $\displaystyle
a_{S}\left|n_{S},n_{I},n_{R},n_{D}\right>=$ $\displaystyle
n_{S}\left|n_{S}-1,n_{I},n_{R},n_{D}\right>~{};$ $\displaystyle
a_{I}\left|n_{S},n_{I},n_{R},n_{D}\right>=$ $\displaystyle
n_{I}\left|n_{S},n_{I}-1,n_{R},n_{D}\right>~{};$ $\displaystyle
a_{R}\left|n_{S},n_{I},n_{R},n_{D}\right>=$ $\displaystyle
n_{R}\left|n_{S},n_{I},n_{R}-1,n_{D}\right>~{};$ $\displaystyle
a_{D}\left|n_{S},n_{I},n_{R},n_{D}\right>=$ $\displaystyle
n_{D}\left|n_{S},n_{I},n_{R},n_{D}-1\right>~{}.$

Since

$\displaystyle a^{\dagger}_{S}a_{S}\,\left|n_{S},n_{I},n_{R},n_{D}\right>$
$\displaystyle=$ $\displaystyle n_{S}\left|n_{S},n_{I},n_{R},n_{D}\right>~{};$
$\displaystyle a^{\dagger}_{I}a_{I}\,\left|n_{S},n_{I},n_{R},n_{D}\right>$
$\displaystyle=$ $\displaystyle n_{I}\left|n_{S},n_{I},n_{R},n_{D}\right>~{};$
$\displaystyle a^{\dagger}_{R}a_{R}\,\left|n_{S},n_{I},n_{R},n_{D}\right>$
$\displaystyle=$ $\displaystyle n_{R}\left|n_{S},n_{I},n_{R},n_{D}\right>~{};$
$\displaystyle a^{\dagger}_{D}a_{D}\,\left|n_{S},n_{I},n_{R},n_{D}\right>$
$\displaystyle=$ $\displaystyle n_{D}\left|n_{S},n_{I},n_{R},n_{D}\right>~{},$
(6)

we may conclude that the operators $n_{S}=a^{\dagger}_{S}a_{S}$,
$n_{I}=a^{\dagger}_{I}a_{I}$, $n_{R}=a^{\dagger}_{R}a_{R}$ and
$n_{D}=a^{\dagger}_{D}a_{D}$ just count the number of individuals in a
definite state, reason why they are called number operators. Note that the
vector state of our system may be then rewritten in terms of creation and
annihilation operators as

$\displaystyle\left|\nu_{S},\nu_{I},\nu_{R},\nu_{D}\right>_{t}$
$\displaystyle=$
$\displaystyle\sum_{n_{S},n_{I},n_{R},n_{D}}p_{(n_{S},n_{I},n_{R},n_{D})}(t)$
(7)
$\displaystyle\times\left(a^{\dagger}_{S}\right)^{n_{S}}\left(a^{\dagger}_{I}\right)^{n_{I}}\left(a^{\dagger}_{R}\right)^{n_{R}}\left(a^{\dagger}_{D}\right)^{n_{D}}\left|0,0,0,0\right>~{}.$

The time evolution of our system will then be generated by a linear operator
$\mathcal{H}$ (Hamiltonian) which may be constructed directly from the
transition rates present in our model according to Table 1 (cf. Ref.
[MondainiBMST2015], Table 1). In Figure 1 we present a compartmental diagram
for the structure of our SIR-type model.

The terms presented in Table 1 can be summed up and rearranged, so that the
Hamiltonian can be written down as

$\displaystyle\mathcal{H}=\beta n_{S}+\gamma n_{I}-\left[\beta
a^{\dagger}_{I}a_{S}+(1-\sigma)\gamma a^{\dagger}_{R}a_{I}+\sigma\gamma
a^{\dagger}_{D}a_{I}\right]~{}.$ (8)

Table 1: Transition rates presented in our model and corresponding terms in
the Hamiltonian $\mathcal{H}$.
Transition | Description | Contribution to $\mathcal{H}$
---|---|---
$\mathcal{S}\xrightarrow{\beta}\mathcal{I}$ | infection (rate $\beta$) | $\beta(a^{\dagger}_{S}a_{S}-a^{\dagger}_{I}a_{S})$
$\mathcal{I}\xrightarrow{(1-\sigma)\gamma}\mathcal{R}$ | change infectious $\rightarrow$ recovered (rate $(1-\sigma)\gamma$) | $(1-\sigma)\gamma(a^{\dagger}_{I}a_{I}-a^{\dagger}_{R}a_{I})$
$\mathcal{I}\xrightarrow{\sigma\gamma}\mathcal{D}$ | change infectious $\rightarrow$ deceased (rate $\sigma\gamma$) | $\sigma\gamma(a^{\dagger}_{I}a_{I}-a^{\dagger}_{D}a_{I})$
Figure 1: The structure of our SIR-type model. The rate $\beta$ is related to
the infection rate of the virus and $\gamma$ is related to the time required
to infected people to become recovered or deceased. Both parameters have the
dimension of time-1. Indeed, $\gamma=1/D$, where $D$ is the average time
period during which an individual remains infectious, before one gets
recovered or deceased. The parameter $\sigma$ (the case fatality rate) is
dimensionless and represents the percentage of infected people who become
deceased.

The notational advantage of this field theoretical description is made clear
at this point if we observe that, from the above definitions, the equation
which represents the flux of probability between states at rates defined by
our model (the so-called master equation or forward Kolmogorov differential
equation [allen2010introduction]) takes the very compact form of an imaginary-
time Schrödinger-type linear equation, namely

$\displaystyle\frac{d}{dt}\left|\nu_{S},\nu_{I},\nu_{R},\nu_{D}\right>_{t}=-\mathcal{H}\left|\nu_{S},\nu_{I},\nu_{R},\nu_{D}\right>_{t}~{}.$
(9)

We can get, after some algebra, a more common representation for the master
equation by substituting the expressions for the Hamiltonian, (8), and the
vector state, (3), into (9)

$\displaystyle\frac{d}{dt}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)$ $\displaystyle=$
$\displaystyle-(\beta n_{s}+\gamma n_{I})\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)$
(10) $\displaystyle+\beta(n_{S}+1)\,p_{(n_{S}+1,n_{I}-1,n_{R},n_{D})}(t)$
$\displaystyle+(1-\sigma)\gamma(n_{I}+1)\,p_{(n_{S},n_{I}+1,n_{R}-1,n_{D})}(t)$
$\displaystyle+\sigma\gamma(n_{I}+1)\,p_{(n_{S},n_{I}+1,n_{R},n_{D}-1)}(t)\,~{}.$

In order to compute the time-dependent expectation values of the observables
$\mathcal{S}(t)$, $\mathcal{I}(t)$, $\mathcal{R}(t)$ and $\mathcal{D}(t)$
through the above master equation, we follow the well-established methodology
presented in Ref. [allen2010introduction] and introduce the following moment
generating function (mgf)

$\displaystyle M(\phi_{S},\phi_{I},\phi_{R},\phi_{D};t)$ $\displaystyle=$
$\displaystyle\left<e^{\phi_{S}\mathcal{S}(t)}e^{\phi_{I}\mathcal{I}(t)}e^{\phi_{R}\mathcal{R}(t)}e^{\phi_{D}\mathcal{D}(t)}\right>$
(11) $\displaystyle=$
$\displaystyle\sum_{n_{S},n_{I},n_{R},n_{D}}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)e^{n_{S}\phi_{S}+n_{I}\phi_{I}+n_{R}\phi_{R}+n_{D}\phi_{D}}\,~{}.$

Note that from the above equation we have

$\displaystyle\left[\frac{\partial
M}{\partial\phi_{S}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\sum_{n_{S},n_{I},n_{R},n_{D}}n_{S}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)=\left<\mathcal{S}(t)\right>~{};$
$\displaystyle\left[\frac{\partial
M}{\partial\phi_{I}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\sum_{n_{S},n_{I},n_{R},n_{D}}n_{I}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)=\left<\mathcal{I}(t)\right>~{};$
$\displaystyle\left[\frac{\partial
M}{\partial\phi_{R}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\sum_{n_{S},n_{I},n_{R},n_{D}}n_{R}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)=\left<\mathcal{R}(t)\right>~{};$
$\displaystyle\left[\frac{\partial
M}{\partial\phi_{D}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\sum_{n_{S},n_{I},n_{R},n_{D}}n_{D}\,p_{(n_{S},n_{I},n_{R},n_{D})}(t)=\left<\mathcal{D}(t)\right>~{}.$
(12)

Multiplying (10) by
$\exp(n_{S}\phi_{S}+n_{I}\phi_{I}+n_{R}\phi_{R}+n_{D}\phi_{D})$ and summing on
($n_{S},n_{I},n_{R},n_{D}$), leads, after some algebra, to

$\displaystyle\frac{\partial M}{\partial
t}=\left[\beta\left(e^{-\phi_{S}+\phi_{I}}-1\right)\right]\frac{\partial
M}{\partial\phi_{S}}+\left[\gamma\left(e^{-\phi_{I}+\phi_{R}}-1\right)+\sigma\gamma\left(e^{-\phi_{I}+\phi_{D}}-e^{-\phi_{I}+\phi_{R}}\right)\right]\frac{\partial
M}{\partial\phi_{I}}~{}.$ (13)

Hence, by differentiating the above equation with respect to
$\phi_{S},\phi_{I},\phi_{R},\phi_{D}$ and evaluating the result at
$\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0$ we get the following differential
equations for $\left<\mathcal{S}(t)\right>$, $\left<\mathcal{I}(t)\right>$,
$\left<\mathcal{R}(t)\right>$ and $\left<\mathcal{D}(t)\right>$:

$\displaystyle\left[\frac{\partial^{2}M}{\partial
t\,\partial\phi_{S}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\frac{d}{dt}\left<\mathcal{S}(t)\right>=-\beta\left<\mathcal{S}(t)\right>~{};$
$\displaystyle\left[\frac{\partial^{2}M}{\partial
t\,\partial\phi_{I}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\frac{d}{dt}\left<\mathcal{I}(t)\right>=\beta\left<\mathcal{S}(t)\right>-\gamma\left<\mathcal{I}(t)\right>~{};$
$\displaystyle\left[\frac{\partial^{2}M}{\partial
t\,\partial\phi_{R}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\frac{d}{dt}\left<\mathcal{R}(t)\right>=(1-\sigma)\gamma\left<\mathcal{I}(t)\right>~{};$
$\displaystyle\left[\frac{\partial^{2}M}{\partial
t\,\partial\phi_{D}}\right]_{\phi_{S},\phi_{I},\phi_{R},\phi_{D}=0}$
$\displaystyle=$
$\displaystyle\frac{d}{dt}\left<\mathcal{D}(t)\right>=\sigma\gamma\left<\mathcal{I}(t)\right>~{}.$
(14)

Last but not least, we must observe at this point that the above differential
equations would be the same if we have considered time-dependent rates
($\beta(t)$,$\gamma(t)$) from the very beginning. This turns out to be an
important aspect of our model as we shall see in the next section.

##  3 Comparison to COVID-19 Epidemic Data - South Korea

We now compare our model against public epidemic data from South Korea, which
was taken from the “Worldometer” reference website [worldometers]. The modeled
data for the number of deaths and number of active cases was compared to the
data reported in South Korea during 96 days of the epidemic. In order to take
into account the slowdown of the spread due the effects of quarantine and
social distancing we use a similar approach to Ref. [palladino2020modelling],
where the infection rate $\beta$ is a function of time
($\beta\equiv\beta(t)$). Specifically, $\beta(t)=\beta_{0}$ before a time
$t_{th}$, while it exponentially decays for $t>t_{th}$. In our approach,
however, we treat the period prior to lockdown and social actions with a a
linearly growing $\beta(t)$ function, instead of a constant, i.e.:

$\displaystyle\beta(t)=\beta_{0}t\ \ \ \mbox{for }0\leq t<t_{th}~{},$ (15)

and due to function continuity the expression for $\beta(t)$ after lockdown
measures becomes:

$\displaystyle\beta(t)=t_{th}\beta_{0}\ e^{-(t-t_{th})/\tau}\ \ \ \mbox{for
}t\geq t_{th}~{},$ (16)

where the time $t_{th}$ represents the starting time of the quarantine
actions, while $\tau$ refers to the decaying period and it has the dimension
of time.

##  4 Number of active cases

For the number of active cases, we perform a weighted fit of the model
allowing four free parameters to vary (the number of active cases is
independent of $\sigma$, the case fatality rate, discussed in Section 5) using
least squares minimization with the LMFIT [newville2016lmfit] package in
Python 3.6.8. We use a population size of $N=5.164\times 10^{7}$ people and as
initial condition we use $19$ active cases, the number of registered cases on
the 18th of February in South Korea according to the Worldometer database. The
end date of this first wave of infections is set to May 23, 2020.

The resulting fit, obtained by solving (14) for $\left<\mathcal{I}(t)\right>$,
overlaid with the real data points is shown in Figure 2. We also obtain the
mean absolute percentage error ($M$) of the real data and the fitted data
points and achieve $M=14.7\%$, which indicates acceptable agreement.

Figure 2: Number of active cases of COVID-19 in South Korea during 96 days of
the epidemic corresponding to first wave of infections in the country. Modeled
data is shown by the continuous red line, while measured data is depicted by
the blue points.

###  4.1 Fitted parameters

The best fit is given by the following set of parameters:

$\displaystyle\beta_{0}=1.898\times 10^{-6}\pm 1.202\times 10^{-7}\quad\quad$
, $\displaystyle\quad\quad\gamma=0.0710\pm 0.0127~{},$ $\displaystyle
t_{th}=10.008\pm 0.868\quad\quad$ , $\displaystyle\quad\quad\tau=19.010\pm
3.029~{},$

where the estimated uncertainties, calculated by LMFIT, are done by inverting
the Hessian matrix, which represents the second derivative of fit quality for
each variable parameter. It is interesting to note that 1/$\gamma$, the time
required to infected people to become recovered or deceased, is estimated to
be 1/0.0710 $\sim$ 14 days. The uncertainty on this value is, however, large
(18%), but even the shortest (12 days) and longest (17 days) values are within
acceptable ranges for the average duration of the disease
[fernandez2020estimating]. The $\beta_{0}$ parameter, which is the initial
constant for the infection rate, is a small number, which is unsurprising
given that we have considered the infection rate to grow linearly with time
until lockdown measures are in place. The value of $t_{th}$ is interpreted as
the day lockdown and social distancing measures started in South Korea during
the first wave of infections.

##  5 Case fatality rate

Figure 3: Number of deaths in South Korea during 96 days of the epidemic.
Modeled data is shown by the continuous red line, while measured data is
depicted by the blue points.

The case fatality rate ($\sigma$) in South Korea between February 18 and May
23, 2020 was obtained from the number of registered deaths. Figure 3 shows the
modeled data (red), the solution for $\left<\mathcal{D}(t)\right>$ in (14),
overlaid with the measured data (blue), and reasonable agreement can be seen
between them. The best fit of the measured data to our model is given by:

$\sigma=0.0109\pm 0.0011~{},$

which corresponds to a rate of $\sim 1.1\%$.

##  6 Second and third outbreaks in South Korea

Figure 4: Number of active cases of COVID-19 in South Korea between February
18, 2020 and March 16, 2021. Three distinct waves of epidemic activity can be
clearly seen from the above data.

South Korea has experienced three waves of the COVID-19 pandemic, as can be
seen on Figure 4, which shows the total number of active cases as a function
of time (days) using data taken from the Worldometer reference website. As our
model was able to reasonably model the first wave of infections, it is worth
checking its explanatory power during the second and third COVID-19 waves for
the same country. There is no strict definition for what an epidemic wave is,
with the term simply implying natural pattern of peaks and valleys of
contaminated individuals. Thus, we define the second wave between August 10,
20120 and October 1st 2020 and fit the data to our model using the same
approach as for the first wave. The result for the second wave can be seen on
Figure 5, with the best fit given by the following set of parameters:

$\displaystyle\beta_{0}=7.666\times 10^{-7}\pm 2.086\times 10^{-8}\quad\quad$
, $\displaystyle\quad\quad\gamma=0.0741\pm 4.110\times 10^{-5}~{},$
$\displaystyle t_{th}=15.653\pm 0.388\quad\quad$ ,
$\displaystyle\quad\quad\tau=12.589\pm 0.643~{}.$

Figure 5: Number of active cases of COVID-19 in South Korea during 53 days of
the epidemic corresponding to the so-called second wave of infections. Modeled
data is shown by the continuous red line, while measured data is depicted by
the blue points.

The mean absolute percentage error for the fit is of 5%, which points to an
even better agreement compared to the first wave, a fact that is partially
explained by the better available statistics. It must also be noted that the
end of the second wave of infections quite rapidly became the start of the
third wave, reason the data was truncated already on October 1st, i.e.,
equation 16 no longer holds, after around that date. The good agreement during
the expanding phase, followed by a worse agreement during the decay phase seen
in Figure 5 can be partially attributed to it.

Figure 6 shows the result for the third epidemic wave in South Korea, between
November 06 2020 and February 06 2021. The $\beta_{0}$ parameter is of the
same order of magnitude as for the second wave and there is also similar
agreement, with a mean absolute percentage error of 5.7%.

Figure 6: Number of active cases of COVID-19 in South Korea during 93 days of
the epidemic corresponding to the so-called third wave of infections. Modeled
data is shown by the continuous red line, while measured data is depicted by
the blue points.

##  7 Discussion of the model approach and results

The second quantization approach described in the present article is generic
and can be applied to a variety of problems for which stochastic modeling is
appropriate, being the main constraint of its application to classical
systems, as opposed to quantum systems. In, e.g., Ref. [MondainiBMST2015] this
approach was applied to the case of hepatitis C virus epidemic, but no
confrontation with real data was possible, and as such, it served only as an
initial presentation of the theoretical framework and a glimpse of its full
potential. To further reinforce the generality of the second quantization
approach to stochastic models it is also worth mentioning Ref.
[Mondaini_2017], which presented a simple two-variable (healthy/cancer cells)
system to the description of tumour growth, but again without yet successful
confrontation with real data. The present article presents therefore, not only
an application of the theoretical framework established in Ref.
[MondainiBMST2015] to a new problem, but also its first successful description
of real data.

Moreover, on the specific problem studied in the current article, namely the
COVID-19 epidemic, the model was successfully expanded by the use of a time-
dependent rate of infection, with a linear growing function for the initial
phase of the infection rate representing the most appropriate way to describe
the data. This is interesting from a public policy perspective, since it
reinforces the social distancing measures in the control of the COVID-19
epidemic. According to our model, on the absence of controls, the rate of
infection itself is accelerating.

Finally, the tension observed between model and data during the second and
third waves of the epidemic can, to a reasonable degree, be attributed to the
fact that the rate of infection decrease was not truly exponential, as
observed in the first wave. Again this might also have some interest for
policy makers since a quantifiable disagreement with an exponential decrease,
when using our model, indicates a potential change in the mathematical
behavior of the infection rate, which can quickly shift the epidemic phase
from infection decline to infection grow.

##  8 Concluding Remarks

Using an SIR-type stochastic model built on field theory techniques, we have
investigated the time evolution of the mean number of infectious (active
cases) and deceased individuals for COVID-19 epidemic in South Korea. Our
model was able to describe the available data successfully. To the best of our
knowledge, it is the first COVID-19 model built from standard field
theoretical language. Despite its simplicity, the model was able to describe
three waves of infections of the COVID-19 pandemic in South Korea, which
indicates its robustness and versatility. Studies including a non-susceptible
subset of the population due to vaccinations and their effects on the
infection rates, as well as comparisons with data from other countries, are
possible future directions to be explored with our model.

## References

## References

  * [1]
  * [2]
  * [3]