Mathematical theory of epidemics and pandemics based on the laws of chain and thermal explosions

Full Text

Как показывает история, эпидемии ежегодно приходят на Землю, а пандемия являлась редко, но есть опасения, что все поменяется к худшему. Поэтому развитие эпидемии и пандемии необходимо изучать не только на основе ряда приблизительных эмпирических теорий SIS, SIRS, SEIR, MSEIR, но опираясь на хорошо изученные процессы в точных науках. В работе мы попытались создать теорию развития вирусных и бактериальных заболеваний, в том числе эпидемию и пандемию, на основе хорошо развитого математического аппарата теорий цепного и теплового взрывов в горючем газе.

1. Цепной взрыв. Теория цепных реакций и цепного взрыва была развита Н.Н. Семеновым [1], за что он получил Нобелевскую премию. Цепные реакции идут посредством активных частиц n, скорость изменения которых описывается уравнением:

dn/dt = W₀ + fn – gn, (1)

где W₀ – скорость образования активных центров, f и g – константы скоростей процессов разветвления и обрыва цепей. Обозначив f – g = ϕ, получим:

dn/dt = W₀ + ϕn. (2)

Люди находятся в постоянном движении, вступая в различные социальные контакты, что похоже на броуновское движение молекул в газе, где концентрация компонентов на единицу объема меняется из-за их диффузии в соседние зоны и химических реакций с другими реагентами. В физической химии этот процесс описывается уравнением диффузии, в котором химическая реакция может резко ускориться, приводя к взрыву по механизму самовоспламенения или зажигания [2–5]. Все это похоже на переход обычного заражения в эпидемию и пандемию, а представленная теория позволит найти требуемые параметры развития болезни в относительно стабильных промежутках времени ее течения.

Рассмотрим общину (город, область и др.) с количеством жителей N, где распространяется один вирус, с общим количеством инфицированных людей равным N_i. Через время (t_f) с момента заражения инфицированный человек перестает быть переносчиком вируса, поэтому их число N_t выразим через безразмерный параметр ჰ:

N_t = ჰN_i. (3)

Рассмотрим распространение вируса в общине, полагая, что: 1) активным центром является инфицированный человек, являющийся переносчиком вируса; 2) W₀ = 0, так как точно неизвестно, где, как и когда возник вирус; 3) n₀ – количество первых инфицированных людей, прибывших в общину из другого региона; 4) количество неинфицированных людей неограничено.

С учетом вышеизложенных предположений условие (2) примет вид

dN_i /dt = φn. (4)

Проинтегрируем выражение (4) при условии n = n₀ при t = 0:

N_i = n₀ e^φt, (5)

где N_i и n₀ – текущая и начальная концентрация инфицированных людей.

При φ < 0 число N_i уменьшается со временем, поэтому нет риска возникновения эпидемии. Но при φ > 0 число N_i экспоненциально растет, поэтому φ = 0 является критическим значением. Под обрывом цепи будем считать смерть инфицированного или утрату им способности заражать не инфицированного человека, что происходит через время t_f с момента заражения, оцениваемого по разным источникам от восьми дней до месяца. Так как доля умерших людей за время t_f мала по отношению к числу инфицированных, то в дальнейшем не будем ее учитывать, а в период времени от начала заражения до времени t_f будем считать величину g равной нулю, т.е. до времени t_f рост инфицированных людей N_i вычисляется по формуле:

N_i = n₀ e^ft. (6)

При t > t_f общее количество инфицированных N_i экспоненциально растет, но число людей N_t, передающих вирус, уменьшается на число реинфицированных (не способных заражать), равное числу инфицированных людей до момента времени t – t_f, т.е.:

N_t = n₀ e^ft – n₀ e^f(t–tf) = n₀[1 – e^–ftf] e^ft. (7)

Используя формулы (3), (6), найдем:

ჰ = 1–e^–ftf. (8)

Используя формулу (8), возьмем логарифм уравнения (7):

lnN_t = lnn₀ჰ + ft. (9)

Уравнение (9) показывает, что в логарифмических координатах скорость роста N_t от времени остается постоянной, даже если инфицированные теряют способность заражать здоровых людей через время t_f. В табл. 1 даны значения параметра ჰ, вычисленного по формуле (8), при разных t_f и f. В статистике часто используется дискретное время, равное одним суткам, так будем считать и мы в табл. 1 и во всех формулах работы.

Таблица 1. Параметр ჰ (2–5 строки) в зависимости от безразмерного параметра f и t_f (дни)

