Ovil Aviation High Technologies
Vol. 20, No. 02, 2017
УДК 519.626.2
PROBLEM OF OPTIMAL CONTROL OF EPIDEMIC IN VIEW OF LATENT PERIOD
N.I. OVSYANNIKOVA1
1Moscow State Technical University of Civil Aviation, Moscow, Russia
The problem of optimal control of epidemic through vaccination and isolation, taking into account latent period is considered. The target function is minimized - functionality summarizing costs on epidemic prevention and treatment and also considering expenses on infected people left at the end of control T who may be a new source of epidemic. On the left endpoint of the integration segment initial data is given - quantity of infected and confirmed people at the moment t, the right endpoint is free. The dynamic constraints are written by way of a system of simple differential equations describing the speed of changes of number of subjected to infection and number of already infected. Besides the inhomogeneous community is considered, consisting of four age groups (babies, preschool children, school children and adults). The speed of vaccination (number of vaccinated per a time unit) and isolation speed are used as the control functions. There are some restrictions on control above and below. The latent period is described by the constant h and is part of the equation describing the contamination speed of people as a retarding in argument t, i.e. a person being in a latent period infects others not being aware of his disease. For problem solving Pontryagin maximum principle is used where it can be seen that the control is piecewise constant. The result of numerical implementation of discrete problem of optimal control is given. The conclusions are made that the latent period significantly influence the incidence rate and as consequence the costs on epidemic suppression. The programme based on the programming language Delphi gives an opportunity to estimate the scale of epidemic at different initial data and restrictions on control as well as to find an optimal control minimizing costs on epimedic suppression.
Key words: optimal control of epidemic, latent period, vaccination and isolation, minimization of epidemic elimination costs.
INTRODUCTION
Problems of optimal control of systems with delayed argument have extensive applications in engineering, economics, medicine, automatic control theory, theory of self-oscillatory systems and other sciences [1, p.78]. The development of approximate and numerical methods for solving such problems is a pressing issue. In this paper, taking into account latent period, we consider the problem of optimal control of the epidemic, which is described by a system of ordinary differential equations with delayed argument. The author constructed a model of the epidemic with a constant delay on the basis of Cermak-McKendrick model, formulated necessary conditions for optimality of the Pontryagin maximum principle, developed a numerical method for solving the problem for an inhomogeneous community, consisting of four groups, and carried out an analysis of the impact of the delay on optimal control.
RESEARCH METHODS AND METHODOLOGY
Let us consider the problem of optimal control of epidemic through vaccination and isolation in an inhomogeneous community, consisting of n age groups, taking into account the latent period h. The dynamics of the epidemic is described by the following system of differential equations with delay:
X (t) = - x (t )£ ß^ (t) - (t) - v, (t) + Л г,
j=i
i = 1, n, t e [0;T],
у, (t) = x, (t - hj£ ßjy} (t - h) - yt (t) - yt (t) - yy (t) - щ (t), i = 1, n, t e [0;T],
j=i
with the initial conditions on the interval of delay T0
<
Том 20, № 02, 2017_Научный Вестник МГТУ ГА
Vol. 20, No. 02, 2017 Civil Aviation High Technologies
X(t) = xi0(t), yt(t) = yi0(t), i = \~n, t e To = [—h;0], (2)
where x (t) - number of people exposed to infection in the i-group at time t (t = 1, n ), yi (t) - number of people infected in the i-group, Xt (t) - rate of change in the number of people exposed to infection in the
n
i-group at time t, уг (t) - rate of change of the number of infected people in the i-group, xt (t)^fitjyj (t) -
j=i
rate of infection of healthy people out of the i-group from the infected people of the j-group at time t, given the fact that the infection could occur from an infected person from any of the j-group (j = 1, n), yiy (t) -number of patients in the i-group, regained their health without the action of external agents: quarantine, vaccination and so forth - average time of natural recovery for the given infectious disease), - growth coefficient, characterizing the frequency of meetings of healthy people out of the i-group with infected people out of the j-group and the probability of infection at the meeting, jut - natural mortality coefficient of people in the i-group, Д. - mortality coefficient from the given infection in the i-group, Лг. - average birth rate in the i-th group, h - constant quantity, characterizing the latent time of the disease. Let us denote vt (t) - population vaccination rate among those who are exposed to infection in i-group at time t, ut (t) - isolation rate of the infected in i-group at time t. Restrictions on the control functions are set in the following form:
0 < v < M, t = in (3)
0 < u < N, t = in (4)
where M , N - maximum rates of vaccination and isolation (are limited by technical and material
possibilities) in the i-group. Functional:
T n n
J(u,v) = (Ay, (t) + Divi (t) ^Cu(t))dt + ^Biyi (T) ^ inf (5)
0 i=1 i=1
characterizes the aim of control, which is to minimize the cost of infection elimination, where A - the average cost per patient for the society per unit of time (known quantity, for Russia it is about $ 50 [4, p. 7]), D - the cost of vaccination per person in the i-group, C - the cost of isolation per person in the i-group. The latter summand refers to the cost of residual patients, which may cause secondary infection, so their cost Bt (i = 1; n) should be large (as a penalty for undertreated patients). Unfortunately, in the latter summand we can not take into account the ones who are in incubation period, as they have not moved in the patients' group yet.
If we take A, the average cost of a patient for the society per a unit of time, equal to one standard monetary unit, then (5) can be rewritten as following:
T n n
J(u, v) = (y (t) + dv (t) +CiUi (t))dt + £ by (T) ^ inf (6)
0 i=1 i=1
where d - relative cost of vaccination per person in the i-group, C - relative cost of isolation of the patient in the i-group, b - relative cost of one undertreated patient
Ovil Aviation High Technologies
Vol. 20, No. 02, 2017
in i-group at time T (and b > 1, since each undertreated patient can infect more people in the future).
Necessary conditions of optimality: let's construct Pontryagin's function:
H(x,y,g,z,v,u,p,q) = -Л0£(yt (t) +dtvt(t) + c,u,(t)) + £(-x, (t)£ßl}y} (t) - M,x,(t) + Лг - v,(t)) +
i=l
i=1
j=1
£qt (t )(g, (t )£ ßjZj (t ) -Mtyt (t ) - My (t ) - y y (t ) - u (t)),
i=1
j=1
whiere gt (t) = x, (t - h), zt (t) = y, (t - h), i = 1; n . Switching functions:
(p.(t) = -A0dt -pt (t\ (t) = -A0ct - qt (t\ i = \;n.
On the basis of the theorem about the necessary conditions of optimality in systems with after effects, formulated and proved in [2], [3], [6], we can write the necessary conditions of optimality for the previously constructed model.
If permissible process w = (X,y,v,u) for each t e [0;T] is optimal in the problem (1)-(4), (6),
then there exist not all at a time equal to zero number A0 and vector functions p (t), q (t) such that
the optimal control is determined by the following conditions:
V = 0, if (-¿0 dt - pt (t)) < 0, i =
v = Mi, if (-Äo d - p, (t)) > 0, i = l;n_
0 < v, < Ml, if (-¿0 d, - p, (t)) = 0, i = l; n
v = 0, if (-¿0 d% - p, (t)) < 0, i = 1n
v% = M,, if (-¿0d% -p, (t)) > 0, i = Vn_
0 < v < Ml, if (-¿0 d, - p, (t)) = 0, i = V in
(7)
(8)
and conjugate functions satisfy the system of differential equations (here and below A0 = 1):
да
dx,
SH
dgi
pi (£ j + Mi)
t +h
j=1
q £ ß
vz1
j=1
t+h
pt (t)(£ßl]yi (t) + M ) +q (t + h)£ßjz3 (t + h), i = 1, n
l=l
j=1
q, =
ж
fyl
<9H
ôz,
1+(£ ßnxiPi +4i (Pi + ßi+ÏI ))
t +h
£ßug14r
' \t+h
= l, n
= \ + £p, (t)ßüXi (0 - £ Чг (' + h)ßax,. (0 + q, (tXu +fi,+r,), 7 = 1.
i=1
p, (T) = 0; qi (T) = -b, ; p, (t) = 0, qt (t) = 0, t > T.
(9)
(10) - transversality conditions at the right end of integration.
Vol. 20, No. 02, 2017
Civil Aviation High Technologies
Numerical realization of the task: We divide the interval of integration [0; T] into q equal subintervals by the points 0 = t0 < tt <... < t = T so that the length of each subinterval is At = ti+l - tt, i = 0, q -1.
Let's divide delay interval [-h; 0] with a step At, rounding the obtained result to an integer m > 0.
Integral I (v, u) may be found by a numerical method of rectangles, the error is estimated by the formula:
R < max I f( ,
[ab]
where f(x) - integrand, [a; b] - the integration interval, n is the number of elementary partition segments of [a; b] [9, p. 99], and differential equations - by the Euler method, the error is estimated by the formula |;y -y(x)| < Cxh (eCl(Xi-X0) -1), where \y -y(x)| - deviation of the approximate solution from the exact one, Ci, C2 - positive constants determined by the right side of the equation y = f (x, y) and its derivatives in the neighborhood of the point (x, y), h = const - step of the numerical method, from h ^ 0 should be uniform aspiration [8 , p. 37].
Then the problem (1)-(5) is approximated by the following discrete problem of optimal control with accuracy O(At):
I (v, u) =
n q-1
(yj + d,v', + c,u', )At + Yb,yq, ^inf
Z—i Z—i v-r j
j=1 i=0
J J J J
jj
j=1
(11)
xj = xj xj-
JM^ßjky'k-At(vj + j-Aj), j = 1;n, i = 0; q -1,
к=1
yf = yj + xjm A^ß^T-Ai((/j +pj ) yJ + uJ ), j = 1; n, i = 0; q -1,
к=1
(12)
= ÇJ, y,0 = Ç), j = 1, n, i = -m,0
(13)
in the delay interval To
0 < v\ < A ., 0 < u'. < B , j = 1,n, i = 0, q -1.
(14)
The Lagrangian function:
L = Л
n q-1
H( yj + J + CjuJ )At + X b
J=1i=0 J =1
_J=1'=0
/ V \ n q_1 n , V
At (v+Jj-AJ ))+££ ; [ у;1 - y - S-mAtx ß]kyk-m+At (rj ) yJ +
" j=11=0 v к=1
y
jj
+ 11 j fXT - X1; + X1; At± ß]kyk +
j =1 i =0
J J J
V k=1
q-1
j=11=0 v
u
(15)
We calculate the derivatives of the Lagrangian function, write down stationary conditions in the problem (11)—(14):
n
n
<
n
Научный Вестник МГТУ ГА_Том 20, № 02, 2017
Ovil Aviation High Technologies Vol. 20, No. 02, 2017
ÔT n n
— = pj - j + PГ* £р]кук+- qrm A Xj = 0
rxj k=1 k=1
r)T " "
= t^pfß^ +ч'-чТ- +(rj+fij+ßj) = o
ryj k=1 k=1
or
Pj = P,T - Pl!l^t£ß]kyk - Pfbtp, T j+m At^ßjky'k, j = 1,n, i = о, q -1,
k=1 k=\
qj = j -At -At X P^ßjXk + At X qk+l+mßjXk -Mjifj T ^ ), j = 1, n, i = 0, q -1.
(16)
= Lnj =
' rvj rv J V, j . j • J ' ' " 77
k=1 k=1
q ÔL n q ÔL . 1
Pq =^ =0' qq =ryq = -Wj if l>q Pj= 0qj=0 (17)
j ^j
Partial derivatives of the Lagrangian function with respect to the control variable are the following:
p) = ^T = (^ dJ+pj1 )&t =^T = (Ao cj+q,)1 )&t
dv) v J J ) fa) \ J J ) (18)
j = I,n, i = 0,q-I.
To build the solution to the discrete problem of optimal control (11)—(18) we use the gradient projection method [10].
RESULTS OF THE RESEARCH
To carry out the experiment on a real model there were used statistical data on influenza epidemic in the city of Arkhangelsk. The data were provided by the Territorial Administration of the epi-demiological surveillance in the city of Arkhangelsk over the past twenty years. On the basis of these data there were allocated four age groups:
I group - children from 0 to 2 years old,
II group - children from 3 to 6 years old,
III group - children from 7 to 14 years old,
IV group - people over 15 years old.
Mortality coefficients (independent of the flu) for each group, calculated in accordance with statistical data for Arkhangelsk: p = 0.008, p2 = 0.001, p = 0.001, pA = 0.028. p. = 0 (j = 14) -
take death coefficients, where death occurred because of the flu, zero, as there were no registered deaths caused by the flu over the past twenty years in Arkhangelsk. Take one week per a unit of time. The average birth rate in the 1st group is 35 people a week, in the rest groups A. = 0 (j = 2;4). y. -
natural recovery coefficient, it can be taken equal to one, as y- - the average time of natural recovery coefficient, as for flu can be taken equal to one week. The coefficients ¡3jk are found by solving the inverse problem:
Vol. 20, No. 02, 2017
Civil Aviation High Technologies
#4 = 4,22 • 10 5; #u * # * P13 * 4,22 • 10 6; f322 = 6,50 • 10 5; #21 * #23 * #24 * 6,50 • 10 6; #33 = 2,41 • 10 5; #31 * #32 * #34 * 2,41 • 10 6; #44 = 3,07 • 10-6; #41 * #42 * #43 * 3,07 • 10-7;
The latent period of the flu disease is from 1 to 10 days (usually 3-5 days) [4, p. 2]. Let us consider the time period of 10 weeks (T = 10). The relative vaccination cost is d = 0,01 in all the groups. The relative isolation cost is c = 3 in all the groups. The calculations were made in Delphi. In case when b1 = b2 = b3 = b4 = 0 the solution of the problem is shown on the following graphs (the epidemic dynamics and optimal control for the first, second and third groups are similar, therefore the graphics are shown only for the third and fourth groups).
Fig. 1. Dynamics of the infected y3 (a), isolation control u3 (b), vaccination control v3 (c) depending on h
Fig. 2. Dynamics of the infected y4 (a), isolation control u4 (b), vaccination control v4 (c) depending on h
Analysis of the results shows that with the growth of h increases the number of infected people yj {j = 1,4) on the interval [0; T], decreases duration of vaccination control as well as duration of isola-
Ovil Aviation High Technologies
Vol. 20, No. 02, 2017
tion control, although not all the patients are isolated and cured. With the growth of h increase the total costs of infection elimination.
Now take bj ^ 0 (j = 1,2,3,4), i.e. introduce a fine for undertreated patients. Let bj = 5
(j = 1, 2, 3, 4).
Fig. 3. Dynamics of the infected y3(a) and y4 (b) under constant isolation control u3 = u4 = 100 (pers./week) and dynamics of the infected y4 (c) under constant isolation control u4 = 500 (pers./week) depending on h
So now it is obvious that it is more profitable to heal the patients than to leave them undertreat-ed, as they may cause secondary infection. Therefore, control should be of maximum rate and duration. Control in the first three groups is enough to eliminate the infection in the period under review (fig. 3a), in the fourth group - not enough (fig. 3b). To eliminate the infection in the fourth group, you have to either strengthen control or increase its duration. Increase in the isolation rate of patients up to five times gives the desired result - an epidemic in the fourth group will be eliminated within the required time frame (fig. 3c).
DISCUSSION OF OBTAINED RESULTS AND CONCLUSION
We have considered the problem of optimal control of epidemic through vaccination and isolation in an inhomogeneous community, consisting of four age groups, taking into account latent period. The aim of the control is to minimize the costs to fight the epidemic at existing control restrictions. The objective of numerical research is to reveal the impact of the latent period duration on program control. To carry out numerical experiment on a real model there were used statistical data on influenza epidemic in the city of Arkhangelsk. The obtained result shows that with the growth of latent period increases the number of the infected in all the groups, decreases duration of vaccination control as well as duration of isolation control, although not all the patients are isolated and cured. When introducing a fine for undertreated patients, control maximizes with respect to rate and duration, as it is more profitable to cure the patients rather than fight off the secondary infection.
REFERENCES
1. Andreeva E.A. Optimal'noe upravlenie dinamicheskimi sistemami [Optimal control of dynamic systems]. Tver, TvGU Publ., 1999, pp. 72-120. (in Russian)
Vol. 20, No. 02, 2017
Ovil Aviation High Technologies
2. Andreeva E.A., Kolmanovskij V.B., Shejhet L.E. Upravlenie sistemami s posledejstviem [Control systems with aftereffect]. Moscow, Nauka Publ, 1992. (in Russian)
3. Haratishvili G.L. Princip maksimuma v teorii optimal'nyh processov s zapazdyvaninem [The maximum principle in the theory of optimal processes with retardation]. Report of the USSR Academy of Sciences. 1961, vol. 136, no. 1, pp. 39-41. (in Russian)
4. Informacionnyj bjulleten' «Vakcinacija. Novosti vakcinoprofilaktiki» [Information bulletin "Vaccination. News of vaccinal prevention»]. Moscow, May/June 2003, vol. 3 (27).
5. Andreeva E.A., Semykina N.A. Optimal'noe upravlenie [Optimal control]. Tver, Tver branch of MESI, 2006, pp. 184-211. (in Russian)
6. Kolmanovskij V.B. Ob approksimacii linejnyh upravljaemyh sistem s posledejstviem [About approximation of linear control systems with delay]. Problemy upravlenija i teorii informacii [Problems of control and information theory]. 1974, vol. 3, no.1. (in Russian)
7. Andreeva E.A., Ciruleva V.M. Variacionnoe ischislenie i metody optimizacii [The calculus of variations and optimization methods]. Orenburg-Tver, 2005. (in Russian)
8. Mihlin S.G., Smolickij H.L. Priblizhjonnye metody reshenija differencial'nyh uravnenij [Approximate methods of solving differential equations]. Reference mathematical library edited by Lyusternik L.A. and Yanpolskii A.R. Moscow, Science Publ., Main edition of Physical and Mathematical Literature, 1965. (in Russian)
9. Bahvalov N.S. Chislennye metody, t. 1 [Numerical methods, vol. 1]. Moscow, Nauka Publ., Main edition of Physical and Mathematical Literature, 1975.
10. Evtushenko Ju.G. Metody reshenija jekstremal'nyh zadach i ih primenenie v sistemah op-timizacii [Methods of solution of extremal problems and their application in optimization systems]. Moscow, Nauka Publ., 1982. (in Russian)
INFORMATION ABOUT THE AUTHOR
Natalia I. Ovsyannikova, PhD in Physical and Mathematical Sciences (Candidate of Physical and Mathematical Sciences), Associate Professor of the Department of Applied Mathematics of Moscow State Technical University of Civil Aviation (MSTU CA), natmat68@mail.ru.
ЗАДАЧА ОПТИМАЛЬНОГО УПРАВЛЕНИЯ ЭПИДЕМИЕЙ С УЧЕТОМ ЛАТЕНТНОГО ПЕРИОДА
Н.И. Овсянникова1
1 Московский государственный технический университет гражданской авиации,
г. Москва, Россия
Рассматривается задача оптимального управления эпидемией путём вакцинации и изоляции с учётом латентного периода. Минимизируется целевая функция - функционал, суммирующий затраты на лечение и профилактику эпидемии, а также учитывающий стоимость инфицированных людей, оставшихся на момент окончания управления Т, которые могут явиться источником новой эпидемии. На левом конце отрезка интегрирования заданы начальные условия - количество инфицированных и подверженных заражению в момент времени t, правый конец - свободный. Динамические ограничения записаны в виде системы обыкновенных дифференциальных уравнений, описывающих скорость изменения числа подверженных заражению и числа уже зараженных. Причем рассматривается неоднородное общество, состоящее из четырех возрастных групп (младенцы, дошкольники, школьники, взрослые). В качестве управляющих функций взяты скорость вакцинации (число вакцинированных в единицу времени) и скорость изоляции. Имеются ограничения на управление сверху и снизу. Латентный период описывается константой h, и входит в уравнение, описывающее скорость инфицирования людей как запаздывание в аргументе t, то есть человек, находящийся в латентном периоде, заражает окружающих, не зная, что он уже болен. Для решения задачи записывается Принцип максимума Понтрягина, откуда видно, что управление является кусочно-постоянным. В работе приводится результат численной реализации дискретной задачи оптимального управления, сделаны выводы о том, что латентный период существенно влияет на рост заболеваемости и, как следствие, расхо-
Civil Aviation High Technologies
Vol. 20, No. 02, 2017
дов на погашение эпидемии. Программа, написанная на языке программирования Delphi, дает возможность оценить масштабы эпидемии при различных начальных данных и ограничениях на управление, а также найти оптимальное управление, минимизирующее расходы на погашение эпидемии.
Ключевые слова: оптимальное управление эпидемией, латентный период, вакцинация и изоляция, минимизация затрат на погашение эпидемии.
СПИСОК ЛИТЕРАТУРЫ
1. Андреева E.A. Оптимальное управление динамическими системами. Тверь: ТвГУ, 1999. С. 72-120.
2. Андреева E.A., Кялмановский В.Б., Шейхет Л.Е. Управление системами с последействием. Москва: Наука, 1992.
3. Харатишвили Г.Л. Принцип максимума в теории оптимальных процессов с запаздыванием // Доклады АН СССР. 1961. Т. 136, № 1. С. 39-41.
4. Информационный бюллетень «Вакцинация. Новости вакцинопрофилактики». № 3 (27). май/июнь 2003. Москва.
5. Андреева Е.А., Semykina N.A. Optimal'noe upravlenie [Optimal control]. Tver, Tver branch of MESI, 2006. pp.184-211. (in Russian)
6. Колмановский В.Б. Об аппроксимации линейных управляемых систем с последействием// Проблемы управления и теории информации. 1974. Т. 3, № 1.
7. Aндреева E.A., Семыкина Н.А. Оптимальное управление. Тверь: ТФ МЭСИ, 2006. С. 184-211.
8. Михлин С.Г., Смолицкий Х.Л. Приближённые методы решения дифференциальных уравнений // Справочная математическая библиотека / под общ. ред. Л.А. Люстерника и А.Р. Янпольского. М.: Наука, Главная редакция физико-математической литературы, 1965.
9. Бахвалов Н.С. Численные методы. Т. 1. М.: Наука, Главная редакция физико-математической литературы, 1975.
10. Евтушенко Ю.Г. Методы решения экстремальных задач и их применение в системах оптимизации. М.: Наука, 1982.
СВЕДЕНИЯ ОБ АВТОРЕ
Овсянникова Наталья Игоревна, кандидат физико-математических наук, доцент, доцент кафедры прикладной математики Московского государственного технического университета гражданской авиации (МГТУ ГА), natmat68@mail.ru.