# Timescale dynamics of COVID-19 pandemic waves: The case of Greece

Dimitris M. Manias1, ,  Dimitris G. Patsatzis2,  Dimitris A. Goussis3,4
(1) Department of Mathematics, _Khalifa University of Science and Technology_
,
Abu Dhabi 127788, United Arab Emirates
(2) Modelling Engineering Risk and Complexity, _Scuola Superiore Meridionale_
, Naples 80138, Italy
(3) Department of Mechanical Engineering, _Khalifa University of Science and
Technology_ ,
Abu Dhabi 127788, United Arab Emirates
(4) Research and Innovation Center on CO2 and H2, _Khalifa University of
Science and Technology_ ,
Abu Dhabi 127788, United Arab Emirates Corresponding Author:
<EMAIL_ADDRESS>
###### Abstract

The results of an alternative methodology for making predictions about the
COVID-19 pandemic in Greece are presented. Instead of focusing on the various
population profiles (subjected to instabilities introduced by the fitting
process), this methodology focuses on the time scale that characterises the
intensity and duration of the _outbreak phase_. Therefore, instead of
predicting the peak of active cases, here their inflection point is predicted
(the point where the increase of active cases stops accelerating and starts
decelerating). Since the inflection point precedes the peak, this methodology
can serve as an early warning of the peak. In addition, the paths between the
various populations (healthy, exposed, infected, etc) that contribute the most
to the _outbreak phase_ are identified.

_Keywords_ Population dynamics $\cdot$ COVID-19 $\cdot$ Computational Singular
Perturbation $\cdot$ time scale analysis

##  1 Introduction

The major effort in making predictions on the spread of COVID-19 focuses on
the temporal evolution of various populations, notably active cases and
deaths. Here, the employment of Computational Singular Perturbation (CSP) [1,
2, 3] is proposed as an alternative methodological approach for the analysis
of COVID-19 pandemics compartmental models. The analysis is based on the time
scale that characterises the _outbreak phase_ (the period in which the number
of active cases increases) and employs tools used for the multi-scale analysis
of models in chemical kinetics, biology and pharmacokinetics.CSP is an
algorithmic method of Geometrical Singular Perturbation Analysis (GSPA), which
is employed in multi-scale dynamical systems; i.e., systems that are driven by
processes evolving in a wide range of fast/slow time scales. In such systems,
the processes characterized by the fast time scales become quickly exhausted
and form low dimensional surfaces in the tangent space, well-known as Slow
Invariant Manifolds (SIM) [4, 5, 6]. The trajectories of the system are then
confined to evolve on the SIM, governed by the processes which are
characterized by the remaining slow time scales. Recently, a physics-informed
machine learning approach was introduced for the computation of SIMs, in the
spirit of GSPA [7]. CSP offers a complete set of diagnostic tools for
identifying the variables and physical processes related to (i) the fast time
scales, which contribute to the formation of the SIM and (ii) the slow time
scales that govern the dominant slow dynamics of the system. The CSP tools
have been employed for gaining physical understanding of chemical kinetic
mechanisms [8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18] and other multi-scale
physical problems encountered in systems biology [19, 20] and population
dynamics [21].

The utilization of CSP offers significant advantages for analyzing
epidemiological models, particularly in compartmental contexts. CSP’s
relevance stems from its compatibility with multi-scale epidemiological
models, as demonstrated in various studies [22, 23, 24, 25, 26]. Unlike
solution-based approaches, CSP focuses on system dynamics, avoiding issues
like non-identifiability. It also provides a systems-level understanding,
surpassing the capabilities of sensitivity analysis in distinguishing fast and
slow dynamics. CSP’s effectiveness in exploring dynamical systems with
explosive eigenvalues, as seen in the initial exponential growth phase of
COVID-19 [27, 28, 29]. The presence of a positive or explosive eigenvalue,
which is a key characteristic during the exponential growth phase of COVID-19
pandemics, is a well-established concept [30, 31]. CSP methodology, known for
its effectiveness in studying dynamical systems with such eigenvalues [32, 33,
34, 21], is therefore a logical choice for analyzing the rapid expansion of
COVID-19.

In this study, we apply a novel approach to analyze COVID-19 pandemics in
Greece using the SEInsRD model, shown in Fig. 1. We focus on the outbreak
phases of the 4th, 5th, and 6th waves, identifying key populations and
processes driving the pandemic dynamics using CSP diagnostic tools.
Additionally, we explore the relationship between the system’s time scales and
the inflection point of the infected population. This approach allows early
estimation of inflection points, aiding in predicting the plateau phase of the
pandemic, as supported by previous studies [35, 36].

Figure 1: The SEInsRD population model employed; SP: susceptible population
(not infected), EP: exposed population (infected but not yet able to infect),
INSP: infected normally susceptible population (able to infect), ISSP:
infected severely susceptible population (able to infect), IDP: infected dead
population, RP: recovered population. $R^{i}$ denotes the transmission rate of
the $i$-th path.

##  2 Mathematical Model

In this study, the SEInsRD population model was employed, a description of
which is depicted in Fig. 1. Since research indicates that COVID-19 symptoms
vary from mild to severe, the SEInsRD model is an adaptation of the SEIRD
model, which expands the later by splitting the infected compartment into two
sub-compartments: ’normally infected’ (IN) and ’severely infected’ (IS)
individuals. This division reflects the immediate healthcare needs of severely
infected individuals and leads to the modified SEInsRD model:

$\dfrac{d}{dt}\begin{bmatrix}S\\\ E\\\ IN\\\ IS\\\ R\\\
D\end{bmatrix}=\begin{bmatrix}-\beta_{N}S.IN-\beta_{S}S.IS-\mu_{TP}S\\\
\beta_{N}S.IN+\beta_{S}S.IS-\sigma E-\mu_{TP}E\\\ (1-ss)\sigma E-\gamma
IN-\mu_{N}IN\\\ ss\sigma E-\gamma IS-\mu_{S}IS\\\ \gamma(IN+IS)-\mu_{TP}R\\\
\mu_{N}IN+\mu_{S}IS\end{bmatrix}$ (1)

where the subscripts $N$ and $S$ indicate the normally and severely infected
transmission $\beta$ and fatality $\mu$ ratios, $ss$ denotes the fraction of
severely over normally infected individuals and $\mu_{TP}$ is the
physiological death ratio.

The SEInsRD model can be rewritten to its - equal to density dependent -
frequency dependent formulation, as:

$\displaystyle\dot{SP}$ $\displaystyle=-R^{1}-R^{2}-R^{3}$
$\displaystyle\dot{EP}$ $\displaystyle=R^{1}+R^{2}-R^{4}-R^{5}$
$\displaystyle\dot{INSP}$ $\displaystyle=(1-ss)R^{4}-R^{6}-R^{8}$
$\displaystyle\dot{ISSP}$ $\displaystyle=ssR^{4}-R^{7}-R^{9}$ (2)
$\displaystyle\dot{RP}$ $\displaystyle=R^{6}+R^{7}-R^{10}$
$\displaystyle\dot{IDP}$ $\displaystyle=R^{8}+R^{9}$