Из табл. 1 видно, что параметр ჰ растет при увеличении f и t_f . Функция N_t (7) не имеет точки перегиба, так как ее вторая производная по времени нигде не равна нулю, поэтому нельзя найти условие, когда обычное распространение вируса перейдет в эпидемию и пандемию. В теории цепных реакций [1] вычисляют период индукции t_i, по истечении которого появившиеся изменения становятся заметны. Это время можно высчитать на начальной стадии распространения вируса, используя уравнения (7) и (8):

t_i = [ln(N_ti /ჰn₀)] / f, (10)

где N_ti – первое обнаруженное врачами количество людей, зараженных вирусом в общине. Если N_ti равно 5 человек, при f = 0.2 и n₀ = 1, то t_i ≈ 8 дней.

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

2. Тепловой взрыв. В теории теплового взрыва рассматривают уравнения теплопроводности для определения критического профиля температур в зоне химической реакции [2–5], по достижении которого может произойти взрыв. Так как уравнения теплопроводности и диффузии в газах подобны, то для нашего случая следует рассмотреть последнее:

dN_tg / dt = D∆N_tg + K₁N_tg (N_g – N_tg), (11)

где Δ – оператор Лапласа, D – коэффициент диффузии, K₁ – константа автокаталитической реакции, N_tg и N_g – концентрации активного и исходного реагентов.

Так как распространение вируса происходит не в объеме, а на поверхности общины, то для использования формулы (11) в процессе заражения предложим несколько аналогов для замены объемной концентрации.

2.1. Будем считать N_i, N_t и N аналогами концентраций инфицированных людей, переносчиков вируса и исходного количества людей в общине соответственно.

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

dN_t / dt = D∆N_t + K₂N_t (N – N_t), (12)

где Δ – оператор Лапласа, K₂ – аналог константы скорости передачи вируса от переносчика вируса к неинфицированному человеку, D – аналог коэффициента диффузии, отвечающий за социальную активность переносчиков вируса в общине, который не одинаков у различных групп населения.

В работах [2, 3] показано, что решить нестационарное уравнение теплопроводности, а по аналогии и уравнение диффузии, возможно лишь в некоторых случаях, используя ряд упрощений. Рассмотрим три возможных варианта.

2.1.1. Стационарное уравнение диффузии (dN_t /dt = 0). В теории теплового взрыва известны условия самовоспламенения горючих газов в симметричных сосудах [2, 3], полученные при решении стационарных уравнений теплопроводности, описывающих предельный профиль температуры. В общине людей, вероятно, существует предельный профиль переносчиков вируса, по достижении которого он начнет распространяться взрывным образом. Но в настоящее время нет возможности найти значение D и другие параметры, чтобы проверить эту гипотезу.

2.1.2. Медленное распространение вируса. В этом случае вторым членом в правой части уравнения (12) можно пренебречь, а оставшееся уравнение будет описывать перемещение носителей вируса практически без заражения людей, но к пандемии это не имеет отношения.

2.1.3. Быстрое распространение вируса среди людей. В этом случае первым членом в правой части уравнения (12) можно пренебречь:

dN_t / dt = K₂N_t (N – N_i). (13)

Будем считать, что община похожа на многослойную луковицу, где процесс инфицирования идет слоями, по малым частям N_a: в первую очередь заболеют люди, не выполняющие санитарные меры или по профессиональным обязанностям находящиеся в постоянном контакте с людьми; затем заболеют те, кто соблюдают минимальные санитарные требования, и т.д. В центре луковицы будут люди, ведущие изолированную жизнь, которые смогут избежать заболевания. Выразим N_a через безразмерный параметр ẟ_a, зависящий от типа общины и других факторов.

N_a ≈ ẟ_aN. (14)

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

dN_t / dt = K₂N_t (ჰẟ_aN – N_t) / ჰ. (15)

Проинтегрируем выражение (15):

ln[N_t / (ჰẟ_aN – N_t)] = K₂ẟ_aNt + C₁, (16)

где С₁ – константа интегрирования. Полагая N_0t << ẟ_aN, при t = 0 найдем:

ჰჰ (17)

где N_0t и N_0i – количество переносчиков вируса и инфицированных людей при t = 0.

Используя условие (17), уравнение (16) примет вид

ln[N_t / (l – N_t /ჰẟ_aN)N_0t] = K₂ẟ_aN t. (18)

Полагая, что при t = 0 верно условие (3), т.е. N_0t = ჰN_0i, найдем количество переносчиков вируса в зависимости от времени:

