Оценка уровня смертности, обусловленной SARS-CoV-2 в регионах России, на основе эконометрической модели

Полный текст

Прогнозирование смертности от COVID-19 на территориях стран и их регионов является актуальной проблемой. Это важно при разработке оперативных профилактических мер в регионах, а также для принятия мер по снижению смертности. Последнее — стратегическая цель Правительства РФ.

Вопросы математического моделирования в эпидемиологии сжато освещаются в [15]; небольшой обзор подходов к моделированию заболеваемости и смертности от COVID-19 приводится в [12]. Дополняя его, отметим, что эконометрическое моделирование этой смертности выполнено в [18]; к временным рядам также применялся подход Бокса и Дженкинса [19]. Сопоставление точности моделей, созданных известными научными центрами, сделано в [16]. Из весьма успешных работ по моделированию смертности в регионах можно выделить [2, 8, 9, 17].

Ряд исследователей в качестве инструментария выбирают регрессионный анализ [3, 4, 6, 12, 18]. Так, в [4] рассматривалась «избыточная смертность» по данным ЗАГС и строились линейные модели; в результате на тех или иных отрезках пандемии наиболее информативными оказались различные факторы. Высокая информативность одного из показателей, связанных со смертностью от COVID-19, была выявлена в [3]. Этот показатель рассчитывается как доля летальных случаев среди случаев с известным исходом (умерших и выздоровевших), которые ранее входили в список зарегистрированных лиц. Он был назван Э.М. Гольдштейном «коэффициентом летальности», поэтому ниже он именуется нами как «индикатор Гольдштейна» (ИГ) (использование «летальности» является здесь нестандартным, согласно О.М. Драпкиной и соавт. [2020]). Также важными оказались распространенность активных случаев и обеспеченность населения врачами-специалистами [12]. Влияние загрязненности атмосферного воздуха на заболеваемость изучалось в [6].

Целью исследования является построение (с оценкой по точности) модели линейной регрессии, которая связывает смертность от COVID-19 в регионах с комплексом объясняющих факторных переменных. На основе этой модели можно выполнять краткосрочный прогноз смертности и более эффективно планировать медицинские ресурсы.

Материалы и методы

Регионами, попавшими в таблицу с обучающими данными, были: 14 областей ЦФО, 2 региона из СЗФО и 4 — из ПФО; в состав южных округов входило 11, а 9 оставшихся находятся на Урале, Дальнем Востоке и в Сибири. В выбранных регионах проживает почти ½ населения страны; в итоге таблица включала n = 40 строк. В таблицу не вошли города федерального уровня, ряд промышленных регионов (Липецкая, Вологодская области и многие другие), а также — территории Арктики, ряд регионов СКФО, Сибири, Дальнего Востока, Урала.

Таблица. Результаты оценивания, коэффициенты эластичности, статистики / Table. Least square estimates for all parameters, coefficients of elasticity, statistics

Некоторые объясняющие факторы были качественными (ординальными) и номинальными по своей природе, поэтому их градации кодируются фиктивными переменными (dummy). В качестве факторных переменных были выбраны (индекс t здесь для простоты не пишется):

• Х₁ — распространенность активных случаев COVID-19 на сутки t — 10, относительно суток с отчетом t, к которым в [14] был привязан показатель смертности (т. е. общее количество подтвержденных случаев COVID-19 за вычетом выздоровевших и умерших к суткам t–10 в соотношении со среднегодовой численностью населения в 2020 г.) (Active Cases);
• Х₂ — показатель обеспеченности населения врачами-специалистами (для трех групп: инфекционисты (И), врачи санитарно-противоэпидемической группы (СПГ) и анестезиологи-реаниматологи (АР), но без учета токсикологов; здесь бралась их линейная комбинация: Х₂ = 0,423 ∙ И + 0,267 ∙ СПГ + + 0,118 ∙ АР, где И – количество врачей-инфекционистов на 10 тыс. жителей за 2020 г. и т. д.) [5, 7] (Doctors);
• Х₃ — осредненное значение для ИГ, измеренное в логарифмической шкале, с использованием 7-дневного окна (при обучении это 23.х-29.х), Log(Mortality) [14];
• d₁ — dummy переменная, равная 1 для регионов, в которых наблюдается более низкий уровень популяционного здоровья (Type 1);
• d₂ — то же для варианта, когда уровень здоровья был сравнительно высок (Type 2);
• d₃ — то же для южных регионов, т. е. регионов из ЮФО, СКФО (South Region);
• d₄ — то же, когда наблюдался повышенный уровень загрязнения атмосферного воздуха в населенных пунктах региона за ряд последних лет (Air Pollution).

