Минимизация нелинейных наименьших квадратов для оценки параметров моделей SEIR и SEIRD

В этом посте мы определим модели SEIR и SEIRD и минимизируем нелинейный метод наименьших квадратов, чтобы оценить их параметры на основе наблюдаемых данных о случаях коронавируса. Для хорошей работы методы нелинейной оптимизации, такие как алгоритм Левенберга – Марквардта (LM), нуждаются в хорошем начальном предположении для решения нелинейной оптимизации наименьших квадратов. Мы будем использовать моделирование SEIR и SEIRD, чтобы угадать оценки параметров и использовать эти значения в оптимизации.

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

Прежде чем мы углубимся в суть дела, давайте познакомимся с типом моделей, к которым относятся SEIR и SEIRD.

Однако, если вы хотите сразу перейти к коду, записные книжки Jupyter с полными кодами доступны в разделах Моделирование и оценка параметров модели SEIR и Моделирование и оценка параметров модели SEIRD.

Механистические модели

Механистические модели используются для моделирования и прогнозирования распространения болезней и инфекций, таких как грипп, лихорадка Эбола или даже коронавирус. Отличительной чертой механистических моделей по сравнению с моделями машинного обучения является их способность кодировать предполагаемую причинную природу передачи болезней. Обычно это включает построение математических формулировок причинного механизма и использование аналитических инструментов для проверки этих гипотез с помощью наблюдаемых данных [R. Baker et al., 2018].

Короче говоря, механистические модели поддаются интерпретации, в то время как модели машинного обучения, как правило, нет. Из-за этого механистическая модель, которая хорошо объясняет это явление, может дать проницательную информацию о том, как распространяется болезнь и, возможно, как с ней бороться. Различия в типах данных и информации, которые предоставляют эти две модели, перечислены в Таблице 1 ниже.

Компартментные модели

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

SEIR & SEIRD

Компартменты и прогрессия между эпидемиологическими состояниями в моделях SEIR и SEIRD показаны на Рисунке 1 и Рисунке 2 ниже.

Чтобы узнать больше о динамических системах этих моделей и о том, как их моделировать на Python, обратитесь к Моделирование компартментных моделей в эпидемиологии с использованием виджетов Python и Jupyter.

Оценка параметров

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

Данные

Мы будем использовать данные о случаях коронавируса в Индии, доступные в COVID-19 REST API для Индии. Структура JSON, возвращаемого из API, выглядит следующим образом:

{
'success': boolean, 
'data': [
         {'day': 'YYYY-MM-DD', 
          'summary': {
                     'total': int,
                     'confirmedCasesIndian': int,
                     'confirmedCasesForeign': int,
                     'discharged': int,
                     'deaths': int,
                     'confirmedButLocationUnidentified': int
                     }, 
          'regional': [
                       {
                       'loc': string,
                       'confirmedCasesIndian': int,
                       'discharged': int,
                       'deaths': int,
                       'confirmedCasesForeign': int,
                       'totalConfirmed': int
                       }                     
                      ]
          }
         ], 
'lastRefreshed': timestamp, 
'lastOriginUpdate': timestamp
}

Пример данных выглядит следующим образом:

{
 'success': True,
 'data': [
          {'day': '2020-03-10',
           'summary': {
                      'total': 47,
                      'confirmedCasesIndian': 31,
                      'confirmedCasesForeign': 16,
                      'discharged': 0,
                      'deaths': 0,
                      'confirmedButLocationUnidentified': 0
                      },
           'regional': [
                        {
                        'loc': 'Delhi',
                        'confirmedCasesIndian': 4,
                        'confirmedCasesForeign': 0,
                        'discharged': 0,
                        'deaths': 0,
                        'totalConfirmed': 4},
                        {
                        'loc': 'Haryana',
                        'confirmedCasesIndian': 0,
                        'confirmedCasesForeign': 14,
                        'discharged': 0,
                        ...
                        }
                        ...
                        ]
            }
            ...
           ],
 'lastRefreshed': '2020-04-27T04:30:02.090Z', 
 'lastOriginUpdate': '2020-04-27T02:30:00.000Z'
}

Эти данные можно легко преобразовать во фрейм данных pandas, мы знакомы с использованием следующего фрагмента кода

Давайте посмотрим первые 5 строк фрейма данных pandas в таблице 2 ниже,

Общее количество положительных случаев коронавируса в Индии продолжает расти из-за ежедневных новых ненулевых положительных случаев. Тенденция случаев показана на Рисунке 3 ниже.

Теперь давайте определим модели на Python.

Определение моделей SEIR и SEIRD

Модели SEIR и SEIRD определяются системой обыкновенных дифференциальных уравнений (ОДУ). Давайте посмотрим, что это за ODE, и как их кодировать и решать на Python,

SEIR

Математическая модель

Код Python

