Метка: графики

Коронавирус: предсказание на Python

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

КДПВ: вирус + Китай + график

Нам понадобятся библиотеки:

  • pandas – для загрузки данных
  • matplotlib – для построения графика
  • scipy – для построения предсказания
  • numpy – для работы с массивами

Установите их, если они еще не установлены, набрав в терминале:

pip install pandas matplotlib scipy numpy

Подключим модули:

import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit
from datetime import timedelta

Я собрал данные по заболевшим в последние дни в CSV файл и закачал его на свой сайт по адресу: https://tirinox.ru/wp-content/uploads/2020/01/corona.csv

Pandas умеет читать данные прямо из интернета по URL:

df = pd.read_csv('https://tirinox.ru/wp-content/uploads/2020/01/corona.csv', parse_dates=['day'])

print(df.head())

""" Вывод:
         day  infected  dead
0 2020-01-20       278     4
1 2020-01-21       326     6
2 2020-01-22       547     8
3 2020-01-23       639    14
4 2020-01-24       916    25
"""

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

df.plot(kind='bar', x='day', y=['infected', 'dead'])
plt.show()
График роста вируса по дням
infected = зараженные люди, dead = погибшие

Очевидно, что прирост с каждым днем все больше и больше. Явно нелинейная зависимость! Это и понятно. Представим, что 1 человек в среднем заражает двоих. Те двое, каждый заразив еще 2, добавят 4 новых человека. А те 4 заразят 8. Затем 16, 32, 64. К этому нужно не забыть добавить уже больных.

Такой закон называется экспоненциальным. Опишем экспоненциальную функцию в параметрическом виде. Параметры a, b, c еще предстоит подобрать, чтобы максимально точно удовлетворить существующим данным:

def fit_func_exp(x, a, b, c):
    return a * np.exp(c * x) + b

Теперь вычислим оптимальные параметры, обеспечивающие наименьшую ошибку:

# зависимая переменная - ее будем предсказывать
infected = df['infected']

# дни - преобразуем их в целые числа от 0 до максимального
days = range(len(df['day']))

# у нас 3 параметра в функции: a, b, c – начальное приближение - единицы
p0 = np.ones(3)

# впишем кривую в наши данные
(a, b, c), _ = curve_fit(fit_func_exp, days, infected, p0=p0)

Теперь, зная параметры, рассчитаем функцию, скажем, до 60 дней с момента начала эпидемии:

# предскажем динамику вируса на 60 дней (начало соответствует 20 января)
MAX_DAY = 60

x_days = np.linspace(0, MAX_DAY - 1, MAX_DAY)
y_infect = fit_func_exp(x_days, a, b, c)

Построим график и убедимся, что он хорошо описывает экспериментальные данные:

plt.xlabel('Дни')
plt.ylabel('Больные')

# график в log шкале
plt.yscale('log')

# это данные текущей статистики - нарисуем их синими точками
plt.scatter(days, infected, marker='D', label='Реальные')

# это красная линия – предсказание (первые 22 дня)
plt.plot(x_days[:22], y_infect[:22], 'r', label='Предсказание')
plt.legend()

plt.show()

Вот, что у нас получится:

График предсказания.

Ого! Рост поражает: сотни, тысячи, десятки тысяч! Узнаем, на какой день число зараженных людей достигнет населения всей Земли:

# население Земли
EARTH_POPULATION = 7_530_000_000

# найдем номер дня, когда количество зараженных достигнет популяции Земли
doom_index = np.argmax(y_infect >= EARTH_POPULATION)
doom_day = x_days[doom_index]

# вычислим дату
day0 = df['day'][0]
doom_date = day0 + timedelta(days=int(doom_day))

# дата конца!
print(f'Doom date: {doom_date:%d.%m.%Y}')

Doom date: 13.03.2020

13 марта…

Полный код загрузил в gist.

UPD: добавил табличку:

