と。

Github: https://github.com/8-u8

COVID-19の日本の感染者数の推移を曜日単位で眺めてみた話

01_統計モデリング

編集履歴

2020/04/22現在で、下記で行っている分析に用いたデータは削除されている模様です。
(gitリンクは更新しました)
現在のデータで再現する記事はきぬいとの気分次第で書きますので、気長にお待ちください。

今日は2本上げます

昨今国内でも猛威を振るう新型コロナウイルス。

東洋経済では国内の感染者データなど、細かい情報を可視化し、データソースも公開するなど、
かなり挑戦的でかつ社会貢献度の高い試みを行っています。

toyokeizai.net

データがオープンになっていたら触るしかない。触ってみました。

github.com

ご注意

この記事はあくまで公開されているデータに基づいて、
きぬいとが「勝手に眺めてみた」ものになります。
あくまで「オープンデータを触ってみた」という範疇でご覧ください。

新型コロナウイルスの対策については厚生労働省、
及び各自治体の方針に従ってください。

厚生労働省のサイト

注意書きは読んだな？

データは下記からDLします。リスペクトをもってStarしましょう(しました)。

github.com

目的は「日次での感染報告者数の推移を理解する」ことにあります。理由は知的好奇心。

時系列については2020年1月15日〜2020年3月24日までです。
元データは感染者ユニークとなっています。

2月12日までは不定期なデータとなっているので、ここまでのデータは一旦除外します。
また、データにはyyyymmdd形式での時系列がないのでさくっと作ります。
日次でみたいので、上記日付でgroup_byし、人数などをカウントしましょう。 その他、「2週間のラグ」とか、累積度数とか、見たいものを好き勝手作りましょう*1

そんな感じで書いたのが以下です。

library(tidyverse)

UseData_JP <- readr::read_csv("covid19/data/individuals.csv") %>% 
  dplyr::mutate(yyyymmdd = lubridate::ymd(paste0(.$`確定年`,"-",.$`確定月`,"-",.$`確定日`))) 

UseData_JP_2 <- UseData_JP %>% 
  dplyr::filter(yyyymmdd >= lubridate::ymd("2020-02-13")) %>% 
  group_by(yyyymmdd) %>% 
  summarise(amounts_by_day = n()) %>% 
  ungroup() %>% 
  dplyr::mutate(lags = dplyr::lag(amounts_by_day),
                lag_2wk = dplyr::lag(amounts_by_day, n=14),
                diff = amounts_by_day - lags,
                diff_2wk = amounts_by_day - lag_2wk,
                cum  = cumsum(amounts_by_day),
                ratio_by_lag = amounts_by_day/lags,
                ratio_by_lag2wk = amounts_by_day / lag_2wk,
                inc_ratio = amounts_by_day/sum(amounts_by_day),
                pareto = cum / sum(amounts_by_day),
                weekday = lubridate::wday(yyyymmdd))

感染者数の推移

早速可視化していきましょう。

plot(UseData_JP_2$yyyymmdd, UseData_JP_2$diff_wkmean, type="b")

f:id:kinuit:20200331192317p:plain
日別感染報告者数

昨今報道でも指摘されているように、増加傾向にはあります。
ただ、細かく見ると若干上下があり、「指数関数的上昇」という感じでもないようにも見えます。
「休日とか報告者少なくなりそう？」とか思ったので、まずは土日(青、赤)、
週明けの動きも気になるので、月火(それぞれ緑、オレンジ)で切ってみましょう。

上のplotに以下を加筆します。

abline(v=UseData_JP_2$yyyymmdd[UseData_JP_2$weekday==7], col ="blue") # sat
abline(v=UseData_JP_2$yyyymmdd[UseData_JP_2$weekday==1], col ="red")  # sun
abline(v=UseData_JP_2$yyyymmdd[UseData_JP_2$weekday==2], col ="green") # mon
abline(v=UseData_JP_2$yyyymmdd[UseData_JP_2$weekday==3], col ="orange")# tue

f:id:kinuit:20200331192356p:plain
曜日の挙動

なんとなく、「オレンジの火曜日」にハネる傾向がありそうです。
この傾向は他の曜日の規則性よりもなかなか顕著に見えます。
とはいえあくまで「なんとなく」なので、もうちょっと詳しく見たいなあと思っていた矢先にTJOさんからこんな感じでアドバイスが。

ほう……と思いつつ。 つまりはこの変化の挙動を、「1週間の規則性」と「トレンド」に分解してみればよいんじゃない？
という話です。これはRではstl関数で簡単に分解ができてしまいます。

早速実装します。

trendline <- stl(ts(as.numeric(UseData_JP_2$amounts_by_day), frequency = 7), s.window = 'per')$time.series
UseData_JP_2 <- UseData_JP_2 %>% 
  dplyr::mutate(seasonal_shift = trendline[,1],
                trend_shift = trendline[,2],
                # amount_diff_seasonal = amounts_by_day - seasonal,
                # amount_diff_trend = amounts_by_day - trend,
  )


plot(UseData_JP_2$yyyymmdd, UseData_JP_2$seasonal_shift, type="b")

1週間の規則性の結果はこう。

f:id:kinuit:20200331192324p:plain

さらに印をつけるとこう。

f:id:kinuit:20200331193630p:plain

「火曜日にハネる」というよりは、「日、月が異常に低くなる」というトレンドのようです。

ちなみにトレンドはこんな感じで、やはり上昇トレンド(しかも波がある？)にあるようです。

f:id:kinuit:20200331193958p:plain

「ふーん、そんな感じか」と思っていました。
今日(2020年3月31日)は都内で78人の感染が確認されました。
一方、昨日は13人と、非常に少なかったです。

……ずいぶん一致するなあ。

そう思っているとこんなご指摘がありました。

それもそうですね、検査から報告までのタイムラグ、というのはあるかもしれないなと。

で？

とりあえず我々ができることは、「正しく恐れてできることをやる」ということに尽きると思っています。
ただ、その数字がある規則性を持っていることを知らないと、単純な増減で一喜一憂してしまい、疲れちゃいます。
今回の検討からは「月曜日の結果は参考までに、より現実的な結果は火曜日まで待とう」という学びが(あくまで個人的に)得られました。
難局ではありますが、できることからやって、日常を取り戻していきましょう。

*1:好き勝手作った結果あんまり見ないんですけど。