ჰ (19)

В начале распространения вируса, когда N_i << ẟ_aN, формула (19) переходит в выражение (7) при f = K₂ẟ_aN и n₀ = N_0i, т.е. наша модель согласуется с теорией цепного взрыва.

В.М. Гольдберг в замечательных работах [6, 7] рассматривал модель бимолекулярной химической реакции в двухкомпонентной общине людей, состоящей из инфицированных и здоровых людей, т.е. без учета реинфицированных людей (ჰ = 1), и получил выражения

dN_i / dt = K₂N_i (N_max – N_i), (20)

1n[N_i / (N_max – N_i)] = K₂N_maxt + C₂. (21)

При t = 0 найдем

C₂ = ln(N_0i/N_max). (22)

Используя выражения (21), (22), найдем:

 (23)

При N_i << N_max формула (23) примет вид уравнения (7), сочетаясь с решением, вытекающим из теории цепного взрыва при f = K₂N_max и ჰ = 1.

В.И. Гольдберг доказал, что его модель развития вируса на некоторых отрезках времени очень хорошо подчиняется экспоненциальному закону. Поэтому при обработке статистических данных в те же периоды времени наша математическая модель (16)–(19) должна иметь в экспоненте тот же показатель степени, и, сравнивая уравнения (16) и (21), найдем

N_max ≈ N_a = ẟ_aN. (24)

Если уравнения (16) и (19) выразить через N_i и N_a, то получим формулы для обработки статистических данных (25) и расчета числа больных в новой волне распространения вируса (26):

ln[N_i / (N_a – N_i)] = K₂N_at + ln(N_0t / ჰ N_a), (25)

ჰ (26)

Уравнения (25), (26) удобны для обработки статистических данных, так как нужно считать N_i, а не число N_t, имеющее большую статистическую ошибку из-за малого времени вычислений. Формулы (25), (26) и (21), (23) отличаются лишь наличием параметра ჰ, поэтому линейная зависимость, обнаруженная ранее Гольдбергом [6, 7], будет наблюдаться и при использовании формул (25), но с учетом ჰ. Подробное изучение влияния параметра ჰ на начальные и переходные моменты распространения вируса в общине не позволяет формат сообщения.

В.М. Гольдберг ввел труднопонимаемый параметр N_max, при котором S-образная функция, представляющая уравнения (20)–(23), при малых N_i / N_max дает экспоненциальный рост N_i (23), а с ростом N_i приобретает S-образный вид, что подтвердилось в реальности [6, 7]. Так, в работе [7] при обсчете статистических данных числа заболевших в Москве в период с 09.10.2020 по 10.01.2021 он нашел N_max≈ 800 000 чел., при котором статистические данные ложатся на одну прямую, с линейной зависимостью: y ≈ 0.034x – 9.185. В этих расчетах Гольдберг вычислял N_i с начала волны, протекавшей в Москве с 22.07.2020 по 08.02.2021. Варьируя параметром ẟ_aN в уравнении (25), мы нашли значение N_a ≈ 790 000 чел. (близкое к N_max Гольдберга), при котором статистические данные хорошо ложатся на прямую с линейной зависимостью: y ≈ 0.0353t – 9.747, где f = K₂N_а = 0.0353 и C₂ ≈ –9.747, рис. 1. Зная f, из табл. 1 при t_f = 8 дней найдем ჰ ≈ 0.24, при котором выражение (26) показывает рост числа инфицированных в 4 раза меньше, чем дает формула (23). По формуле (17) найдем число переносчиков вируса в начале возникновения пандемии: N_0t ≈ 11 чел., а по формуле (22) получим результат: N_0i ≈ 46 чел. В Москве проживает 12*10⁶ чел., отсюда найдем ẟ_a ≈ 0.067. Заметим, что параметры ẟ_a и N_a взаимосвязаны (14), поэтому независим только один из них, но смысл имеют оба параметра.

Рис. 1. Рост числа заболевших в Москве в период новой волны распространения вируса 09.10.2020–10.01.2021, в логарифмических координатах от времени, в соответствии с уравнением (25), при котором статистические данные (точки) ложатся на теоретическую прямую, с линейной зависимостью: y ≈ 0.0353x – 9.747. В расчетах было найдено оптимальное число Na ≈ 790 000 чел., а Ni отсчитывали от начала возникновения волны с 09.10.2020 г.