where $SP$, $EP$, $INSP$, $ISSP$, $RP$ and $IDP$ denotes the number of
susceptible, exposed, infected normally symptomatic, infected severely
symptomatic, recovered and infected deceased individuals, respectively. The
rates $R^{i}$ are defined as:

$\displaystyle R^{1}$
$\displaystyle=\dfrac{\beta_{INSP}}{TP}INSP\leavevmode\nobreak\ SP$
$\displaystyle R^{2}$
$\displaystyle=\dfrac{\beta_{ISSP}}{TP}ISSP\leavevmode\nobreak\ SP$
$\displaystyle R^{3}$ $\displaystyle=\mu_{TP}SP$ $\displaystyle R^{4}$
$\displaystyle=aincp^{-1}EP$ $\displaystyle R^{5}$ $\displaystyle=\mu_{TP}EP$
$\displaystyle R^{6}$ $\displaystyle=aip^{-1}INSP$ $\displaystyle R^{7}$
$\displaystyle=aip^{-1}ISSP$ $\displaystyle R^{8}$
$\displaystyle=\mu_{INSP}INSP$ $\displaystyle R^{9}$
$\displaystyle=\mu_{ISSP}ISSP$ $\displaystyle R^{10}$
$\displaystyle=\mu_{TP}RP$

where the transmission ratios are expressed in a frequency dependent
representation as $\beta_{N}=\beta_{INSP}/TP$ and $\beta_{S}=\beta_{ISSP}/TP$,
$\mu_{N}=\mu_{INSP}$, $\mu_{S}=\mu_{ISSP}$, $\gamma=aip^{-1}$ and
$\sigma=aincp^{-1}$. In addition, Eq. (2) yields the forming equilibration:

$\dfrac{d(SP+EP+INSP+ISSP+RP+IDP)}{dt}=R^{3}+R^{5}+R^{10}=-\mu_{TP}(SP+EP+RP)$
(3)

More details on the inflection points through analytic calculations are
included in Appendix A.

##  3 Computational Singular Perturbation

The SEInsRD compartmental model of Eq. (2) can be formed as an Ordinary
Differential Equation (ODE) system:

$\dfrac{d\mathbf{y}}{dt}=\mathbf{g(y)}=\sum_{k=1}^{K}\mathbf{S}_{k}R^{k}(\mathbf{y})$
(4)

where $\mathbf{y}$ is the $N$-dim. column state vector, the elements of which
are the number of individuals in the population group or their fraction over
the population and $\mathbf{g(y)}$ is the $N$-dim. column vector field that
incorporates the transition laws from a compartmental group to another. In
particular, the vector field $\mathbf{g(y)}$ is described by the $N$-dim.
stoichiometric vector $\bm{S}_{k}$, which indicates the direction of the $K$
transitions from a group to another and $R_{k}$, which is the related
transition rate/interaction; e.g., transmission, recovery, death rate.
According to CSP, the system in Eq. (4) can be cast in the following form [8,
37]:

$\dfrac{d\mathbf{y}}{dt}=\mathbf{g}(\mathbf{y})=\sum_{n=1}^{N}\mathbf{a}_{n}(\mathbf{y})f^{n}(\mathbf{y}),\qquad
f^{n}(\mathbf{y})=\mathbf{b}^{n}(\mathbf{y})\cdot\mathbf{g(y)}=\sum_{k=1}^{K}\left(\mathbf{b}^{n}(\mathbf{y})\cdot\mathbf{S}_{k}\right)R^{k}(\mathbf{y}),$
(5)

where $\mathbf{a}_{n}(\mathbf{y})$ and $\mathbf{b}^{n}(\mathbf{y})$ are the
$N$-dim. CSP column and row, respectively, basis vectors of the $n$-th mode,
which satisfy the orthogonality conditions
$\mathbf{b}^{i}(\mathbf{y})\cdot\mathbf{a}_{j}(\mathbf{y})=\delta_{j}^{i}$ [8,
1]. In this way, the vector field $\mathbf{g(y)}$ is decomposed in $n$ CSP
modes, each having a different impact in the vector field and acting in a
different timescale. The impact of the $n$-th CSP mode is measured by the
amplitude $f^{n}(\mathbf{y})$ (proper adjustment of the signs of the CSP basis
vectors, set $f^{n}(\mathbf{y})$ always positive), which provides a measure of
the projection of the vector field $\mathbf{g(y)}$ on the CSP vector
$\mathbf{a}_{n}$. When the system in Eq. (5) exhibits $M$ timescales that are
[i] of dissipative nature (i.e., when the components of the system that
generate them tend to drive the system towards equilibrium) and [ii] much
faster than the rest, the model can be reduced in:

$f^{r}(\mathbf{y})\approx
0\quad(r=1,\ldots,M),\qquad\dfrac{d\mathbf{y}}{dt}\approx\sum_{s=M+1}^{N}\mathbf{a}_{s}(\mathbf{y})f^{s}(\mathbf{y}),$
(6)

the first relation of which is an $M$-dim. system of algebraic equations
defining the manifold $\mathscr{M}$ (a low dimensional surface in phase-space,
where the system is confined to evolve), while the second relation is an
$N$-dim. system of ODEs governing the slow evolution of the system on this
manifold. The objective of asymptotic analysis techniques is to deliver the
expressions in Eq. (6), which are provided, order by order, by CSP in
algorithmic fashion [38, 39]. This is accomplished through the CSP vectors,
$\mathbf{a}_{i}(\mathbf{y})$ and $\mathbf{b}^{i}(\mathbf{y})$
($i=1,\ldots,N$), which can be approximated in leading order accuracy by the
right and left, respectively, eigenvectors of the $N\times N$-dim. Jacobian
$\mathbf{J(y)}$ of $\mathbf{g(y)}$; i.e.,
$\mathbf{a}_{i}(\mathbf{y})=\bm{\alpha}_{i}(\mathbf{y})$ and
$\mathbf{b}^{i}(\mathbf{y})=\bm{\beta}^{i}(\mathbf{y})$ [37, 8, 40]. In the
following, for simplicity reasons, the dependency of $R^{k}$, $\mathbf{g}$,
$\bm{\alpha}$, etc, to $\mathbf{y}$ will be removed, and be mentioned only
when explicitly needed.

The interactions among the population groups have a different impact to the
evolution of the system, either to the constraints along which the system is
confined to evolve or to its slow evolution along them. In particular, the $M$
constraints in Eq. (6)
$f^{r}=(\bm{\beta}^{r}\cdot\mathbf{S}_{1})R^{1}+\ldots+(\bm{\beta}^{r}\cdot\mathbf{S}_{K})R^{K}\approx
0$ ($r=1,\ldots,M$) result from significant cancellations among some of the
additive terms $(\bm{\beta}^{r}\cdot\mathbf{S}_{k})R^{k}$ ($k=1,\ldots,K$). In
addition, the slow evolution of the system is characterized by the slow CSP
modes, the impact of which is measured by the slow amplitudes $f^{s}$ in Eq.
(6). The interactions that contribute significantly to (i) the formation of
each of the $M$ constraints and (ii) the governing components of the slow
system, are identified by the Amplitude Participation Index (API):

