реализация модели SIRD на Python

Вопрос или проблема

import numpy as np
import matplotlib.pyplot as plt

#sirdv
#a-alpha b-beta g-gamma d-detla m-mu e=эффективность вакцины
#N=население s(t)=подверженные d(t)=умершие i(t)=инфицированные r(t)=выздоровевшие v(t)=вакцинированные
#n=s+i+r+d+v
# ds/dt=-b/N(si) dr/dt=gi dd/dt=mi di/dt=b/N(si-gi-mi)-a*s) dv/dt=a*s-(e*b*i*v)/N
N=1000 #общее население
b=0.5 # коэффициент передачи от 0.4 до 0.6
a=0.05 #5 % подверженного населения
e=0.7 #на данный момент 
g=1/7 # коэффициент выздоровления
m=0.01 #коэффициент смертности от 0.01 до 0.1 %
M=(m)/(m+g) #=смертность/выздоровление
R0=b/g
print(M)
print(R0)

steps=np.arange(0,365,dtype=int)
time=1*steps
def SIRDV(s,i,r,d,b,g,m,a,v,e):
    sdash= -s*i*(b/(s+i+r+d+v))-(a*s)
    rdash=g*i
    idash=((s*i*b)/(s+i+r+d+v))-(g*i)-(m*i)
    ddash=m*i
    vdash=(a*s)-((e*b*i*v)/(s+i+r+d+v))

    return sdash, idash, rdash, ddash, vdash

tlen=len(time)
s=np.empty([tlen])
i=np.empty([tlen])
r=np.empty([tlen])
d=np.empty([tlen])
v=np.empty([tlen])

r[0]=0
i[0]=1
d[0]=0
v[0]=0
s[0]=N-i[0]-r[0]-d[0]-v[0]

dt=0.05

for t in steps:
    if t<len(steps)-1:
        sdash, idash, rdash, ddash, vdash= SIRDV(s[t],i[t],r[t],d[t],v[t],b,g,m,e,a)
        s[t+1]=s[t]+dt*sdash
        i[t+1]=i[t]+dt*idash
        r[t+1]=r[t]+dt*rdash
        d[t+1]=d[t]+dt*ddash
        v[t+1]=v[t]+dt*vdash

plt.plot(time,s,label="подверженные")
plt.plot(time,i,label="инфицированные")
plt.plot(time,r,label="выздоровевшие")
plt.plot(time,d,label="умершие")
plt.plot(time,v,label="вакцинированные")
plt.xlabel('время в днях')
plt.ylabel('число людей')
plt.legend()
plt.pause(0.01)
plt.show()

Ответ или решение

Чтобы решить вашу проблему с реализацией модели SIRD (Susceptible-Infected-Recovered-Deceased) с вакцинацией в Python, давайте подробно рассмотрим ваш код и выявим возможные недочёты. Опираясь на описанные вами графики и уравнения, мы постараемся улучшить вашу модель, чтобы она более точно отображала динамику вовлечения вакцинированного населения.

Анализ кода

Основные параметры и уравнения

Вы определили основные параметры и уравнения для модели, что является правильным шагом. Тем не менее, есть несколько моментов, которые могут вызывать проблему с динамикой графиков:

1. Уравнение для вакцинированных (vdash): Ваше уравнение выглядит следующим образом:
   [
   vdash = a \cdot s – \frac{e \cdot b \cdot i \cdot v}{s+i+r+d+v}
   ]
   Здесь нас интересует, что часть (- \frac{e \cdot b \cdot i \cdot v}{s+i+r+d+v}) может ограничивать рост вакцинированных при низких значениях (v). Тем не менее, если количество вакцинированных остаётся низким, это может привести к отсутствию их значимого увеличения.

Потенциальные исправления кода

1. Измените уравнение для вакцинированных

Попробуйте упростить уравнение для (vdash):

vdash = a * s

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

2. Корректировка расчётов в цикле

Имеет смысл пересмотреть ваше уравнение для изменяемого населения в цикле:

if t<len(steps)-1:
    # Изменяем максимум и адаптируем формулы

Исправьте здесь лишние проверки для увеличения читаемости кода.

3. Точность временного шага

Ваш временной шаг (dt) равен 0.05. Убедитесь, что он достаточно мал, чтобы изменения были гладкими и реалистичными, но не слишком мал, чтобы код не стал ненужным образом медлительным.

Обновлённый код

Ниже приведённый код с изменениями, которые помогут улучшить вашу реализацию:

import numpy as np
import matplotlib.pyplot as plt

# Константы
N = 1000  # общее население
b = 0.5   # скорость передачи
a = 0.05  # скорость вакцинации
e = 0.7   # эффективность вакцины
g = 1/7   # скорость выздоровления
m = 0.01  # смертность

# Временные параметры
steps = np.arange(0, 365, dtype=int)
time = 1 * steps

def SIRDV(s, i, r, d, v, b, g, m, a, e):
    sdash = -b/(N) * (s * i) - a * s
    rdash = g * i
    idash = b/(N) * (s * i) - g * i - m * i
    ddash = m * i
    vdash = a * s  # Упрощенное уравнение вакцинации
    return sdash, idash, rdash, ddash, vdash

# Инициализация массивов
tlen = len(time)
s = np.zeros(tlen)
i = np.zeros(tlen)
r = np.zeros(tlen)
d = np.zeros(tlen)
v = np.zeros(tlen)

# Начальные условия
s[0] = N - 1
i[0] = 1
r[0] = 0
d[0] = 0
v[0] = 0

# Расчет
dt = 0.05
for t in steps[:-1]:
    sdash, idash, rdash, ddash, vdash = SIRDV(s[t], i[t], r[t], d[t], v[t], b, g, m, a, e)
    s[t + 1] = s[t] + dt * sdash
    i[t + 1] = i[t] + dt * idash
    r[t + 1] = r[t] + dt * rdash
    d[t + 1] = d[t] + dt * ddash
    v[t + 1] = v[t] + dt * vdash

# Визуализация
plt.plot(time, s, label="Susceptible")
plt.plot(time, i, label="Infected")
plt.plot(time, r, label="Recovered")
plt.plot(time, d, label="Deceased")
plt.plot(time, v, label="Vaccinated")
plt.xlabel('Time in Days')
plt.ylabel('Number of People')
plt.legend()
plt.show()

Заключение

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