На наш взгляд, параметр N_a ≈ N_max имеет ясную статистическую природу – число здоровых людей, имеющих близкую степень уязвимости к инфекции, в составе общины N. В.М. Гольдберг в работе [6], рассматривая пандемию в России с ее начала (03.03.2020), обнаружил четыре волны роста заболевших, каждая из которых сначала развивалась по экспоненте со своими параметрами, а затем внезапно начиналась новая волна, с новыми показателями. То есть представление об общине как о многослойной луковице близко к реальности, в которой процесс инфицирования людей идет “рывками”, по слоям, по малым группам. Так как показатель степени K₂N_а остается постоянным в каждой волне, то, вероятно, количество заболевших людей из других слоев общины невелико и слабо влияет на параметры развития заболевания части общины, числом N_a.

Отметим, что в работе можно было ввести аналоги площадных концентраций, равных тем же, что и в п. 2.1, но деленных на площадь S общины. Но этот аналог концентраций не точен, с учетом разницы плотности населения в общине. Если же принять это определение, то правая часть уравнения (15) будет поделена на S², что изменит показатель степени экспоненциальной функции уравнений (16)–(19), который отражает способность вируса распространяться в общине. Однако для оценки опасности разных вирусов важны относительные значения этих показателей друг к другу. Можно было ввести пропорциональные концентрации, равные тем же, что и в пункте 2.1, но деленные на количество людей N в общине. Но если принять это определение, то правая часть уравнения (15) будет поделена на постоянный множитель N², что, как было показано выше, не даст новой информации для относительной оценки опасности вирусов.

Параметр K₂N_a входит в формулы (16)–(19), (25), (26) и является важным показателем опасности вируса, характеризующим скорость заражения людей в части общины, числом N_a. Используя формулы (25), (26), можно предсказать число инфициро- ванных людей в новой волне и ее скорость, высчитав параметры K₂N_a и N_a по первым дням выявленной волны. На рис. 1, а также на рисунках, представленных в работах [6, 7], видно – точки ложатся кучно, так как большое количество статистических данных уменьшают статистическую ошибку, что и позволяет получать гладкую зависимость требуемых функций от времени.

Несмотря на то, что на стабильных отрезках времени теория хорошо описывает основные особенности распространения волн вируса в общине, нет ответа на многие вопросы: каков механизм перехода от одной волны к другой? что произойдет при появлении мутации старого вируса или появлении нового штамма? как оценить действие вакцин? Ответы на эти и многие другие вопросы можно найти при анализе расчетов, полученных при обработке по формулам (25), (26) статистических данных по разным общинам, в разных условиях распространения вируса.

В работе представлена теория распространения вируса, в том числе эпидемии и пандемии, в изолированной общине людей, полученная с помощью теорий цепного и теплового взрывов, куда введены параметры ჰ и ẟ_a, адаптирующие картину развития вируса к вышеупомянутым теориям. Предложена модель общины в виде многослойной луковицы, в которой процесс инфицирования идет по слоям. Полученные формулы с большой точностью описывают статистические данные по распространению вируса на стабильных этапах развития пандемии. При ряде упрощений наши формулы переходят в формулы, вытекающие из теории цепного взрыва [1] и в выражения, данные в работах В.М. Гольдберга [6, 7]. Наша теория находится в качественном согласии с эмпирическими теориями SIS, SIRS, SEIR, MSEIR. Решения описывают распространение вируса среди усредненных людей, в усредненной общине.

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

БЛАГОДАРНОСТИ

Автор благодарит академика РАН А.А. Берлина за интерес к работе.

ИСТОЧНИК ФИНАНСИРОВАНИЯ

Данная работа финансировалась за счет средств бюджета Федерального исследовательского центра химической физики им. Н.Н. Семенова Российской академии наук. Никаких дополнительных грантов на проведение или руководство данным конкретным исследованием получено не было.

КОНФЛИКТ ИНТЕРЕСОВ

Автор данной работы заявляет, что у него нет конфликта интересов.

About the authors

A. A. Philippov

Author for correspondence.
Email: 7987961@mail.ru
Russian Federation, Moscow

References

Supplementary files

Supplementary Files

Action

1. JATS XML

Download

2. Fig. 1. The increase in the number of cases in Moscow during the new wave of virus spread 09.10.2020–10.01.2021, in logarithmic coordinates from time, in accordance with equation (25), in which statistical data (points) fall on a theoretical straight line, with linear dependence: y ≈ 0.0353x – 9.747. The optimal number of Na ≈ 790,000 people was found in the calculations, and Ni was counted from the beginning of the wave occurrence from 09.10.2020.

Download (77KB)

Indexing metadata