Commit 2c00d2ad authored by Luca Lista's avatar Luca Lista
Browse files

added vax vs incidence chart

parent 09c5dc35
......@@ -16,7 +16,7 @@
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": false
"scrolled": true
},
"outputs": [],
"source": [
......@@ -27,6 +27,8 @@
"import matplotlib.colors as mcolors\n",
"from scipy import optimize\n",
"from scipy import stats\n",
"from scipy.stats import chi2\n",
"from scipy.stats import norm\n",
"import plotly\n",
"import plotly.graph_objects as go\n",
"import plotly.express as px\n",
......@@ -579,11 +581,78 @@
"# Plot results"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"vaccines = pd.read_csv('vaccine/data/csv/somministrazioni-vaccini-summary-latest.csv',\n",
" index_col=['region','data_somministrazione'])"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def plot_vax_inc(cfg):\n",
" fig = go.Figure()\n",
" col = 0\n",
" regs = sorted(italia_regioni)\n",
" vaxs, incs = [], []\n",
" for region in regs:\n",
" dat = regions[region][['data', 'totale_casi_giornalieri_somma1w_rel', 'totale_casi_giornalieri_somma2w_rel', 'rt']].set_index('data')\n",
" dat.index = dat.index.str.slice(stop=10)\n",
" vax = vaccines.loc[region.replace('uli Ven', 'uli-Ven').replace('ino Alt', 'ino-Alt')].sort_index()[['prima_dose', 'seconda_dose']].rename_axis('data')\n",
" df = dat.merge(vax, how='inner', on='data')\n",
" df = df[df.prima_dose>0]\n",
" pops = population.loc[population[\"Regione/Provincia\"]==region.replace('ia-Ro', 'ia Ro')]['Popolazione'].tolist()\n",
" if len(pops)>0: pop = pops[0]\n",
" else:\n",
" print(region, \" not found!\")\n",
" pop = 1\n",
" df['vax'] = df.seconda_dose.cumsum()/pop\n",
" fig.add_trace(go.Scatter(x=df.vax, y=df.totale_casi_giornalieri_somma2w_rel*100000, customdata=df.index,\n",
" hovertemplate = '%{customdata}<br>percentuale seconda dose: %{x:.2%}<br>positivi, 14gg/100k ab.: %{y:.2f}',\n",
" name=region, line=dict(color=rainbow_colors[col])))\n",
" col += 1\n",
" sel = df[~df.totale_casi_giornalieri_somma2w_rel.isnull() & ~df.seconda_dose.isnull()]\n",
" sel = sel[sel.index == sel.index.max()]\n",
" vaxs.append(sel.vax.values[0])\n",
" incs.append(sel.totale_casi_giornalieri_somma2w_rel.values[0]*100000)\n",
" fig.add_trace(go.Scatter(x=vaxs, y=incs, mode = 'markers+text', textposition=\"middle right\", text=regs,\n",
" hovertemplate = '%{text}<br>percentuale seconda dose: %{x:.2%}<br>positivi, 14gg/100k ab.: %{y:.2f}',\n",
" marker=dict(size=10, color=rainbow_colors), name='Attuale'))\n",
" add_logo(fig)\n",
" fig.update_layout(\n",
" autosize = True, paper_bgcolor=\"white\",\n",
" title_text = 'Vaccini vs incidenza',\n",
" xaxis_title = \"vaccini: percentuale seconda dose\",\n",
" yaxis_title = \"positivi, 14gg / 100.000 abitanti\",\n",
" xaxis_tickformat='.0%', xaxis_range=[min(vaxs)*0.98,max(vaxs)*1.02], yaxis_range=[0,max(incs)*1.3],\n",
" template='plotly_white',\n",
" title_x=0.5,\n",
" legend=(dict(bgcolor=\"rgba(255,255,255,0)\"))\n",
" )\n",
" filename = ('italia_vaccini_vs_incidenza')\n",
" file = path + filename\n",
" if not batch_mode: fig.show(renderer='iframe')\n",
" plot_div = py.plot(fig, output_type = 'div', include_plotlyjs = False, config=img_button_cfg(filename))\n",
" with open(file + \".div\", \"w\") as text_file: text_file.write(plot_div)\n",
" fig.write_image(file + \".png\", width=800, height=600, scale=1)\n",
" if plot_pdf:\n",
" fig.write_image(file + \".pdf\", width=pdf_width, height=pdf_height)\n",
" return (\"plot_vax_inc\", \"success\")\n",
"batch_mode = True\n",
"cfg = {}\n",
"plot_vax_inc(cfg)"
]
},
{
"cell_type": "markdown",
"metadata": {
"scrolled": true
},
"metadata": {},
"source": [
"Make summary plots"
]
......@@ -596,7 +665,6 @@
},
"outputs": [],
"source": [
"vaccines = pd.read_csv('vaccine/data/csv/somministrazioni-vaccini-summary-latest.csv', index_col=['region','data_somministrazione'])\n",
"def summary(cfg):\n",
" region = cfg['region']\n",
" ratios = 'ratios' in cfg.keys()\n",
......@@ -949,9 +1017,7 @@
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": true
},
"metadata": {},
"outputs": [],
"source": [
"def prov_map(cfg):\n",
......@@ -1023,6 +1089,106 @@
" prov_map(dict(var='totale_casi_giornalieri_1w'))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"with open('italy-regions.json') as json_file: italy_regions_json = json.load(json_file)\n",
"with open('italy-provinces.json') as json_file: italy_provinces_json = json.load(json_file)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def benford_map(cfg):\n",
" zone = cfg['zone']\n",
" var = cfg['var']\n",
" ben_data = []\n",
" ord0 = ord('0')\n",
" db = np.array([math.log10(1 + 1/d) for d in range(1,10)])\n",
" def zscore(data):\n",
" dd = np.zeros(9)\n",
" for n in data:\n",
" if n.is_integer():\n",
" d1 = ord(str(n)[0])-ord0\n",
" dd[d1-1] +=1\n",
" dbe = db*dd.sum()\n",
" c2 = (((dd-dbe)/dbe)**2).sum()\n",
" pval=chi2.sf(c2, 8)\n",
" zsc = -norm.ppf(pval)\n",
" return pval, zsc\n",
" if zone == 'regions':\n",
" for region in regions:\n",
" data = regions[region][var]\n",
" reg_name = region\n",
" if reg_name == 'Trentino Alto Adige': reg_name = 'Trentino-Alto Adige/Südtirol'\n",
" if reg_name == 'Friuli Venezia Giulia': reg_name = 'Friuli-Venezia Giulia'\n",
" pval, zsc = zscore(data)\n",
" ben_data.append([reg_name, pval, zsc])\n",
" ben_df = pd.DataFrame(ben_data, columns=['region', 'pval', 'zsc'])\n",
" reg_names=regions_only\n",
" else:\n",
" provs=frame_p['denominazione_provincia'].unique()\n",
" for prov in provs:\n",
" data = frame_p[frame_p['denominazione_provincia']==prov][var]\n",
" prov_name = prov.replace('Bolzano','Bolzano/Bozen')\n",
" pval, zsc = zscore(data)\n",
" ben_data.append([prov_name, pval, zsc])\n",
" ben_df = pd.DataFrame(ben_data, columns=['province', 'pval', 'zsc'])\n",
" zdata = ben_df.zsc\n",
" zmin, zmax = zdata.min(), zdata.max()\n",
" zrange = zmax - zmin\n",
" ref = 0\n",
" zref = (ref - zmin)/zrange\n",
" if zref > 0 and zref < 1:\n",
" colorscale = [[0, 'green'], [zref, 'yellow'], [1, 'orangered']]\n",
" elif zref < 0:\n",
" colorscale = [[0, 'yellow'], [1, 'orangered']]\n",
" else:\n",
" colorscale = [[0, 'green'], [1, 'yellow']]\n",
" if zone == 'regions':\n",
" fig = go.Figure(data=go.Choropleth(geojson=italy_regions_json, locations=ben_df.region,\n",
" featureidkey = \"properties.reg_name\", colorbar_title='Z score', colorscale=colorscale,\n",
" hovertemplate='', z=zdata), layout=dict(geo=dict(projection = {'type':'transverse mercator'})))\n",
" else:\n",
" fig = go.Figure(data=go.Choropleth(geojson=italy_provinces_json, locations=ben_df.province,\n",
" featureidkey = \"properties.prov_name\", colorbar_title='Z score', colorscale=colorscale,\n",
" hovertemplate='', z=zdata), layout=dict(geo=dict(projection = {'type':'transverse mercator'})))\n",
" fig.update_geos(fitbounds=\"locations\", visible=False)\n",
" fig.update_layout(margin={\"r\":0,\"l\":0,\"b\":0})\n",
" fig.update_layout(autosize = True, paper_bgcolor=\"white\", showlegend=True,\n",
" title_text = \"Italia - Bendford - \" + labels[var] + ' - '+ last_day, legend=(dict(bgcolor=\"rgba(255,255,255,0)\")),\n",
" template='plotly_white', title_x=0.5)\n",
" fig.update_coloraxes(colorbar_tickformat='.2f')\n",
" add_logo(fig, 1.1)\n",
" filename = ('italia_regions_bendford_'+var+'_map').lower()\n",
" file = path + filename\n",
" if not batch_mode: fig.show(renderer='iframe')\n",
" plot_div = py.plot(fig, output_type = 'div', include_plotlyjs = False, config=img_button_cfg(filename))\n",
" with open(file + \".div\", \"w\") as text_file: text_file.write(plot_div)\n",
" fig.write_image(file + \".png\", width=800, height=600, scale=1)\n",
" if plot_pdf: fig.write_image(file + \".pdf\", width=pdf_width, height=pdf_height)\n",
" return (\"benford map: \"+var, \"success\")\n",
"\n",
"batch_mode = True\n",
"if batch_mode:\n",
" pass\n",
"# cfgs = []\n",
"# for var in ['totale_casi_giornalieri', 'deceduti_giornalieri', 'terapia_intensiva', 'totale_ospedalizzati']:\n",
"# cfgs.append(dict(zone='regions', var=var))\n",
"# with concurrent.futures.ProcessPoolExecutor(max_workers=32) as executor:\n",
"# res = executor.map(prov_map, cfgs)\n",
"# result = list(res)\n",
"# print_job_results(result)\n",
"else:\n",
" benford_map(dict(zone='provinces', var='totale_casi_giornalieri'))"
]
},
{
"cell_type": "markdown",
"metadata": {
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment