2020-04-12

新型肺炎COVID-19の感染者数をgooglevis を使って表示する

こんな記事を観測した。
www.ft.com
新型肺炎の感染者数および死亡者数を、ECDCという機関のデータから取得してプロットしている。
Homepage | European Centre for Disease Prevention and Control
この機関、いいところがcsv データを整備してくれていて、かつ、Rでデータを使えるようにデータ取得のスクリプトまで書いてくれているという優しさ。

library(utils)
#read the Dataset sheet into “R”. The dataset will be called "data".
data <- read.csv("https://opendata.ecdc.europa.eu/covid19/casedistribution/csv", na.strings = "", fileEncoding = "UTF-8-BOM")

前述のサイトでは7日間の移動平均を取っている。Rでやるなら

filter(d, rep(1, 7))/7

でできる。そのままパクるのは芸がないので時系列でもやろうと思ったが200カ国近くやるのは時間がかかるのでspline で適当に補完した。
前述のサイトにならって最初に30人以上の感染者が報告された日付を0にそろえてプロットしている。
googlevisならインタラクティブに、見たい国のデータを触れば見れるのに、と思ってやった。
googlevis のgvisLineChart関数のオプションはここを見ながら頑張る。
Line Chart | Charts | Google Developers

LineChartID1192a1e30d4

Data from European Centre for Disease Prevention and Control (ECDC) on 2020-04-12

最初の感染報告から今日までを見たいと思ったのでベタにやってみた。見にくい。

LineChartID11922e2c0b29

Data from European Centre for Disease Prevention and Control (ECDC) on 2020-04-12

library(utils)
#read the Dataset sheet into “R”. The dataset will be called "data".
data <- read.csv("https://opendata.ecdc.europa.eu/covid19/casedistribution/csv", na.strings = "", fileEncoding = "UTF-8-BOM")

# spline
library(stringr)
library(googleVis)

daily_case <- 30
p <- vector("list", length(s))
for(i in seq(s)){
  y <- rev(s[[i]]$case)
  if(length(y) > 8){
    sp <- smooth.spline(y, spar=0.5)
    y1 <- round(replace(sp$y, sp$y < 0, 0), 2)
    p[[i]] <- cbind(x=sp$x - which(y > daily_case)[1] + 1, y=y1)
  }
}
names(p) <- names(s)

lm <- apply(do.call(rbind, p), 2, range, na.rm=TRUE)
X <- seq(lm[1,1], lm[2,1])
Y0 <- mapply(function(z) replace(rep(NA, length(X)), na.omit(match(z[,1], X)), z[,2]), p)
#Y0 <- log10(Y0)
#Y0 <- replace(Y0, Y0 == -Inf, NA)
Y <- cbind.data.frame(X, Y0)
Y <- Y[-c(1:min(unlist(sapply(apply(!is.na(Y[,-1]), 2, which), head, 1)))), ]
Y <- subset(Y, X >= 0)

url <- "https://opendata.ecdc.europa.eu/covid19/casedistribution/csv"
webpage <- "European Centre for Disease Prevention and Control (ECDC)"
headertxt <- sprintf("Data from <span><a href='%s' target='_blank'>%s</a></span> on %s", url, webpage, Sys.Date())
# googlevis 用のプロットオプション
ops <- list(
  fontSize=12, height=800, width=800, lineWidth=1,
  vAxis="{title:'Number [log10]', fontSize: 20, titleTextStyle: {fontSize: 20, bold: 'true', italic: 'false'}, viewWindow: {min: 0, max: 40000}, gridlines: {color: '#00F', minSpacing: 100}, logScale: 'true'}",
  hAxis="{title:'Days since average daily cases passed 30', fontSize: 5, titleTextStyle: {fontSize: 10, bold: 'true', italic: 'false'}, viewWindow: {min: 0, max: 85}}",
  explorer="{actions: ['dragToZoom', 'rightClickToReset']}",
  focusTarget="datum",
  selectionMode="multiple",
  legend="{position: 'right', textStyle: {color: 'black', fontSize: 10}, pageIndex: 3}"
)
gmap <- gvisLineChart(Y, options=ops)

gmap$html$caption <- NULL
#gmap$html$footer <- NULL
gmap$html$footer <- str_replace(gmap$html$footer,
                                "<span> \n  .+\n  .+\n</span>",
                                headertxt)
plot(gmap)


# 日付すべてプロットする
daily_case <- 30
XX <- as.Date(sprintf("%d-%02d-%02d", data$year, data$month, data$day))
X <- seq.Date(min(XX), max(XX), by=1)
p <- vector("list", length(s))
for(i in seq(s)){
  y <- rev(s[[i]]$case)
  if(length(y) > 8){
    x <- as.Date(sprintf("%d-%02d-%02d", s[[i]]$year, s[[i]]$month, s[[i]]$day))
    idx <- match(rev(x), X)
    sp <- smooth.spline(idx, y, spar=0.5)
    px <- min(idx):max(idx)
    py <- predict(sp, px)
    y1 <- round(replace(py$y, py$y < 0, 0), 2)
    p[[i]] <- cbind.data.frame(x=seq.Date(min(x), max(x), by=1), y=y1)
  }
}
names(p) <- names(s)
lm <- do.call(rbind, p)
Y0 <- mapply(function(z) replace(rep(NA, length(X)), na.omit(match(z[,1], X)), z[,2]), p)
#Y0 <- log10(Y0)
#Y0 <- replace(Y0, Y0 == -Inf, NA)
Y <- cbind.data.frame(X, Y0)
Y <- Y[-c(1:min(unlist(sapply(apply(!is.na(Y[,-1]), 2, which), head, 1)))), ]
colnames(Y)[1] <- "Date"

Y$Date <- as.character(Y$Date)

url <- "https://opendata.ecdc.europa.eu/covid19/casedistribution/csv"
webpage <- "European Centre for Disease Prevention and Control (ECDC)"
headertxt <- sprintf("Data from <span><a href='%s' target='_blank'>%s</a></span> on %s", url, webpage, Sys.Date())
# googlevis 用のプロットオプション
ops <- list(
  fontSize=12, height=800, width=800, lineWidth=1,
  vAxis="{title:'Number [log10]', fontSize: 20, titleTextStyle: {fontSize: 20, bold: 'true', italic: 'false'}, viewWindow: {min: 0, max: 40000}, gridlines: {color: '#00F', minSpacing: 100}, logScale: 'true'}",
  sizeAxis="{minValue: -80, maxValue: 75, maxSize: 15}",
  hAxis="{title:'Date', fontSize: 5, titleTextStyle: {fontSize: 10, bold: 'true', italic: 'false'}, viewWindow: {min: 0, max: 105}}",
  explorer="{actions: ['dragToZoom', 'rightClickToReset']}",
  focusTarget="datum",
  selectionMode="multiple",
  legend="{position: 'right', textStyle: {color: 'black', fontSize: 10}, pageIndex: 3}"
)
gmap <- gvisLineChart(Y, options=ops)

gmap$html$caption <- NULL
#gmap$html$footer <- NULL
gmap$html$footer <- str_replace(gmap$html$footer,
                                "<span> \n  .+\n  .+\n</span>",
                                headertxt)
plot(gmap)