Играем в эпидемиологов

Для того чтобы создать простейшую модель по распространению инфекции, нам потребуется немного математики и небольшой скрипт на python. Модель будет упрощенная, т.е. мы пропустим ряд параметров, влияющих на распространение болезни. Поэтому нельзя сказать, что реальная инфекция распространяется именно так, как на графиках ниже. Хотя некоторые тенденции все-таки будут прослеживаться.

Будем наблюдать за группой мышей, в которой доля привитых (вакциной против некого патогена) равна vасс_ratio. Допустим нам известно, что патоген распространяется по воздуху не дальше, чем на расстояние distination. И на расстоянии менее или равным distination заражает некоторый процент среди привитых мышек и некоторый процент среди непривитых. Особь продолжает выделять патоген еще некоторое время после заражения – ill_durat, далее она становится безопасной для окружающих.

Импортирую необходимые библиотеки:

import numpy as np import random import pandas as pd import matplotlib.pyplot as plt import matplotlib.animation

Задаю начальные условия:

population = 100# Количество точек vacc_ratio = 0.9 # Доля вакцинированных vacc_infect_prb = 0.05 # Вероятность заражения вакцинированной точки not_vacc_infect_prb = 0.4 # Вероятность заражения невакцинированной точки distination = 0.08 # Расстояние для заражения ill_durat = 12 # Длительность выделения патогена точкой timer = 100 # Время наблюдений

Одна мышка — это точка на плоскости, имеющая координаты (x, y). Привитых мышек обозначим зеленым цветом, а непривитых - синим. Мышки время от времени перемещаются (меняют свои координаты). Первую точку мы заразим искусственно, и отметим красным цветом.

Создаем Dataframe наблюдений, содержащий идентификатор точки, временную метку, цвет, координаты точки и время заболевания.

# Задаем координаты, начальное время, ID, цвет (состояние), время заражения X = pd.Series(np.random.rand(population)) Y= pd.Series(np.random.rand(population)) time = [0.0]*population id = pd.Series(range(0,population)) color = ['green']*int(population*vacc_ratio)+['blue']*(population-int(population*vacc_ratio)) # Первое заражение first_infect_id = random.randint(0,population-1) color[first_infect_id] = 'red' infect_date = pd.Series([-1]*population) infect_date[first_infect_id] = 0 # Сохраняем в DateFrame observ_df = pd.DataFrame({ 'ID': id, 'time':time, 'color': color, 'X': X, 'Y':Y,'infect_date':infect_date })

Далее запустим цикл, который обновляет состояния на каждой итерации (одна итерация - одна единица времени) и сохраняет эти состояния в общей таблице.

for i in range(1, timer): # Cоединяем таблицу наблюдений саму с собой, чтобы получить все пары точек (сочетание без повторений). cross_observ = observ_df.merge(observ_df, how='cross',suffixes=('_1', '_2')) cross_observ = cross_observ[cross_observ['ID_1'] < cross_observ['ID_2']] # Считаем евклидово расстояние между каждой парой точек и записываем его в столбец destinaton cross_observ['destinaton']=cross_observ.apply(lambda obs: np.sqrt((obs['X_1'] - obs['X_2']) ** 2 + (obs['Y_1'] - obs['Y_2']) ** 2), axis=1) # Рассчет новых состояний точек a = cross_observ.apply(func=change_color, axis=1) # Обновляем DataFrame observ df в связи с изменением состояний некоторых точек после применения функции change_color df_new_1 = pd.DataFrame() df_new_2 = pd.DataFrame() df_new = pd.DataFrame() df_new_1[['ID','color','infect_date']] = a[a['color_1']=='red'][['ID_1', 'color_1','infect_date_1']].drop_duplicates() df_new_2[['ID','color','infect_date']] = a[a['color_2'] == 'red'][['ID_2','color_2','infect_date_2']].drop_duplicates() df_new[['ID', 'color','infect_date']] = pd.concat([df_new_1, df_new_2],ignore_index=True).drop_duplicates() observ_df = observ_df.merge(df_new, how='left', on=['ID'], suffixes=('', '_new')) observ_df['color'] = np.where(pd.notnull(observ_df['color_new']), observ_df['color_new'], observ_df['color']) observ_df['infect_date'] = np.where(pd.notnull(observ_df['infect_date_new']), observ_df['infect_date_new'], observ_df['infect_date']) observ_df.drop(['color_new','infect_date_new'], axis=1, inplace=True) # Добавляем observ_df в общий DataFrame observ_all = observ_all.append(observ_df) # Обновляем время и координаты точек observ_df['time'] = i if i%5 == 0: observ_df['X'] = pd.Series(np.random.rand(population)) observ_df['Y'] = pd.Series(np.random.rand(population))

