CC BY
50
Статистический анализ параметров модели эпидемической ситуации
Ю.Б. Гришунина1 (juliagri@list.ru), Н.А. Контаров 2 3, Г.В. Архарова3, Н.В. Юминова2 (yuminova@mail.ru)
Московский институт электроники и математики Национального исследовательского университета «Высшая школа экономики» 2ФГБУ «НИИ вакцин и сывороток им. И.И. Мечникова» РАН, Москва 3ГБОУ ВПО «Первый МГМУ им. И.М. Сеченова» Минздрава России
Резюме
Предлагается способ вычисления и метод вариации параметров модели эпидемической ситуации, учитывающей внешние риски. Приведен доступный в приложении Exel подробный поэтапный алгоритм расчета и анализа параметров модели по результатам наблюдений. Построен прогноз развития эпидемической ситуации с оценкой длительности вспышки и оценкой общего числа вовлекаемых в нее членов популяции. Исследована эпидемическая ситуация на конкретном примере группы заболеваний ОРВИ-гриппом в одном из регионов Московской области.
Ключевые слова: математическая модель эпидемии, внешние риски, оценка параметров, заболеваемость, популяционный иммунитет
Statistical Analysis of the Model parameters of the Epidemic Situation
Yu.B. Grishunina1 (juliagri@list.ru), N.A. Kontarov2,3, G.V. Arkharova3, N.V. Yuminova2 (yuminova@mail.ru)
1 Moscow Institute of Electronics and Mathematics of National Research University «Higher School of Economics»
2 I.I. Mechnikov Research Institute of Vaccines and Sera of Russian Academy of Sciences, Moscow
3 I.M. Sechenov First Moscow State Medical University, State Educational Institution of Higher Professional Training of the Ministry of Healthcare of the Russian Federation
abstract
A calculation method and a variation method of parameters of the epidemic situation model taking into account external risks is proposed. The detailed step-by-step algorithm available in the application Exel for the calculation and analysis of the model parameters based on the results of observations is given. The prognosis of the epidemic situation with the estimation of the disease outbreak duration and the total number of population members involved in the outbreak is constructed. We investigated the epidemic situation by an example of reports of a diseases group acute viral respiratory infections-flu in one of the districts of the Moscow region. Key words: mathematical model of the epidemic, external risks, parameter estimation, morbidity, population immunity
Введение
Работы по моделированию важных в медицинском отношении и социально значимых процессов предпринимаются давно [1], а в последние годы активизировались в связи с усложнением экономической ситуации, что потребовало более рационального использования финансовых средств в экономике и социальной сфере, в том числе и на профилактические мероприятия [2 - 8].
Для анализа и прогнозирования развития эпидемической ситуации в заданном регионе утвердилась в практическом использовании модель SIR [1 - 5], основанная на разделении населения на три группы: S - восприимчивые (Susceptible), I -инфицированные (Infectious) и R - имеющие иммунитет (Removed): N = S + I + R, где N - общая численность населения региона.
Эпидемическую ситуацию в регионе характеризует распределение людей по указанным группам в зависимости от их статуса к рассматриваемому
заболеванию. Динамика изменения численности групп описывается нелинейной системой дифференциальных уравнений, основными параметрами которых являются частота контактов (интенсивность заражения от одного больного) и интенсивность выздоровления - величина, обратная средней продолжительности заболевания.
В современных работах по моделированию эпидемической ситуации рассматриваются различные модификации SIR-модели, учитывающие, например, латентную заболеваемость, рождаемость, смертность, вакцинацию и т.д., для чего в уравнения вводятся новые переменные и параметры, описывающие соответствующие процессы. В большинстве работ значения параметров модели предполагаются известными или взяты из литературы. При этом не учитывается тот факт, что в различных регионах эти параметры могут существенно отличаться, а каждая вспышка заболеваемости может иметь разные причины и разную интенсивность,
поэтому для построения адекватной модели развития эпидемической ситуации эти параметры каждый раз должны вычисляться заново на основании статистических данных конкретного региона и конкретной вспышки заболевания. Следует также отметить, что параметры моделей могут меняться на фоне развития эпидемической ситуации или зависеть от времени года, погодных условий и других факторов.
Цель данной работы - предложить способ вычисления и метод вариации параметров модели эпидемической ситуации на основании статистических данных по заболеваемости.
Материалы и методы
В основе работы лежит модель эпидемической ситуации, описанная в [9]. Эта модель, построенная на основе классической детерминированной SIR-модели, позволяет учесть влияние внешних источников заражения: для этого вводится новый параметр - А - агрессивность внешних рисков, равный среднему числу зараженных от внешних источников в единицу времени (за неделю).
Следует отметить, что развитие эпидемического процесса в значительной степени определяется уровнем популяционного иммунитета, в формировании которого важную роль играет разъяснительная работа, проводимая медицинскими работниками среди населения, а также в средствах массовой информации, направленная на предупреждение заражения с помощью индивидуальных профилактических мер. К таким мерам относятся прием витаминов, применение защитных масок, правильный режим и т.д., а также вакцинация. При построении модели эпидемической ситуации влияние перечисленных факторов учитывается с помощью еще одного нового параметра - г, равного среднему числу людей, еженедельно приобретающих иммунитет в результате профилактических мер. Динамика численности групп больных I и невосприимчивых R описывается рекуррентными соотношениями:
vt = (l - C/tA + Л)
It+l = /t(l - /л)+ vt flt+1 = at(l-y) + v£ + r
где I > 1 - номер недели, г/г - среднее количество заболевших в единицу времени (за неделю), Я1 -количество невосприимчивых к данной инфекции, N - общее число жителей в наблюдаемом регионе, 1г - число больных в начале текущей недели, Л - интенсивность заражения, ц - интенсивность выздоровления одного больного, величина, обратная средней продолжительности заболевания; Y - скорость потери иммунитета одним человеком, величина, обратная средней продолжительности сохранения иммунитета к данному заболеванию.
Вычисления оценок параметров модели проводятся с помощью метода наименьших квадратов.
Для расчетов использована сводка понедельной заболеваемости ОРВИ - гриппом в одном из регионов Московской области с населением N = 272 390 человек. Наблюдения проводились в течение 23 недель, полученные данные приведены в таблице 1.
Результаты и обсуждение
Графическое представление данных о заболеваемости ОРВИ и гриппом реализовано в Excel (Меню-Вставка-Диаграммы) и представлено на рисунке 1.
Визуально на построенном графике (см. рис. 1) можно выделить три этапа развития эпидемической ситуации. Первый этап продолжается с 1-й по 14-ю неделю наблюдений и характеризуется отсутствием резких скачков заболеваемости, на втором этапе. с 15-й по 17-ю неделю, происходит ее резкий рост, а затем, на третьем этапе, который длится с 18-й по 23-ю неделю, заболеваемость постепенно снижается. Таблица 1.
Данные о заболеваемости ОРВИ и гриппом
№ недели Количество заболевших (чел./нед)
1 1525
2 1338
3 1435
4 1498
5 1262
6 1305
7 1696
8 1685
9 1480
10 1443
11 1418
12 1401
13 1464
14 1663
15 2700
16 2509
17 3380
18 3029
19 2871
20 3046
21 3098
22 2971
23 2292
Рисунок 1.
Заболеваемость ОРВИ и гриппом
При моделировании эпидемической ситуации необходимо иметь в виду, что значения некоторых параметров, входящих в рекуррентные соотношения (*), могут меняться на разных этапах развития эпидемической ситуации, а другие могут быть постоянными величинами.
К параметрам, значения которых не изменяются, относятся прежде всего ц и у, поскольку эти параметры характеризуют заболевание в целом и определяются усредненными величинами; ц -интенсивность выздоровления, величина, обратная средней продолжительности заболевания, у -скорость потери иммунитета, величина, обратная средней продолжительности сохранения иммунитета к данному заболеванию. Значения этих параметров могут быть взяты из доступных медицинских источников информации. Для ОРВИ - гриппа средняя продолжительность заболевания составляет 10 дней, то есть 1,5 недели, ц = 2/3 нед, а продолжительность сохранения иммунитета после перенесенной болезни или вакцинации, по разным данным, составляет в среднем от 6 недель до года и даже больше. В данной работе будем рассматривать вариант со сроком сохранения иммунитета 6 месяцев, это примерно 26 недель, то есть у = 1/26 нед-1. К постоянным параметрам относится также и интенсивность заражения Л, т.е. среднее число людей, которых может заразить в течение недели один больной. Параметр Л определяется прежде всего, частотой контактов, а также зави-
сит от свойств возбудителя, особенно контагиоз-ности: если она близка к 1, как при ОРВИ - гриппе, то Л можно интерпретировать как среднее число контактов в единицу времени (под контактом понимается взаимодействие между двумя людьми, продолжительность которого достаточна для заражения). Частота контактов определяется плотностью населения, а также рядом социальных и экономических факторов, характеризующих уровень жизни в регионе: структурой рынка труда, уровнем развития транспорта и т.д.
Параметры А и г относятся к другой группе параметров, их можно назвать вариативными, поскольку их значения зависят от этапа эпидемической ситуации, более того, можно предположить, что именно они и определяют стадию эпидемического процесса. Суть предлагаемого метода вариации параметров модели состоит в том, что для построения адекватной модели для каждого этапа развития эпидемической ситуации необходимо строить такие оценки этих параметров, чтобы при их подстановке в соотношения (*) отклонение результатов вычислений от реальной заболеваемости было минимальным.
Рассмотрим далее отличительные особенности каждого этапа
Первый этап наблюдений, с 1-й по 14-ю неделю, отличается тем, что в этот период не происходит значительных изменений заболеваемости, и мож-
но предположить, что она находится в пределах естественного фона, когда заражение обусловлено только непосредственными контактами больных и восприимчивых; агрессивность внешних рисков на этом этапе - А = 0. Поскольку при этом эпидемическая ситуация остается спокойной, нет оснований для проведения профилактических мероприятий, поэтому параметр r, характеризующий их эффективность, также можно считать равным нулю. Таким образом, развитие эпидемической ситуации на данном этапе определяется количеством невосприимчивых в начале наблюдений и интенсивностью заражения Л. Для оценки этих параметров с помощью метода наименьших квадратов можно воспользоваться надстройкой «Поиск решения» в Меню «Данные»; в качестве целевой функции выбирается сумма квадратов отклонений расчетных значений заболеваемости от результатов наблюдений:
il>t-vt)2
где - v результаты вычислений по соотношениям (*), vt- данные наблюдений, а ее минимизация проводится по двум переменным: Л и R1. При решении задачи оптимизации получены следующие значения параметров: R1 = 0 чел., Л = 0,72126 чел./нед. Отметим, что интенсивность заражения Л является постоянным параметром и для данного региона характеризует возможность передачи заболеваний с контагиозностью, близкой к 1, как при ОРВИ - гриппе, от больных к восприимчивым, поэтому полученное значение может быть использовано в дальнейшем для моделирования и прогнозирования эпидемической ситуации по таким заболеваниям в этом регионе. Кроме того, по данным наблюдений на этом этапе можно оценить границы значений естественного фона заболеваемости ОРВИ - гриппом в данном регионе. Для этого воспользуемся правилом «трех сигм», смысл которого состоит в том, что вероятность отклонения случайной величины от среднего значения больше, чем на три среднеквадратических отклонения, не превосходит 0,1. В Excel среднее значение вычисляется с помощью функции «СРЗНАЧ», а средне-квадратическое отклонение - с помощью функции «СТАНДОТКЛОН», где в качестве аргументов используются данные о заболеваемости. Исходя из этих данных получаем: среднее значение а = 1472 чел./нед, среднеквадратическое отклонение о = 135 чел./нед, тогда границы а ± 3 о равны 1068 - 1876 чел./нед. На практике правило трех сигм означает, что при неизменных внешних условиях фоновая заболеваемость должна находиться в указанных пределах. Используя полученные значения параметров модели, можно построить прогноз заболеваемости в предположении, что внешние условия не изменятся. Результаты моделирования, прогнозирования, а также диапазон фоновых значений заболеваемости представлены на рисунке 2.
Сравнительный анализ прогноза и статистических данных показывает, что фактический рост заболеваемости значительно превышает значения, полученные из рекуррентных соотношений (*) в предположении, что агрессивность внешних рисков сохраняется на нулевом уровне. Отсюда следует, что на втором этапе наблюдений (15 - 17-я недели), который характеризуется резким ростом заболеваемости, значение параметра А изменится и для построения адекватной модели его необходимо заново оценивать по данным наблюдений на данном этапе. Для этого снова можно воспользоваться методом наименьших квадратов и надстройкой «Поиск решения», взяв в качестве целевой функции:
Т7 V
При проведении вычислений следует учесть, что период роста заболеваемости достаточно короткий, в связи с чем профилактические и санитарно-гигиенические мероприятия на данном этапе только начинают проводиться и либо реализуются в недостаточном объеме, либо эти меры не успевают дать заметных результатов, поэтому параметр г, как и на первом этапе, можно считать равным 0. По фактическим данным заболеваемости получена оценка агрессивности внешних рисков на втором этапе наблюдений: А = 1134 чел./нед. Повышение уровня агрессивности внешних рисков может быть связано с появлением новых источников заражения, помимо больных, инфицированных на первом этапе. Появление этих внешних источников может быть,в частности обусловлено тем, что складываются некоторые благоприятные для определенного вируса условия (например, погодные), при которых он способен сохранять жизнеспособность и высокую концентрацию во внешней среде в течение некоторого времени, достаточного для массового заражения людей, то есть для начала вспышки заболеваемости. Кроме того, рост заболеваемости может быть обусловлен и недостаточным уровнем популяционного иммунитета, который сформировался на первом этапе только у перенесших заболевание. Для полученного значения агрессивности внешних рисков можно построить прогноз заболеваемости в предположении, что она останется на этом уровне и профилактические мероприятия не будут эффективны, то есть при г = 0. Результаты моделирования и прогнозирования для второго этапа представлены на рисунке 3.
Построенный прогноз показывает, что при сохранении значения агрессивности внешних рисков на прежнем уровне рост заболеваемости должен продолжаться, что противоречит данным наблюдений: начиная с 18-й недели фактическая заболеваемость постепенно снижается. Снижение заболеваемости на третьем этапе наблюдений может быть обусловлено двумя причинами. Во-первых, вероятно уменьшение интенсивности заражения от внешних источников: это может произойти по естественным причинам, например в связи с изменением погоды,
Рисунок 2.
Результаты моделирования, прогнозирования и диапазон фоновых значений заболеваемости
Рисунок 3.
Результаты моделирования и прогнозирования для второго этапа развития эпидемической ситуации
когда внешние условия будут уже не столь благоприятны для сохранения высокой концентрации вируса в местах массового скопления людей, или скажутся результаты проводимых санитарно-гигиенических мероприятий. Во-вторых, проходит достаточное время для формирования соответствующего уровня по-пуляционного иммунитета. На этом этапе рост числа невосприимчивых происходит уже не только за счет переболевших, но и за счет лиц, предпринимающих профилактические меры. Кроме того, формируется иммунитет у лиц, прошедших вакцинацию, а у некоторых это происходит естественным путем на фоне воздействия внешних факторов риска. Ослабление влияния внешних источников заражения приводит к уменьшению значения параметра А, характеризующего агрессивность внешних рисков, а накопление эффекта от лечебно-профилактических мероприятий - к увеличению параметра г. Для оценки значений этих параметров модели на третьем этапе эпидемического процесса снова можно воспользоваться методом наименьших квадратов и надстройкой «Поиск решения», где целевая функция будет иметь вид:
23 (v -v): М8 t V
а ее минимизация проводится по переменным А и г. Получены следующие значения вариативных параметров: А = 674 чел./нед, г = 12 299 чел./ нед. Таким образом, на третьем этапе наблюдений агрессивность внешних рисков снизилась в среднем в 1,7 раза, а еженедельное увеличение
количества невосприимчивых в результате проведения лечебно-профилактических мероприятий составило примерно 4,5% от общей численности населения региона. По полученным на третьем этапе оценкам параметров А и г можно построить прогноз дальнейшего развития эпидемической ситуации в предположении, что значения этих параметров не изменятся. Результаты моделирования и прогнозирования представлены на рисунке 4.
Построенный прогноз показывает, что при сохранении значений вариативных параметров продолжится снижение заболеваемости и ее значения достаточно быстро, в течение 2 - 3-х недель, вернутся в границы фонового диапазона. Этот результат может быть использован для прогнозирования длительности рассматриваемой вспышки заболеваемости: если не произойдет резких изменений внешних условий, то продолжительность данной вспышки составит примерно 12 недель, что соответствует ежегодным сезонным вспышкам заболеваемости ОРВИ - гриппом.
Можно также провести сравнительный анализ между количеством переболевших за время вспышки (по фактическим данным за 15 - 23-ю недели) и оценочными данными (по результатам моделирования). Для этого надо просуммировать данные наблюдений и расчетные значения заболеваемости. Сложив результаты наблюдений, получаем, что за время вспышки реальное количество заболевших составило 25 896 человек. Вычисление смоделиро-
Рисунок 4.
Прогноз дальнейшего развития эпидемической ситуации
Рисунок 5.
Результаты моделирования эпидемической ситуации при раннем начале лечебно-профилактических мероприятий
ванных показателей дает близкий результат: 26 058 человек.
Визуальное сравнение графиков фактической и расчетной заболеваемости (см. рис. 4), а также небольшая величина отклонения при вычислении количества заболевших по фактическим данным и по результатам моделирования позволяют сделать вывод, что построенная модель и найденные оценки ее параметров на различных этапах эпидемического процесса обеспечивают достаточно высокий уровень соответствия реальной эпидемической ситуации в конкретном регионе. Вывод - описанная модель может быть использована для исследования возможности управления эпидемической ситуацией в регионе.
В качестве примера можно рассмотреть сценарий, когда профилактические мероприятия начинают проводиться в более ранние сроки и дают эффект уже на втором этапе эпидемического процесса, когда происходит рост заболеваемости в связи с появлением внешних источников заражения. Для моделирования такого сценария будем использовать значение агрессивности внешних рисков для второго этапа наблюдений: А = 1134 чел./нед, и оценку эффективности профилактических мероприятий: г = 12 299 чел./нед, которая ранее была полу-
чена для третьего этапа эпидемического процесса, но теперь это значение параметра г будет фигурировать также и в рекуррентных соотношениях для 15 - 17-й недель наблюдений. Результаты моделирования представлены на рисунке 5.
Из полученного прогноза следует, что при раннем начале лечебно-профилактических мероприятий происходит, во-первых, уменьшение длительности вспышки заболеваемости - она составляет теперь 8 недель, так как на 23-й неделе расчетное значение заболеваемости уже находится в пределах фонового диапазона, а во-вторых, значительно, примерно на 5,5 тыс., сокращается количество перенесших заболевание - 20 410 человек, то есть ранняя профилактика дает и ощутимый экономический эффект.
Таким образом, предлагаемый метод позволяет построить оценки параметров модели, не зависящих от времени, и проследить изменение остальных параметров по статистическим данным, во-вторых, спрогнозировать предполагаемую длительность вспышки и оценить возможное число вовлекаемых в нее членов популяции, в-третьих, провести сравнительный анализ возможной и реальной заболеваемости, а также оценить уровень естественно накопленного или сформированного
лечебно-профилактическими мероприятиями попу-ляционного иммунитета к заболеванию в заданном регионе. Расчеты по приведенным рекуррентным соотношениям доступны любому уверенному пользователю приложения Excel.
Используя полученные оценки параметров модели на разных этапах эпидемического процесса, можно рассмотреть и другие сценарии, направленные на сокращение негативных последствий вспышек заболеваемости, например связанные со снижением агрессивности внешних рисков путем выявления и устранения внешних источников заражения или с повышением эффективности профилактических мероприятий путем увеличения числа вакцинируемых и т.д. Результаты таких исследований могут стать основой для планирования сроков и объемов лечебно-профилактических и санитарно-гигиенических мероприятий в регионе.
При этом возникает еще одна задача - определение сроков начала роста заболеваемости. Поскольку начало вспышки связано с повышением агрессивности внешних рисков, необходимо включение в систему мониторинга эпидемической ситуации контроля показателей состояния окружающей среды, характеризующих возможность заражения от внешних источников. Например, для ОРВИ - гриппа такими показателями могли бы стать концентрация РНК вирусных антигенов и инфекционный титр вирусного материала в пробах воздуха и/или смывах с поверхности в местах массового скопления людей. Ухудшение этих показателей может свидетельствовать о возрастающем риске возникновения вспышки заболеваемости и должно послужить сигналом к началу массовых профилактических мероприятий.
Интересно было бы проследить зависимость между значениями этих показателей и полученными при моделировании оценками агрессивности внешних рисков, которые характеризуют интенсивность вспышки заболеваемости. Поскольку для этого необходимы данные наблюдений хотя бы за 4 - 5 лет, а в сложившейся на сегодняшний день системе мониторинга эпидемической ситуации такие наблюдения не проводятся, изучение влияния
характеристик внешней среды на показатели эпидемической ситуации могло бы стать предметом
дальнейших исследований.
Выводы
1. Параметры модели эпидемической ситуации можно разделить на две группы: постоянные, значения которых не изменяются, и вариативные, значения которых зависят от стадии эпидемического процесса. Для построения адекватной модели эпидемической ситуации на каждом этапе ее развития необходимо заново вычислять оценки вариативных параметров.
2. Частота контактов является постоянным параметром, характеризующим конкретный регион с точки зрения возможности распространения заболеваний от больных к восприимчивым. Однажды построенная оценка этого параметра может быть в дальнейшем использована для моделирования и прогнозирования эпидемической ситуации.
3. К вариативным параметрам относятся агрессивность внешних рисков, определяющая интенсивность заражения от внешних источников, и эффективность профилактических мер, характеризующая уровень популяционного иммунитета, сформированного в результате профилактики.
4. Управление эпидемической ситуацией может осуществляться путем воздействия на вариативные параметры с помощью лечебно-профилактических и санитарно-гигиенических мероприятий.
5. Предметом дальнейших исследований может стать изучение зависимости значений агрессивности внешних рисков от показателей состояния окружающей среды, связанных с конкретным заболеванием.
6. Контроль показателей состояния окружающей среды по социально значимым для конкретного региона заболеваниям должен быть включен в систему мониторинга эпидемической ситуации в данном регионе. ш
литература
1. Бейли Н. Математика в биологии и медицине. Пер. с англ. Е.Г. Коваленко. Москва: Мир; 1970.
2. Мастихин А.В. Финальные вероятности для марковских процессов эпидемии. Москва: МГТУ им. Н.Э. Баумана; 2011.
3. Авилов К. Математическое моделирование в эпидемиологии как задача анализа сложных данных. Доступно на: http://download.yandex.ru/company/ experience/seminars/KAvilovmatmodelirovanie.pdf
4. Allen L.J.S. An Introduction to stochastic epidemic models. Доступно на: http://eaton.math.rpi.edu/csums/papers/epidemic/allenstochasticepidemic.pdf
5. Gray A., Greenhalgh D., Mao X., Pan J. The SIS epidemic model with markovian switching. Доступно на: http://strathprints.strath.ac.uk/41322/
6. Асатрян М.Н., Салман Э.Р., Боев Б.В., Кузин С.Н., Глиненко В.М., Ефимов М.В. Моделирование и прогнозирование эпидемического процесса гепатита В. Эпидемиология и Вакцинопрофилактика. 2012: 1 (62): 49 - 54.
7. Боев Б.В. Моделирование развития эпидемии гриппа A(H1N1) в России в сезон 2009 - 2010 годов. Эпидемиология и Вакцинопрофилактика. 2010; 1 (50): 52 - 58.
8. Боев Б.В., Салман Э.Р., Асатрян М.Н. Применение компьютерного инструментария для прогнозирования водных вспышек гепатита А техногенного характера c оценкой эффективности мер противодействия. Эпидемиология и Вакцинопрофилактика. 2010; 3 (52): 57 - 62.
9. Гришунина Ю.Б., Контаров Н.А., Архарова Г.В., Юминова Н.В. Моделирование эпидемической ситуации с учетом внешних рисков. Эпидемиология и Вакцинопрофилактика. 2014; 5 (78): 61 - 66.
References
1. Bailey N. Mathematics in biology and medicine. Trans. from engl. E.G. Kovalenko. Moscow: Mir; 1970 (in Russian).
2. Mastihin A.V. Final probabilities for markov processes epidemic: Diss. of Phd. of math. sci. Moscow: N.E. Bauman Moscow State Technical University; 2011 (in Russian).