Эти четыре переменные и X₂ влияют лишь на изменение свободного члена в (1).

Для материала обучения были взяты значения фактора Х₁ на дату 21.x.2021, а значения Y ^t включались здесь на 10 суток позднее, т. е. на 1.xi.2021. При этом на этапе «обучение» использовалось значение Y ^t, осредненное за 5 суток: с 30.x по 3.xi.2021. При обработке данных по ИГ использовалась операция осреднения по окну шириной 7 суток, после чего брался натуральный логарифм, в итоге получались значения Х₃. Факторы Х₂, d₁,…,d₄ предполагались постоянными по t.

Затем оценивались параметры модели регрессии с переменной структурой:

Y ^t = b ₀ + b ₁X ₁^t–10 + b ₂X ₂ + b ₃X ₃^t–7 + b ₄d ₁ + + b ₅d ₂ + b ₆d ₃ + b ₇d ₄ + e, (1)

где Y ^t — смертность от COVID-19 (отношение количества умерших от этого заболевания за все время пандемии к суткам t к среднегодовой численности населения в 2020 г., умноженное на 100 тыс. [11, 14]). Здесь X₁, X₂, X₃, d₁, d₂, d₃, d₄ — объясняющие факторы, β₀ отражает влияние на Y ^t неучтенных переменных, ε — случайная величина с E(ε) = 0 и дисперсией σ ². Неизвестные параметры β₀, β₁,…,β₇ и σ ² оцениваются по таблице c обучающими данными, причем оценки b₀,…,b₇ для β-параметров и σ ² находятся методом наименьших квадратов (МНК) средствами MS Excel. Каждая строка этой таблицы содержит набор Y ^t, X₁,…., X₃, d₁,…,d₄ для того или иного региона; факторы X₁, X₃ входили здесь с разными лагами относительно даты привязки показателя смертности в [14]. Отчасти близкий подход применялся ранее при анализе данных по онкологии [13].

Нами выделялись регионы типа 1, в которых в 2015–2020 гг. наблюдалась более высокая заболеваемость населения по некоторым основным классам болезней, зарегистрированных у пациентов с диагнозом, установленным впервые в жизни (на 1 тыс. человек населения). В частности, здесь были рассмотрены болезни системы кровообращения; болезни эндокринной системы, расстройства питания и нарушения обмена веществ; врожденные аномалии (пороки развития), деформации и хромосомные нарушения; болезни крови, кроветворных органов и отдельные нарушения, вовлекающие иммунный механизм; болезни органов дыхания; болезни органов пищеварения; болезни мочеполовой системы и новообразования. Также было рассмотрено среднее значение контингентов пациентов со злокачественными новообразованиями за 2016–2020 гг. (численность пациентов, состоящих на учете в лечебно-профилактических организациях, на конец года). Дополнительно учитывались общие коэффициенты смертности за 2020 г. Часть таких регионов сильно пострадала после катастрофы на ЧАЭС; в других регионах или вблизи них действуют АЭС (с не самыми передовыми технологиями). На территории некоторых, или сравнительно недалеко от них, еще во времена СССР был ряд технологических аварий или взрывов [1]. Для 17-ти регионов этого типа выбирались значения переменных d₁ = 1, d₂ = 0. К этому типу можно отнести почти половину областей из ЦФО: Брянскую, Владимирскую, Ивановскую, Курскую, Орловскую, Рязанскую, Смоленскую; 7 регионов из СЗФО, ПФО, ЮФО: Архангельскую, Астраханскую, Волгоградскую, Оренбургскую, Ростовскую, Саратовскую области и Республику Удмуртию, а также некоторые регионы из УФО, СФО, ДФО: Тюменскую область, Алтайский край и Еврейскую автономную область.

Для регионов типа 2, наоборот, характерны довольно низкие значения заболеваемости вышеупомянутых восьми классов болезней, а также — относительно малые радиационные загрязнения и/или загрязнения из-за техногенных аварий. Для них принималось d₁ = 0, d₂ = 1; к этому типу были отнесены пять регионов: три южных (республики Адыгея, Кабардино-Балкария и Ставрополье), Хабаровский край и Московская область; здесь проживает больше молодого населения, а также в них сравнительно мала доля лиц старших возрастных групп. Все прочие регионы, вошедшие в материал обучения, кроме регионов типов 1 и 2, имели нулевые значения dummy переменных d₁ и d₂.