Реализация функции change_color():

#Если 2 точки находятся на расстоянии не более destination # и только одна из точек помечена красным(заражена), # и с момента ее заражения прошло не более ill_durat, # тогда вторая точка окрашивается красным (заражается) с вероятностью #vacc_infect_prb, если она была привита и с вероятностью not_vacc_infect_prb, если непривита. def change_color(row): if row['destinaton'] <= distination and row['color_1']== 'red' and (row['time_1'] - row['infect_date_1']) < ill_durat: if row['color_2']=='blue': row['color_2'] = np.random.choice(['blue', 'red'], p=[1 - not_vacc_infect_prb, not_vacc_infect_prb]) if row['color_2']=='red': row['infect_date_2'] = row['time_1'] if row['color_2'] == 'green': row['color_2'] = np.random.choice(['green', 'red'], p=[1 - vacc_infect_prb, vacc_infect_prb]) if row['color_2'] == 'red': row['infect_date_2'] = row['time_1'] if row['destinaton'] <= distination and row['color_2']== 'red' and (row['time_1'] - row['infect_date_2']) < ill_durat: if row['color_1']=='blue': row['color_1'] = np.random.choice(['blue', 'red'], p=[1 - not_vacc_infect_prb, not_vacc_infect_prb]) if row['color_1'] == 'red': row['infect_date_1'] = row['time_1'] if row['color_1'] == 'green': row['color_1'] = np.random.choice(['green', 'red'], p=[1 - vacc_infect_prb, vacc_infect_prb]) if row['color_1'] == 'red': row['infect_date_1'] = row['time_1'] return row

Осталось создать анимированные графики и сохранить их в формате gif.

fig, ax = plt.subplots() x, y = [], [] sc = ax.scatter(x, y) sc.set_cmap("brg") plt.xlim(0, 1) plt.ylim(0, 1) ttl = ax.text(-0.05, 1.05, "", transform=ax.transAxes, fontsize=14) def animate(i): ill_prc = int(observ_all[(observ_all['time']==i) & (observ_all['color']==0.5)]['ID'].count()/population *100) ttl.set_text('vacc_ratio = {}, time = {}, ill_prc = {}'.format(vacc_ratio, i,ill_prc)) x = list(observ_all[observ_all['time']==i]['X']) y = list(observ_all[observ_all['time'] == i]['Y']) colors = list(observ_all[observ_all['time'] == i]['color']) sc.set_offsets(np.c_[x,y]) sc.set_array(colors) ani = matplotlib.animation.FuncAnimation(fig, animate,frames=timer-1, interval=1, repeat=True) # Сохраняем анимацию в виде gif файла: my_writer=matplotlib.animation.PillowWriter(fps = 8) ani.save('vacc_ratio='+str(vacc_ratio)+'.gif',writer=my_writer)

Теперь, изменяя параметры vacc_ratio, vacc_infect_prb, not_vacc_infect_prb, distination, ill_durat можно наблюдать за охватом и скоростью распространения инфекции. Изменяя параметр vacc_ratio c 30% до 90%, я построила следующие графики:

Те же результаты, но в табличном виде (в абсолютных и относительных значениях):

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