Created
July 30, 2020 20:54
-
-
Save dgleich/c697c41e48b0f7a74bcf0a2ab6b3e768 to your computer and use it in GitHub Desktop.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
using Plots, DataFrames, CSV, Statistics | |
## Get positions for each state | |
# From https=> //github.com/jkbren/covid19-cases-grid-gartogram/blob/master/covid19-cases-grid-gartogram.ipynb | |
state_posx = Dict("ak"=> (0, 0), "me"=> (0, 10),"gu"=> (7, 0), "vi"=> (7, 9), "pr"=> (7, 8), "mp"=> (7, 1), | |
"vt"=> (1, 9), "nh"=> (1, 10),"wa"=> (2, 0), "id"=> (2, 1), "mt"=> (2, 2), "nd"=> (2, 3), "mn"=> (2, 4), | |
"il"=> (2, 5), "wi"=> (2, 6), "mi"=> (2, 7), "ny"=> (2, 8), "ri"=> (2, 9), "ma"=> (2, 10),"or"=> (3, 0), | |
"nv"=> (3, 1), "wy"=> (3, 2), "sd"=> (3, 3), "ia"=> (3, 4), "in"=> (3, 5), "oh"=> (3, 6), "pa"=> (3, 7), | |
"nj"=> (3, 8), "ct"=> (3, 9), "ca"=> (4, 0), "ut"=> (4, 1), "co"=> (4, 2), "ne"=> (4, 3), "mo"=> (4, 4), | |
"ky"=> (4, 5), "wv"=> (4, 6), "va"=> (4, 7), "md"=> (4, 8), "de"=> (4, 9), "az"=> (5, 1), "nm"=> (5, 2), | |
"ks"=> (5, 3), "ar"=> (5, 4), "tn"=> (5, 5), "nc"=> (5, 6), "sc"=> (5, 7), "dc"=> (5, 8), "ok"=> (6, 3), | |
"la"=> (6, 4), "ms"=> (6, 5), "al"=> (6, 6), "ga"=> (6, 7), "hi"=> (7, 0), "tx"=> (7, 3), "fl"=> (7, 7)) | |
## Downloaded from | |
# https://d14wlfuexuxgcm.cloudfront.net/covid/rt.csv | |
# on 2020-07-30 | |
# MD5 (rt.csv) = 5bf496e12311b32a3b2979ba403fcfff | |
data = CSV.read("rt.csv") | |
## Check that we can map all states in the data | |
states = unique(data.region) | |
map(x -> state_posx[lowercase(x)], states) | |
## Look at data as of May 1 | |
df = filter(row -> row.date >= Dates.Date("2020-05-01"), data) | |
## Extract state-by-date data as a matrix. | |
# There is probably a better way of doing this... | |
function build_state_by_date_matrix(df, col::Symbol; states=nothing) | |
if states == nothing | |
states = sort(unique(df.region)) | |
end | |
dates = sort(unique(df.date)) | |
state2index = Dict( states .=> 1:length(states) ) | |
date2index = Dict( dates .=> 1:length(dates) ) | |
X = Array{Union{Missing, Float64}}(undef, length(state2index), length(date2index)) | |
for row in eachrow(df) | |
if haskey(state2index, row.region) | |
X[ state2index[row.region], date2index[row.date] ] = row[col] | |
end | |
end | |
return X, states, dates | |
end | |
X, Xs, Xd = build_state_by_date_matrix(df, :median) | |
## Animate the results with the media Rt value | |
using Measures | |
function make_anim(X, states, dates; kwargs...) | |
state_coords = map(st -> state_posx[lowercase(st)], states) | |
ycoords = -map(first, state_coords) | |
xcoords = map(x -> x[2], state_coords) | |
anim = @animate for (i,date) in enumerate(dates) | |
scatter(xcoords, ycoords, marker=:square, | |
markerstrokewidth=0, size=(800,550), | |
markersize=25, | |
margin=10mm, | |
framestyle=:none, | |
aspect_ratio=:equal, | |
xlims=(-0.5,11), | |
ylims=(-7.5,1), | |
marker_z = X[:,i],legend=false; kwargs...) | |
annotate!(collect(zip(xcoords, ycoords, Plots.text.(states, :white)))) | |
title!("$date") | |
end | |
return anim | |
end | |
# This centers the Rt values at 1 and shows higher ones in purple | |
# and lower ones (good ones) in green. | |
anim = make_anim(X, Xs, Xd; color=cgrad(:PRGn_11, rev=true), clim=(0.7, 1.3)) | |
gif(anim, "rt-vals.gif", fps=3) | |
## Now shows positive test fractions | |
P = build_state_by_date_matrix(df, :new_cases)[1] | |
T, S, D = build_state_by_date_matrix(df, :new_tests) | |
Pf =clamp!(P ./ (T .+ 1), 0, 1) | |
anim = make_anim(Pf, S, D; color=cgrad(:thermal, rev=true), clim=(0.05, 0.15)) | |
gif(anim, "posf-vals.gif", fps=3) | |
# Too noisy! | |
## Need moving averages! | |
# https://stackoverflow.com/questions/59562325/moving-average-in-julia | |
moving_average(vs,n) = [sum(@view vs[i:(i+n-1)])/n for i in 1:(length(vs)-(n-1))] | |
nday = 7 | |
Pa = copy(hcat(moving_average.(eachcol(P'), nday)...)') | |
Ta = copy(hcat(moving_average.(eachcol(T'), nday)...)') | |
Pf = clamp!(Pa ./ Ta, 0, 1) | |
Da = D[nday:end] | |
anim = make_anim(Pf, S, Da; colorbar=true,color=cgrad(:thermal, rev=true), clim=(0.05, 0.15), ) | |
gif(anim, "posf-mov-vals.gif", fps=3) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment