1. Introduction
Infectious diseases are some of the modern world’s most important problems, e.g., influenza, the periodic outbreaks of coronavirus infections caused by COVID-19, and many other respiratory diseases have significantly affected the economies of many countries (lockdowns, masks, and reduced production efficiency) [
1,
2]. It is therefore important to minimize the economic damage, which includes both loss of work and treatment costs. Recent studies have presented many different models describing the dynamics of virus spread, which help to analyze and explain epidemic outbreaks [
3,
4]. Previously, it was shown that several infections can circulate in the population at the same time. In this case, the controlled Susceptible–Infected–Recovered model with multiple infected subgroups can be used [
5]. The solution of this model is usually found by using classical optimality criteria, e.g., the Pontryagin’s maximum principle [
6]. In this study, we focus on finding solutions that are robust to noise and take into account the dynamics of future changes in the process. We extend previous results by using a nonlinear model-predictive-control (MPC) controller to find effective controls [
7]. MPC is a computational mathematical method used in dynamically controlled systems with observations to find effective controls.
It has been shown that several virus species, or several different virus strains of the same species, can circulate simultaneously in a population. In this case, the subpopulation of infected individuals can be divided into several subgroups. In our study, two types of virus were considered to simplify the calculations, but the results obtained can be generalized to a larger number of circulating strains.
Experience with influenza epidemics and the SARS-CoV-2 pandemic has shown that there can be several waves of illness in a population during a single epidemic season. To reduce the morbidity of an epidemic wave, the government and health authorities apply various measures to prevent and contain the epidemic. The use of different preventive and curative measures during a single wave reduces the number of cases and the rate at which the virus spreads in the population. Thus, the rise of the disease in the next wave of the epidemic begins with new initial states of the system. A review of the literature led to the idea of using predictive control models to develop effective control that takes into account changes in the environment. MPC controllers allow for the resulting control to be adapted by dividing the entire time interval into smaller sub-intervals. At each interval, the control system is checked to see if the constraints are satisfied, and if not, the control action is adjusted.
The main applications of this method are in the petrochemical [
8], woodworking [
9], and energy [
10,
11] industries, where, among other things, the stability of the obtained controllers and the possibility of their automatic adjustment are important. The MPC controller finds a solution in the form of piecemeal functions on the control interval, but the optimization problem is solved over a larger time interval—the forecast horizon. The MPC controller finds the solution in the form of piecewise constant functions on the control interval, but the optimization problem is solved on the forecast horizon.
As the main contribution of this work, we formulated a control problem and solved it using an MPC controller, a series of numerical experiments are carried out to confirm the effectiveness of the solutions obtained by the MPC controller for the SIIR model.
The paper is organized as follows.
Section 2 presents the epidemic models and modified infection rate parameter.
Section 3 shows the application of MPC-controllers to a SIIR epidemic model. In
Section 4, we use numerical simulation to illustrate our results. The paper is concluded in
Section 5.
2. Epidemic Process in Total Urban Population
Our model is based on the extended Susceptible–Infected–Recovered model with two types of viruses, which circulate simultaneously in an urban population [
6]. A size of population is equal to
N. According to these assumptions, we can consider four groups: the Susceptible, the Infected by virus
, the Infected by virus
, and the Recovered group.
Susceptible is a group of people who are not infected by viruses, but can be infected by one or both types of viruses and who have no immunity to the viruses. Depending on the viral strength, we observe that the number of people infected by virus 1 or by virus 2 can be different, and people infected by virus
or
belong to the subgroup
Infected. The
Recovered subgroup consists of people who have recovered from the infection. The mixing of populations allows for viruses to spread quickly, and each person in the population is assumed to have an equal probability of coming into contact with others. Therefore, when an infected person interacts with a susceptible person, the virus can spread. A virus with higher virulence, we assume, will have a higher probability of success in spreading when an encounter occurs between the infected and the susceptible. The scheme of process is represented in the
Figure 1.
We model the spread of the virus in the population based on the classical SIR epidemiological model, where a system of differential equations is used to describe the proportion of the population as a function of time. Then, at time t, , , , and correspond to the fractions of the population that are Susceptible, Infected by virus , Infected by virus , and Recovered, respectively, and for all t the condition is justified. Define as proportions of Susceptible, Infected, and Recovered. At the beginning of the epidemic , most people in the population belong to the susceptible subpopulation, a small group in the total population is infected and other people are in the recovered subpopulation. Therefore, initial states are , , , .
We have extended the simple SIR model introduced by [
12] to describe the situation with two virus types:
where
are infection rates for virus
,
,
are recovery rates (see
Figure 1). Infection rate is defined as a product of the contact rate
l and transmissibility of infection, i.e., probability of transmission infection during the contact,
Infection rate is integrated into the evolution of mutation process to epidemics in the urban population.
Medical treatment or quarantine reduces the number of the infected individuals in the urban population. These prevention measures can be interpreted as control parameters in the system defined as ; here, is a fraction of the infected which are quarantined or under intensive medical treatment. Recovered rates are inversely proportional to disease duration , hence .
The objective function: We will minimize the overall cost in time interval
[
13,
14,
15,
16]. At any given
t following costs exist in the system:
is a function measuring the direct costs of infection, with a non-decreasing, twice differentiable function, such that , for . These costs can be interpreted as the price the government or medical organization pays for the health system to treat the infected and increase the probability of recovery for the treated individuals. Infected individuals who are not detected because they are asymptomatic are not treated, and therefore there are no direct costs for them.
The treatment or medical cure of those infected produces costs for society. The functions represents the costs of treating the infected to keep them out of the labor market until recovery. The function is increasing on the argument , is twice differentiable, and is such that .
The function is a benefit rate. These are the bonuses that society receives from Recovered individuals who have returned to work and started contributing to society. The function is non-decreasing and differentiable and .
Constants , , represent the costs and benefit for Infected and Recovered in the end of the epidemic.
Aggregated system costs is given by
and the optimal control problem is to minimize the cost, i.e.,
To simplify the analysis, we consider the case where .
Here, we assume that as a control parameters which help to reduce numbers of infected in human population, medical treatment or isolation can be considered, then define variable as control , for all t.
2.1. SIIRS Models
The SIIR models discussed in the previous section are best suited to describe epidemic processes in which viruses mutate slowly. Because of the small number of mutations, individuals can develop lifelong immunity to the virus. However, there are viruses that mutate rapidly and for these, SIIRS (Susceptible–Infected–Recovered–Susceptible) models have been developed (see
Figure 2) to account for the decline of immunity in the population over time.
The system of differential equations describing the dynamics of state change is equal to
Here, as in the SIIR model, the parameters
are responsible for the intensity of infection and cure,
—the proportion of infected people on medication, and
—the intensity of immune decline in the population. The objective function for such tasks is no different from (
2).
Note that due to the gradual decline of immunity, the epidemic process in the SIIRS problem can continue indefinitely. For this reason, when solving this problems we should choose large values for T. However, in such a case there is a risk (due to the high mutation rate of viruses) that the model parameters will no longer be relevant. Therefore, methods that can adapt to changes in the process by predicting its behavior can be effective in solving such problems.
2.2. SWIIRS Model
In 2020, the SWIIRS (Susceptible–Warned–Infected–Recovered–Susceptible) model was proposed, which takes into account citizens’ awareness of the epidemic [
17]. This model divides the population into five groups instead of four. An additional group of individuals is the group of warned
W individuals who are aware of the virus (see
Figure 3). This approach allows a more accurate description of epidemic processes in large countries, as the speed of spread of the virus often depends not only on the characteristics of the virus itself, but also on the behavior of members of the population.
The changes in population groups for the controlled SWIIRS model are described by a system of five differential equations:
In the system (
4), the constants
and
describe the infection rates for groups
S and
W, respectively, and
– the cure rate
. The value
is responsible for the probability that a knowledgeable individual will acquire immunity to viruses (e.g., go to a health facility to get a vaccine). The intensity of the decline in immunity is described by
. Two parameters are responsible for informing (warning) members of the population:
Parameter —percentage of susceptible individuals informed by the media (become Warned);
Parameter —probability of informing the susceptible individual when communicating with the warned individual.
2.3. Controllable Infection Rate Parameter
Consider a model (
1) with one virus
and let
be the probability of virus transmission
from a diseased individual to a susceptible one, and
l be the average number of contacts of an individual (with other members of the population) per unit of time. Consider a random contact of a susceptible individual. Infection with the virus will only occur if two events occur simultaneously: the contact was with the infected individual, and the virus was transmitted on contact.
Define two incompatible events:
Event A— a susceptible individual becomes infected through one of an l contacts;
Event B—with a similar number of contacts, infection will not occur.
Obviously that events
A and
B form a complete group of events (see
Figure 4), then
The event
B will only occur if there is no infection at any of the
l contacts. In other words
From Equations (
5) and (
6) we see that the probability of transmission to a susceptible individual by
l contact with infected individuals can be defined as
The same formula can be obtained in another way. The probability of the event
A is given by
This probability is called the modified infection rate parameter and it is denoted by .
Let us redefine
(intensive virus treatment) as
. Let us introduce two types of control into the resulting system: isolation of infected individuals and reduction in the number of contacts in the whole system. Let
be the fraction of infected individuals isolated, and
be the fraction of contacts prevented by reducing the mobility of individuals in the whole population. In this case we obtain
where
. Then the system of differential equations of the corresponding single-variable SIR problem is
If we perform similar transformations to (
8) for the SIIR dual virus model, the calculated formulae for the controllable infection parameters are as follows
where
. If we substitute these parameters into the system (
1), the differential equations take the following form:
For the introduced controls we also need to redefine the cost functions. Let
,
,
be non-reducing functions corresponding to the costs of the controls
,
, and
. Then the function (
2) takes the form
4. Numerical Simulations
This section presents the results of a numerical simulation that corroborates our theoretical results. The numerical experiments are aimed at obtaining effective treatment strategies using nonlinear MPC controllers. We selected MATLAB as our platform for programming and numerical calculations for a series of experiments. MATLAB offers an extensive library of the model-predictive-control method, including various examples. Nonlinear-MPC controllers were employed in our study. These controllers utilize the linearization algorithm described in the article to find solutions. Nonlinear MPC is adequate to address these issues. However, larger systems may demand more precise predictive models based on machine learning and neural networks. For first three experiments, two variants of the SARS-CoV-2 virus were selected for study: Alpha [
23] and Omicron [
24]. The epidemic situation was considered for a period of 180 days (6 months). During this period, six MPC controllers were used sequentially with a forecast horizon of 30 days (1 month) and a control horizon of 21 days (3 weeks). Time dependence of MPC controller accuracy and performance on
k parameter are shown in
Table 1. The values of the main system parameters are shown in
Table 2.
The serial application of multiple nonlinear MPC controllers allows the constraints, system parameters and objective function coefficients to be adjusted step by step. To simplify the calculations, we will use the same parameters in each iteration of the process.
Experiment 1. In the first experiment, we considered the epidemic process with a small number of contacts, which is possible in cases where most of the infected comply with the quarantine. It can be seen that applying control is having a positive impact on stopping the virus spreading (
Figure 7). On the effective control graph (
Figure 8), you can notice “jumps” due to the fact that a new MPC controller is connected each new month (30th day). However, the resulting control greatly reduces the damage from the epidemic compared to a situation where we do not control the epidemic process. So the cost in the uncontrolled case is
monetary units (m.u.) and with the effective strategy
m.u. The benefit of using control is therefore 2143 m.u.
Experiment 2. A second experiment (
Figure 9 and
Figure 10) was run with the average number of contacts per day
. This experiment is the most different from the others. The main wave of the epidemic does not occur in the first few days, but only at the beginning of the third month. The highest costs are in the same period. After
, the effective control is similar to that of experiments 1 and 3. With this strategy, it is possible to reduce the cost of the epidemic from
m.u. to
m.u., which saves 2568 monetary units.
It can also be seen that as the average number of contacts increased, so did the values of the objective functions. Just as importantly, the benefits of using MPC controllers have also increased.
Experiment 3. In the third simulation the number of contacts is
, which approximates the real situation.
Figure 11 shows that with a large number of contacts, almost all agents are infected by
. However, when efficient control is used, the infection peaks are significantly reduced (from
to
and from
to
). As in the previous experiments, the use of efficient control is economically beneficial, reducing costs from
m.u. to
m.u. The benefit of using an effective strategy is 3422 m.u. (see
Figure 12).
It should also be noted that despite the largest value of k the costs of using effective management in Experiment 3 turned out to be the smallest. This may be due to the fact that when solving the problem, they are linearized and the accuracy of the solution can vary from the parameters of the problem. At the same time, all three experiments showed that effective control can reduce the economic damage from the epidemic, which shows the practical importance of using the model predictive control in SIIR problem.
Experiment 4. In the fourth experiment we examined the SIIR model with the data taken from the paper “Optimal control of heterogeneous mutating viruses” [
5].
Table 3 shows the parameter values used in the simulation.
Figure 13 shows the population states in the uncontrolled and controlled cases. It can be seen that the use of the control found with MPC significantly reduces the number of infections. In the uncontrolled case, the maximum proportion of the number of infections is
,
, and in the controlled case,
,
, corresponding to the initial proportion of infected in the population.
After the computer simulation, the solution of the formulated problem was obtained and the values of the functional were found (see
Figure 14).
The strategies that increased to
and
during the first three days of the epidemic and then decreased to zero were found to be effective. The resulting control signals are similar in appearance to the continuous controls presented in the experiments of this paper. The aggregated costs are equal to
m.u. in the controlled case and
m.u. in the uncontrolled case. Thus, the benefit of using the control is
monetary units, which is slightly less than using the optimal control
m.u. In
Figure 15, we can see the structure of the optimal control for the model and the data from this experiment. It can be observed that the control structure in the case of effective control is similar to the structure in the case of optimal control.
Experiment 4 showed that the developed software package finds near-optimal solutions to classical SIIR problems and can be used as an alternative to classical optimal control methods.
Experiment 5. As the software implementation of the MPC method showed good accuracy for standard tasks, Experiment 5 was conducted using a modified SIIRS model (with controlled infection rate parameter). Data for the experiment were obtained from sources such as the World Health Organization, Our World in Data, and Worldometer. The parameters of the experiment are shown in
Table 4.
Using effective strategies (see
Figure 16), it was possible to reduce the maximum proportion infected with
virus from
to
, but the proportion infected with
virus increased from
to
, which is due to the fact that the proportion of people (30% of the population) who were not infected with
virus (due to the measures introduced) managed to become infected with
virus. It can also be seen that from
, the proportion infected with both viruses does not exceed 1%.
By running the computer simulation program, effective controls for the formulated problem were obtained (see
Figure 17). It can be observed that the treatment of individuals for both viruses takes a maximum value starting from
and decreases for virus
after
and for virus
after
. For
, a small fraction (
) of those infected during
are isolated, and after
the control
turns on again, reaching a maximum value, and turns off after
. For
, maximum isolation of the infected occurs at
,
,
, and partial isolation at
,
. Reducing the total number of contacts was most effective in periods
.
For problems with a modified infection rate parameter, no optimality conditions have yet been found, so the only way to assess the quality of the solution is the cost savings provided by the constructed efficient control. For this reason, the values of the target functional and the total control for the two viruses were calculated for this experiment (see
Figure 18). With the combined effect of all controls it was possible to significantly reduce the value of the target function from 33,286 to 11,180 m.u. The benefit of using effective controls was therefore 22,106 m.u.
Experiment 6. The SWIIRS model was used for Experiment 6, and the main objective was to construct long term strategies (
) using multiple MPC controllers. The experimental parameter values are available in
Table 5.
Figure 19 and
Figure 20 show the results of Experiment 6. The effective control plots show that when several MPC controllers are connected, there are pulses (switching points) in the control. They are related to the fact that the prediction model is recalculated to a new prediction interval. You may notice that from the third controller (
) the constructed controllers have a similar structure. This is reflected both in the behavior of the process and in the behavior of the total cost function (they become periodic).
By applying efficient management, the value of the target function was significantly reduced from to m.u. Thus, the benefit of applying management was 2748 monetary units. It can be seen that model-predictive-control solutions remain effective over long time intervals, but require the connection of several MPC controllers.