ГИБРИДНАЯ SIRS-МОДЕЛЬ РАСПРОСТРАНЕНИЯ ИНФЕКЦИЙ

Известия высших учебных заведений. Прикладная нелинейная динамика. 2022. Т. 30, № 6 Izvestiya Vysshikh Uchebnykh Zavedeniy. Applied Nonlinear Dynamics. 2022;30

Научная статья DOI: 10.18500/0869-6632-003014

УДК 517.9, 621.372 EDN: EJQVCA

Гибридная SIRS-модель распространения инфекций

А. В. Шабунин

Саратовский национальный исследовательский государственный университет имени Н. Г. Чернышевского, Россия E-mail: shabuninav@info.sgu.ru Поступила в редакцию 29.05.2022, принята к публикации 25.07.2022, опубликована онлайн 31.10.2022, опубликована 30.11.2022

Аннотация. Цель — построение модели распространения инфекции в виде системы дифференциальных уравнений, учитывающей инерционный характер передачи инфекции между особями. Методы. В работе проводится теоретическое и численное исследование устройства фазового пространства системы обыкновенных дифференциальных уравнений модели среднего поля. Результаты. Построена модифицированная SIRS-модель распространения эпидемий в виде системы обыкновенных дифференциальных уравнений третьего порядка. От стандартных моделей она отличается учетом инерционного характера процесса передачи инфекции между особями популяции, что реализуется посредством введения в модель «агента-переносчика». В модели не учитывается влияние заболевания на численность популяции, при этом плотность населения рассматривается как параметр, влияющий на ход эпидемии. Динамика модели демонстрирует хорошее качественное соответствие с рядом наблюдаемых при развитии заболеваний явлений. Обсуждение. Предложенное усложнение стандартной SIRS-модели посредством добавления в него уравнения для динамики возбудителя инфекции предоставляет перспективы для ее уточнения посредством более точной настройки на конкретные заболевания, а также для учета неоднородности в распределении особей и возбудителя в пространстве. Модификация модели может идти по пути усложнения вида функций, регулирующих вероятность заражения, генерации и инактивации возбудителя, влияния климатических факторов и т. п., а также по пути перехода к пространственно распределенным системам, например решеткам вероятностных клеточных автоматов.

Ключевые слова: популяционная динамика, SIRS-модель, динамические системы.

Благодарности. Исследование выполнено при финансовой поддержке РФФИ и DFG в рамках научного проекта № 20-51-45001.

Для цитирования: Шабунин А. В. Гибридная SIRS-модель распространения инфекций // Известия вузов. ПНД. 2022. T. 30, № 6. С. 717-731. DOI: 10.18500/0869-6632-003014. EDN: EJQVCA

Статья опубликована на условиях Creative Commons Attribution License (CC-BY 4.0).

Article

DOI: 10.18500/0869-6632-003014

Hybrid SIRS model of infection spread

A. V. Shabunin

Saratov State University, Russia E-mail: shabuninav@info.sgu.ru Received 29.05.2022, accepted 25.07.2022, available online 31.10.2022, published 30.11.2022

Abstract. Purpose of this work is to build a model of the infection spread in the form of a system of differential equations that takes into account the inertial nature of the transfer of infection between individuals. Methods. The paper presents a theoretical and numerical study of the structure of the phase space of the system of ordinary differential equations of the mean field model. Results. A modified SIRS model of epidemic spread is constructed in the form of a system of ordinary differential equations of the third order. It differs from standard models by considering the inertial nature of the infection transmission process between individuals of the population, which is realized by introducing a «carrier agent» into the model. The model does not take into account the influence of the disease on the population size, while population density is regarded as a parameter influencing the course of the epidemic. The dynamics of the model shows a good qualitative correspondence with a variety of phenomena observed in the evolution of diseases. Discussion. The suggested complication of the standard SIRS model by adding to it an equation for the dynamics of the pathogen of infection presents prospects for its specification via more precise adjustment to specific diseases, as well as taking into account the heterogeneity in the distribution of individuals and the pathogen in space. Further modification of the model can go through complicating the function which defines the probability of infection, generation and inactivation of the pathogen, the influence of climatic factors, as well as by means of transition to spatially distributed systems, for example, networks of probabilistic cellular automata.

Keywords: population dynamics, SIRS model, dynamical systems.

Acknowledgements. The reported study was funded by RFBR and DFG according to the research project No. 20-51-45001.

For citation: Shabunin AV. Hybrid SIRS model of infection spread. Izvestiya VUZ. Applied Nonlinear Dynamics. 2022;30(6): 717-731. DOI: 10.18500/0869-6632-003014

This is an open access article distributed under the terms of Creative Commons Attribution License (CC-BY 4.0).

Введение

Методы нелинейной динамики позволяют изучать явления природы путем построения и анализа простых (так называемых качественных) математических моделей. Качественные модели, в отличие от имитационных, не претендуют на предсказание детального поведения моделируемой системы в конкретных условиях, однако, позволяют понять характер наблюдаемых явлений, выявлять их закономерности и внутренние механизмы, а также определять, за счет каких параметров или воздействий можно изменять их характеристики. Одной из областей применения качественных моделей является исследование распространения эпидемий инфекционных заболеваний в биологических популяциях [1-4].

Математическое моделирование эпидемий является одним из разделов популяционной динамики [5], используемой как в биологии, так и в других науках, например, социологии. Для предсказания хода эпидемий используются различные математические методы: анализ временных рядов [6], регрессионные [7] и авторегрессионные [8] модели, решетки клеточных автоматов [9-12], искусственные нейронные сети [13-15] и др. Классические модели распространения инфекций представляют собой системы обыкновенных дифференциальных уравнений (ОДУ). Наиболее известной из них является модель SIRS, предложенная в 1920-х годах Кермаком и МакКендриком [16]. В SIRS-модели популяцию разбивают на группы здоровых и восприимчивых (S — Susceptible), больных (I — Infectious) и выздоровевших (R — Recovered) особей и строят системы уравнений, определяющие закон изменения относительного числа особей в каждой из

групп, опираясь на предположение о случайном и равномерном распределении особей в популяции. Такую систему называют моделью среднего поля. При этом процессы инфицирования особей описывают подобно столкновениям частиц идеального газа в статистической физике. Основой этого подхода является предположение, сделанное В. Хамером [17], о том, что скорость распространения эпидемии зависит от частоты контактов между восприимчивыми и инфицированными особями; частота контактов определяется произведением плотностей населения восприимчивых и инфицированных особей в популяции. Такой подход, несмотря на очевидную простоту и наглядность, не всегда адекватно описывает реальные процессы заражения, которые могут характеризоваться определенной нелокальностью и инерционностью. Инерционность биологических процессов может быть учтена посредством введения в модель времени задержки, то есть путем использования уравнений с отклоняющимся аргументом. Одним из первых такой подход был применен в 1948 г. в работе Хатчирсона [18]. В дальнейшем уравнения с запаздыванием использовались во множестве работ, обзор которых можно найти в [19,20]. Однако модели с запаздыванием представляют собой системы с бесконечным числом степеней свободы, что затрудняет их анализ. Другой возможный подход — введение дополнительного уравнения (и соответственно дополнительной переменной), описывающего механизм инерционности. Именно такой метод и предлагается в работе.

В данной работе предлагается модификация SIRS-модели, в которой передача инфекции происходит опосредованно, за счет взаимодействия с агентом-переносчиком, в качестве которого могут выступать вирусы, бактерии или паразиты1. Такой подход особенно актуален для прогнозирования распространения респираторных вирусных инфекций, при которых агент, вызывающий заражение, является чрезвычайно подвижным и относительно долгоживущим, поэтому акт заражения может происходить в отрыве (во времени и пространстве) от инфицированных особей. Такая модель является моделью взаимодействия двух систем: популяции особей и популяции вирусов, каждая из которых живет по своим законам. В частности, под особями мы будем понимать обособленных индивидуумов, состояние которых меняется дискретным образом, а под вирусами — внешнее поле, воздействующее на особей и приводящее к изменению их состояния, то есть заражению.

1. Модифицированная двухкомпонентная SIRS-модель распространения инфекционных заболеваний

Несмотря на различия, обусловленные спецификой того или иного инфекционного заболевания, большинство из них имеют общие черты:

заболевание распространяется через заражение, источником которого является ранее инфицированная особь (/), а реципиентом — восприимчивая особь (S). В результате восприимчивая особь заболевает и сама становится источником последующих заражений: S ^ I;

• заболевшая особь со временем излечивается, приобретая при этом иммунитет к последующим заражениям (R): I ^ R;

• наконец, с течением времени приобретенный иммунитет утрачивается, возвращая особь к исходному состоянию: R ^ S.

Таким образом, в эволюции каждой особи мы наблюдаем циклическую цепочку превращений между дискретным и конечным набором состояний S ^ I ^ R ^ S. Отсюда и название данной модели — SIRS.

В стандартной SIRS-модели акт заражения описывается как результат локального контакта особей S и I: S + I ^ 21. Однако на практике заражение может происходить и опосредованно,

1В дальнейшем всех таких агентов-посредников будем условно называть вирусами.

без непосредственного взаимодействия особей. В настоящей работе предлагается именно такая схема заражения, основанная на «обменном взаимодействии» особей вирусными частицами. В этой схеме зараженная особь (!) выступает генератором вирусов (V), которые затем вследствие диффузии или перемешивания распространяются по ареалу обитания, заражая восприимчивых

особей (5): в I. Вирусы в данной схеме являются посредниками между больной особью и восприимчивой. Таким образом, вместо стандартной SIRS-модели нами предлагается двухкомпо-

нентная (особи + вирусы) модель: в I — К — Б.

В предлагаемой модели динамика вирусных частиц принципиально отличается от поведения особей популяции. Последние представляют собой объекты с дискретным набором состояний (5,I, К], переходы между которыми являются случайными событиями и характеризуются своими значениями вероятностей (Р^). Кроме того, каждой из особей требуется некоторый ареал обитания, вследствие чего число особей на заданной площади всегда ограничено некоторым максимальным числом N. В противоположность этому, вирусные частицы способны неограниченно накапливаться в каждой точке пространства, поэтому их количество может принимать произвольные положительные значения. Они также способны перемещаться в процессе диффузии. Таким образом, взаимодействие особей с вирусами является некоторой аналогией взаимодействия частиц с полем. Описанные превращения, произошедшие в каждой элементарной ячейке пространства за время ДЬ, могут быть представлены в виде следующей схемы (см. рис. 1):

I — I + Б + ^ I + пУ,

I % К,

к 5, У — (1 - ц) У,

где буквы над стрелками обозначают вероятности соответствующих переходов. Кратко опишем последовательность операций схемы (1).

1. Зараженная особь (!) генерирует V вирусных частиц (V).

2. Восприимчивая особь (5) заражается при контакте с п вирусами с вероятностью Р(п).