Еще одна переменная кодировала южные регионы на карте страны: d₃ = 1. Наконец, высокий уровень загрязнения воздуха регулярно выявлялся в Архангельской, Белгородской, Воронежской, Курской, Новосибирской, Ростовской, Саратовской и Смоленской областях, а также в Алтайском, Забайкальском, Приморском и Хабаровском краях, республиках Кабардино-Балкария, Крым и Северная Осетия — Алания [10]. Для всех этих 15-ти регионов, принималось, что d₄ = 1.

Математические методы включали корреляционный и регрессионный анализ, методы проверки статистических гипотез. Использовались различные тесты (через t-статистику Стьюдента для гипотез H₀: β_j = 0 и DW-статистику Дарбина–Уотсона для выявления автокорреляции) и традиционный дисперсионный анализ Р. Фишера. Для корректного использования тестов также применялись критерии согласия распределения остатков с нормальным законом. Так как дисперсия σ ² неизвестна, то были взяты тесты χ² Пирсона и Лилиефорса после анализа эмпирической функции распределения остатков. Условие гомоскедастичности проверялось на диаграммах рассеяния остатков регрессии. Весь статистический анализ и обработка данных выполнялись в MS Excel.

Результаты и обсуждение

В верхней части таблицы приводятся b-коэффициенты, или МНК-оценки β.

В 3-й строке в скобках приводятся модули t-статиcтик, из которых после проверки соответствующих гипотез следует, что β-коэффициенты отличны от 0. В последней строке приводится F-статистика Р. Фишера и число степеней свободы ν. Все достигнутые уровни значимости (p-values) в тестах Стьюдента и Фишера были ниже 0,001. Также в 5-й строке приводятся стандартная ошибка s, значение DW-статистики, несмещенная оценка коэффициента детерминации R2, а также две статистики для критериев согласия: χ² и статистика, содержащая расстояние Колмогорова D_n. Наконец, в 4-й строке приводятся средние коэффициенты эластичности. Каждый из них показывает, на сколько процентов изменится смертность Y ^t в среднем по таблице обучения, если Х_j увеличить на 1%, и при этом зафиксировать значения прочих факторов. Сильнее всего на переменную Y ^t влияют обеспеченность врачами и распространенность активных случаев, а несколько слабее — d₁ и X₃ (логарифм скользящего среднего для ИГ). Так, если значение фактора X₂ с обеспеченностью врачами увеличить на 1% (соблюдая нормативы по балансу между специализациями), то смертность снизится на 0,5%.

Модель получилась довольно точной для многих регионов, попавших в таблицу с обучающими данными и сформировавших материал обучения, когда отсчет t привязан в [14] к суткам 1.xi.2021. А именно для 30 регионов ошибка аппроксимации MAPE = 1,94% для модели (1) с b-коэффициентами и со стандартной ошибкой 1,5%. Иногда она была сравнима по точности с моделью из [2]. Для всего материала обучения MAPE (mean absolute percentage error) была 3%; при этом модель объяснила 97,8% разброса переменной Y ^t вокруг ее среднего значения. На рисунке показано, как модель (1) с b-коэффициентами работает на «экзамене» при изменении X₁ и X₃ по оси 0-t. Модель можно улучшить, если уточнить данные по врачам-специалистам; другой вариант ее улучшения — регулярно проводить новое обучение: через 5–7 дней.

Рисунок. Прогноз смертности Y на 1,...,17.xi на 100 тыс. жителей / Figure. Predicting the mortality Y on the 1st… 17^th Nov on 10⁵ population

Каждый точечный график связан с прогнозом Y ^t — распространенности смертности от COVID-19 за сутки t из промежутка [1.xi,17.xi] — для одного из семи регионов. Здесь дате 1.xi.2021 для показателя Y ^t будет соответствовать отсчет t; дате 2.xi — t+1,… Непрерывными и более контрастными линиями на рисунке изображается изменение смертности по факту.

Например, чтобы получить прогнозное значение Y ^t на 5 ноября для Подмосковья, надо взять его значения dummy-переменных d₁,…, d₄, а также — переменной с врачами X₂. Все эти пять переменных не изменяются по t ; они лишь скорректируют для региона значение b₀ (см. табл.). От суток к суткам изменяется количество активных случаев; оно легко находится из [14]; для прогноза на 5.xi.2021 надо брать накопленные данные на 25.x, когда значение фактора X₁ = 8,267 на 1 тыс. жителей (с учетом того, что здесь проживало около 7700 тыс. граждан). Аналогично изменяется и фактор X₃; он определяется по набору из нескольких значений ИГ. Для его расчета здесь берутся данные за 22.х,…, 28.x по количеству выздоровевших r пациентов в регионе и умерших d за каждые сутки [12], затем находятся значения ИГ = 100 ∙ d/(d+r) за 7 суток: с 22.x по 28.x (Предполагается, что за эти дни имелся хотя бы 1 летальный случай). Потом берется их среднее значение, а от него — натуральный логарифм, равный 0,3863; это и есть значение X₃ на 25.х, т. е. середина 7-дневного окна. В результате находится сумма: 9,141 ∙ X₁ + 12,76 ∙ X₃ – 96,99 ∙ X₂ + (сумма четырех произведений с dummy…) + 106,462, которая считается через b-коэффициенты из таблицы. При этом выражение в скобках равно –21,82, так как d₂ = 1, а все d₁,…, d₃ нулевые (см. текст ниже). Итак, для региона будет Y ^t= 9,14X₁ + 12,76X₃ + 37,9186, так как здесь X₂ = 0,48792. В итоге на дату 5.xi.2021 получаем прогноз показателя Y ^t, равный 118. Для сравнения: смертность из-за новообразований в регионе к концу 2020 г. была 172 [7].