+------------+------------------+----------------+
|    Дата    | Число заболевших | Число погибших |
+------------+------------------+----------------+
| 20.01.2020 |        54        |       7        |
| 21.01.2020 |       239        |       9        |
| 22.01.2020 |       492        |       11       |
| 23.01.2020 |       839        |       15       |
| 24.01.2020 |       1314       |       21       |
| 25.01.2020 |       1963       |       30       |
| 26.01.2020 |       2853       |       45       |
| 27.01.2020 |       4072       |       68       |
| 28.01.2020 |       5740       |      104       |
| 29.01.2020 |       8024       |      161       |
| 30.01.2020 |      11150       |      252       |
| 31.01.2020 |      15431       |      394       |
| 01.02.2020 |      21293       |      618       |
| 02.02.2020 |      29317       |      971       |
| 03.02.2020 |      40304       |      1528      |
| 04.02.2020 |      55347       |      2405      |
| 05.02.2020 |      75942       |      3787      |
| 06.02.2020 |      104139      |      5965      |
| 07.02.2020 |      142745      |      9398      |
| 08.02.2020 |      195601      |     14807      |
| 09.02.2020 |      267968      |     23331      |
| 10.02.2020 |      367047      |     36763      |
| 11.02.2020 |      502698      |     57930      |
| 12.02.2020 |      688423      |     91288      |
| 13.02.2020 |      942703      |     143854     |
| 14.02.2020 |     1290844      |     226690     |
| 15.02.2020 |     1767494      |     357228     |
| 16.02.2020 |     2420088      |     562938     |
| 17.02.2020 |     3313572      |     887108     |
| 18.02.2020 |     4536865      |    1397953     |
| 19.02.2020 |     6211706      |    2202971     |
| 20.02.2020 |     8504777      |    3471566     |
| 21.02.2020 |     11644280     |    5470689     |
| 22.02.2020 |     15942657     |    8621023     |
| 23.02.2020 |     21827679     |    13585498    |
| 24.02.2020 |     29885019     |    21408803    |
| 25.02.2020 |     40916535     |    33737214    |
| 26.02.2020 |     56020077     |    53165030    |
| 27.02.2020 |     76698735     |    83780495    |
| 28.02.2020 |    105010434     |   132026095    |
| 29.02.2020 |    143772730     |   208054274    |
| 01.03.2020 |    196843216     |   327863828    |
| 02.03.2020 |    269503423     |   516666580    |
| 03.03.2020 |    368984436     |   814192761    |
| 04.03.2020 |    505186524     |   1283051543   |
| 05.03.2020 |    691664408     |   2021906043   |
| 06.03.2020 |    946976216     |   3186235247   |
| 07.03.2020 |    1296530371    |   5021051835   |
| 08.03.2020 |    1775114218    |   7912460814   |
| 09.03.2020 |    2430356032    |  12468908548   |
| 10.03.2020 |    3327464945    |  19649219634   |
| 11.03.2020 |    4555720506    |  30964364746   |
| 12.03.2020 |    6237357709    |  48795417935   |
| 13.03.2020 |    8539731720    |  76894611953   |
| 14.03.2020 |   11691972927    |  121174929895  |
| 15.03.2020 |   16007789810    |  190954388899  |
| 16.03.2020 |   21916688952    |  300916853606  |
| 17.03.2020 |   30006719188    |  474201998218  |
| 18.03.2020 |   41082993743    |  747274645537  |
| 19.03.2020 |   56247814448    | 1177598150075  |
+------------+------------------+----------------+

Предостережение

Данные расчеты чисто теоретические. Как будет в реальности, мне неизвестно. Но я уверен, что власти Китая и других стран справятся с эпидемией быстрее, чем она поразит весь мир.

Узнавать оперативно о статистике вируса можно через телеграм бот @NovelCoronaVirusBot.

Желаю вам здоровья и не подхватить даже обычной простуды! Меньше бывайте в людных местах, часто мойте руки с мылом, носите маски, при любых симптомах обращайтесь к врачу!

Если заболели, не пытайтесь скрыть болезнь, сбив температуру!

Специально для канала @pyway. Подписывайтесь на мой канал в Телеграм @pyway 👈 

​​Анимация Jupyter Notebook

Сегодня мы будем анимировать график прямо внутри Jupyter Notebook. Сперва сделаем плавную отрисовку графика. Переключим режим отображения графиков в notebook:

%matplotlib notebook

Импортируем все, что нужно:

import matplotlib.pyplot as plt
from matplotlib import animation
import numpy as np

Сгенерируем наши данные:

# время (200 точек)
t = np.linspace(0, 2 * np.pi, 200)
x = np.sin(t)  # синусоида

Создадим пустой график:

fig, ax = plt.subplots()
# пределы отображения
ax.axis([0, 2 * np.pi, -2, 2])
l, = ax.plot([], [])

Функция animate будет вызываться при отрисовка каждого кадра, аргумент i – номер кадра:

def animate(i):
    # рисуем данные только от 0 до i
    # на первом кадре будет 0 точек, 
    # а на последнем - все
    l.set_data(t[:i], x[:i])

Запускаем анимацию:

fps = 30  # карды в сек
# frames - число кадров анимации
ani = animation.FuncAnimation(fig, animate, frames=len(t), interval=1000.0 / fps)

Если мы хотим анимировать сами данные, например, заставить синусоиду «плясать», то на каждом шаге перегенерируем данные заново, используя переменную i:

def animate(i):
    x = np.sin(t - i / len(t) * np.pi * 2) * np.sin(t * 15)
    l.set_data(t, x)

Можно сохранить в GIF:

ani.save('myAnimation.gif', writer='imagemagick', fps=30)

Сам ноутбук я загрузил на GitHub, но поиграться онлайн с ним не получится, надо скачать себе и запустить локально. Анимированные графики отрисовываются в реальном времени, поэтому требуют достаточно много ресурсов. Пример 3D анимации:

Специально для канала @pyway. Подписывайтесь на мой канал в Телеграм @pyway 👈