$P^{n}_{k}=\frac{\displaystyle(\bm{\beta}^{n}\mathbf{S}_{k})R^{k}}{\displaystyle\sum\nolimits_{i=1}^{K}|(\bm{\beta}^{n}\mathbf{S}_{i})R^{i}|}\qquad(k=1,\ldots,K),$
(7)

where by definition $\sum_{k=1}^{K}|P_{k}^{n}|=1$ [1, 41, 42]. When employed
to the fast $r=1,\ldots,M$ CSP modes, $P^{r}_{k}$ measures the relative
contribution of the $k$-th interaction to the cancellations among the additive
terms in $f^{r}\approx 0$, while when employed to the slow $s=M+1,\ldots,N$
CSP modes, $P^{s}_{k}$ measures the relative contribution of the $k$-th
interaction to governing slow system. Both $P^{r}_{k}$ and $P^{s}_{k}$ can be
either positive or negative; the sum of positive and negative $P^{r}_{k}$
terms equaling 0.5, by definition [43, 44, 45, 46].

The formation of the $M$ constraints is characterized by the $M$ fastest
timescales, which by definition are of dissipative nature; see requirement [i]
immediately before Eq. (6). On the contrary, the dynamics of the slow system
in Eq. (6) is characterized either by the fastest of the $N-M$ dissipative
slow timescales or by explosive timescales, when the later are present, since
they tend to drive the system away from equilibrium. The timescales are
approximated by the inverse of the eigenvalues of the Jacobian $\mathbf{J}$,
$\tau_{n}=|\lambda_{n}|^{-1}$ ($n=1,\ldots,N$). The CSP diagnostic tool,
Timescale Participation Index (TPI), identifies the reactions significantly
contributing to the generation of the timescales:

$J^{n}_{k}=\dfrac{c_{k}^{n}}{\sum\nolimits_{i=1}^{K}|c^{n}_{i}|}\qquad(k=1,\ldots,K),$
(8)

where ${\lambda}_{n}=c_{1}^{n}+\ldots+c_{K}^{n}$ and by definition
$\sum_{k=1}^{K}|J_{k}^{n}|=1$ [47, 41, 48]. $c_{k}^{n}$ denotes the
contribution of the $k$-th interaction to the $n$-th eigenvalue and can be
calculated as
$c_{k}^{n}=\bm{\beta}^{n}\nabla\left(\mathbf{S}_{k}R^{k}\right)\bm{\alpha}_{n}$,
where $\sum_{k=1}^{K}\nabla\left(\mathbf{S}_{k}R^{k}\right)$ is the Jacobian
$\mathbf{J}$. $c_{k}^{n}$ can be either positive or negative and therefore, a
negative (positive) $J^{n}_{k}$ implies that the $k$-th interaction
contributes to a dissipative (explosive) character of the $n$-th timescale
$\tau_{n}$ [49, 50]. By definition, dissipative (explosive) timescales relate
to the components of the system that tend to drive it towards (away from)
equilibrium [8, 1]. TPI has been succesfuly used

Each population group/compartment is associated differently to each CSP mode;
e.g., a fast CSP mode is expected to be much more related to the infected
compartment than the recovered one, since the former are playing a major role
as opposed to the latter at the beginning of an epidemics. The relation of the
$m$-th CSP mode ($m=1,\dots,M$) to the various population groups is identified
by the Pointer (Po):

$\mathbf{D}^{m}=diag\left[\bm{\alpha}_{m}\bm{\beta}^{m}\right]=\left[\alpha^{1}_{m}\beta^{m}_{1},\alpha^{2}_{m}\beta^{m}_{2},\ldots,\alpha^{N}_{m}\beta^{m}_{N}\right]\qquad(m=1,\ldots,M),$
(9)

where, due to the orthogonality condition
$\bm{\beta}^{i}\cdot\bm{\alpha}_{j}=\delta^{i}_{j}$, the sum of all $N$
elements of $\mathbf{D}^{m}$ equals unity, i.e.
$\sum_{i=1}^{N}\alpha^{i}_{m}\beta^{m}_{i}=1$ [3, 1, 42, 51]. A relatively
large value of $\alpha^{i}_{m}\beta^{m}_{i}$ indicates that the $i$-th
population group is strongly associated to $m$-th CSP mode and the $m$-th
timescale. A value of $D^{m}_{i}$ close to unity suggests that the $i$-th
population group is in Quasi Steady-State (QSS) [51].

##  4 Results

The fastest of the _explosive time scales_ in the dynamics of the process
characterises the _outbreak phase_ from the start of the process, up to the
point where these time scales cease to exist (this is the period in which the
system accelerates at the start of the outbreak phase) and that the remaining
period of this phase is characterised by the fastest of the _dissipative time
scales_ (this is the period in which the system decelerates as it approaches
the end of the outbreak phase, known as the _peak_). The transition from the
accelerating to the decelerating period, defines the _inflection point_ in the
profile of the active cases, during the _outbreak phase_. The development of
the _inflection point_ is an early warning that the _peak_ will be reached
soon.

The explosive and dissipative time scales that characterize the outbreak phase
were computed by fitting the available data (from Worldometer database, number
of exposed, infected and dead) to the SEInsRD model, shown in Fig. 1.

###  4.1 The 4th wave

In this section, results for the 4th wave will be presented. For validation
purposes, the analysis was done by considering solutions obtained by fitting
data from two periods in 2021:

  1. i)

Period A, July 1 - July 8, 8 days long,

  2. ii)

Period B, July 1 - July 18, 18 days long,

where July 1 marks the start of the 4th wave.

The evolution of the active cases from the Worldometer database are shown in
Fig. 2 with circles (up to 12/8/21). In addition, Fig. 2 displays the number
of active cases, as it was computed from the SEInsRD model, when fitted to
Period A (8 days long) and Period B (18 days long) data.

Figure 2: The active cases (INSP+ISSP, see Fig. 1) according to the available
data (circles) and the solution of the SEInsRD model (curves), fitted to
Period A (July 1-8) and Period B (July 1-18) data. The two periods are shown
with horizontal arrows.

The _explosive_ and _dissipative_ (attenuating) time scales that characterise
the dynamics in the outbreak phase are shown in Fig. 3. The _fast explosive_
time scale characterises the evolution of the process, from July 1st until the
day when this time scale ceases to exist; i.e., July 14 according to Period A
and July 12 according to period B.

Figure 3: The explosive (in red) and dissipative (attenuating, in gray) time
scales; profiles based by fitting data from Period A (8 days) and Period B (18
days). The vertical blue dashed line denotes the point where the two explosive
time scales coalesce, indicating that their disappearance is imminent. The
vertical green dashed line denotes the point where the explosive time scales
disappear. The fastest explosive time scale characterises the initial part of
the outbreak phase and the fastest dissipative (attenuating) time scale
characterises the remaining part.

According to the profiles shown in Fig. 2, the _inflection point_ in the
profile of the active cases was manifested on July 15 and July 13, according
to the solutions based on data spanning Period A and Period B, respectively.
These two days are approximated very well by the days in which the _explosive_
time scales disappear; i.e., July 14 and July 12 for the two periods, as shown
in Fig. 3. The disappearance of the _explosive_ time scales is preceded by the
their coalescence, as shown in Fig. 3; July 9 for Period A and July 8 for
Period B. According to these findings, by July 8-9 it was possible to predict
that the development of a peak (plateau in this case) was imminent.

