В предыдущей статье мы обсуждали, почему мы моделируем болезнь и какие методы мы можем использовать. В этой статье мы рассмотрим, как разработать модель системы дифференциальных уравнений.
Мы обсудили 3 дифференциальных уравнения, использованные в предыдущей статье, чтобы дать вам краткое описание того, что каждое из них делает. Первый показывает, что градиент графика восприимчивости отрицательно пропорционален скорости инфицирования людей. Второй показывает, что градиент инфицированного графика равен скорости, с которой люди заражаются, минус скорость, с которой люди удаляются из модели. Последний показывает, что градиент удаленного графика пропорционален количеству инфицированных людей.
Построение модели
Python - идеальный выбор языка для разработки этой модели благодаря обширной библиотеке графических инструментов и быстрому времени разработки.
Сначала нам нужно решить, какую библиотеку построения графиков мы собираемся использовать, matplotlib - очевидный выбор из-за простоты использования. Если у вас не установлен matplotlib, вам придется установить его с помощью pip.
Окна
py -m pip install matplotlib
Unix / OSX
python3 -m pip install matplotlib
При добавлении matplotlib в наш код python нам нужно указать, какую часть matplotlib мы хотим добавить, в этом случае мы хотим добавить pyplot. Мы также собираемся установить его как plt, чтобы не писать повсюду matplotlib.pyplot.
install matplotlib.pyplot as plt
Для нашей модели нам также нужно импортировать два других модуля: numpy и scipy. Вы можете сделать это с помощью pip. Затем нам нужно добавить их в наш код Python. Для scipy, как и для matplotlib, нам нужно указать, какую часть мы хотим, нам нужен решатель дифференциального уравнения, который можно найти как интегрированную часть scipy. Итак, теперь код должен выглядеть так
install matplotlib.pyplot as plt import numpy as np from scipy.integrate import odeint
В нашей основной функции мы хотим получить наши значения для: всего населения, начала заражения, начала удаления, вероятности заражения, вероятности удаления и времени выполнения модели. Вы можете сделать это как хотите, я решил, что пользователь вводит их
def main():
N = int(input("Total Population: "))
I_0 = int(input("Total Number of Infected at time = 0: "))
R_0 = int(input("Total Number of Removed at time = 0: "))
S_0 = N - I_0 - R_0 #S_0 is starting susceptible
alpha = float(input("Chance of Infection: "))
beta = float(input("Fraction of People who recover per day: "))
MaxTime = int(input("Time for which the model should run: "))
RunTime = np.linspace(0, MaxTime, MaxTime)
Несмотря на то, что можно решить систему дифференциальных уравнений вручную, будет быстрее выбрать функцию Python, которая сделает это за нас. Вот тут-то и появляется scipy.
def Find_Deriv(y, RunTime, N, alpha, beta):
S, I, R = y
dsdt = -alpha * S * I/N
didt = (alpha * S * I/N) - (beta * I)
drdt = beta * I
return dsdt, didt, drdt
Эта функция сообщает модели, какие дифференциальные уравнения находятся в данной точке модели. Теперь нам нужно добавить строку в нашу основную функцию
y_0 = S_0, I_0, R_0 coord = odeint(Find_Deriv, y_0, RunTime, args = (N, alpha, beta)) S, I, R = coord.T
Эти 3 строки генерируют координаты для 3 графиков для каждого подсчета. Теперь самое простое - построить график.
fig = plt.pyplot.figure(facecolor ='w')
ax = fig.add_subplot(111, facecolor = '#dddddd', axisbelow = True)
ax.plot(RunTime, S/1000, 'b', alpha = 0.5, lw = 2, label = "Susceptible")
ax.plot(RunTime, I/1000, 'r', alpha = 0.5, lw = 2, label = "Infected")
ax.plot(RunTime, R/1000, 'g', alpha = 0.5, lw = 2, label = "Removed")
ax.set_xlabel('Time/count')
ax.set_ylabel('Number Of people / 1000')
ax.yaxis.set_tick_params(length = 0)
ax.xaxis.set_tick_params(length = 0)
ax.grid(b=True, which = 'major', c= 'w', lw=2, ls='-')
legend = ax.legend()
legend.get_frame().set_alpha(0.5)
for x in ('top', 'right', 'bottom', 'left'):
ax.spines[x].set_visible(False)
plt.pyplot.show()
Как только это добавлено к основной функции, у нас есть наша модель. Итак, сложив все это вместе, мы получаем модель SIR.
import matplotlib.pyplot as plt
import numpy as np
from scipy.integrate import odeint
def Find_Deriv(y, RunTime, N, alpha, beta):
S, I, R = y
dsdt = -alpha * S * I/N
didt = (alpha * S * I/N) - (beta * I)
drdt = beta * I
return dsdt, didt, drdt
def main():
N = int(input("Total Population: ")) #N = TotalPopulation
I_0 = int(input("Total Number of Infected at time = 0: "))
R_0 = int(input("Total Number of Removed at time = 0: "))
S_0 = N - I_0 - R_0
alpha = float(input("Chance of Infection: ")) #Chance of Infection
beta = float(input("Fraction of People who recover per day: ")) #Fraction Of people who recover per count
MaxTime = int(input("Time for which the model should run: "))
RunTime = np.linspace(0, MaxTime, MaxTime)
y_0 = S_0, I_0, R_0
coord = odeint(Find_Deriv, y_0, RunTime, args = (N, alpha, beta))
S, I, R = coord.T
fig = plt.pyplot.figure(facecolor ='w')
ax = fig.add_subplot(111, facecolor = '#dddddd', axisbelow = True)
ax.plot(RunTime, S/1000, 'b', alpha = 0.5, lw = 2, label = "Susceptible")
ax.plot(RunTime, I/1000, 'r', alpha = 0.5, lw = 2, label = "Infected")
ax.plot(RunTime, R/1000, 'g', alpha = 0.5, lw = 2, label = "Removed")
ax.set_xlabel('Time/count')
ax.set_ylabel('Number Of people / 1000')
ax.yaxis.set_tick_params(length = 0)
ax.xaxis.set_tick_params(length = 0)
ax.grid(b=True, which = 'major', c= 'w', lw=2, ls='-')
legend = ax.legend()
legend.get_frame().set_alpha(0.5)
for x in ('top', 'right', 'bottom', 'left'):
ax.spines[x].set_visible(False)
plt.pyplot.show()
main()
Когда вы запускаете код, пока вы не сделали никаких ошибок, вы получите график, который выглядит следующим образом
В этом примере:
- Общая численность населения = 500000 человек.
- Начальный зараженный = 1
- Пуск удален = 0
- Вероятность заражения = 0,08
- Шанс быть удаленным = 0,03
Надеюсь, эта статья была полезной для понимания того, как вы будете строить базовую модель дифференциальных уравнений. В следующих парах статей мы рассмотрим разработку агентных моделей.