Заметим, что распространенность активных случаев COVID-19 в последней декаде октября была примерно одинаковой в Ленинградской, Московской, Тверской и Костромской областях [14]. Среди них в первой декаде ноября в Тверской области была самая высокая смертность от COVID-19, а также одно из наихудших в РФ мест по коэффициенту общей смертности в 2019–2020 гг. Значения Y ^t были здесь близки к значениям в Костромской и Ленинградской областях, а в Подмосковье умирали реже.

Сравнивая средние значения для этой четверки регионов по набору более чем 10 медико-демографических показателей, взятых за 2015–2020 гг. [7], можно предположить, почему смертность из-за COVID-19 в Подмосковье оказалась в это время ниже. Вблизи Москвы проживает более молодое население и меньше людей старше 70 лет (по отношению к тем, кто младше 25), также здесь стабильно высокие душевые расходы бюджета на здравоохранение. Из-за перечисленного здесь невелик фактор X₃; его снизила и более широкая вакцинация взрослых. Также X₃ бывает меньше в регионах, где активнее проводится тестирование на антитела к SARS-CoV-2, отчего возрастает доля нетяжелых случаев, и поэтому больше выздоравливающих [3]. Еще в Московской области в 2019 и 2020 гг. была наименьшая общая смертность в ЦФО (исключая Москву), сравнительно низкие уровни заболеваемости вышеупомянутых восьми классов болезней (исключая болезни органов дыхания и пищеварения). Также здесь относительно низка онкологическая смертность, нет острых проблем с радиоактивным загрязнением территории и загрязнением атмосферного воздуха; поэтому для Подмосковья принимались следующие значения dummy-переменных: d₁ = d₃ = d₄ = 0, d₂ = 1.

Для читателя с более «продвинутой» математической подготовкой может быть интересно следующее. Нами проверялось выполнение предпосылок теоремы Гаусса–Маркова, в которой обоснованы хорошие свойства МНК-оценок. Для этого сначала рассматривались коэффициенты корреляций Y ^t с каждым из факторов, а также коэффициенты корреляции между факторами. Явление мультиколлинеарности не обнаружено. Гипотеза о нормальности по критерию χ² на 5%-м уровне не отвергалась. График эмпирической функции распределения для остатков сравнивался с графиком распределения Гаусса с нулевым средним и стандартом s, где s — несмещенная оценка σ. Наибольшее расхождение графиков D_n = 0,100, так что из-за малости n^½ ∙ D_n, аналогичный вывод сделан тестом Лиллиефорса [13]. Итак, распределение остатков не противоречит закону Гаусса. Анализ диаграмм рассеяния остатков относительно каждого Х_j, а также относительно оценки Y ^t, не выявил явного нарушения условия гомоскедастичности. Тест об отсутствии автокорреляции остатков не отвергался. В итоге все предпосылки теоремы были выполнены.

В результате была построена статистическая модель для краткосрочного прогноза смертности в регионах от COVID-19 на текущие сутки, и для нее найдены коэффициенты эластичности. Модель может оказаться полезной при разработке мер борьбы против пандемии в ряде регионов, а также при создании улучшенных секторных моделей (типа SEIRF, SIRD и др.).

Об авторах

Владимир Сергеевич Степанов

Автор, ответственный за переписку.
Email: _stepanov@cemi.rssi.ru
ORCID iD: 0000-0002-4478-376X
http://istina.msu.ru/profile/VSStepanov

к.ф.-м.н., старший научный сотрудник

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

Дополнительные файлы

Доп. файлы

Действие

1. JATS XML

Скачать

2. Рисунок. Прогноз смертности Y на 1,...,17.xi на 100 тыс. жителей

Скачать (151KB)

Метаданные