3. Зараженная особь (!) выздоравливает с вероятностью Р2 и становится иммунной (Д).

4. Иммунная особь (К) теряет иммунитет с вероятностью Р3 и возвращается к восприимчивому состоянию (5).

5. Инактивация ц вирусов.

Рис. 1. Схема переходов между состояниями SIRS+V-модели: сплошные линии показывают направления переходов; штриховая линия указывает на влияние вирусов на переход между состояниями особей

Fig. 1. Diagram of transitions between states of the SIRS+V model: solid lines indicate the direction of transitions; a dashed line indicates the influence of viruses on the transition between states of individuals

Все трансформации однонаправленны и порядок смены состояний особей жестко детерминирован, в то время как сами превращения, а следовательно, и моменты времени, в которые они происходят, являются случайными. Динамика данной модели определяется как параметрами, так и выбором функции Р(п).

При рассмотрении процессов распространения эпидемии основной интерес представляет не динамика индивидуальных особей, а изменение их количеств ^ (к € {5,I, Щ). Как уже говорилось выше, общая численность популяции ^+N2ограничена по величине некоторым максимально допустимым значением N, которое назовем максимальной емкостью популяции. Величина N обусловлена некоторым общим ресурсом (количеством пищи, площадью территории и т. д.), сдерживающим неограниченный рост численности. Отношение полной численности популяции к максимально возможной — есть относительная населенность

^ + N1 + ° = N ,

принимающая значение от нуля до единицы. Поскольку полное число частиц в схеме (1) не меняется, С является параметром. Величина параметра С — важный фактор, определяющий ход инфекционных процессов в популяции.

2. Модель среднего поля

В приближении среднего поля для схемы (1) можно составить систему уравнений, управляющих изменением числа особей в каждом из состояний (к € {5,I, Щ) и числа вирусов Иу за малый интервал Д1 Вводя в качестве переменных относительные плотности населения к = Ык /Ы и переходя к пределу Дí ^ 0, запишем систему обыкновенных дифференциальных уравнений

I = Р

г = Р2г - Р3г,

V = VI - ^у

(здесь учтено, что 8 + г + г = С). Для определения вида функции Р(у) естественно предположить, что она должна монотонно увеличиваться с концентрацией вирусов от Р(0) =0 до Р(ж) = 1. В качестве такой функции выберем Р(у) = (1 - (1 - Р\)и), где Р\ € [0 : 1[ — вероятность заразиться «единичной» порцией вирусов2. Обозначив а = - 1п(1 - Р\), можно записать функцию вероятности в более удобной форме: Р(у) = 1 - ехр (-ау).

В уравнении (2) можно уменьшить число независимых параметров, введя новую переменную: г = ау, и новый параметр: о = ау. В итоге получим

; =

г = Р2г - Р3г, г = ог - ^г.

Уравнение (3) будет являться основой для анализа поведения системы (1). Из условий задачи следует, что все переменные и параметры (3) являются неотрицательными числами; кроме того, переменные г и г суммарно ограничены сверху: г + г ^ С ^ 1.

2Выбор данной зависимости определяется тем, что вероятность не заразиться V порциями вирусов составит , где ^ — вероятность остаться здоровым после контакта с одной порцией вирусов.

3. Стационарные решения

Проанализируем стационарные решения системы (3) в области допустимых значений параметров и переменных. Предельными траекториями системы (3) являются два состояния равновесия: тривиальное Ео = (0, 0, 0) и нетривиальное

Ег =

где введены обозначения: А = Р2/Р3, В = о/ц, а значение го определяется как корень трансцендентного уравнения

20 = ТГА

Точка Ео соответствует случаю полного выздоровления популяции. Ее устойчивость определяется корнями характеристического уравнения

= — Р3,

— ( Р2 + ц) ( Р2 — ц)2 + 4оС

2,3 =

2

которые в силу неотрицательности всех входящих в формулу величин являются вещественными. В зависимости от соотношения между ( = ВС и Р2 она может быть либо устойчивым узлом (при Р2/(1> 1), либо седлоузлом (при Р2/(1 < 1). Бифуркационное условие для Ео

Р2 =

одновременно соответствует прохождению точки Е1 через начало координат, при котором она приобретает устойчивость. Таким образом, при любых значениях параметров в фазовом пространстве (3) существует только одна устойчивая неподвижная точка. При Р2 > ( это точка Ео и при Р2 < ( — точка Е1. В первом случае популяция не подвержена заражению и любое случайное проникновение инфекции затухает во времени; во втором случае при наличии первоначального инфицирования возникает эпидемия.

Из выражения (6) видно, что определяющую роль в эпидемических процессах играет отношение (1/?2, являющееся произведением нескольких параметров: avСт,0Т2, где х,0 = ц-1 — среднее время жизнеспособности вирусов, Т2 = Р2-1 — средняя продолжительность заболевания (так называемый «период инфицирования»). Каков биологический смысл данного фактора? Если учесть, что средняя концентрация вирусов, генерируемых одной инфицированной особью составляет vт^, а доля восприимчивых особей близка к 100%, то величина ( определяет среднюю скорость заражений, производимых одной инфицированной особью на начальном этапе эпидемии. Умножив ее на средний период инфицирования Т2, получим, что (/Р2 представляет собой хорошо известную в моделировании эпидемий характеристику — базовый индекс репродукции, обозначаемый обычно как Ко. Таким образом, условие потери устойчивости состоянием равновесия Ео в уравнении (3) полностью согласуется с известным в эпидемиологии условием возникновения эпидемии Ко > 1.

4. Установившийся уровень заболевания в популяции

Далее рассмотрим систему (3) при Р2 < то есть в условиях развития заболевания. В этом случае координаты устойчивой точки Е1 определяют установление динамического равновесия между числом больных и здоровых индивидуумов. Обычно представляет интерес относительное

число больных особей в популяции Хо = го/С, которое характеризует уровень заболевания. Оно определяется трансцендентным уравнением

Хо =

1 + А

1 -

Р2Х0

1 — exp

Проведем качественную оценку формулы (7) для типичных значений входящих в нее параметров. Прежде всего, заметим, что величина в скобках меньше единицы, соответственно, уровень

заболевания ограничен сверху множителем (1 + А)-1. Для большинства инфекций средняя длительность иммунитета Т3 = Р-1 во много раз превосходит длительность самого заболевания Т2 [4]. Поэтому А — достаточно большое число3. Отсюда следует, что: (а) единицей в знаменателе (7) можно пренебречь и (б) величина Хо не превосходит нескольких процентов. Поэтому при не слишком больших значениях фактора (I экспоненту в знаменателе можно линеаризовать: ехр (-(1Х0) — 1 — (1Х0. Тогда формула для Х0 существенно упростится

-о -1 (i—£)•

Здесь 1 — Л-1 — это тот установившийся уровень заболевания в популяции, который существовал бы там при отсутствии иммунитета. Как видно из формулы (8), наличие иммунитета снижает средний уровень заболевших в А раз.

Нужно отметить, что величины А и Ко не являются независимыми, поскольку обе содержат параметр Р2. Для исключения этой неоднозначности формула (8) может быть переписана в виде

Хо - -Тз

(Т2 — £> •

—т3=300 —т,=100

В этом случае, как видно из (9), установившийся уровень заболевания в популяции определяется действием двух сомножителей, одним из которых является т-1, а вторым — превышение длительности инфекционного периода над его критическим значением В целом, наличие длительного иммунитета (тз ^ 1) предохраняет популяцию от высоких уровней заражения. Другими сдерживающими Хо факторами в равной степени являются уменьшение периода инфекции Т2 и увеличение обратного фактора й-1. ¡^

На рис. 2 приведены графики Хо(т2) при длительных интервалах иммунитета тз ~ 100, построенные как по полной формуле (7), так и по ее упрощенной версии (9); последние изображены пунктирными линиями. Мы видим, что даже при больших й — 4, если средний период инфекции не продолжается слишком долго (т2 ^ 10), обе формулы дают близкие величины. Поэтому, для многих заболеваний соотношение (9) может быть хорошим приближением.

0.25

0.20 -

0.15 -

0.10 -

0.05 -

20

Рис. 2. Зависимость уровня инфицирования популяции от периода инфицирования т2 при длительном иммунитете; d = 4.16

Fig. 2. Dependence of the population infection level on infection period т2 with long-term immunity; d = 4.16

3Например, для большинства респираторных вирусных инфекций средний период инфицирования составляет одну-две недели, а иммунитет может сохраняться в течение года и более.

1

5. Первоначальное развитие заболевания

Координаты точки Е\ определяют уровень заражения в популяции на больших временах, то есть при £ —> те. Однако на практике часто бывает важно знать, как будет развиваться заболевание на начальном этапе после первоначального заражения. Поэтому нас также интересует переходный процесс от произвольного начального состояния в окрестности начала координат к устойчивой точке Е1. Переходный процесс определяется устройством фазового пространства в окрестности Е1, то есть типом данного состояния равновесия. Его можно определить, исходя из собственных значений А,1_з якобиана системы (3):

exp

(--¿0) - Р2 - 1 exp -г0) - 1 (с -

Р2 о

--

-Рз 0

— I

Аналитический вид собственных значений является достаточно громоздким, поэтому воспользуемся численными расчетами. Выберем значения параметров, которые использовались в предыдущем разделе: С = 0.6, о = 0.7, ц = 0.3, т2 ~ 10, т3 ~ 100 и построим графики зависимости собственных значений от среднего периода инфекции Т2 (рис. 3, а) и средней длительности иммунитета т3 (рис. 3, Ь) Как видно из графиков, в рассматриваемой области параметров точка Е1 имеет два комплексно-сопряженных собственных значения Х1 и Х2 и одно вещественное — Хз, причем вещественные части Х12 очень малы и монотонно уменьшаются с ростом длительности иммунитета. Таким образом, Е1 представляет собой устойчивый фокус, и при подходе к этой точке траектория будет совершать множество витков, демонстрируя вначале значительные осцилляции числа заболевших, постепенно уменьшающиеся по амплитуде. Временная реализация данного переходного процесса определяется начальными условиями. Выберем в качестве таковых точку вблизи начала координат: г о = 0.001, го = 0 и го = 0, которая соответствует ситуации проникновения нескольких зараженных особей в изначально здоровую популяцию. Типичный

<и

-0.9 -1.0

0

-0.1 -0.2 - Im Х2 Re A,1j2 Х2 << 0 -0.1

-0.3 - - а

-0.4 - — Хз : <<

-0.5 - - <u -0.3

-0.6 - [ -

-0.7 - - -0.4

-0.8 - I - -0.5

0 4

1 < 1 <

^ReX^ Im X2

- 1

, 1,1,1,1,1,1,1,1,1,

12 16 20 24 28 32 36 40 t2

b

50 100 150 200 250 300 350 400 450 500

Рис. 3. Зависимость собственных чисел состояния равновесия Е1: а — от т2 при хз = 200 и b — от хз при т2 = 10; оставшиеся параметры зафиксированы: о = 0.7, - = 0.3, С = 0.6. Сплошные кривые соответствуют вещественным величинам, штриховые — мнимым; серая область на рис. a выделяет зону, где Е1 является фокусом

Fig. 3. Dependence of the eigenvalues of the equilibrium state Ei: a — on т2 at т3 = 200 and b — from т3 at т2 = 10; the other parameters are fixed: о = 0.7, - = 0.3, С = 0.6. The solid curves correspond to real values, dashed lines correspond to imaginary ones; gray region on the fig. a marks the zone where Ei is the focus

8

a

и

0.40 0.35 0.30 0.25 0.20 0.15 0.10 0.05 0

i0 =0.001, r0= =0, Z0= 0

¿0= =0.001, r0= =0.2, Z0= 0

¿0= =0.001, r0= =0.4, Z0= 0

0 50 100 150 200 250 300 350 400 450 500 t

Рис. 4. Уровни заражения в популяции X(t) на начальном этапе эпидемии при разных начальных условиях; значения параметров: т2 = 5, хз = 300, а = 0.7, ц = 0.3, С = 0.6

Fig. 4. Population infection levels X (t) at initial stage of the epidemic under different initial conditions; the parameter values are: x2 =5, тз = 300, а = 0.7, ц = 0.3, С = 0.6

вариант реализующегося при этом переходного процесса изображена на рис. 4 (линия, отмеченная кружками). Здесь мы видим чрезвычайно быстрый рост числа заболевших непосредственно после инфицирования, за которым следует такой же быстрый спад почти до нуля. Далее следует второй пик существенно меньшей величины, после чего траектория выходит на уровень, близкий к . С точки зрения популяционной динамики, данное поведение можно интерпретировать как затухающую последовательность волн заражения, при которой уровень заражения популяции демонстрирует последовательность резких пиков, прежде чем эпидемия достигнет динамического равновесия, характеризуемого малым уровнем заражения. В данном случае, как следует из рисунка, пик заражения в самом начале заболевания затрагивает больше трети численности популяции. Большая величина этого пика определяется тем, что на данном этапе иммунитет в популяции еще не сформирован.

Таким образом, несмотря на сравнительно небольшой устанавливающийся стационарный уровень заражения, число зараженных в отдельные моменты времени может достигать значительных величин. При тяжелом протекании заболевания, сопровождающимся вероятными летальными исходами, наличие высокого пика первой волны заражения может приводить для популяции к серьезным последствиям.

При других начальных условиях динамика заражения может выглядеть иначе. Если на момент первоначального заражения в популяции уже присутствует достаточное число особей, обладающих иммунитетом (например, вследствие предварительной вакцинации), переходный процесс становится более плавным. При этом нарастание числа зараженных происходит медленнее, а максимальный уровень заражения оказывается существенно меньшим. Так, например, при го = 0.2 (что составляет 33% от численности популяции) на пике заражения число инфицированных составляет около 0.2, а при го = 0.4 (67% от численности популяции) — менее 0.075 (см. рис. 4). В целом, зависимость максимального числа инфицированных от уровня иммунных особей, имеющегося в наличии на момент заражения, представлена на рис. 5. Как видно из графика, при уровне иммунных особей от нуля примерно до 60% величина пика заражения снижается с ростом го почти линейно, после чего становится

0.40 0.35 0.30 0.25

S

^f 0.20 0.15 0.10 0.05 0

0 0.1 0.2 0.3

0.4 0.5

rJC

0.6 0.7 0.8 0.9

Рис. 5. График зависимости максимального значения относительного числа инфицированных от уровня сформированного иммунитета, имеющегося в популяции на момент заражения. Штриховая линия отмечает уровень заражения, соответствующий установившемуся состоянию

Fig. 5. Dependence of maximum value of infected relative number on the level of formed immunity available in the population at the time of infection. Dashed line marks the infection level corresponding to the steady state

более пологой и, начиная с уровня 80%, происходит стабилизация около величины Xmax ~ 0.075. При других значениях параметров количественные значения будут другими, но качественный вид данной зависимости сохраняется.

5.1. Влияние сезонных факторов на ход заболевания. В предыдущих разделах была рассмотрена автономная модель распространения инфекции, где все параметры — постоянные числа. Однако процессы, связанные с инфицированием, как правило, зависят от условий окружающей среды, которые для большинства климатических зон периодически меняются из-за сезонных факторов. Поэтому автономная модель для них является недостаточной.

Попробуем модифицировать рассматриваемую систему (3), чтобы учесть влияние изменяющихся условий. Смена времен года будет приводить к периодической модуляции параметров уравнений вокруг некоторых средних значений, что может отражаться на сезонных колебаниях числа заболевших. Из всех используемых в (3) параметров наиболее чувствительным к влиянию внешних факторов является скорость инактивации вирусов ц, обратная величина к которой (т^) определяет время, в течение которого вирусные частицы могут оставаться активными вне организма зараженной особи. Этот интервал существенным образом зависит от температуры и влажности окружающей среды, а также от ряда других факторов, таких как количество солнечной радиации, средняя интенсивность ветра и т. п. Для учета сезонных факторов полагаем, что параметр ц в уравнении (3) модулируется гармонической функцией времени вокруг среднего значения цо по закону: ц(£ ) = ц0 (1 + mF (cos (2n/i))), где f — частота, а т — индекс модуляции; функция F задает влияние внешних факторов на параметр ц, она выбирается так, что |F(х)| ^ 1 при |х| ^ 1.

Рассмотрим сначала случай линейной зависимости: F(х) = х. Выберем среднее значение параметра ц0 = 0.3, используемое ранее. При наличии модуляции параметра мы наблюдаем периодические колебания уровня заражения X (t) в виде регулярной последовательности острых пиков заражения, вид которых показан на рис. 6, a. Период колебаний равен периоду функции модуляции Т = /-1. В результате заболевание периодически достигает больших значений уровня заражений на пиках и спадает почти до нуля в интервалах между ними. Амплитуда установившихся колебаний Ах очевидным образом растет с ростом коэффициента модуляции. Что касается ее зависимости от частоты модуляции /, то она носит выраженный резонансный характер, что видно по резонансным кривым, построенным на рис. 6, b. При малых индексах модуляции резонансная частота близка к собственной частоте осциллятора, определяемой как

И

0.40 0.35 0.30 0.25 0.20 0.15 0.10 0.05

500

1000 1500 2000 2500 3000 t

0.005

b

0.005 0.010 0.015 0.020 0.025 0.030 0.035 0.040 f

Рис. 6. Параметрические колебания в системе

Fig. 6. Parametric oscillations in the system

0

0

a

0.40 г 0.35

—m=0.25 —m=0.40

0 30 ' " л л Л Л л л Л Л Л л л л л л л Л 030 л; \/ V/ V V V \i V/ V V V м м V V V м \

и

0.25 0.20 0.15 0.10 0.05 0

a

0

1000

2000

3000

b

0.005 0.010 0.015 0.020 i

0.025

Рис. 7. Удвоение периода колебаний: а — временная реализация X

Fig. 7. Doubling of the oscillation period: а — time-series X

0

t

0.25 0.20 ^ 0.15 0.10 0.05 0

—m=0.3 —m=0.5

Хш (А,1,2)/(2л) (последняя отображена на рисунке пунктирной линией). С увеличением глубины модуляции наблюдается небольшое смещение резонансной частоты в сторону нижних частот.

Помимо колебаний с частотой внешнего параметрического воздействия при некоторых значениях параметров в осцилляторе наблюдается удвоение периода колебаний. В этом случае период колебаний удваивается по сравнению с периодом воздействия (рис. 7, a), а фазовый портрет приобретает вид двухоборотного предельного цикла (рис. 7, Ь). Такой тип колебаний означает на практике наличие вспышек заболевания с двухгодичным интервалом. Подобная динамика действительно наблюдается для некоторых заболеваний. Например, она отмечалась для заболевания корью до введения массовой вакцинации [21].

Зависимость скорости инактивации вирусов ц от условий внешней среды может выражаться и нелинейной функцией Р (х). В этом случае в спектре колебаний численности заболевших можно ожидать пики на гармониках, кратных базовой частоте. Например, если мы выберем зависимость в виде квадратичного полинома Р(х) = ах2 + Ъх + с, то максимумы заболевания появляются дважды за период. Такой вид колебаний показан на рис. 8. Он соответствует хорошо известному явлению сезонных вспышек заболевания респираторными вирусными инфекциями.

Поскольку неавтономная модель не обладает, как известно, инвариантностью по отношению ко времени, существенное влияние на процессы в ней может оказывать выбор момента начального инфицирования. Очевидно, что развитие эпидемии может происходить

0

1000

2000

3000

4000

5000

t

Рис. 8. Временная реализация параметрических колебаний при нелинейной зависимости параметра инактивации:

F

Fig. 8. Time-series of parametric oscillations with nonlinear dependence of the inactivation parameter: F

Рис. 9. Зависимость максимального заражения от начальной фазы функции модуляции параметра

Fig. 9. Dependence of maximum infection on initial phase of the parameter modulation function

по-разному, в зависимости от того, в какое время года произошло первоначальное заражение. Чтобы учесть данный фактор, введем в закон модуляции параметра начальную фазу (ф): = ц0 (1 + mF (cos (2п ft + ф))) и будем измерять уровень заболевания на пике эпидемии Xmax в зависимости от ф. Результаты расчетов подтверждают предположение о важности момента начального заражения на ход эпидемии. Типичный вид зависимости приведен на рис. 9. Как видно из графика, величина пика заражения может меняться в несколько раз в зависимости от того, насколько благоприятным для распространения эпидемии оказался момент первоначального проникновения инфекции в популяцию.

Заключение

Предложенная модифицированная 81Я8-модель демонстрирует динамику, характерную для развития эпидемических процессов, наблюдаемых в природных популяциях. Если не учитывать периодическое изменение условий, вызываемое климатическими факторами, то при наличии первоначального заражения в сообществе возникает последовательность из нескольких затухающих «волн» заражения, сходящаяся к небольшому среднему установившемуся уровню. При этом амплитуда первой волны может достигать очень больших значений, сравнимых с численностью всей популяции. Учет климатических факторов ведет к периодически повторяющимся вспышкам заболевания, которые могут происходить как с удвоенной частотой (сезонные вспышки респираторных инфекций), так и с половинной частотой воздействия (двухгодичные вспышки заболевания некоторых инфекций). Период этих процессов определяется видом функции, описывающей влияние климатических факторов на скорость инактивации вирусов в окружающей среде. Конкретизация вида этой функции возможна за счет учета биологических факторов, характеризующих поведение возбудителя инфекции вне организма зараженной особи.

Рассматривая система описывает многие характерные явления, наблюдаемые при развитии эпидемий. В то же время она является достаточно грубой, поскольку не учитывает конкретные особенности заболеваний, а также возможное перемещение особей и возбудителя в пространстве. Модификация модели может идти по пути уточнения вида функции Р(у), например с учетом порогового характера заражения, наблюдаемого для некоторых заболеваний; использования более реалистичного (нелинейного) уравнения, описывающего генерацию вирусов заболевшей особью; уточнения вида функции модуляции Р(х) и т. п. Такие уточнения представляют собой подстройку модели под конкретный вид заболевания и учитывают его биологические особенности, то есть соответствуют переходу от качественного моделирования к имитационному. Другим направлением модификации предложенной модели может быть переход от системы ОДУ к решеткам вероятностных клеточных автоматов. Такой подход даст возможность естественным образом учесть стохастический характер эпидемиологических процессов, а также рассмотреть влияние на эти процессы неравномерности в распределении особей и вирусов по ареалу обитания и влияние процессов их перемещения.

Список литературы

1. Бейли Н. Математика в биологии и медицине. М.: Мир, 1970. 326 с.

2. Марчук Г. И. Математические модели в иммунологии. Вычислительные методы и эксперименты. М.: Наука, 1991. 304 c.

3. Hethcote H. W. The mathematics of infectious diseases // SIAM Review. 2000. Vol. 42, no. 4. P. 599-653. DOI: 10.1137/S0036144500371907.

4. Андерсон Р., Мэй Р. Инфекционные болезни человека. Динамика и контроль. М.: Мир, 2004. 784 c.

5. Базыкин А. Д. Нелинейная динамика взаимодействующих популяций. Москва - Ижевск: Институт компьютерных исследований, 2003. 368 c.

6. Serfling R. E. Methods for current statistical analysis of excess pneumonia-influenza deaths // Public Health Reports. 1963. Vol. 78, no. 6. P. 494-506. DOI: 10.2307/4591848.

7. Burkom H. S., Murphy S. P., Shmueli G. Automated time series forecasting for biosurveillance // Statistics in Medicine. 2007. Vol. 26, no. 22. P. 4202-4218. DOI: 10.1002/sim.2835.

8. Pelat C., Boeile P.-Y., Cowling B. J., Carrat F., Flahault A., Ansart S., Valleron A.-J. Online detection and quantification of epidemics // BMC Medical Informatics and Decision Making. 2007. Vol. 7. P. 29. DOI: 10.1186/1472-6947-7-29.

9. Boccara N, Cheong K. Automata network SIR models for the spread of infectious diseases in populations of moving individuals // Journal of Physics A: Mathematical and General. 1992. Vol. 25, no. 9. P. 2447-2461. DOI: 10.1088/0305-4470/25/9/018.

10. Sirakoulis G. C., Karafyllidis I., Thanailakis A. A cellular automaton model for the effects of population movement and vaccination on epidemic propagation // Ecological Modelling. 2000. Vol. 133, no. 3. P. 209-223. DOI: 10.1016/S0304-3800(00)00294-5.

11. Шабунин А. В. SIRS-модель распространения инфекций с динамическим регулированием численности популяции: Исследование методом вероятностных клеточных автоматов // Известия вузов. ПНД. 2019. T. 27, № 2. C. 5-20. DOI: 10.18500/0869-6632-2019-27-2-5-20.

12. Шабунин А. В. Синхронизация процессов распространения инфекций во взаимодействующих популяциях: Моделирование решетками клеточных автоматов // Известия вузов. ПНД. 2020. T. 28, № 4. С. 383-396. DOI: 10.18500/0869-6632-2020-28-4-383-396.

13. Фирсов О. В. Гибридное прогнозирование заболеваемости раком почки и смертности от него на основе нейросетевых и статистических технологий // Врач-аспирант. 2006. Т. 10, № 1. C. 15-32.

14. Ефимова Н. В., Горнов А. Ю., Зароднюк Т. C. Опыт использования искусственных нейронных сетей при прогнозировании заболеваемости населения (на примере г. Братска) // Экология человека. 2010. № 3. C. 3-7.

15. Белецкая C. Ю., Коровин В. Н., Родионов О. В. Разработка прогностических моделей развития заболеваемости детей в городском административном районе на основе нейросетевых технологий // Вестник Воронежского государственного технического университета. 2010. Т. 6, № 12. C. 201-205.

16. Kermack W. O., McKendrick A. G. A contribution to the mathematical theory of epidemics // Proc. R. Soc. Lond. A. 1927. Vol. 115, no. 772. P. 700-721. DOI: 10.1098/rspa.1927.0118.

17. Hamer W. H. The Milroy lectures on epidemic disease in England: The evidence of variability and persistence of type // The Lancet. 1906. Vol. 1. P. 733-739.

18. Hutchinson G.E. Circular casual systems in ecology // Annals of the New York Academy of Sciences. 1948. Vol. 50, no. 4. P. 221-246. DOI: 10.1111/j.1749-6632.1948.tb39854.x.

19. Gopalsamy K. Stability and Oscillations in Delay Differential Equations of Population Dynamics. Netherlands: Springer, 1992. 502 p. DOI: 10.1007/978-94-015-7920-9.

20. Переварюха А.Ю. Непрерывная модель трех сценариев инфекционного процесса при

факторах запаздывания иммунного ответа // Биофизика. 2021. Т. 66, № 2. С. 384-407. DOI: 10.31857/S0006302921020204.

21. Anderson R. M., May R. M. Spatial, temporal, and genetic heterogeneity in host populations and the design of immunization programmes // Mathematical Medicine and Biology: A Journal of the IMA. 1984. Vol. 1, no. 3. P. 233-266. DOI: 10.1093/imammb/1.3.233.

References

1. Bailey N. Mathematical Approach to Biology and Medicine. London: John Wiley and Sons; 1967. 296 p. DOI: 10.2307/2982529.

2. Marchuk GI. Mathematical Models in the Immunology. Simulation Methods and Experiments. Moscow: Nauka; 1991. 304 p. (in Russian).

3. Hethcote HW. The mathematics of infectious diseases. SIAM Review. 2000;42(4):599-653. DOI: 10.1137/S0036144500371907.

4. Anderson RM, May RM. Infectious Diseases of Humans: Dynamics and Control. Oxford: Oxford University Press; 1991. 768 p.

5. Bazykin AD. Nonlinear Dynamics of Interacting Populations. Singapore: World Scientific; 1998. 216 p. DOI: 10.1142/2284.

6. Serfling RE. Methods for current statistical analysis of excess pneumonia-influenza deaths. Public Health Reports. 1963;78(6):494-506. DOI: 10.2307/4591848.

7. Burkom HS, Murphy SP, Shmueli G. Automated time series forecasting for biosurveillance. Statistics in Medicine. 2007;26(22):4202-4218. DOI: 10.1002/sim.2835.

8. Pelat C, Boelle PY, Cowling BJ, Carrat F, Flahault A, Ansart S, Valleron AJ. Online detection and quantification of epidemics. BMC Medical Informatics and Decision Making. 2007;7:29. DOI: 10.1186/1472-6947-7-29.

9. Boccara N, Cheong K. Automata network SIR models for the spread of infectious diseases in populations of moving individuals. Journal of Physics A: Mathematical and General. 1992;25(9): 2447-2461. DOI: 10.1088/0305-4470/25/9/018.

10. Sirakoulis GC, Karafyllidis I, Thanailakis A. A cellular automaton model for the effects of population movement and vaccination on epidemic propagation. Ecological Modelling. 2000;133(3): 209-223. DOI: 10.1016/S0304-3800(00)00294-5.

11. Shabunin AV. SIRS-model with dynamic regulation of the population: Probabilistic cellular automata approach. Izvestiya VUZ. Applied Nonlinear Dynamics. 2019;27(2):5-20 (in Russian). DOI: 10.18500/0869-6632-2019-27-2-5-20.

12. Shabunin AV. Synchronization of infections spread processes in populations interacting: Modeling by lattices of cellular automata. Izvestiya VUZ. Applied Nonlinear Dynamics. 2020;28(4):383-396 (in Russian). DOI: 10.18500/0869-6632-2020-28-4-383-396.

13. Firsov OV. Hybrid prediction of kidney cancer incidence and mortality based on neural network and statistical technologies. Postgraduate Doctor. 2006;10(1):15-32 (in Russian).

14. Efimova NV, Goryunov AY, Zarodnyuk TS. Experience in using the artificial neural networks in prognosing the population morbidity (on an example of Bratsk). Human Ecology. 2010;(3):3-7 (in Russian).

15. Beletzkaya SY, Korovin VN, Rodionov OV. Development of the prognostic models of development deseases of children in city administrative area on the basis of neural network technologies. Bulletin of the Voronezh State Technical University. 2010;6(12):201-205 (in Russian).

16. Kermack WO, McKendrick AG. A contribution to the mathematical theory of epidemics. Proc. R. Soc. Lond. A. 1927;115(772):700-721. DOI: 10.1098/rspa.1927.0118.

17. Hamer WH. The Milroy lectures on epidemic disease in England: The evidence of variability and persistence of type. The Lancet. 1906;1:733-739.

18. Hutchinson GE. Circular casual systems in ecology. Annals of the New York Academy of Sciences. 1948;50(4):221-246. DOI: 10.1111/j.1749-6632.1948.tb39854.x.

19. Gopalsamy K. Stability and Oscillations in Delay Differential Equations of Population Dynamics. Netherlands: Springer; 1992. 502 p. DOI: 10.1007/978-94-015-7920-9.

20. Perevaryukha AY. A continuous model of three scenarios of the infection process with delayed immune response factors. Biophysics. 2021;66(2):327-348.

DOI: 10.1134/S0006350921020160.

21. Anderson RM, May RM. Spatial, temporal, and genetic heterogeneity in host populations and the design of immunization programmes. Mathematical Medicine and Biology: A Journal of the IMA. 1984;1(3):233-266. DOI: 10.1093/imammb/1.3.233.

Шабунин Алексей Владимирович — родился в Саратове (1966). Окончил с отличием физический факультет Саратовского государственного университета по направлению «Радиофизика и электроника» (1990). Защитил диссертацию на соискание учёной степени кандидата физико-математических наук по специальности «Радиофизика» (1998, СГУ) и доктора физико-математических наук по специальности «Радиофизика» (2007, СГУ). С 1990 года работает на кафедре радиофизики и нелинейной динамики Саратовского государственного университета, в настоящее время — в должности профессора. Научные интересы — нелинейная динамика, синхронизация, мультистабильность, клеточные автоматы, искусственные нейронные сети. Опубликовал свыше 80 научных статей по указанным направлениям.

Россия, 410012 Саратов, Астраханская, 83

Саратовский государственный университет имени Н. Г. Чернышевского E-mail: shabuninav@info.sgu.ru ORCID: 0000-0002-3495-9418 AuthorID (eLibrary.Ru): 34839