SEIRD

Математическая модель

Код Python

Мы будем использовать модуль odeint в scipy для решения этих систем уравнений. Как и в системе ODE, функция решения ODE очень похожа. Давайте посмотрим, как его можно закодировать на Python для модели SEIRD.

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

Далее мы переходим к оптимизации параметров ОДУ с помощью нелинейной минимизации методом наименьших квадратов.

Минимизация нелинейных наименьших квадратов

Что такое нелинейный метод наименьших квадратов?

Нелинейная модель - это модель, в которой отношения между выходными и входными данными не являются линейными по параметрам модели. Например.

Модель Михаэлиса – Ментен не может быть выражена как линейная комбинация двух бета-версий.

Теперь нелинейный метод наименьших квадратов можно определить как: «Нелинейный метод наименьших квадратов - это форма анализа наименьших квадратов, используемая для согласования набора из m наблюдений с моделью, которая является нелинейной по n неизвестным параметрам (m ≥ n ). » - Википедия

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

Моделирование

Мы будем использовать симулятор модели SEIR и SEIRD, построенный в посте Моделирование компартментных моделей в эпидемиологии с использованием виджетов Python и Jupyter с некоторыми модификациями для этой цели.

При испытании различных комбинаций параметров бета (скорость заражения) = 1,14, сигма (скорость инкубации) = 0,02, гамма (скорость восстановления) = 0,02, mu (уровень смертности) = 0,01, кажется, разумными начальными оценками для модели SEIRD. Настройки симулятора можно увидеть на Рисунке 7 ниже.

Аналогичный метод можно использовать для модели SEIR. Теперь давайте воспользуемся первоначальным предположением для дальнейшей оптимизации параметров с помощью алгоритма LM.

Оптимизация

Не вдаваясь в подробности, в алгоритме Левенберга-Марквардта, также известном как метод наименьших квадратов с затуханием (DLS), изменение весов в последовательные эпохи определяется выражением

Мы будем использовать модуль минимизировать из библиотеки lmfit в Python для минимизации нелинейных наименьших квадратов с использованием алгоритма LM. Для минимизировать требуются остатки в качестве входных данных вместе с другими параметрами модели для оптимизации. Давайте закодируем остатки и передадим их минимизировать, как показано ниже,

Теперь минимизация нелинейных наименьших квадратов с помощью алгоритма LM (leastsq) является простой задачей, как показано в коде ниже,

После завершения оптимизации вы можете распечатать оптимизированные параметры модели и визуализировать соответствие.

Показатели

Мы также можем вычислять такие метрики, как MAE и RMSE, используя оптимизированные параметры и вызывая ode_solver, определенный выше. На момент написания этой статьи я обнаружил, что эти значения следующие:

Fitted MAE
Infected:  587.5587793520101
Recovered: 216.5347633701369
Dead:      21.532532235433454

Fitted RMSE
Infected:  702.1643902204037
Recovered: 284.9841911847741
Dead:      26.4463792416788

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

Полный код

Блокноты jupyter с полным кодом доступны в разделах Моделирование и оценка параметров модели SEIR и Моделирование и оценка параметров модели SEIRD.

Заключение

В этом посте мы узнали о том,

  1. Механические, секционные, модели SEIR и SEIRD
  2. Кодирование ODE на Python и их решение
  3. Нелинейные модели и нелинейный метод наименьших квадратов
  4. Оптимизация параметров ODE с помощью симулятора SEIR / SEIRD и модуля минимизации в библиотеке lmfit
  5. Использование оптимизированных параметров в симуляторе для визуализации кривых эпидемиологического состояния

* Обратите внимание, что тот же подход можно использовать для кодирования и оптимизации любой системы ODE.

То, чему мы не научились (сознательный выбор)

  1. Как интерпретировать результаты этих эпидемиологических моделей.
  2. Когда кривая станет плоской? Сколько из них заразятся, выздоровеют, умрут и т. Д.? Когда мы можем перезапустить экономику?

Оставим эти важные вопросы специалистам.

использованная литература

  1. Бейкер Р.Э., Пенья Дж.М., Джаямохан Дж. И Джерусалем А., 2018, Механистические модели против машинного обучения, борьба, за которую стоит бороться за биологическое сообщество? Биол. Lett.1420170660
    http://doi.org/10.1098/rsbl.2017.0660
  2. Човелл Г., Подгонка динамических моделей к эпидемическим вспышкам с количественной неопределенностью: учебник по неопределенности, идентифицируемости и прогнозам параметров (2017 г.), Infectious Disease Modeling Volume 2, Issue 3, August 2017, Pages 379–398
  3. Нелинейный метод наименьших квадратов (Википедия)
  4. Нелинейная регрессия (Википедия)
  5. Самарасингхе С., 2007, Нейронные сети для прикладных наук и инженерии