Ebola epidemiology data scraper

Note: The following post is current as of WHO report of 09-18-14, which includes data to 09-14-14. [I’ve altered this code a number of times because of the nearly constantly changing format and location of the WHO data.]
#######

I wanted to find certain statistics for the West African Ebolavirus (EBV) outbreak from its inception, e.g. recent case and death rates. By “recent” I mean at WHO GAR reporting intervals, typically a few days. But the sites where I’d expect to find such (WHO, CDC etc) didn’t have them, at least in a synthesized form. So, I wrote an R script to scrape, compile, rearrange, and plot the data, for any country or all together, as taken from two sources. Starting with the WHO July 1 GAR report, tables of new and cumulative cases and deaths, categorized by degree of certainty, for each of the three countries in the outbreak, are given. Wikipedia has a less detailed table for all dates from March 25 on. I used it to obtain total cases and deaths up to July 1.

Below are graphs of the recent, per-day (1) case rates, and (2) death rates, from March 22. Each shows the raw data (thinner line), and a loess-fitted trend (thicker line). Note that reporting issues are partly, perhaps largely, responsible for extreme fluctuations in the raw data.

library(XML);library(chron);library(gdata)
rm(list=ls())

## set working directory containing data table (as "ebola data.csv") 
setwd("")

## 1. Read data table:
options(stringsAsFactors = FALSE)
EBV = read.csv("ebola data.csv",header=T)
EBV[,4:7] = as.integer(unlist(EBV[,4:7]))
colnames(EBV) = c("Country", "Date", "Type", "C", "P", "S", "CPS")
for (i in 1:3) EBV[,i] = trim(EBV[,i])
EBV = EBV[order(EBV$Country, EBV$Type, EBV$Date),]
EBV$Date = as.Date(EBV$Date)
str(EBV)

## 2. Compute new cases from differences, compute time in days from first report (March 22)
dates = unique(EBV$Date); n.dates = length(dates); t1 = dates[1]
EBV$T.sum = as.integer(EBV$Date - t1)
EBV$dCPS = c(NA, diff(EBV$CPS)); EBV$dC = c(NA, diff(EBV$C)); EBV$dT = c(NA, diff(EBV$T.sum))
temp2 = rep(NA, nrow(EBV))
for (i in 2:nrow(EBV)) {if (EBV$Type[i] != EBV$Type[i-1]) temp2[i] = 1}
EBV$dCPS[temp2==1] = NA; EBV$dC[temp2==1] = NA; EBV$dT[temp2==1] = NA 
EBV$rCPS = EBV$dCPS/EBV$dT; EBV$rC = EBV$dC/EBV$dT

## 3. Compute and plot recent (since last update report) daily case and mortality rates for given country and case type 
# a) compute
strata = expand.grid(sort(unique(EBV$Country)), sort(unique(EBV$Type)), stringsAsFactors=F); strata
pick = which(strata[,1]=="All countries")			# choose desired country; from names given in strata
stiff = 1; flex = 1/2
i = pick[1]

 subset1 = EBV[which(EBV$Country==strata[i,1] & EBV$Type==strata[i,2]),][-1,]
 subset2 = subset1[-(1:21),c(8,12)]; rownames(subset2) = 1:nrow(subset2)
 smooth1 = loess(subset1$rCPS ~ subset1$T.sum, span=stiff)[[2]]
 smooth2 = loess(subset2$rCPS ~ subset2$T.sum, span=stiff)[[2]]

 mn = smooth2[which(smooth2==min(smooth2))];  mx = smooth2[which(smooth2==max(smooth2))] 
 (diff = mx-mn); subset2[,1] = subset2[,1] - min(subset2[,1]) + 1
 (days = max(subset2[,1]))
 (x = diff^(1/days)); x^days; log(x)
 (tests =  exp(seq(0.5*log(x), 1.5*log(x), length.out=201))) 
 testing = matrix(NA, nrow=nrow(subset2), ncol = length(tests))
 for (j in 1:ncol(testing)) testing[,j] = tests[j]^subset2[,1]
 
 a = (testing - subset2[,2])^2
 b = apply(a,2,sum)
 (daily.spread = tests[which(b==min(b))])
 (weekly.spread = daily.spread^(6:12))

# b) plot, with loess smoothing, reported cases and deaths, (except WHO-reported "New" cases)
 for (i in pick){
 subset3 = EBV[which(EBV$Country==strata[i,1] & EBV$Type==strata[i,2]),][-1,]
 smooth3 = loess(subset3$rCPS ~ subset3$T.sum, span=stiff)[[2]]
 windows()
 plot(subset3$T.sum, subset3$rCPS, type="l",
  xlab = "Days from March 22 2014", ylab = "Recent Daily Rate",
  main = paste(strata[i,1], tolower(strata[i,2]), "all reported", sep=", "))
  lines(subset3$T.sum, smooth3, lwd=4)
  print(smooth3)
}