Given the significance of the _fast explosive_ time scale (the faster this
scale the more intense the outbreak), it is of interest to identify the paths
of the model in Fig. 1, which either promote or oppose its appearance. The
results displayed in Table 1 suggest that the intensity of the initiation of
the outbreak is mainly promoted by Paths 1 and 4 (susceptible becoming exposed
and then infected), while it is opposed by Path 6 (infected becoming
recovered). In addition, it is shown that - as the end of the explosive period
is reached - the promoting influence of Paths 1 and 4 diminishes, while the
opposing one of Path 6 increases (doubles).

Table 1: The percentage contribution of the various paths to the development
of the fastest explosive time scale, estimated on July 1 and 7.

| Period A solution
---|---
Path | July 1 | July 7
$R^{1}$ |  51.8% |  43.8%
$R^{4}$ |  32.6% |  28.9%
$R^{6}$ | -14.8% | -26.3%

| Period B solution
---|---
Path | July 1 | July 7
$R^{1}$ |  47.6% |  35.1%
$R^{4}$ |  38.3% |  33.4%
$R^{6}$ | -13.6% | -30.3%

The analysis of available data extending in the period July 1-8 produced
similar results with the analysis of data extending in the period July 1-18.
By concentrating on the inflection point in the profile of active cases, it
was shown that we could know on July 8 or 9 - the latest - that the peak of
active cases was within reach (eventually it was recorded around July 22).
Finally, it was shown that it is possible to quantify, the degree to which the
intensity and duration of the _outbreak phase_ can be controlled by decreasing
the rate by which susceptible become infected and by increasing the rate by
which infected become recovered.

###  4.2 The 5th wave

In this section, results for the 5th wave will be presented. For validation
purposes, the analysis was done by considering solutions obtained by fitting
data from two periods in 2021:

  1. i)

Period C, October 10 - October 23, 14 days long,

  2. ii)

Period D, October 10 - November 6, 28 days long,

where October 10 marks the start of the 5th wave.

The evolution of the active cases from the Worldometer database are shown in
Fig. 4 with circles (up to 23/11/21). In addition, Fig. 4 displays the number
of active cases, as it was computed from the SEInsRD model, when fitted to
Period C (14 days long) and Period D (28 days long) data; the extend of the
two periods shown in Fig. 4.

Figure 4: The active cases (INSP+ISSP, see Fig. 1) according to the available
data (circles) and the solution of the SEInsRD model (curves), fitted to
Period C (October 10 - 23) and Period D (October 10 - November 6) data. The
two periods are shown with horizontal arrows. For comparison, the solution of
the model fitted to Periods A and B are displayed.

The _explosive_ and _dissipative_ (attenuating) time scales that characterise
the dynamics in the outbreak phase are shown in Fig. 7. The _fast explosive_
time scale is present from the assumed start of the 5th wave in October 10,
until 2 January 2022 according to the Period C data and 26 December 2021
according to Period D data. Since the influence of the fast explosive time
scale weakens after it meets the slow one (21-22 November), it is expected
that the inflection point of the infected population will manifest early in
December.

Figure 5: The explosive (in red) and dissipative (attenuating, in gray) time
scales; profiles based by fitting data from Period C (14 days) and Period D
(28 days). The vertical blue dashed line denotes the point where the two
explosive time scales coalesce, indicating that their disappearance is
imminent. The vertical green dashed line denotes the point where the explosive
time scales disappear. The fastest explosive time scale characterises the
initial part of the outbreak phase and the fastest dissipative (attenuating)
time scale characterises the remaining part.

Similarly to 4th wave, it is of interest to identify the paths of the model in
Fig. 1, which either promote or oppose the generation of the _fast explosive_
time scale. The results displayed in Table 2 suggest that the intensity of the
initiation of this outbreak is mainly promoted by Paths 1 and 4 (susceptible
becoming exposed and then infected), while it is opposed by Path 6 (infected
becoming recovered). It is shown that over a 2-weeks period from the start of
the 5th wave, there are no significant changes in the major contributions from
paths 1, 4 and 6. This is in contrast to the 4th wave (see Table 1), where
over a 1-week period from the start there was a significant decrease in the
promoting to the explosive stage contribution of Paths 1 and 4 and a
significant increase in the opposing contribution of Path 6. This feature is
related to the expectation that the duration of the 5th wave will be longer
than that of the 4th one; see Section 5.

In addition, a comparison of Tables 1 and 2 reveals that the opposing to the
explosive time scale influence of Path 6 (infected becoming recovered) is much
stronger in wave 5 (about -33% during the first two weeks of the outbreak)
than in wave 4 (about -15% at the start of the wave to about -30% a week
later).

Table 2: The percentage contribution of the various paths to the development
of the fastest explosive time scale, estimated on October 10 and 24.

| Period C solution
---|---
Path | October 10 | October 24
$R^{1}$ |  48.9% |  48.4%
$R^{6}$ | -33.4% | -34.5%
$R^{4}$ |  17.1% |  16.5%

| Period D solution
---|---
Path | October 10 | October 24
$R^{1}$ |  48.9% |  48.5%
$R^{6}$ | -32.2% | -33.2%
$R^{4}$ |  18.3% |  17.7%

The analysis of available data extending in the period October 10-23 produced
qualitatively similar results with the analysis of data extending in the
period October 10 - November 6. Both periods predict the development of an
inflection point (the point where the increase of active cases stops
accelerating and starts decelerating) early in December.

_Although the explosive period during the 5th wave is longer than that during
the 4th wave, the opposing Path 6 (infected population becoming recovered one)
is stronger in the 5th wave from its start, most likely due to the larger
portion of vaccinated people in the infected population. Consequently, if this
opposing factor had not being so strong, the promoting influence of Paths 1
and 4 (susceptible population becoming exposed and then infected one) would
have made the _outbreak phase_ even stronger. As in the 4th wave, this
opposing influence of Path 6 is expected to grow further towards the end of
the explosive period._

###  4.3 The 6th wave

In this section, results for the 6th wave will be presented. The analysis was
done by considering a solution obtained by fitting data from a period of ten
days, say Period E, from 26 December 2021 to 4 January 2022, where 26 December
2021 marks the start of the 6th wave.

Figure 6: The active cases (INSP+ISSP, see Fig. 1) according to the available
data (circles) and the solution of the SEInsRD model (curves), fitted to
Period E (December 26 ’21- January 4 ’22) data. For comparison, the solution
of the model fitted to Periods C and D are displayed. The various periods are
shown with horizontal arrows.

The evolution of the active cases from the Worldometer database are shown in
Fig. 6 with circles (up to 08/01/22), along with the number of active cases,
as it was computed from the SEInsRD model, when fitted to the Period E data.

The _explosive_ and _dissipative_ (attenuating) time scales that characterise
the dynamics in the outbreak phase are shown in Fig. 7. The _fast explosive_
time scale is present from the assumed start of the 6th wave in December 26
until 31 December 2021. Since the influence of the fast explosive time scale
weakens after it meets the slow one (30 December), it is expected that the
inflection point of the infected population will manifest the first days of
January 2022.

Figure 7: The explosive (in red) and dissipative (attenuating, in gray) time
scales; profiles based by fitting data from Period E (10 days). The vertical
blue dashed line denotes the point where the two explosive time scales
coalesce, indicating that their disappearance is imminent. The vertical green
dashed line denotes the point where the explosive time scales disappear. The
fastest explosive time scale characterises the initial part of the outbreak
phase and the fastest implosive (dissipative) time scale characterises the
remaining part.

The use of the CSP tools reveals the paths that affect the generation of the
_fast explosive_ time scale. The results displayed in Table 3 suggest that the
intensity of the initiation of this outbreak is mainly promoted by Path 4
(exposed becoming effected) and in a lesser degree by Path 1 (susceptible
becoming exposed), while it is opposed by Path 6 (infected becoming
recovered).

Table 3: The percentage contribution of the various paths to the development of the fastest explosive time scale, estimated on December 26 and 28.  | Period E solution
---|---
Path | December 26 | December 28
$R^{4}$ |  56.8% |  59.1%
$R^{1}$ |  35.6% |  28.6%
$R^{6}$ | -7.4% | -11.9%

According the Table 3, during the outbreak of the 6th wave, the relevant
contributions of the main paths do not change much. Furthermore, it is shown
that the opposing to the explosive time scale influence of Path 6 (infected
becoming recovered) is much more weaker in the 6th wave, compared to the one
of both waves 4 and 5; see Tables 1 and 2.

The analysis of available data extending in the period 26 December 2021 to 4
January 2022 predicted the development of an inflection point of active cases
(the point where the increase of cases stops accelerating and starts
decelerating) early in January. This implies that the relaxation period
(decrease of active cases) will initiate soon.

In contrast to the case of the 4th and 5th waves, the promoting the outbreak
phase influence of Path 4 (exposed population becoming infected) was more
significant than that of Path 1 (susceptible population becoming exposed),
confirming the wider dispersion of the omicron variant to the general
population

It was also shown that the opposing influence of Path 6 (infected population
becoming recovered one) is very small; smaller than in the 4th and 5th waves,
most likely due to the weaker effects of the vaccines when dealing with the
omicron variant.

##  5 Comparison between waves 4, 5 and 6

Tables 4 and 5 allow for a direct comparison of the dynamics of the three
outbreak phases and the origin of their differences.

Table 4 lists the explosive time scale $\tau$ at the start of the outbreak
phase and at the point where this time scale disappears; the latter coinciding
with the inflection point of the active cases. It is shown that $\tau$ is
about 8 d, 28 d and 4 d at the start of the outbreak phase in the case of the
4th, 5th and 6th waves, respectively. This implies that the process is twice
as fast at the start of the outbreak phase of the 6th wave in comparison to
that of the 4th wave and seven times as fast in comparison to that of the 5th
wave. Regarding the point where the fast time scale disappears, the 4th and
6th waves are equally fast and about twice as fast in comparison to the 5th
wave.

Table 4: Comparison between the developing explosive time scales of the 4th, 5th and 6th waves. $\tau_{e}^{0}$ and $\tau_{e}^{end}$ represent the value of the fast explosive time scale at the start of each wave and at the time when the explosive modes disappear, respectively. Wave 4 | Wave 5 | Wave 6
---|---|---
Period A | Period B | Period C | Period D | Period E
$\tau_{e,A}^{0}$ = 8.33 d |  $\tau_{e,B}^{0}$ = 7.24 d |  $\tau_{e,C}^{0}$ = 29.5 d |  $\tau_{e,D}^{0}$ = 25.89 d |  $\tau_{e,E}^{0}$ = 3.8 d
$\tau_{e,A}^{end}$ = 30.53 d |  $\tau_{e,B}^{end}$ = 29.14 d |  $\tau_{e,C}^{end}$ = 78.7 d |  $\tau_{e,D}^{end}$ = 71.1 d |  $\tau_{e,E}^{end}$ = 34.9 d
Table 5: Comparison between the major contributors to the explosive time scale at the initiation of the 4th, 5th and 6th waves. Wave 4 | Wave 5 | Wave 6
---|---|---
$R^{1}$: $\sim$ 50% |  $R^{1}$: $\sim$ 49% |  $R^{4}$: $\sim$ 59%
$R^{4}$: $\sim$ 36% |  $R^{6}$: $\sim$-33% |  $R^{1}$: $\sim$ 29%
$R^{6}$: $\sim$-14% |  $R^{4}$: $\sim$ 18% |  $R^{6}$: $\sim$-12%

Table 5 lists the contributions of the various paths of the model employed to
the explosive time scale at the start of the outbreak in each wave considered.
The two major findings are (i) the switch from the non-infected population to
the exposed one, as the population contributing the most to the outbreak
phase, which confirms the wider dispersion of the omicron variant to the
general population and (ii) the weaker opposition to the outbreak phase in the
6th wave relative to the 5th one (the first dominated by the omicron variant
and the latter by the delta one), which is in accordance to the higher
infection rate due to the omicron variant and the decreased efficiency of the
vaccines when dealing with this variant..

##  6 Conclusions

Results based on an alternative and robust methodology for predicting the
evolution of each COVID-19 wave were presented. The methodology focuses on the
_fast explosive_ time scale that characterises the intensity and duration of
the _outbreak phase_. The point were this time scale ceases to characterize
the process coincides with the inflection point of active cases; the point
where the increase of active cases stops accelerating and starts decelerating.
Since the inflection point precedes the peak, this methodology serves as an
early warning of the peak. In addition, this methodology allows for the
identification of the factors (paths) that can influence the intensity and
length of the _outbreak phase_.

The last major three epidemic waves were analyzed, as they were recorded in
the periods 1.7.2021-9.10.2021 (4th wave), 10.10.2021-19.12.2021 (5th wave)
and 20.12.2021-today (6th wave, still developing); the first two periods
dominated by the _delta_ variant of the virus and the last one by the
_omicron_ variant. By examining the dynamics of their _outbreak phase_ , it is
concluded that:

  1. i)

The most intense _outbreak phase_ is that of the 6th wave, followed by that of
the 4th wave; the one of the 5th wave being by far the least intense. This
finding confirms (i) _the dominance of the omicron variant, relative to the
delta one_ , (ii) _the relaxing influence of vaccinations_ (much higher
percentage of vaccinated population during the 5th wave relative to the 4th
and not significant difference between the 5th and the 6th waves) and (iii)
the weak influence of vaccines when dealing with the omicron variant.

  2. ii)

The major contribution to the _outbreak phase_ of the 4th and 5th waves
originates from the non-infected population becoming exposed (infected but not
yet able to transmit the virus), followed by the exposed population becoming
infected (able to infect). In the 6th wave, this ordering is reversed, so that
the population group dominating the _outbreak phase_ is the exposed population
and not the non-infected ones of the 4th and 5th waves. This switch _confirms
the wider dispersion of the omicron variant in the general population_.

  3. iii)

The major opposing influence to the intensity of the _outbreak phase_ of all
three waves originates from the infected population becoming recovered. This
opposition is strongest in the 5th wave (least intense wave) and weakest in
the 6th wave (most intense wave).

    1. a)

The stronger opposition to the _outbreak phase_ in the 5th wave relative to
the 4th one (both dominated by the delta variant) is _in accordance to the
higher vaccinated percentage of the general population during the 5th wave_.

    2. b)

The weaker opposition to the _outbreak phase_ in the 6th wave relative to the
5th one (the first dominated by the omicron variant and the latter by the
delta one) is _in accordance to the higher infection rate due to the omicron
variant and the decreased efficiency of the vaccines when dealing with this
variant_.

  4. iv)

The inflection point of active cases in the current 6th wave occurred in the
first days of January 2022, so the outbreak is about to start relaxing.

NOTE 1: _The predictions reported here for the 4th wave were based on the
prevailing conditions in Periods A (July 1-8) and B (July 1-18). It was
demonstrated that it was possible to predict by July 8 or 9 - the latest -
that the peak of active cases was within reach; eventually, it was recorded
around July 22._

_The predictions reported for the 5th wave are valid under the assumption that
the prevailing conditions in Periods C (October 10-23) and D (October 10 -
November 6) keep prevailing in the period that follows._

_Similarly, the predictions reported for the 6th wave are valid under the
assumption that the prevailing conditions in Periods E (26 December 2021 to 4
January 2022) keep prevailing in the period that follows._

NOTE 2: _In accordance to recent studies, it is assumed here that the
probability of transmission between vaccinated or not-vaccinated does not
differ significantly_.

##  Appendix A Appendix

According to the frequency dependent representation of SEInsRD model in Eq.
(2), the following conclusions can be drawn, regarding to the inflection
points $SP$, $EP$, $INSP$ and $ISSP$, through analytic calculations.

###  A.1 Inflection point of SP close to max R1 and R2

Differentiation of the SP differential equation in Eq. (2) yields:

$\ddot{SP}=-\dot{R}^{1}-\dot{R}^{2}-\dot{R}^{3}=-(\dot{R}^{1}+\dot{R}^{2})+\mu_{TP}(R^{1}+R^{2}+R^{3})$
(10)

Neglecting the terms $\mu_{TP}(R^{1}+R^{2}+R^{3})$, since $\mu_{TP}\ll 1$,
yields:

$\ddot{SP}\approx-(\dot{R}^{1}+\dot{R}^{2})$ (11)

During the explosive stage both $R^{1}$ and R2 increase with time; i.e.
initially the rate by which $SP$ gets infected increases with time
($\dot{SP}<0$ and $\ddot{SP}<0$). According to Eq. (11), an inflection point
of $SP$ ($\ddot{SP}=0$) might occur when R1 and R2 attain their maximum values
($\dot{R}^{1}=\dot{R}^{2}=0$). The inflection point of SP and the points in
time where R1 and R2 attain their maximum values are indeed very close. The
differences are related to the magnitude of the term
$\mu_{TP}(R^{1}+R^{2}+R^{3})$.

###  A.2 Inflection points of SP and EP

Differentiation of the sum of SP and EP equations in Eq. (2) yields:

$\ddot{SP}+\ddot{EP}=-\dot{R}^{3}-\dot{R}^{4}-\dot{R}^{5}=-\dfrac{\dot{EP}}{aincp}-\mu_{TP}(\dot{SP}+\dot{EP})$
(12)

Neglecting the terms $\mu_{TP}(\dot{SP}+\dot{EP})$, since $\mu_{TP}\ll 1$,
yields:

$\ddot{SP}+\ddot{EP}\approx-\dfrac{1}{aincp}\dot{EP}$ (13)

Given that $\dot{EP}>0$ during the explosive stage, Eq. (13) indicates that
$\ddot{SP}+\ddot{EP}<0$ there. In addition, with the exception of the very
early and very late parts of the explosive stage, $\ddot{SP}\leq 0$ and
$\ddot{EP}\geq 0$ throughout this stage; i.e., the rate of increase of $EP$
and the rate of decrease of $SP$ intensify. Moving towards the end of the
explosive stage, $\dot{EP}$ decreases, being reflected with a drop of the
increasing trend of EP. At the end of the explosive stage, $\dot{EP}$ attains
small positive values, so that Eq. (13) yields that the inflection point of EP
($\ddot{EP}=0$) comes shortly before the end of the explosive stage, followed
by the inflection point of SP ($\ddot{SP}=0$).

###  A.3 Inflection points of INSP and ISSP coinciding with max R4

Differentiation of the sum of INSP and ISSP equations in Eq. (2) yields:

$\displaystyle\ddot{INSP}+\ddot{ISSP}$
$\displaystyle=\dot{R}^{4}-\dot{R}^{6}-\dot{R}^{7}-\dot{R}^{8}-\dot{R}^{9}$
$\displaystyle=\dfrac{\dot{EP}}{aincp}-\left(\dfrac{1}{aip}+\mu_{INSP}\right)\dot{INSP}-\left(\dfrac{1}{aip}+\mu_{ISSP}\right)\dot{ISSP}$
(14)

The summation of all the differential equations in Eq. (2) yields:

$\dot{SP}+\dot{EP}+\dot{INSP}+\dot{ISSP}+\dot{RP}+\dot{IDP}=-R^{3}-R^{5}-R^{10}=$
(15)

in which the terms $\dot{RP}$, $\dot{IDP}$ and $-R^{3}-R^{5}-R^{10}$ are very
small during the explosive stage, such that:

$\dot{SP}+\dot{EP}+\dot{INSP}+\dot{ISSP}=0$ (16)

Differentiation of Eq. (16) yields:

$\ddot{SP}+\ddot{EP}+\ddot{INSP}+\ddot{ISSP}=0$ (17)

which, after algebraic calculations, takes the form:

$\mu_{TP}(\dot{SP}+\dot{EP})+\left(\dfrac{1}{aip}+\mu_{INSP}\right)\dot{INSP}+\left(\dfrac{1}{aip}+\mu_{ISSP}\right)\dot{ISSP}=0$
(18)

Differentiation of Eq. (18) implies:

$\mu_{TP}(\ddot{SP}+\ddot{EP})+\left(\dfrac{1}{aip}+\mu_{INSP}\right)\ddot{INSP}+\left(\dfrac{1}{aip}+\mu_{ISSP}\right)\ddot{ISSP}=0$
(19)

in which the substitution of Eq. (13) yields:

$-\left(\mu_{TP}^{2}+\dfrac{\mu_{TP}}{aincp}\right)\dot{EP}+\left(\dfrac{1}{aip}+\mu_{INSP}\right)\ddot{INSP}+\left(\dfrac{1}{aip}+\mu_{ISSP}\right)\ddot{ISSP}=0$
(20)

Given that $\dot{EP}>0$, $\ddot{INSP}>0$ and $\ddot{ISSP}>0$ during the
explosive stage and that $\dot{INSP}$ and $\dot{ISSP}$ attain small positive
values, when moving towards to the end of it, Eq. (20) denotes that the
inflection points of INSP and ISSP ($\ddot{INSP}\approx\ddot{ISSP}=0$)
coincide with the time point at which EP attains its maximum value
($\dot{EP}=0$).

