covid.md

   1 ---
   2 jupyter:
   3   jupytext:
   4     formats: ipynb,md
   5     text_representation:
   6       extension: .md
   7       format_name: markdown
   8       format_version: '1.2'
   9       jupytext_version: 1.3.4
  10   kernelspec:
  11     display_name: Python 3
  12     language: python
  13     name: python3
  14 ---
  15
  16 Data from [European Centre for Disease Prevention and Control](https://www.ecdc.europa.eu/en/publications-data/download-todays-data-geographic-distribution-covid-19-cases-worldwide)
  17
  18 ```python
  19 import itertools
  20 import collections
  21 import json
  22 import pandas as pd
  23 import numpy as np
  24 from scipy.stats import gmean
  25 import datetime
  26
  27 import matplotlib as mpl
  28 import matplotlib.pyplot as plt
  29 %matplotlib inline
  30 ```
  31
  32 ```python
  33 DEATH_COUNT_THRESHOLD = 10
  34 COUNTRIES_CORE = 'IT DE UK ES IE FR'.split()
  35 COUNTRIES_NORDIC = 'SE NO DK FI UK'.split()
  36 # COUNTRIES_FRIENDS = 'IT UK ES BE SI MX'.split()
  37 COUNTRIES_FRIENDS = 'IT UK ES BE SI PT'.split()
  38
  39 COUNTRIES_AMERICAS = ['AG', 'AR', 'AW', 'BS', 'BB', 'BZ', 'BM', 'BO', 'BR', 'VG', 'KY', # excluding Canada and USA
  40        'CL', 'CO', 'CR', 'CU', 'CW', 'DM', 'DO', 'EC', 'SV', 'GL', 'GD', 'GT',
  41        'GY', 'HT', 'HN', 'JM', 'MX', 'MS', 'NI', 'PA', 'PY', 'PE', 'PR', 'KN',
  42        'LC', 'VC', 'SX', 'SR', 'TT', 'TC', 'VI', 'UY', 'VE']
  43 COUNTRIES_OF_INTEREST = list(set(COUNTRIES_CORE + COUNTRIES_FRIENDS))
  44 COUNTRIES_ALL = list(set(COUNTRIES_CORE + COUNTRIES_FRIENDS + COUNTRIES_NORDIC + COUNTRIES_AMERICAS))
  45 ```
  46
  47 ```python
  48 !curl https://opendata.ecdc.europa.eu/covid19/casedistribution/csv/ > covid.csv
  49 ```
  50
  51 ```python
  52 # First col is a date, treat geoId of NA as 'Namibia', not "NA" value
  53 raw_data = pd.read_csv('covid.csv',
  54                        parse_dates=[0], dayfirst=True,
  55                        keep_default_na=False, na_values = [''],
  56 #                        dtype = {'day': np.int64,
  57 #                                 'month': np.int64,
  58 #                                 'year': np.int64,
  59 #                                 'cases': np.int64,
  60 #                                 'deaths': np.int64,
  61 #                                 'countriesAndTerritories': str,
  62 #                                 'geoId': str,
  63 #                                 'countryterritoryCode': str,
  64 #                                 'popData2019': np.int64,
  65 #                                 'continentExp': str,
  66 #                                 }
  67                       )
  68 ```
  69
  70 ```python
  71 raw_data.size
  72 ```
  73
  74 ```python
  75 raw_data.fillna(0, inplace=True)
  76 ```
  77
  78 ```python
  79 raw_data.head()
  80 ```
  81
  82 ```python
  83 raw_data.dtypes
  84 ```
  85
  86 ```python
  87 raw_data = raw_data.astype({'dateRep': np.datetime64,
  88     'day': np.int64,
  89     'month': np.int64,
  90     'year': np.int64,
  91     'cases': np.int64,
  92     'deaths': np.int64,
  93     'countriesAndTerritories': str,
  94     'geoId': str,
  95     'countryterritoryCode': str,
  96     'popData2019': np.int64,
  97     'continentExp': str })
  98 ```
  99
 100 ```python
 101 raw_data.dtypes
 102 ```
 103
 104 ```python
 105 raw_data[((raw_data.geoId == 'UK') & (raw_data.dateRep >= '2020-07-10'))]
 106 ```
 107
 108 ```python
 109 # raw_data = raw_data[~ ((raw_data.geoId == 'ES') & (raw_data.dateRep >= '2020-05-22'))]
 110 ```
 111
 112 ```python
 113 base_data = raw_data.set_index(['geoId', 'dateRep'])
 114 base_data.sort_index(inplace=True)
 115 base_data
 116 ```
 117
 118 ```python
 119 base_data.loc['ES'].loc['2020-05-10':]
 120 ```
 121
 122 ```python
 123 countries = raw_data[['geoId', 'countriesAndTerritories', 'popData2019', 'continentExp']]
 124 countries = countries[countries['popData2019'] != '']
 125 countries = countries.drop_duplicates()
 126 countries.set_index('geoId', inplace=True)
 127 countries = countries.astype({'popData2019': 'int64'})
 128 countries.head()
 129 ```
 130
 131 ```python
 132 countries.shape
 133 ```
 134
 135 ```python
 136 countries[countries.countriesAndTerritories == 'Finland']
 137 ```
 138
 139 ```python
 140 countries.loc[COUNTRIES_OF_INTEREST]
 141 ```
 142
 143 ```python
 144 countries[countries.continentExp == 'America'].index
 145 ```
 146
 147 ```python
 148 data_by_date = base_data[['cases', 'deaths']]
 149 data_by_date.head()
 150 ```
 151
 152 ```python
 153 data_by_date.loc['UK']
 154 ```
 155
 156 ```python
 157 # data_by_date.deaths.drop_duplicates().sort_values().to_csv('dth.csv', header=True)
 158 ```
 159
 160 ```python
 161 data_by_date.groupby(level=0).cumsum()
 162 ```
 163
 164 ```python
 165 data_by_date = data_by_date.merge(
 166     data_by_date.groupby(level=0).cumsum(),
 167     suffixes=('', '_culm'),
 168     left_index=True, right_index=True)
 169 data_by_date
 170 ```
 171
 172 ```python
 173 data_by_date = data_by_date.merge(
 174     data_by_date[['cases', 'deaths']].groupby(level=0).diff(),
 175     suffixes=('', '_diff'),
 176     left_index=True, right_index=True)
 177 data_by_date
 178 ```
 179
 180 ```python
 181 data_by_date.loc['UK', '2020-04-17']
 182 ```
 183
 184 ```python
 185 data_by_date.loc['UK']
 186 ```
 187
 188 ```python
 189 # data_by_date[data_by_date.deaths_culm > DEATH_COUNT_THRESHOLD]
 190 ```
 191
 192 ```python
 193 # days_since_threshold = data_by_date[data_by_date.deaths_culm > DEATH_COUNT_THRESHOLD].groupby(level=0).cumcount()
 194 # days_since_threshold.rename('since_threshold', inplace=True)
 195 ```
 196
 197 ```python
 198 dbd = data_by_date[data_by_date.deaths_culm > DEATH_COUNT_THRESHOLD].reset_index(level=1)
 199 dbd['since_threshold'] = dbd.dateRep
 200 dbd.set_index('dateRep', append=True, inplace=True)
 201 dbd.sort_index(inplace=True)
 202 days_since_threshold = dbd.groupby(level=0).diff().since_threshold.dt.days.fillna(0).astype(int).groupby(level=0).cumsum()
 203 # days_since_threshold.groupby(level=0).cumsum()
 204
 205 # days_since_threshold = dbd.rename('since_threshold')
 206 days_since_threshold
 207 ```
 208
 209 ```python
 210 # days_since_threshold = (data_by_date[data_by_date.deaths_culm > DEATH_COUNT_THRESHOLD]
 211 #                         .reset_index(level=1).groupby(level=0)
 212 #                         .diff().dateRep.dt.days
 213 #                         .groupby(level=0).cumcount()
 214 #                        )
 215 # days_since_threshold.rename('since_threshold', inplace=True)
 216 # days_since_threshold
 217 ```
 218
 219 ```python
 220 data_since_threshold = data_by_date.merge(days_since_threshold,
 221     left_index=True, right_index=True)
 222 data_since_threshold
 223 ```
 224
 225 ```python
 226 data_since_threshold = data_since_threshold.set_index('since_threshold', append=True
 227                               ).reorder_levels(['since_threshold', 'geoId', 'dateRep']
 228                                               ).reset_index('dateRep')
 229 data_since_threshold.sort_index(inplace=True)
 230 data_since_threshold
 231 ```
 232
 233 ```python
 234 data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT']), :]
 235 ```
 236
 237 ```python
 238 data_since_threshold.loc[(slice(None), ['ES']), :].tail(8)
 239 ```
 240
 241 ```python
 242 data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT']), ['deaths_culm']].unstack().plot(logy=True)
 243 ```
 244
 245 ```python
 246 # deaths = data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT', 'IE']), ['deaths_culm']].unstack().xs('deaths_culm', axis=1, drop_level=True)
 247 ```
 248
 249 ```python
 250 deaths = data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['deaths_culm']].unstack().sort_index().xs('deaths_culm', axis=1, drop_level=True)
 251 ```
 252
 253 ```python
 254 cases = data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['cases_culm']].unstack().sort_index().xs('cases_culm', axis=1, drop_level=True)
 255 ```
 256
 257 ```python
 258 COUNTRIES_AMERICAS_DEAD = list(set(deaths.columns) & set(COUNTRIES_AMERICAS))
 259 ```
 260
 261 ```python
 262 data_since_threshold.reset_index().merge(countries, on='geoId').set_index(['since_threshold', 'geoId'])
 263 ```
 264
 265 ```python
 266 data_since_threshold.reset_index().merge(countries, on='geoId').set_index(['since_threshold', 'geoId']).sort_index(inplace=True)
 267 ```
 268
 269 ```python
 270 data_since_threshold_per_capita = data_since_threshold.reset_index().merge(countries, on='geoId').set_index(['since_threshold', 'geoId'])
 271 data_since_threshold_per_capita['cases_culm_pc'] = data_since_threshold_per_capita.cases_culm / data_since_threshold_per_capita.popData2019
 272 data_since_threshold_per_capita['deaths_culm_pc'] = data_since_threshold_per_capita.deaths_culm / data_since_threshold_per_capita.popData2019
 273 data_since_threshold_per_capita
 274 ```
 275
 276 ```python
 277 deaths_pc = data_since_threshold_per_capita.loc[(slice(None), ['UK', 'DE', 'IT', 'IE']), ['deaths_culm_pc']].unstack().sort_index().xs('deaths_culm_pc', axis=1, drop_level=True)
 278 ```
 279
 280 ```python
 281 deaths_pc.index
 282 ```
 283
 284 ```python
 285 deaths_pc = data_since_threshold_per_capita.loc[(slice(None), COUNTRIES_ALL), ['deaths_culm_pc']].unstack().xs('deaths_culm_pc', axis=1, drop_level=True)
 286 ```
 287
 288 ```python
 289 deaths[COUNTRIES_CORE].plot()
 290 ```
 291
 292 ```python
 293 deaths[COUNTRIES_FRIENDS].plot()
 294 ```
 295
 296 ```python
 297 ax = deaths[COUNTRIES_FRIENDS].plot(figsize=(10, 6), title="Total deaths, linear")
 298 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
 299 for c in COUNTRIES_FRIENDS:
 300     lvi = deaths[c].last_valid_index()
 301     ax.text(x = lvi + 1, y = deaths[c][lvi], s = f"{c}: {deaths[c][lvi]:.0f}")
 302 # plt.savefig('covid_deaths_total_linear.png')
 303 ```
 304
 305 ```python
 306 ax = deaths[COUNTRIES_CORE].plot(figsize=(10, 6), title="Total deaths, linear")
 307 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
 308 for c in COUNTRIES_CORE:
 309     lvi = deaths[c].last_valid_index()
 310     ax.text(x = lvi + 1, y = deaths[c][lvi], s = f"{c}: {deaths[c][lvi]:.0f}")
 311 plt.savefig('covid_deaths_total_linear.png')
 312 ```
 313
 314 ```python
 315 deaths_prime = deaths[COUNTRIES_CORE].copy()
 316 deaths_prime.loc[73:, 'ES'] = np.NaN
 317 # deaths_prime['ES'][70:]
 318 ```
 319
 320 ```python
 321 ax = deaths_prime[COUNTRIES_CORE].plot(figsize=(10, 6), title="Total deaths, linear")
 322 for c in COUNTRIES_CORE:
 323     lvi = deaths_prime[c].last_valid_index()
 324     ax.text(x = lvi + 1, y = deaths_prime[c][lvi], s = f"{c}: {deaths_prime[c][lvi]:.0f}")
 325 # plt.savefig('covid_deaths_total_linear.png')
 326 ```
 327
 328 ```python
 329 ax = cases[COUNTRIES_CORE].plot(figsize=(10, 6), title="Total cases, linear")
 330 for c in COUNTRIES_CORE:
 331     lvi = cases[c].last_valid_index()
 332     ax.text(x = lvi + 1, y = cases[c][lvi], s = c)
 333 plt.savefig('covid_cases_total_linear.png')
 334 ```
 335
 336 ```python
 337 ax = deaths[COUNTRIES_AMERICAS_DEAD].plot(figsize=(10, 6), title="Total deaths, linear")
 338 for c in COUNTRIES_AMERICAS_DEAD:
 339     lvi = deaths[c].last_valid_index()
 340     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 341 # plt.savefig('covid_deaths_total_linear.png')
 342 ```
 343
 344 ```python
 345 ax = deaths[COUNTRIES_CORE + ['BR', 'MX']].plot(figsize=(10, 6), title="Total deaths, linear")
 346 for c in COUNTRIES_CORE + ['BR', 'MX']:
 347     lvi = deaths[c].last_valid_index()
 348     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 349 # plt.savefig('covid_deaths_total_linear.png')
 350 ```
 351
 352 ```python
 353 ax = deaths[COUNTRIES_NORDIC].plot(figsize=(10, 6), title="Total deaths, linear")
 354 for c in COUNTRIES_NORDIC:
 355     lvi = deaths[c].last_valid_index()
 356     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 357 # plt.savefig('covid_deaths_total_linear.png')
 358 ```
 359
 360 ```python
 361 ax = deaths[COUNTRIES_OF_INTEREST].plot(figsize=(10, 6), title="Total deaths, linear")
 362 for c in COUNTRIES_OF_INTEREST:
 363     lvi = deaths[c].last_valid_index()
 364     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 365 plt.savefig('covid_deaths_total_linear_of_interest.png')
 366 ```
 367
 368 ```python
 369 ax = deaths[COUNTRIES_CORE].plot(logy=True, figsize=(10, 6), title="Total deaths, log")
 370 for c in COUNTRIES_CORE:
 371     lvi = deaths[c].last_valid_index()
 372     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 373
 374 plt.savefig('covid_deaths_total_log.png')
 375 ```
 376
 377 ```python
 378 ylim = (5*10**3, 5*10**4)
 379 ax = deaths[COUNTRIES_CORE].plot(logy=True, figsize=(10, 6), ylim=ylim, title="Total deaths, log")
 380 for c in COUNTRIES_CORE:
 381     lvi = deaths[c].last_valid_index()
 382     if ylim[0] < deaths[c][lvi] < ylim[1]:
 383         ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 384
 385 # plt.savefig('covid_deaths_total_log.png')
 386 ```
 387
 388 ```python
 389 ax = deaths[COUNTRIES_FRIENDS].plot(logy=True, figsize=(10, 6), title="Total deaths, log")
 390 for c in COUNTRIES_FRIENDS:
 391     lvi = deaths[c].last_valid_index()
 392     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 393
 394 # plt.savefig('covid_deaths_total_log.png')
 395 ```
 396
 397 ```python
 398 ax = deaths[COUNTRIES_NORDIC].plot(logy=True, figsize=(10, 6), title="Total deaths, log")
 399 for c in COUNTRIES_NORDIC:
 400     lvi = deaths[c].last_valid_index()
 401     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 402
 403 # plt.savefig('covid_deaths_total_log.png')
 404 ```
 405
 406 ```python
 407 ax = deaths[COUNTRIES_OF_INTEREST].plot(logy=True, figsize=(10, 6), title="Total deaths, log")
 408 for c in COUNTRIES_OF_INTEREST:
 409     lvi = deaths[c].last_valid_index()
 410     ax.text(x = lvi + 1, y = deaths[c][lvi], s = c)
 411
 412 plt.savefig('covid_deaths_total_log.png')
 413 ```
 414
 415 ```python
 416 ax = deaths_pc.plot(figsize=(10, 6), title="Deaths per capita, linear")
 417 for c in deaths_pc.columns:
 418     lvi = deaths_pc[c].last_valid_index()
 419     ax.text(x = lvi + 1, y = deaths_pc[c][lvi], s = c)
 420 plt.savefig('covid_deaths_per_capita_linear.png')
 421 ```
 422
 423 ```python
 424 ax = deaths_pc.plot(logy=True, figsize=(10, 6), title="Deaths per capita, log")
 425 for c in deaths_pc.columns:
 426     lvi = deaths_pc[c].last_valid_index()
 427     ax.text(x = lvi + 1, y = deaths_pc[c][lvi], s = c)
 428 ```
 429
 430 ```python
 431 deaths_pc[['UK', 'IE']].plot( figsize=(10, 6), title="Deaths per capita, linear")
 432 ```
 433
 434 ```python
 435 deaths_pc[['UK', 'IE']].plot(logy=True, figsize=(10, 6), title="Deaths per capita, log")
 436 ```
 437
 438 ```python
 439 deaths[['UK', 'ES', 'IT']].plot(logy=True, figsize=(10, 6), title="Deaths, log")
 440 plt.savefig('covid_deaths_selected_log.png')
 441 ```
 442
 443 ```python
 444 deaths[['UK', 'ES', 'IT', 'MX']].plot(logy=True, figsize=(10, 6), title="Deaths, log")
 445 ```
 446
 447 ```python
 448 data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT']), :]
 449 ```
 450
 451 ```python
 452 data_since_threshold['deaths_m4'] = data_since_threshold.groupby(level=1)['deaths'].transform(lambda x: x.rolling(4, 1).mean())
 453 data_since_threshold['deaths_m7'] = data_since_threshold.groupby(level=1)['deaths'].transform(lambda x: x.rolling(7, 1).mean())
 454 data_since_threshold['cases_m7'] = data_since_threshold.groupby(level=1)['cases'].transform(lambda x: x.rolling(7, 1).mean())
 455 # data_since_threshold['deaths_diff_m4'] = data_since_threshold.groupby(level=1)['deaths_diff'].transform(lambda x: x.rolling(4, 1).mean())
 456 # data_since_threshold['deaths_diff_m7'] = data_since_threshold.groupby(level=1)['deaths_diff'].transform(lambda x: x.rolling(7, 1).mean())
 457 data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT']), :]
 458 ```
 459
 460 ```python
 461 deaths_m4 = (data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['deaths_m4']]
 462              .unstack().sort_index().xs('deaths_m4', axis=1, drop_level=True))
 463 ```
 464
 465 ```python
 466 deaths_m7 = (data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['deaths_m7']]
 467              .unstack().sort_index().xs('deaths_m7', axis=1, drop_level=True))
 468 ```
 469
 470 ```python
 471 cases_m7 = (data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['cases_m7']]
 472              .unstack().sort_index().xs('cases_m7', axis=1, drop_level=True))
 473 ```
 474
 475 ```python
 476 ax = deaths_m4.plot(figsize=(10, 6), title="Deaths per day, 4 day moving average")
 477 for c in deaths_m4.columns:
 478     lvi = deaths_m4[c].last_valid_index()
 479     ax.text(x = lvi + 1, y = deaths_m4[c][lvi], s = c)
 480 plt.savefig('covid_deaths_per_day.png')
 481 ```
 482
 483 ```python
 484 ax = deaths_m4[COUNTRIES_CORE].plot(figsize=(10, 6), title="Deaths per day, 4 day moving average")
 485 for c in COUNTRIES_CORE:
 486     lvi = deaths_m4[c].last_valid_index()
 487     ax.text(x = lvi + 1, y = deaths_m4[c][lvi], s = c)
 488 plt.savefig('covid_deaths_per_day-core.png')
 489 ```
 490
 491 ```python
 492 ax = deaths_m4[COUNTRIES_FRIENDS].plot(figsize=(10, 6), title="Deaths per day, 4 day moving average")
 493 for c in COUNTRIES_FRIENDS:
 494     lvi = deaths_m4[c].last_valid_index()
 495     ax.text(x = lvi + 1, y = deaths_m4[c][lvi], s = c)
 496 plt.savefig('covid_deaths_per_day-friends.png')
 497 ```
 498
 499 ```python
 500 C7s = 'ES FR IT UK'.split()
 501 ax = deaths_m7[C7s].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average")
 502 for c in C7s:
 503     lvi = deaths_m7[c].last_valid_index()
 504     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = c)
 505 # plt.savefig('covid_deaths_per_day-friends.png')
 506 ```
 507
 508 ```python
 509 ax = deaths_m7[COUNTRIES_CORE].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average")
 510 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
 511 for c in COUNTRIES_CORE:
 512     lvi = deaths_m7[c].last_valid_index()
 513     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = c)
 514 # plt.axhline(0, color='0.7')
 515 plt.savefig('covid_deaths_per_day_7.png')
 516 ```
 517
 518 ```python
 519 ax = deaths_m7[COUNTRIES_FRIENDS].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average")
 520 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
 521 for c in COUNTRIES_FRIENDS:
 522     lvi = deaths_m7[c].last_valid_index()
 523     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = c)
 524 # plt.axhline(0, color='0.7')
 525 # plt.savefig('covid_deaths_per_day_7.png')
 526 ```
 527
 528 ```python
 529 deaths_m7_prime = deaths_m7[COUNTRIES_CORE].copy()
 530 deaths_m7_prime.loc[73:, 'ES'] = np.NaN
 531 deaths_m7_prime['ES'][70:]
 532 ```
 533
 534 ```python
 535 ax = deaths_m7_prime[COUNTRIES_CORE].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average")
 536 for c in COUNTRIES_CORE:
 537     lvi = deaths_m7_prime[c].last_valid_index()
 538     ax.text(x = lvi + 1, y = deaths_m7_prime[c][lvi], s = c)
 539 # plt.savefig('covid_deaths_per_day_7.png') # see below for where this is written, with the projection
 540 ```
 541
 542 ```python
 543 ax = deaths_m7[COUNTRIES_FRIENDS].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average")
 544 for c in COUNTRIES_FRIENDS:
 545     lvi = deaths_m7[c].last_valid_index()
 546     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = c)
 547 plt.savefig('covid_deaths_per_day_friends_7.png')
 548 ```
 549
 550 ```python
 551 ax = deaths_m7[COUNTRIES_CORE + ['BR', 'MX']].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average")
 552 for c in COUNTRIES_CORE + ['BR', 'MX']:
 553     lvi = deaths_m7[c].last_valid_index()
 554     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = c)
 555 # plt.savefig('covid_deaths_per_day_7.png')
 556 ```
 557
 558 ```python
 559 ax = cases_m7[COUNTRIES_CORE].plot(figsize=(10, 6), title="Cases per day, 7 day moving average")
 560 for c in COUNTRIES_CORE:
 561     lvi = cases_m7[c].last_valid_index()
 562     ax.text(x = lvi + 1, y = cases_m7[c][lvi], s = c)
 563 plt.savefig('covid_cases_per_day-core.png')
 564 ```
 565
 566 ```python
 567 ax = cases_m7[COUNTRIES_FRIENDS].plot(figsize=(10, 6), title="Cases per day, 7 day moving average")
 568 for c in COUNTRIES_FRIENDS:
 569     lvi = cases_m7[c].last_valid_index()
 570     ax.text(x = lvi + 1, y = cases_m7[c][lvi], s = c)
 571 # plt.savefig('covid_cases_per_day-core.png')
 572 ```
 573
 574 ```python
 575 def gmean_scale(items):
 576     return gmean(items) / items[-1]
 577 ```
 578
 579 ```python
 580 def doubling_time(df):
 581     return np.log(2) / np.log((df.deaths_culm + df.deaths_g4) / df.deaths_culm)
 582
 583 def doubling_time_7(df):
 584     return np.log(2) / np.log((df.deaths_culm + df.deaths_g7) / df.deaths_culm)
 585 ```
 586
 587 ```python
 588 # data_since_threshold['deaths_g4'] = data_since_threshold.groupby(level=1)['deaths'].transform(lambda x: x.rolling(4, 1).apply(gmean_scale, raw=True))
 589 # data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT']), :]
 590 ```
 591
 592 ```python
 593 data_since_threshold['deaths_g4'] = data_since_threshold.groupby(level=1)['deaths'].transform(lambda x: x.rolling(4, 1).apply(gmean, raw=True))
 594 data_since_threshold['deaths_g7'] = data_since_threshold.groupby(level=1)['deaths'].transform(lambda x: x.rolling(7, 1).apply(gmean, raw=True))
 595 data_since_threshold.loc[(slice(None), ['UK', 'DE', 'IT']), :]
 596 ```
 597
 598 ```python
 599 data_since_threshold['doubling_time'] = data_since_threshold.groupby(level=1).apply(doubling_time).reset_index(level=0, drop=True).sort_index()
 600 data_since_threshold['doubling_time_7'] = data_since_threshold.groupby(level=1).apply(doubling_time_7).reset_index(level=0, drop=True).sort_index()
 601 # data_since_threshold.loc[(slice(None), 'UK'), :]
 602 ```
 603
 604 ```python
 605 doubling_times = (data_since_threshold.loc[(slice(None), COUNTRIES_OF_INTEREST), ['doubling_time']]
 606              .unstack().sort_index().xs('doubling_time', axis=1, drop_level=True))
 607 doubling_times.replace([np.inf, -np.inf], np.nan, inplace=True)
 608 ```
 609
 610 ```python
 611 doubling_times_7 = (data_since_threshold.loc[(slice(None), COUNTRIES_OF_INTEREST), ['doubling_time_7']]
 612              .unstack().sort_index().xs('doubling_time_7', axis=1, drop_level=True))
 613 doubling_times_7.replace([np.inf, -np.inf], np.nan, inplace=True)
 614 ```
 615
 616 ```python
 617 ax = doubling_times.plot(figsize=(10, 6), title="Doubling times, 4 day average")
 618 for c in doubling_times.columns:
 619     lvi = doubling_times[c].last_valid_index()
 620     ax.text(x = lvi + 1, y = doubling_times[c][lvi], s = c)
 621 # plt.savefig('covid_deaths_per_day.png')
 622 ```
 623
 624 ```python
 625 ax = doubling_times_7[COUNTRIES_CORE].plot(figsize=(10, 6), title="Doubling times, 7 day average")
 626 ax.legend(loc="upper left")
 627 for c in COUNTRIES_CORE:
 628     lvi = doubling_times_7[c].last_valid_index()
 629     ax.text(x = lvi + 1, y = doubling_times_7[c][lvi], s = c)
 630 plt.savefig('covid_doubling_times_7.png')
 631 ```
 632
 633 ```python
 634 ax = doubling_times[COUNTRIES_CORE].plot(figsize=(10, 6), title="Doubling times, 4 day average")
 635 for c in COUNTRIES_CORE:
 636     lvi = doubling_times[c].last_valid_index()
 637     ax.text(x = lvi + 1, y = doubling_times[c][lvi], s = c)
 638 plt.savefig('covid_doubling_times.png')
 639 ```
 640
 641 ```python
 642 ax = doubling_times[COUNTRIES_FRIENDS].plot(figsize=(10, 6), title="Doubling times")
 643 for c in COUNTRIES_FRIENDS:
 644     lvi = doubling_times[c].last_valid_index()
 645     ax.text(x = lvi + 1, y = doubling_times[c][lvi], s = c)
 646 plt.savefig('covid_doubling_times_friends.png')
 647 ```
 648
 649 ```python
 650 ax = doubling_times[C7s].plot(figsize=(10, 6), title="Doubling times")
 651 for c in C7s:
 652     lvi = doubling_times[c].last_valid_index()
 653     ax.text(x = lvi + 1, y = doubling_times[c][lvi], s = c)
 654 # plt.savefig('covid_doubling_times_friends.png')
 655 ```
 656
 657 ```python
 658 # deaths_diff_m4 = (data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['deaths_diff_m4']]
 659 #              .unstack().sort_index().xs('deaths_diff_m4', axis=1, drop_level=True))
 660 ```
 661
 662 ```python
 663 # deaths_diff_m7 = (data_since_threshold.loc[(slice(None), COUNTRIES_ALL), ['deaths_diff_m7']]
 664 #              .unstack().sort_index().xs('deaths_diff_m7', axis=1, drop_level=True))
 665 ```
 666
 667 ```python
 668 # deaths_diff_m7
 669 ```
 670
 671 ```python
 672 # data_since_threshold.replace([np.inf, -np.inf], np.nan).groupby(level=1).last().loc[COUNTRIES_ALL]#, [doubling_time]]
 673 ```
 674
 675 ```python
 676 dstl = data_since_threshold.replace([np.inf, -np.inf], np.nan).groupby(level=1).last()
 677 dstl.loc[dstl.index.intersection(COUNTRIES_ALL)]
 678 ```
 679
 680 ```python
 681 # data_since_threshold.replace([np.inf, -np.inf], np.nan).groupby(level=1).last().loc[['UK', 'DE', 'IT']]#, [doubling_time]]
 682 dstl.loc[['UK', 'DE', 'IT', 'FR', 'ES']]
 683 ```
 684
 685 ```python
 686 data_since_threshold.loc[(slice(None), ['UK']), :].tail(20)
 687 ```
 688
 689 ```python
 690 data_since_threshold.loc[(slice(None), ['ES']), :].tail(20)
 691 ```
 692
 693 ## Death projections
 694
 695 ```python
 696 data_since_threshold.loc[(slice(None), ['UK']), :].tail(15)
 697 ```
 698
 699 ```python
 700 it_since_threshold = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['IT']), :]
 701 s_end = it_since_threshold.index.max()[0]
 702 s_end
 703 ```
 704
 705 ```python
 706 uk_projection = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), :]
 707 uk_current_end = uk_projection.index.max()[0] + 1
 708 # s_start = uk_projection.index.max()[0] + 1
 709 uk_current_end
 710 ```
 711
 712 ```python
 713 current_uk_deaths_m7 = uk_projection[uk_projection.deaths_m7 >= 0].iloc[-1].deaths_m7
 714 current_uk_deaths_m7
 715 ```
 716
 717 ```python
 718 it_since_threshold[it_since_threshold.deaths_m7 <= current_uk_deaths_m7].loc[60:].first_valid_index()[0]
 719 ```
 720
 721 ```python
 722 s_start = it_since_threshold[it_since_threshold.deaths_m7 <= current_uk_deaths_m7].loc[60:].first_valid_index()[0]
 723 s_start
 724 ```
 725
 726 ```python
 727 s_start_date = data_since_threshold.loc[(89, 'IT'), 'dateRep']# .iloc[0]
 728 s_start_date
 729 ```
 730
 731 ```python
 732 s_end - s_start
 733 ```
 734
 735 ```python
 736 uk_end = s_end - s_start + uk_current_end
 737 uk_end
 738 ```
 739
 740 ```python
 741 proj = it_since_threshold.loc[(slice(s_start, s_end), slice(None)), ['cases', 'deaths', 'deaths_m7']]
 742 ndiff = uk_current_end - s_start
 743 proj.index = pd.MultiIndex.from_tuples([(n + ndiff, 'UK') for n, _ in proj.index], names=proj.index.names)
 744 proj
 745 ```
 746
 747 ```python
 748 it_since_threshold.loc[(slice(s_start - 8, s_start + 2), slice(None)), ['cases', 'deaths', 'deaths_m7']]
 749 ```
 750
 751 ```python
 752 uk_projection[['cases', 'deaths', 'deaths_m7']].tail()
 753 ```
 754
 755 ```python
 756 # proj['deaths_m7'] = proj['deaths_m7'] + 20
 757 # proj
 758 ```
 759
 760 Projected deaths, UK following IT trend from now.
 761
 762 ```python
 763 uk_projection = uk_projection.append(proj, sort=True)
 764 uk_projection.deaths.sum()
 765 ```
 766
 767 ```python
 768 uk_projection = uk_projection.droplevel(1)
 769 uk_projection
 770 ```
 771
 772 ```python
 773 uk_projection.loc[152, 'deaths']
 774 ```
 775
 776 ## Correction for cumulative deaths correction on 14 August
 777
 778 ```python
 779 # uk_projection.loc[152, 'deaths'] = 50
 780 ```
 781
 782 ```python
 783 uk_projection['deaths_m7'] = uk_projection['deaths'].transform(lambda x: x.rolling(7, 1).mean())
 784 uk_projection.loc[(uk_current_end - 20):(uk_current_end + 5)]
 785 ```
 786
 787 ```python
 788 uk_projection.loc[(uk_current_end - 5):]
 789 ```
 790
 791 ```python
 792 uk_projection.deaths_m7.plot()
 793 ```
 794
 795 ```python
 796 proj.droplevel(level=1)
 797 ```
 798
 799 ```python
 800 ax = deaths_m7[COUNTRIES_CORE].plot()
 801 # uk_projection['deaths_m7'].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average", label="Projection", style='--', ax=ax)
 802 proj.droplevel(level=1)['deaths_m7'].plot(figsize=(10, 6), title="Deaths per day, 7 day moving average", label="Projection", style='--', ax=ax)
 803 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
 804 for c in COUNTRIES_CORE:
 805     lvi = deaths_m7[c].last_valid_index()
 806     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = c)
 807 plt.savefig('covid_deaths_per_day_7.png')
 808 ```
 809
 810 ```python
 811 it_since_threshold.deaths.sum()
 812 ```
 813
 814 # Excess deaths calculation
 815
 816 ```python
 817 with open('excess_deaths.json') as f:
 818     excess_deaths_data = json.load(f)
 819
 820 with open('excess_death_accuracy.json') as f:
 821     excess_death_accuracy = json.load(f)
 822
 823 excess_deaths_data, excess_death_accuracy
 824 ```
 825
 826 ```python
 827 # excess_deaths_upto = '2020-05-08'
 828 # excess_deaths = 54500
 829 ```
 830
 831 ```python
 832 excess_deaths_upto = excess_deaths_data['end_date']
 833 excess_deaths = excess_deaths_data['excess_deaths']
 834 ```
 835
 836 Recorded deaths in period where ONS has reported total deaths
 837
 838 ```python
 839 ons_reported_deaths = base_data.loc['UK'][:excess_deaths_upto]['deaths'].sum()
 840 ons_reported_deaths
 841 ```
 842
 843 ```python
 844 excess_deaths_upto
 845 ```
 846
 847 ## Correction for deaths total correction on 14 August
 848
 849 ```python
 850 ons_unreported_deaths_data = base_data.loc['UK'][excess_deaths_upto:].iloc[1:]['deaths']
 851 # ons_unreported_deaths_data['2020-08-14'] = 50
 852 ```
 853
 854 ```python
 855 ons_unreported_deaths = ons_unreported_deaths_data.sum()
 856 ons_unreported_deaths
 857 ```
 858
 859 ```python
 860 scaled_ons_unreported_deaths = ons_unreported_deaths * excess_death_accuracy
 861 scaled_ons_unreported_deaths
 862 ```
 863
 864 ```python
 865 uk_deaths_to_date = excess_deaths + scaled_ons_unreported_deaths
 866 uk_deaths_to_date
 867 ```
 868
 869 ```python
 870 # data_since_threshold.loc[(slice(None), 'UK'), :][data_since_threshold.dateRep == excess_deaths_data['end_date']]
 871 ```
 872
 873 ```python
 874 data_since_threshold[data_since_threshold.dateRep == excess_deaths_data['end_date']].loc[(slice(None), 'UK'), :]
 875 ```
 876
 877 ```python
 878 ons_unreported_start = data_since_threshold[data_since_threshold.dateRep == excess_deaths_data['end_date']].loc[(slice(None), 'UK'), :].first_valid_index()[0] + 1
 879 ons_unreported_start
 880 ```
 881
 882 ```python
 883 unreported_projected_deaths = uk_projection.loc[ons_unreported_start:].deaths.sum()
 884 unreported_projected_deaths
 885 ```
 886
 887 ```python
 888 unreported_projected_deaths_scaled = unreported_projected_deaths * excess_death_accuracy
 889 unreported_projected_deaths_scaled
 890 ```
 891
 892 ```python
 893 uk_projection.loc[(s_start):].deaths.sum()
 894 ```
 895
 896 ```python
 897 deaths_actual_projected_scaled = uk_deaths_to_date + uk_projection.loc[(s_start):].deaths.sum() * excess_death_accuracy
 898 deaths_actual_projected_scaled
 899 ```
 900
 901 ```python
 902 # excess_deaths / reported_deaths
 903 ```
 904
 905 True deaths to date, if we follow the scaling of excess deaths over reported deaths so far.
 906
 907 ```python
 908 uk_covid_deaths = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), 'deaths_culm'].iloc[-1]
 909 uk_covid_deaths
 910 ```
 911
 912 ```python
 913 # uk_covid_deaths_scaled = excess_deaths + unreported_deaths * excess_death_accuracy
 914 # uk_covid_deaths_scaled
 915 ```
 916
 917 ```python
 918 # data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['IT']), 'dateRep'].iloc[-1] + pd.Timedelta(s_end - s_start, unit='days')
 919 ```
 920
 921 ```python
 922 # data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), 'dateRep'].iloc[-1].strftime("%Y-%m-%d")
 923 ```
 924
 925 ```python
 926 # uk_covid_deaths * excess_deaths / reported_deaths
 927 ```
 928
 929 ```python
 930 # uk_projection.deaths.sum() * excess_deaths / reported_deaths
 931 ```
 932
 933 ```python
 934 # data_since_threshold.loc[(slice(None), 'FR'), :]
 935 # data_since_threshold[data_since_threshold.dateRep == '2020-05-18'].loc[(slice(None), 'FR'), :]
 936 ```
 937
 938 ## School reopenings
 939
 940 ```python
 941 school_reopenings = {
 942     'ES': {'date': '2020-05-18'},
 943     'FR': {'date': '2020-05-18'}, # some areas only
 944 #     'IT': {'date': '2020-09-01'},
 945     # 'IE': {'date': '2020-09-01'},
 946     'DE': {'date': '2020-05-04'},
 947     'UK': {'date': '2020-06-01'}
 948 }
 949 ```
 950
 951 ```python
 952 data_since_threshold[data_since_threshold.dateRep == '2020-05-04'].loc[(slice(None), ['DE']), :].first_valid_index()
 953 ```
 954
 955 ```python
 956 data_since_threshold[data_since_threshold.dateRep == '2020-05-04'].loc[(slice(None), ['DE']), :].iloc[0].deaths_m7
 957 ```
 958
 959 ```python
 960 for cID in school_reopenings:
 961     dst_in = data_since_threshold[data_since_threshold.dateRep == (school_reopenings[cID]['date'])].loc[(slice(None), [cID]), :]
 962     dst_i = dst_in.first_valid_index()
 963     dst_n = dst_in.iloc[0].deaths_m7
 964     school_reopenings[cID]['since_threshold'] = dst_i[0]
 965     school_reopenings[cID]['deaths_m7'] = dst_n
 966 school_reopenings
 967 ```
 968
 969 ```python
 970 ax = deaths_m7[COUNTRIES_CORE].plot(figsize=(15, 9), title="Deaths per day, 7 day moving average")
 971 # uk_projection.deaths_m7.plot(ax=ax)
 972 for c in COUNTRIES_CORE:
 973     lvi = deaths_m7[c].last_valid_index()
 974     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = f"{c}: {deaths_m7[c][lvi]:.0f}")
 975     if c in school_reopenings:
 976         marker_col = [l for l in ax.lines if l.get_label() == c][0].get_color()
 977         ax.plot(school_reopenings[c]['since_threshold'], school_reopenings[c]['deaths_m7'], '*',
 978                 markersize=18, markerfacecolor=marker_col, markeredgecolor=marker_col)
 979         ax.text(x = school_reopenings[c]['since_threshold'] + 1, y = school_reopenings[c]['deaths_m7'],
 980                 s = f"{school_reopenings[c]['date']}: {school_reopenings[c]['deaths_m7']:.0f}")
 981 plt.savefig('school_reopenings.png')
 982 ```
 983
 984 ```python
 985 # ax = deaths_m7[COUNTRIES_CORE].plot(figsize=(15, 9), title="Deaths per day, 7 day moving average",
 986 #                                    xlim=(46, 91), ylim=(0, 400))
 987 # # uk_projection.deaths_m7.plot(ax=ax)
 988 # for c in COUNTRIES_CORE:
 989 #     lvi = deaths_m7[c].last_valid_index()
 990 #     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = f"{c}: {deaths_m7[c][lvi]:.0f}", fontsize=14)
 991 #     if c in school_reopenings:
 992 #         marker_col = [l for l in ax.lines if l.get_label() == c][0].get_color()
 993 #         ax.plot(school_reopenings[c]['since_threshold'], school_reopenings[c]['deaths_m7'], '*',
 994 #                 markersize=18, markerfacecolor=marker_col, markeredgecolor=marker_col)
 995 #         ax.text(x = school_reopenings[c]['since_threshold'] + 1, y = school_reopenings[c]['deaths_m7'],
 996 #                 s = f"{school_reopenings[c]['date']}: {school_reopenings[c]['deaths_m7']:.0f}",
 997 #                 fontsize=14)
 998 # plt.savefig('school_reopenings_detail.png')
 999 ```
1000
1001 # Lockdown graphs
1002
1003 ```python
1004 lockdown_dates = {
1005     'ES': { 'part_start': {'date': '2020-03-14'}
1006           , 'full_start': {'date': '2020-03-15'}
1007           , 'part_finish': {'date': '2020-05-18'}
1008           },
1009     'FR': { 'part_start': {'date': '2020-03-13'}
1010           , 'full_start': {'date': '2020-03-17'}
1011           , 'part_finish': {'date': '2020-05-11'}
1012           },
1013     'IT': { 'part_start': {'date': '2020-03-08'}
1014           , 'full_start': {'date': '2020-03-10'}
1015           , 'part_finish': {'date': '2020-05-04'}
1016           },
1017     'DE': { #'part_start': {'date': '2020-03-13'}
1018           'full_start': {'date': '2020-03-22'}
1019           , 'part_finish': {'date': '2020-05-06'}
1020           },
1021     'UK': { 'part_start': {'date': '2020-03-23'}
1022           , 'full_start': {'date': '2020-03-23'}
1023           , 'part_finish': {'date': '2020-05-31'}
1024           },
1025     'IE': { #'part_start': {'date': '2020-03-12'}
1026           'full_start': {'date': '2020-03-27'}
1027           , 'part_finish': {'date': '2020-05-18'}
1028           },
1029 }
1030 ```
1031
1032 ```python
1033 for cID in lockdown_dates:
1034     for phase in lockdown_dates[cID]:
1035         dst_in = data_since_threshold[data_since_threshold.dateRep == (lockdown_dates[cID][phase]['date'])].loc[(slice(None), [cID]), :]
1036         dst_i = dst_in.first_valid_index()
1037         dst_n = dst_in.iloc[0].deaths_m7
1038         dst_c = dst_in.iloc[0].cases_m7
1039         lockdown_dates[cID][phase]['since_threshold'] = dst_i[0]
1040         lockdown_dates[cID][phase]['deaths_m7'] = dst_n
1041         lockdown_dates[cID][phase]['cases_m7'] = dst_c
1042
1043 lockdown_dates
1044 ```
1045
1046 ```python
1047 ax = deaths_m7[COUNTRIES_CORE].plot(figsize=(15, 9), title="Deaths per day, 7 day moving averagee, with lockdown dates")
1048 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
1049 # uk_projection.deaths_m7.plot(ax=ax)
1050 for c in COUNTRIES_CORE:
1051     lvi = deaths_m7[c].last_valid_index()
1052     if c != 'UK':
1053         ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = f"{c}: {deaths_m7[c][lvi]:.0f}")
1054     if c in lockdown_dates:
1055         for phase in lockdown_dates[c]:
1056             marker_col = [l for l in ax.lines if l.get_label() == c][0].get_color()
1057             ax.plot(lockdown_dates[c][phase]['since_threshold'], lockdown_dates[c][phase]['deaths_m7'], '*',
1058                     markersize=18, markerfacecolor=marker_col, markeredgecolor=marker_col)
1059             if 'start' not in phase:
1060                 ax.text(x = lockdown_dates[c][phase]['since_threshold'] + 1, y = lockdown_dates[c][phase]['deaths_m7'],
1061                         s = f"{lockdown_dates[c][phase]['date']}: {lockdown_dates[c][phase]['deaths_m7']:.0f}")
1062 # plt.savefig('school_reopenings.png')
1063 ```
1064
1065 ```python
1066 ax = cases_m7.iloc[-50:][COUNTRIES_CORE].plot(figsize=(15, 9), title="Cases per day, 7 day moving average, with lockdown dates") #, ylim=(-10, 1500))
1067 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
1068 # uk_projection.deaths_m7.plot(ax=ax)
1069 for c in COUNTRIES_CORE:
1070     lvi = cases_m7[c].last_valid_index()
1071 #     if c != 'UK':
1072     ax.text(x = lvi + 1, y = cases_m7[c][lvi], s = f"{c}: {cases_m7[c][lvi]:.0f}")
1073
1074 ```
1075
1076 ```python
1077 ax = cases_m7[COUNTRIES_CORE].plot(figsize=(15, 9), title="Cases per day, 7 day moving average, with lockdown dates")
1078 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
1079 # uk_projection.deaths_m7.plot(ax=ax)
1080 for c in COUNTRIES_CORE:
1081     lvi = cases_m7[c].last_valid_index()
1082 #     if c != 'UK':
1083     ax.text(x = lvi + 1, y = cases_m7[c][lvi], s = f"{c}: {cases_m7[c][lvi]:.0f}")
1084     if c in lockdown_dates:
1085         for phase in lockdown_dates[c]:
1086             marker_col = [l for l in ax.lines if l.get_label() == c][0].get_color()
1087             if 'start' in phase:
1088                 marker_shape = '^'
1089             else:
1090                 marker_shape = 'v'
1091             ax.plot(lockdown_dates[c][phase]['since_threshold'], lockdown_dates[c][phase]['cases_m7'],
1092                     marker_shape,
1093                     markersize=18, markerfacecolor=marker_col, markeredgecolor=marker_col)
1094             if 'start' not in phase:
1095                 ax.text(x = lockdown_dates[c][phase]['since_threshold'] + 1, y = lockdown_dates[c][phase]['cases_m7'],
1096                         s = f"{lockdown_dates[c][phase]['date']}: {lockdown_dates[c][phase]['cases_m7']:.0f}")
1097 plt.savefig('cases_per_day_with_lockdown.png')
1098 ```
1099
1100 ```python
1101 ax = cases_m7[COUNTRIES_CORE].plot(figsize=(10, 6), title="Cases per day, 7 day moving average")
1102 for c in COUNTRIES_CORE:
1103     lvi = cases_m7[c].last_valid_index()
1104     ax.text(x = lvi + 1, y = cases_m7[c][lvi], s = c)
1105 plt.savefig('covid_cases_per_day-core.png')
1106 ```
1107
1108 ```python
1109 ax = deaths_m7[COUNTRIES_CORE].plot(figsize=(15, 9), title="Deaths per day, 7 day moving average",
1110                                    xlim=(0, 15),
1111                                     ylim=(0, 66)
1112                                    )
1113 # uk_projection.deaths_m7.plot(ax=ax)
1114 for c in COUNTRIES_CORE:
1115     lvi = deaths_m7[c].last_valid_index()
1116     if c in lockdown_dates:
1117         for phase in lockdown_dates[c]:
1118             if 'start' in phase:
1119                 print(c, phase)
1120                 marker_col = [l for l in ax.lines if l.get_label() == c][0].get_color()
1121                 ax.plot(lockdown_dates[c][phase]['since_threshold'], lockdown_dates[c][phase]['deaths_m7'], '*',
1122                         markersize=18, markerfacecolor=marker_col, markeredgecolor=marker_col)
1123                 ax.text(x = lockdown_dates[c][phase]['since_threshold'] + 0.3, y = lockdown_dates[c][phase]['deaths_m7'],
1124                         s = f"{lockdown_dates[c][phase]['date']}: {lockdown_dates[c][phase]['deaths_m7']:.0f}")
1125 # plt.savefig('school_reopenings.png')
1126 ```
1127
1128 ```python
1129
1130 ```
1131
1132 ```python
1133
1134 ```
1135
1136 # Write results to summary file
1137
1138 ```python
1139 with open('covid_summary.md', 'w') as f:
1140     f.write('% Covid death data summary\n')
1141     f.write('% Neil Smith\n')
1142     f.write(f'% Created on {datetime.datetime.now().strftime("%Y-%m-%d")}\n')
1143     f.write('\n')
1144
1145     last_uk_date = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), 'dateRep'].iloc[-1]
1146     f.write(f'> Last UK data from {last_uk_date.strftime("%Y-%m-%d")}\n')
1147 ```
1148
1149 ```python
1150 with open('covid_summary.md', 'a') as f:
1151     f.write('\n')
1152     f.write('## Headlines\n')
1153     f.write('\n')
1154     f.write('| []() | |\n')
1155     f.write('|:---|---:|\n')
1156     f.write(f'| Deaths reported so far | {uk_covid_deaths} | \n')
1157     f.write(f'| Total Covid deaths to date (estimated) | {uk_deaths_to_date:.0f} |\n')
1158     projection_date = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['IT']), 'dateRep'].iloc[-1] + pd.Timedelta(s_end - s_start, unit='days')
1159     f.write(f'| Projected total deaths up to {projection_date.strftime("%Y-%m-%d")} | {deaths_actual_projected_scaled:.0f} | \n')
1160     f.write('\n')
1161 ```
1162
1163 ```python
1164 with open('covid_summary.md', 'a') as f:
1165     f.write('\n')
1166     f.write('## Total deaths\n')
1167     f.write(f'Time based on days since {DEATH_COUNT_THRESHOLD} deaths\n')
1168     f.write('\n')
1169     f.write('![Total deaths](covid_deaths_total_linear.png)\n')
1170     f.write('\n')
1171     f.write('| Country ID | Country name | Total deaths |\n')
1172     f.write('|:-----------|:-------------|-------------:|\n')
1173     for c in sorted(COUNTRIES_CORE):
1174         lvi = deaths[c].last_valid_index()
1175         f.write(f'| {c} | {countries.loc[c].countriesAndTerritories} | {int(deaths[c][lvi])} |\n')
1176     f.write('\n')
1177 ```
1178
1179 ```python
1180 with open('covid_summary.md', 'a') as f:
1181     f.write('\n')
1182     f.write('## All-causes deaths, UK\n')
1183     f.write('\n')
1184     f.write('![All-causes deaths](deaths-radar.png)\n')
1185     f.write('\n')
1186     f.write('### Excess deaths\n')
1187     f.write('\n')
1188     f.write(f'From week ending 20 March 2020 until week ending {pd.to_datetime(excess_deaths_upto).strftime("%d %B %Y")}, ')
1189     f.write(f'there were approximately **{excess_deaths:.0f}** excess deaths, over the average for the previous five years.\n')
1190     f.write('\n')
1191 ```
1192
1193 ```python
1194 with open('covid_summary.md', 'a') as f:
1195     f.write(f'In that period, the UK reported {ons_reported_deaths} Covid deaths.\n')
1196     f.write(f'In the last three weeks for which excess deaths have been reported, the excess deaths have been {excess_death_accuracy:.3f} higher than the Covid-reported deaths.\n')
1197 #     f.write(f'That means the actual number of Covid death is about {excess_deaths / reported_deaths:.2f} times higher than the reported figures.\n')
1198     f.write('\n')
1199     f.write(f'The UK has reported {uk_covid_deaths} deaths so far.\n')
1200     f.write(f'Using the scaling factor above (for Covid-19 deaths after the ONS figures), I infer that there have been **{uk_deaths_to_date:.0f}** total deaths so far.\n')
1201     f.write('\n')
1202 ```
1203
1204 ```python
1205 with open('covid_summary.md', 'a') as f:
1206     f.write('\n')
1207     f.write('## Deaths per day\n')
1208     f.write(f'Based on a 7-day moving average\n')
1209     f.write('\n')
1210     f.write('![Deaths per day](covid_deaths_per_day_7.png)\n')
1211     f.write('\n')
1212 ```
1213
1214 ```python
1215 s_end - s_start - 1
1216 ```
1217
1218 ```python
1219 with open('covid_summary.md', 'a') as f:
1220     f.write('\n')
1221     f.write('## Projected deaths\n')
1222     f.write(f"The UK's daily deaths data is very similar to Italy's.\n")
1223     f.write(f'If I use the Italian data for the next {s_end - s_start - 1} days (from {s_start_date.strftime("%d %B %Y")} onwards),')
1224     f.write(f' the UK will report {uk_projection.deaths.sum()} deaths on day {uk_end} of the epidemic.\n')
1225     f.write('\n')
1226     f.write('Using the excess deaths scaling from above, that will translate into ')
1227     f.write(f'**{deaths_actual_projected_scaled:.0f}** Covid deaths total.\n')
1228     f.write('\n')
1229 ```
1230
1231 ```python
1232 with open('covid_summary.md', 'a') as f:
1233     f.write('\n')
1234     f.write('## Deaths doubling times\n')
1235     f.write(f'Based on a 7-day moving average\n')
1236     f.write('\n')
1237     f.write('![Deaths doubling times](covid_doubling_times_7.png)\n')
1238     f.write('\n')
1239 ```
1240
1241 ```python
1242 with open('covid_summary.md', 'a') as f:
1243     f.write('\n')
1244     f.write('## Cases per day and lockdown dates\n')
1245     f.write(f'Based on a 7-day moving average\n')
1246     f.write('\n')
1247     f.write('![Cases per day](cases_per_day_with_lockdown.png)\n')
1248     f.write('\n')
1249 ```
1250
1251 ```python
1252 with open('covid_summary.md', 'a') as f:
1253     f.write('\n')
1254     f.write('| Country ID | Country name | Most recent daily cases | Most recent daily deaths |\n')
1255     f.write('|:-----------|:-------------|------------------------:|-------------------------:|\n')
1256     for c in sorted(COUNTRIES_CORE):
1257         lvic = cases_m7[c].last_valid_index()
1258         lvid = deaths_m7[c].last_valid_index()
1259         f.write(f'| {c} | {countries.loc[c].countriesAndTerritories} | {cases_m7[c][lvic]:.0f} | {deaths_m7[c][lvid]:.0f} | \n')
1260     f.write('\n')
1261     f.write('(Figures are 7-day averages)\n')
1262     f.write('\n')
1263 ```
1264
1265 ```python
1266 with open('covid_summary.md', 'a') as f:
1267     f.write('# Data sources\n')
1268     f.write('\n')
1269     f.write('> Covid data from [European Centre for Disease Prevention and Control](https://www.ecdc.europa.eu/en/publications-data/download-todays-data-geographic-distribution-covid-19-cases-worldwide)\n')
1270     f.write('\n')
1271     f.write("""> Population data from:
1272
1273 * [Office of National Statistics](https://www.ons.gov.uk/peoplepopulationandcommunity/birthsdeathsandmarriages/deaths/datasets/weeklyprovisionalfiguresondeathsregisteredinenglandandwales) (Endland and Wales) Weeks start on a Saturday.
1274 * [Northern Ireland Statistics and Research Agency](https://www.nisra.gov.uk/publications/weekly-deaths) (Northern Ireland). Weeks start on a Saturday. Note that the week numbers don't match the England and Wales data.
1275 * [National Records of Scotland](https://www.nrscotland.gov.uk/statistics-and-data/statistics/statistics-by-theme/vital-events/general-publications/weekly-and-monthly-data-on-births-and-deaths/weekly-data-on-births-and-deaths) (Scotland). Note that Scotland uses ISO8601 week numbers, which start on a Monday.""")
1276
1277     f.write('\n\n')
1278     f.write('> [Source code available](https://git.njae.me.uk/?p=covid19.git;a=tree)\n')
1279     f.write('\n')
1280
1281 ```
1282
1283 ```python
1284 !pandoc -s covid_summary.md > covid_summary.html
1285 ```
1286
1287 ```python
1288 !scp covid_summary.html neil@ogedei:/var/www/scripts.njae.me.uk/covid/index.html
1289 !scp covid_deaths_total_linear.png neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1290 !scp deaths-radar.png neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1291 !scp covid_deaths_per_day_7.png neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1292 !scp covid_doubling_times_7.png neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1293 !scp cases_per_day_with_lockdown.png neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1294 ```
1295
1296 ```python
1297 with open('uk_covid_deaths.js', 'w') as f:
1298     f.write(f"document.write('{uk_covid_deaths}');")
1299
1300 with open('estimated_total_deaths.js', 'w') as f:
1301     f.write(f"document.write('{uk_deaths_to_date:.0f}');")
1302
1303 with open('projection_date.js', 'w') as f:
1304     f.write(f"document.write(\'{projection_date.strftime('%d %B %Y')}\');")
1305
1306 with open('projected_deaths.js', 'w') as f:
1307     f.write(f"document.write('{uk_projection.deaths.sum():.0f}');")
1308
1309 with open('projected_excess_deaths.js', 'w') as f:
1310     f.write(f"document.write('{deaths_actual_projected_scaled:.0f}');")
1311
1312 with open('excess_deaths_upto.js', 'w') as f:
1313     f.write(f"document.write('{pd.to_datetime(excess_deaths_upto).strftime('%d %B %Y')}');")
1314
1315 with open('excess_deaths.js', 'w') as f:
1316     f.write(f"document.write('{excess_deaths:.0f}');")
1317
1318 with open('reported_deaths.js', 'w') as f:
1319     f.write(f"document.write('{ons_reported_deaths:.0f}');")
1320
1321 with open('scaling_factor.js', 'w') as f:
1322     f.write(f"document.write('{excess_death_accuracy:.2f}');")
1323
1324 with open('projection_length.js', 'w') as f:
1325     f.write(f"document.write('{s_end - s_start - 1}');")
1326
1327 with open('s_end.js', 'w') as f:
1328     f.write(f"document.write('{s_end}');")
1329
1330 s_start_date_str = s_start_date.strftime("%d %B %Y")
1331 with open('s_start_date.js', 'w') as f:
1332     f.write(f"document.write('{s_start_date_str}');")
1333
1334 with open('uk_end.js', 'w') as f:
1335     f.write(f"document.write('{uk_end}');")
1336
1337 with open('last_uk_date.js', 'w') as f:
1338     f.write(f"document.write('{pd.to_datetime(last_uk_date).strftime('%d %B %Y')}');")
1339 ```
1340
1341 ```python
1342 pd.to_datetime(excess_deaths_upto).strftime('%d %B %Y')
1343 ```
1344
1345 ```python
1346 !scp uk_covid_deaths.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1347 !scp estimated_total_deaths.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1348 !scp projection_date.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1349 !scp projected_deaths.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1350 !scp projected_excess_deaths.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1351 !scp excess_deaths_upto.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1352 !scp excess_deaths.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1353 !scp reported_deaths.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1354 !scp scaling_factor.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1355 !scp projection_length.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1356 !scp s_end.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1357 !scp s_start_date.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1358 !scp uk_end.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1359 !scp last_uk_date.js neil@ogedei:/var/www/scripts.njae.me.uk/covid/
1360 ```
1361
1362 ```python
1363 data_by_date.loc['UK'].to_csv('data_by_day_uk.csv', header=True, index=True)
1364 ```
1365
1366 ```python
1367 ukd = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), ['deaths', 'deaths_m7']].droplevel(1)
1368 ax = ukd.deaths.plot.bar(figsize=(12, 8))
1369 ukd.deaths_m7.plot.line(ax=ax, color='red')
1370 # ax = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), 'deaths_m7'].plot.line(figsize=(12, 8), color='red')
1371 # ax = data_since_threshold.replace([np.inf, -np.inf], np.nan).loc[(slice(None), ['UK']), 'deaths'].plot.bar(ax=ax)
1372 ax.set_xticks(range(0, 120, 20))
1373 ```
1374
1375 ```python
1376 np.arange(0, 130, 20)
1377 ```
1378
1379 ```python
1380 data_by_date.loc['UK']
1381 ```
1382
1383 ```python
1384 data_by_date.loc['UK'].plot(x='deaths_culm', y='deaths', logx=True, logy=True)
1385 ```
1386
1387 ```python
1388 data_by_date.loc['UK'].plot(x='cases_culm', y='cases')
1389 ```
1390
1391 ```python
1392 ukdbd = data_by_date.loc['UK'].copy()
1393 ukdbd['deaths_m7'] = ukdbd.deaths.transform(lambda x: x.rolling(7, 1).mean())
1394 ukdbd['cases_m7'] = ukdbd.cases.transform(lambda x: x.rolling(7, 1).mean())
1395 ukdbd
1396 ```
1397
1398 ```python
1399 ukdbd.plot(x='deaths_culm', y='deaths_m7', logx=True, logy=True)
1400 ```
1401
1402 ```python
1403 fig, ax = plt.subplots(figsize=(12, 8))
1404 xmax = 10
1405 for c in COUNTRIES_CORE:
1406     if data_since_threshold.loc[(slice(None), c), 'deaths_culm'].max() > xmax:
1407         xmax = data_since_threshold.loc[(slice(None), c), 'deaths_culm'].max()
1408     data_since_threshold.loc[(slice(None), c), :].plot(x='deaths_culm', y='deaths_m7', logx=True, logy=True, xlim=(10, xmax * 1.1), label=c, ax=ax)
1409 ```
1410
1411 ```python
1412 data_since_threshold.loc[(slice(None), 'UK'), 'deaths_culm'].max()
1413 ```
1414
1415 ```python
1416 countries.continentExp.unique()
1417 ```
1418
1419 ```python
1420 countries.loc['KW']
1421 ```
1422
1423 ```python
1424 data_by_date.groupby(level=0)['deaths'].shift(-25)
1425 ```
1426
1427 ```python
1428 offset_data = data_by_date.loc[:, ['cases']]
1429 offset_data['deaths'] = data_by_date.groupby(level=0)['deaths'].shift(-25)
1430 offset_data['cases_m7'] = offset_data.groupby(level=0)['cases'].transform(lambda x: x.rolling(7, 1).mean())
1431 offset_data['deaths_m7'] = offset_data['deaths'].dropna().groupby(level=0).transform(lambda x: x.rolling(7, 1).mean())
1432 offset_data['deaths_per_case'] = offset_data.deaths_m7 / offset_data.cases_m7
1433 offset_data
1434 ```
1435
1436 ```python
1437 deaths_m7
1438 ```
1439
1440 ```python
1441 offset_deaths_m7 = (offset_data.loc[COUNTRIES_ALL, ['deaths_m7']]
1442              .unstack().sort_index().xs('deaths_m7', axis=1, drop_level=True)).T.sort_index()
1443 offset_deaths_m7
1444 ```
1445
1446 ```python
1447 offset_deaths_m7['UK']
1448 ```
1449
1450 ```python
1451 data_since_threshold.loc[(slice(None), 'UK'), :].tail()
1452 ```
1453
1454 ```python
1455 countries.loc['PT']
1456 ```
1457
1458 ```python
1459 ax = cases_m7.iloc[-50:][COUNTRIES_FRIENDS].plot(figsize=(15, 9), title="Cases per day, 7 day moving average", ylim=(-10, 1500))
1460 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
1461 # uk_projection.deaths_m7.plot(ax=ax)
1462 for c in COUNTRIES_FRIENDS:
1463     lvi = cases_m7[c].last_valid_index()
1464     if c != 'ES':
1465         ax.text(x = lvi + 1, y = cases_m7[c][lvi], s = f"{c}: {cases_m7[c][lvi]:.0f}")
1466
1467 ```
1468
1469 ```python
1470 ax = deaths_m7.iloc[-50:][COUNTRIES_FRIENDS].plot(figsize=(15, 9), title="Cases per day, 7 day moving average", ylim=(-10, 100))
1471 ax.set_xlabel(f"Days since {DEATH_COUNT_THRESHOLD} deaths")
1472 # uk_projection.deaths_m7.plot(ax=ax)
1473 for c in COUNTRIES_FRIENDS:
1474     lvi = deaths_m7[c].last_valid_index()
1475 #     if c != 'ES':
1476     ax.text(x = lvi + 1, y = deaths_m7[c][lvi], s = f"{c}: {deaths_m7[c][lvi]:.0f}")
1477
1478 ```
1479
1480 ```python
1481
1482 ```