## References

  * Lam and Goussis [1994a] S. H. Lam and D. A. Goussis.  The csp method for simplifying kinetics.  _International journal of chemical kinetics_ , 26(4):461–486, 1994a.
  * Lam and Goussis [1994b] SH Lam and DA Goussis.  The csp method for simplifying kinetics.  _International journal of chemical kinetics_ , 26(4):461–486, 1994b.
  * Goussis and Lam [1992] D.A. Goussis and S.H. Lam.  A study of homogeneous methanol oxidation kinetics using csp.  In _Symposium (International) on Combustion_ , volume 24, pages 113–120. Elsevier, 1992.
  * Valorani et al. [2003a] Mauro Valorani, Habib N Najm, and Dimitris A Goussis.  Csp analysis of a transient flame-vortex interaction: time scales and manifolds.  _Combustion and Flame_ , 134(1-2):35–53, 2003a.
  * Goussis [2013] Dimitris A Goussis.  The role of slow system dynamics in predicting the degeneracy of slow invariant manifolds: the case of vdp relaxation–oscillations.  _Physica D: Nonlinear Phenomena_ , 248:16–32, 2013.
  * Maris and Goussis [2015] Dimitris T Maris and Dimitris A Goussis.  The “hidden” dynamics of the rössler attractor.  _Physica D: Nonlinear Phenomena_ , 295:66–90, 2015.
  * Patsatzis et al. [2023] Dimitrios G Patsatzis, Gianluca Fabiani, Lucia Russo, and Constantinos Siettos.  Slow invariant manifolds of singularly perturbed systems via physics-informed machine learning.  _arXiv preprint arXiv:2309.07946_ , 2023.
  * Lam and Goussis [1989] S.H. Lam and D.A. Goussis.  Understanding complex chemical kinetics with computational singular perturbation.  In _Symposium (International) on Combustion_ , volume 22, pages 931–941. Elsevier, 1989.
  * Hadjinicolaou and Goussis [1998] M Hadjinicolaou and D A Goussis.  Asymptotic solution of stiff PDEs with the CSP method: the reaction diffusion equation.  _SIAM J. Sci. Comput._ , 20:781–810, 1998.
  * Manias et al. [2016] Dimitris M Manias, Efstathios Al Tingas, Christos E Frouzakis, Konstantinos Boulouchos, and Dimitris A Goussis.  The mechanism by which ch2o and h2o2 additives affect the autoignition of ch4/air mixtures.  _Combustion and Flame_ , 164:111–125, 2016.
  * Rabbani et al. [2023] Shahid Rabbani, Dimitris M Manias, Dimitrios C Kyritsis, and Dimitris A Goussis.  Comparative review of the chemical dynamics underlying five models of ammonia fuel oxidation.  _Fuel_ , 352:129063, 2023.
  * Rabbani et al. [2022] Shahid Rabbani, Dimitris M Manias, Dimitrios C Kyritsis, and Dimitris A Goussis.  Dominant dynamics of n-butanol/air autoignition and the influence of additives.  _Combustion and Flame_ , 242:112173, 2022.
  * Rabbani et al. [2021] Shahid Rabbani, Dimitris M Manias, Dimitrios C Kyritsis, and Dimitris A Goussis.  Chemical dynamics of the autoignition of near-stoichiometric and rich methanol/air mixtures.  _Combustion Theory and Modelling_ , Accepted:1–31, 2021\.
  * Manias et al. [2021] Dimitris M Manias, Dimitris G Patsatzis, Dimitrios C Kyritsis, and Dimitris A Goussis.  Nh3 vs. ch4 autoignition: A comparison of chemical dynamics.  _Combustion Theory and Modelling_ , 25(6):1110–1131, 2021.
  * Manias et al. [2022a] Dimitris M Manias, Aliou Sow, Efstathios-Al Tingas, Francisco E Hernandez Perez, Hong G Im, and Dimitris A Goussis.  Mechanism of cold-spot autoignition in a hydrogen/air mixture.  _arXiv preprint arXiv:2210.06953_ , 2022a.
  * Khalil et al. [2021] Ahmed T Khalil, Dimitris M Manias, Dimitrios C Kyritsis, and Dimitris A Goussis.  No formation and autoignition dynamics during combustion of h2o-diluted nh3/h2o2 mixtures with air.  _Energies_ , 14(1):84, 2021.
  * Khalil et al. [2019] Ahmed T Khalil, Dimitris M Manias, Efstathios-Al Tingas, Dimitrios C Kyritsis, and Dimitris A Goussis.  Algorithmic analysis of chemical dynamics of the autoignition of nh3–h2o2/air mixtures.  _Energies_ , 12(23):4422, 2019.
  * Tingas et al. [2018] Efstathios Al Tingas, Dimitris M Manias, S Mani Sarathy, and Dimitris A Goussis.  Ch4/air homogeneous autoignition: A comparison of two chemical kinetics mechanisms.  _Fuel_ , 223:74–85, 2018.
  * Patsatzis and Goussis [2019] Dimitris G Patsatzis and Dimitris A Goussis.  A new michaelis-menten equation valid everywhere multi-scale dynamics prevails.  _Mathematical biosciences_ , 315:108220, 2019.
  * Patsatzis et al. [2016] Dimitris G Patsatzis, Dimitris T Maris, and Dimitris A Goussis.  Asymptotic analysis of a target-mediated drug disposition model: algorithmic and traditional approaches.  _Bulletin of mathematical biology_ , 78(6):1121–1161, 2016.
  * Patsatzis [2022] Dimitrios G Patsatzis.  Algorithmic asymptotic analysis: Extending the arsenal of cancer immunology modeling.  _Journal of Theoretical Biology_ , 534:110975, 2022.
  * Gaucel et al. [2009] Sébastien Gaucel, Béatrice Laroche, Pauline Ezanno, Elisabeta Vergu, and Suzanne Touzeau.  Using singular perturbations to reduce an epidemiological model: Application to bovine viral diarrhoea virus within-herd spread.  _Journal of theoretical biology_ , 258(3):426–436, 2009.
  * Brauer [2019] Fred Brauer.  A singular perturbation approach to epidemics of vector-transmitted diseases.  _Infectious Disease Modelling_ , 4:115–123, 2019.
  * Souza [2014] Max O Souza.  Multiscale analysis for a vector-borne epidemic model.  _Journal of mathematical biology_ , 68:1269–1293, 2014\.
  * Zhang et al. [2009] Zhonghua Zhang, Yaohong Suo, Jigen Peng, and Weihua Lin.  Singular perturbation approach to stability of a sirs epidemic system.  _Nonlinear Analysis: Real World Applications_ , 10(5):2688–2699, 2009.
  * Jardón-Kojakhmetov et al. [2021] Hildeberto Jardón-Kojakhmetov, Christian Kuehn, Andrea Pugliese, and Mattia Sensi.  A geometric analysis of the sir, sirs and sirws epidemiological models.  _Nonlinear Analysis: Real World Applications_ , 58:103220, 2021.
  * Chowell et al. [2011] Gerardo Chowell, Santiago Echevarría-Zuno, Cecile Viboud, Lone Simonsen, James Tamerius, Mark A Miller, and Víctor H Borja-Aburto.  Characterizing the epidemiology of the 2009 influenza a/h1n1 pandemic in mexico.  _PLoS Med_ , 8(5):e1000436, 2011.
  * Wallinga and Lipsitch [2007] Jacco Wallinga and Marc Lipsitch.  How generation intervals shape the relationship between growth rates and reproductive numbers.  _Proceedings of the Royal Society B: Biological Sciences_ , 274(1609):599–604, 2007.
  * Park et al. [2019] Sang Woo Park, David Champredon, Joshua S Weitz, and Jonathan Dushoff.  A practical generation-interval-based approach to inferring the strength of epidemics from their speed.  _Epidemics_ , 27:12–18, 2019.
  * Patsatzis [2021] Dimitrios G Patsatzis.  On the relation of the covid-19 reproduction number to the explosive timescales: the case of italy.  _arXiv preprint arXiv:2101.06101_ , 2021.
  * Barbarossa et al. [2020] Maria Vittoria Barbarossa, Jan Fuhrmann, Jan H Meinke, Stefan Krieg, Hridya Vinod Varma, Noemi Castelletti, and Thomas Lippert.  Modeling the spread of covid-19 in germany: Early assessment and possible scenarios.  _Plos one_ , 15(9):e0238559, 2020.
  * Valorani et al. [2020a] M Valorani, F Creta, PP Ciottoli, R Malpica Galassi, DA Goussis, HN Najm, S Paolucci, HG Im, E-A Tingas, DM Manias, et al.  Computational singular perturbation method and tangential stretching rate analysis of large scale simulations of reactive flows: Feature tracking, time scale characterization, and cause/effect identification. part 1, basic concepts.  In _Data Analysis for Direct Numerical Simulations of Turbulent Combustion_ , pages 43–64. Springer, 2020a.
  * Valorani et al. [2020b] M Valorani, F Creta, PP Ciottoli, R Malpica Galassi, DA Goussis, HN Najm, S Paolucci, HG Im, E-A Tingas, DM Manias, et al.  Computational singular perturbation method and tangential stretching rate analysis of large scale simulations of reactive flows: Feature tracking, time scale characterization, and cause/effect identification. part 2, analyses of ignition systems, laminar and turbulent flames.  In _Data Analysis for Direct Numerical Simulations of Turbulent Combustion_ , pages 65–88. Springer, 2020b.
  * Manias et al. [2022b] Dimitris M Manias, Shahid Rabbani, Dimitrios C Kyritsis, and Dimitris A Goussis.  The effect of fuel additives on the autoignition dynamics of rich methanol/air mixtures.  _Fuel_ , 323:124275, 2022b.
  * Fokas et al. [2020a] Athanasios S Fokas, Nikolaos Dikaios, and George A Kastis.  Covid-19: predictive mathematical models for the number of deaths in south korea, italy, spain, france, uk, germany, and usa.  _medRxiv_ , 2020a.
  * Fokas et al. [2020b] Athanasios S Fokas, Nikolaos Dikaios, and George A Kastis.  Predictive mathematical models for the number of individuals infected with covid-19.  _medRxiv_ , 2020b.
  * Lam and Coussis [1991] SH Lam and DA Coussis.  Conventional asymptotics and computational singular perturbation for simplified kinetics modelling.  In _Reduced kinetic mechanisms and asymptotic approximations for methane-air flames_ , pages 227–242. Springer, Berlin, Heidelberg, 1991.
  * Zagaris et al. [2004] Antonios Zagaris, Hans G Kaper, and Tasso J Kaper.  Fast and slow dynamics for the computational singular perturbation method.  _Multiscale Modeling & Simulation_, 2(4):613–638, 2004.
  * Kaper et al. [2015] Hans G Kaper, Tasso J Kaper, and Antonios Zagaris.  Geometry of the computational singular perturbation method.  _Mathematical modelling of natural phenomena_ , 10(3):16–30, 2015.
  * Tingas et al. [2016] Efstathios Al Tingas, Dimitrios C Kyritsis, and Dimitris A Goussis.  Ignition delay control of dme/air and etoh/air homogeneous autoignition with the use of various additives.  _Fuel_ , 169:15–24, 2016.
  * Goussis and Najm [2006] Dimitris A Goussis and Habib N Najm.  Model reduction and physical understanding of slowly oscillating processes: the circadian cycle.  _Multiscale Modeling & Simulation_, 5(4):1297–1332, 2006.
  * Valorani et al. [2003b] Mauro Valorani, Habib N Najm, and Dimitris A Goussis.  Csp analysis of a transient flame-vortex interaction: time scales and manifolds.  _Combustion and Flame_ , 134(1-2):35–53, 2003b.
  * Manias et al. [2019a] Dimitris M Manias, Efstathios-Al Tingas, Yuki Minamoto, and Hong G Im.  Topological and chemical characteristics of turbulent flames at mild conditions.  _Combustion and Flame_ , 208:86–98, 2019a.
  * Manias et al. [2019b] Dimitris M Manias, Efstathios-Al Tingas, Francisco E Hernández Pérez, Riccardo Malpica Galassi, Pietro Paolo Ciottoli, Mauro Valorani, and Hong G Im.  Investigation of the turbulent flame structure and topology at different karlovitz numbers using the tangential stretching rate index.  _Combustion and Flame_ , 200:155–167, 2019b.
  * Manias et al. [2019c] Dimitris M Manias, Alexandros-Efstathios Tingas, Hong G Im, and Yuki Minamoto.  Dynamics analysis of a turbulent methane flame in mild combustion conditions.  In _AIAA Scitech 2019 Forum_ , page 1731, 2019c.
  * Manias et al. [2018] Dimitris M Manias, Alexandros-Efstathios Tingas, Francisco E Hernandez Perez, Hong G Im, Riccardo Malpica Galassi, Pietro Paolo Ciottoli, and Mauro Valorani.  Analysis of hydrogen/air turbulent premixed flames at different karlovitz numbers using computational singular perturbation.  In _2018 AIAA Aerospace Sciences Meeting_ , page 0364, 2018.
  * Goussis and Skevis [2005] D.A. Goussis and G Skevis.  Nitrogen chemistry controlling steps in methane-air premixed flames.  _Computational fluid and solid mechanics_ , 1:650–653, 2005\.
  * Diamantis et al. [2015] Dimitris J Diamantis, Epaminondas Mastorakos, and Dimitris A Goussis.  H2/air autoignition: the nature and interaction of the developing explosive modes.  _Combustion Theory and Modelling_ , 19(3):382–433, 2015.
  * Manias et al. [2015] Dimitris M Manias, Dimitris J Diamantis, and Dimitris A Goussis.  Algorithmic identification of the reactions related to the initial development of the time scale that characterizes ch 4/air autoignition.  _Journal of Energy Engineering_ , 141(2):C4014015, 2015.
  * Manias and Goussis [2015] Dimitris M Manias and Dimitris A Goussis.  The initiation of methane/air autoignition: the important chemical components for various initial conditions.  In _Journal of Physics: Conference Series_ , volume 574, page 012112\. IOP Publishing, 2015.
  * Goussis [2012] Dimitris A Goussis.  Quasi steady state and partial equilibrium approximations: their relation and their validity.  _Combustion Theory and Modelling_ , 16:869–926, 2012.