import pandas as pd # Loading the NYC dataset
import numpy as np # performing fast calculations to build our simulation
import plotly.express as px # Quick visuals for data exploration
from ipywidgets import widgets # sliders, buttons and visual layout of the dashboard
import plotly.graph_objects as go # plotly objects to place on our graph
import math # python mathematical functions
Python 104
Building on the code from last week
Recap - Calculation method to calculate # of people affected
# We'll set some easy defaults for this function so that we immediately have outputs which make some sense when we execute this method.
def calculate_values_SIR(removal_rate = .1
= 1.1
, infection_rate = 1000
, population = 150
, days = 1
, initial_infection
# we've added this additional parameter which allows us to cut the
# infection rate by 2 at a given day to simulate drastic policy action.
= 50
, intervention_day
):
# build out starting positions
= np.linspace(0,days, days)
x_days = [population]
y_susceptible = [initial_infection]
y_infected = [0]
y_removed
for day in x_days[1:days-1]:
# if this day is after the day drastic intervention is made, reduce the infection rate by 2.
if day > intervention_day:
= infection_rate / 2
infection_rate
= int(day)
day
= min((infection_rate * y_infected[day-1] * (y_susceptible[day-1])/population), y_susceptible[day-1])
daily_infected = removal_rate * y_infected[day-1]
daily_removed
-1] - daily_infected )
y_susceptible.append( y_susceptible[day-1] + daily_infected - daily_removed )
y_infected.append( y_infected[day-1] + daily_removed )
y_removed.append( y_removed[day
return {'infected': np.rint(y_infected)
'removed' : np.rint(y_removed)
, 'susceptible': np.rint(y_susceptible)}
,
Recap - Build widgets to control the model inputs
# This slider works with floating point numbers (hence being called Float Slider) and
# allows us to set a variable with this slider. This is going to be the way we set the
# infection rate.
= widgets.FloatSlider(
ir =1.187, # this is the initial value of our slider when it appears
valuemin=0.0, # the minimum value we'll allow
max=5.0, # the maximum value we'll allow
=.001, # by what increments the slider will change when we move it
step='Infection_rate:', # the name of the slider
description=False # Will this slider wait until it stops moving to
continuous_update# call it's update function or does it call the
# update function immediately?
)
= widgets.FloatSlider(
rr =.46,
valuemin=0.1, # this is set to greater than 0 because this is the denominator in the R0 calculation
max=2.0,
=.01,
step='Removal_Rate:',
description=False
continuous_update
)
= widgets.IntSlider(
ii =1,
valuemin=1,
max=50,
=1,
step='Initially Infected:',
description=False
continuous_update
)
= widgets.IntSlider(
ip =1000,
valuemin=500,
max=10_000_000,
=500,
step='Initial Population:',
description=False
continuous_update
)
= widgets.IntSlider(
iday =15,
valuemin=1,
max=500,
=1,
step='Day of intervention (reducing infection rate):',
description=False
continuous_update
)
= widgets.HBox(children=[ir, rr])
first_slider_group = widgets.HBox(children=[ii, ip, iday]) second_slider_group
Recap - Create the first graph with default values
# First, we use the method created above to calculate a model using the initial
# values of the sliders we just created. Given that at this point we haven't
# displayed the sliders yet, their values will be the default values we set above.
= calculate_values_SIR( removal_rate = rr.value
data = ir.value
, infection_rate = ip.value
, population = 150
, days = ii.value
, initial_infection = iday.value
, intervention_day
)
# Next we add all the data traces to the chart
= go.Bar(x = list(range(1,len(data['infected'])))
infected_trace = data['infected']
,y ='Infected'
, name= dict(color='red')
, marker
)
= go.Bar(x = list(range(1, len(data['susceptible'])))
susceptible_trace = data['susceptible']
, y ='Susceptible'
, name= dict(color='rgba(0,0,255,0.5)')
, marker =0.5
, opacity
)
= go.Bar(x =list(range(1, len(data['removed'])))
removed_trace = data['removed']
,y ='Removed'
, name= dict(color='rgba(0,128,0,0.5)')
, marker =0.5)
, opacity
# This trace is interesting as it's basically just drawing a straight line on the
# selected intervention day.
= go.Scatter(x = [iday.value, iday.value]
intervention_day = [0, ip.value]
,y ='Intervention day'
, name= dict(color='darkblue')
, marker = dict(width=5)
, line
)
# We create our figure adding all the traces we created to the data list, and setting some layout values in the layout parameter.
= go.FigureWidget(data=[ infected_trace, removed_trace, susceptible_trace, intervention_day ],
g =go.Layout(
layout={
title'text': f'R0 = {ir.value / rr.value} <br /> Post-Intervention R0: {ir.value/ 2* rr.value} <br />Infection_rate={ir.value} Removal_rate={rr.value}',
'y':.95,
'x':0.5,
'xanchor': 'center',
'yanchor': 'top'}
='stack'
,barmode='x'
,hovermode=900
,height=dict(title='Number of Days')
,xaxis=dict(title='Number of People')
,yaxis
))
# This is to update the x-axis range to show only the days where we have cases
range=[0,np.where(data['infected']==0)[0][0]])
g.update_xaxes(
=1 x
Recap - Connect the widgets to the calculation method to allow the widgets to update the graph
# This method will be called any time one of the sliders is modified. It will re-run our model calculation
# with the new values and update the data for the 4 traces we added to the figure.
def response(change):
=150
num_days
# recalculate the model using the new values defined by the sliders
= calculate_values_SIR(removal_rate = rr.value
pop_values =ir.value
, infection_rate=ii.value
, initial_infection=ip.value
, population=num_days
, days= iday.value)
, intervention_day
# Try to find the first day where we have no more infections,
# if that fails and we get an error, use the maximum number of
# days for which we've calculated the model. We use this later
# to update the x-axis range to keep our curve centered.
try:
= np.where(pop_values['infected']==0)[0][0]
end_infection except IndexError:
= num_days
end_infection
# plotly updates much faster if you send all the updates in one go,
# this construct allows us to open an update session on our chart
# and when we leave it, it will send all the updates to the figure
# at once.
with g.batch_update():
# update the y-axis values from the model
0].y = pop_values['infected']
g.data[1].y = pop_values['removed']
g.data[2].y = pop_values['susceptible']
g.data[
# update the x-axis values
0].x = list(range(1,num_days))
g.data[1].x = list(range(1,num_days))
g.data[2].x = list(range(1,num_days))
g.data[
# Add the intervention day line
3].y = [0 , ip.value]
g.data[3].x = [iday.value, iday.value]
g.data[
# update the title to show the R values as well as the infection rate and the removal rate
= {
g.layout.title 'text': f'R0 = {ir.value / rr.value} <br /> Post-Intervention R0: {ir.value/ 2* rr.value} <br /> Infection_rate={ir.value} Removal_rate={rr.value}',
'x':0.5,
'xanchor': 'center',
'yanchor': 'top'}
# change the x-axis range to keep the infection curve in scope
= dict(range=[0,end_infection])
g.layout.xaxis
# Update each of the widgets and register our update method as the method to
# call when they change.
='value')
ir.observe(response, names='value')
rr.observe(response, names='value')
ii.observe(response, names='value')
ip.observe(response, names='value')
iday.observe(response, names
# put the widgets, and our chart together into our layout
widgets.VBox([first_slider_group, second_slider_group, g])
Add an intervention effectiveness measure
We’re going to add an additional slider to allow us to control how effective the intervention measure is.
# We'll set some easy defaults for this function so that we immediately have outputs which make some sense when we execute this method.
def calculate_values_SIR(removal_rate = .1
= 1.1
, infection_rate = 1000
, population = 150
, days = 1
, initial_infection
# we've added this additional parameter which allows us to cut the
# infection rate by 2 at a given day to simulate drastic policy action.
= 50
, intervention_day = .46
, intervention_effectiveness
):
# build out starting positions
= np.linspace(0,days, days)
x_days = [population]
y_susceptible = [initial_infection]
y_infected = [0]
y_removed
for day in x_days[1:days-1]:
# if this day is after the day drastic intervention is made, reduce the infection rate by 2.
if day > intervention_day:
= infection_rate * max((1-intervention_effectiveness),0)
infection_rate
= int(day)
day
= min((infection_rate * y_infected[day-1] * (y_susceptible[day-1])/population), y_susceptible[day-1])
daily_infected = removal_rate * y_infected[day-1]
daily_removed
-1] - daily_infected )
y_susceptible.append( y_susceptible[day-1] + daily_infected - daily_removed )
y_infected.append( y_infected[day-1] + daily_removed )
y_removed.append( y_removed[day
return {'infected': np.rint(y_infected)
'removed' : np.rint(y_removed)
, 'susceptible': np.rint(y_susceptible)}
,
# This slider works with floating point numbers (hence being called Float Slider) and
# allows us to set a variable with this slider. This is going to be the way we set the
# infection rate.
= widgets.FloatSlider(
ir =1.187, # this is the initial value of our slider when it appears
valuemin=0.0, # the minimum value we'll allow
max=5.0, # the maximum value we'll allow
=.001, # by what increments the slider will change when we move it
step='Infection_rate:', # the name of the slider
description=False # Will this slider wait until it stops moving to
continuous_update# call it's update function or does it call the
# update function immediately?
)
= widgets.FloatSlider(
rr =.46,
valuemin=0.1, # this is set to greater than 0 because this is the denominator in the R0 calculation
max=2.0,
=.01,
step='Removal_Rate:',
description=False
continuous_update
)
= widgets.IntSlider(
ii =1,
valuemin=1,
max=50,
=1,
step='Initially Infected:',
description=False
continuous_update
)
= widgets.IntSlider(
ip =1000,
valuemin=500,
max=10_000_000,
=500,
step='Initial Population:',
description=False
continuous_update
)
= widgets.IntSlider(
iday =15,
valuemin=1,
max=500,
=1,
step='Day of intervention (reducing infection rate):',
description=False
continuous_update
)
= widgets.FloatSlider(
ie =.46,
valuemin=0.0,
max=1.0,
=.01,
step='Intervention effectiveness:',
description=False
continuous_update
)
= widgets.HBox(children=[ir, rr, ie])
first_slider_group = widgets.HBox(children=[ii, ip, iday])
second_slider_group
# First, we use the method created above to calculate a model using the initial
# values of the sliders we just created. Given that at this point we haven't
# displayed the sliders yet, their values will be the default values we set above.
= calculate_values_SIR( removal_rate = rr.value
data = ir.value
, infection_rate = ip.value
, population = 150
, days = ii.value
, initial_infection = iday.value
, intervention_day = ie.value
, intervention_effectiveness
)
# Next we add all the data traces to the chart
= go.Bar(x = list(range(1,len(data['infected'])))
infected_trace = data['infected']
,y ='Infected'
, name= dict(color='red')
, marker
)
= go.Bar(x = list(range(1, len(data['susceptible'])))
susceptible_trace = data['susceptible']
, y ='Susceptible'
, name= dict(color='rgba(0,0,255,0.5)')
, marker =0.5
, opacity
)
= go.Bar(x =list(range(1, len(data['removed'])))
removed_trace = data['removed']
,y ='Removed'
, name= dict(color='rgba(0,128,0,0.5)')
, marker =0.5)
, opacity
# This trace is interesting as it's basically just drawing a straight line on the
# selected intervention day.
= go.Scatter(x = [iday.value, iday.value]
intervention_day = [0, ip.value]
,y ='Intervention day'
, name= dict(color='darkblue')
, marker = dict(width=5)
, line
)
# We create our figure adding all the traces we created to the data list, and setting some layout values in the layout parameter.
= go.FigureWidget(data=[ infected_trace, removed_trace, susceptible_trace, intervention_day ],
g =go.Layout(
layout={
title'text': f'R0 = {ir.value / rr.value} <br /> Post-Intervention R0: {(ir.value * .46) / rr.value} <br />Infection_rate={ir.value} Removal_rate={rr.value}',
'y':.95,
'x':0.5,
'xanchor': 'center',
'yanchor': 'top'}
='stack'
,barmode='x'
,hovermode=900
,height=dict(title='Number of Days')
,xaxis=dict(title='Number of People')
,yaxis
))
# This is to update the x-axis range to show only the days where we have cases
range=[0,np.where(data['infected']==0)[0][0]])
g.update_xaxes(
=1
x
# This method will be called any time one of the sliders is modified. It will re-run our model calculation
# with the new values and update the data for the 4 traces we added to the figure.
def response(change):
=150
num_days
# recalculate the model using the new values defined by the sliders
= calculate_values_SIR(removal_rate = rr.value
pop_values =ir.value
, infection_rate=ii.value
, initial_infection=ip.value
, population=num_days
, days= iday.value
, intervention_day = ie.value
, intervention_effectiveness
)
# Try to find the first day where we have no more infections,
# if that fails and we get an error, use the maximum number of
# days for which we've calculated the model. We use this later
# to update the x-axis range to keep our curve centered.
try:
= np.where(pop_values['infected']==0)[0][0]
end_infection except IndexError:
= num_days
end_infection
# plotly updates much faster if you send all the updates in one go,
# this construct allows us to open an update session on our chart
# and when we leave it, it will send all the updates to the figure
# at once.
with g.batch_update():
# update the y-axis values from the model
0].y = pop_values['infected']
g.data[1].y = pop_values['removed']
g.data[2].y = pop_values['susceptible']
g.data[
# update the x-axis values
0].x = list(range(1,num_days))
g.data[1].x = list(range(1,num_days))
g.data[2].x = list(range(1,num_days))
g.data[
# Add the intervention day line
3].y = [0 , ip.value]
g.data[3].x = [iday.value, iday.value]
g.data[
# update the title to show the R values as well as the infection rate and the removal rate
= ir.value * max((1-ie.value),0)
ie_rate = {
g.layout.title 'text': f'R0 = {ir.value / rr.value} <br /> Post-Intervention R0: {ie_rate / rr.value} <br /> Infection_rate={ir.value} Removal_rate={rr.value}',
'x':0.5,
'xanchor': 'center',
'yanchor': 'top'}
# change the x-axis range to keep the infection curve in scope
= dict(range=[0,end_infection])
g.layout.xaxis
# Update each of the widgets and register our update method as the method to
# call when they change.
='value')
ir.observe(response, names='value')
rr.observe(response, names='value')
ii.observe(response, names='value')
ip.observe(response, names='value')
iday.observe(response, names='value')
ie.observe(response, names
# put the widgets, and our chart together into our layout
widgets.VBox([first_slider_group, second_slider_group, g])
Add display the effectiveness vs the rate had there been no change
Update the calculation method
The first thing we need to do if we want to show the effectiveness of the intervention at the same time as what happens after an intervention is to update the calculation method to calculate both the intervened infectious rates and the un-intervened infectious rates at the same time. This sounds more complex than it actually is, all we need to do during this is keep track of how many people would have been infected/removed if there had been no intervention and subtract the intervened numbers from that.
To that end, we need to keep track of 2 additional numbers for every day: - The number of people we need to add to the infected group to count how many would be infected if there were no intervention. - The number of people we need to add to the removed group if there were no intervention.
To denote such non -intervention figures in the code we’re going to append values with ‘_noi’ so that we can keep track of the changes.
# We'll set some easy defaults for this function so that we immediately have outputs which make some sense when we execute this method.
def calculate_values_SIR(removal_rate = .1
= 1.1
, infection_rate = 1000
, population = 150
, days = 1
, initial_infection
= 10
, incubation_period = .03
, fatality_rate = 30
, immunity_duration
# we've added this additional parameter which allows us to cut the
# infection rate by 2 at a given day to simulate drastic policy action.
= 50
, intervention_day = .46
, intervention_effectiveness
):
# build out starting positions
= np.linspace(0,days, days)
x_days = [population]
y_susceptible = [initial_infection]
y_infected = [0]
y_removed
= [0]
y_infected_noi = [0]
y_removed_noi
for day in x_days[1:days-1]:
# if this day is after the day drastic intervention is made, reduce the infection rate by 2.
if day >= intervention_day:
= infection_rate * max((1-intervention_effectiveness),0)
inter_infection_rate else:
= infection_rate
inter_infection_rate
= int(day)
day
# Calculate the number of people susceptible from the previous day
= y_susceptible[day-1] + y_infected_noi[day-1]
pday_susceptible
= min((inter_infection_rate * y_infected[day-1] * (pday_susceptible)/population), pday_susceptible)
daily_infected = removal_rate * y_infected[day-1]
daily_removed
# If we're after the intevention, caluclate what would have happened and update the
# hypothetical values
if day >= intervention_day:
# calculate the hypothetical number of infected
= y_infected[day-1] + y_infected_noi[day-1]
pday_infected_noi
# Calculate the hypothetical daily infected rates
= min((infection_rate * (pday_infected_noi) * ( y_susceptible[day-1])/population), y_susceptible[day-1]) - daily_infected
daily_infected_noi
# Calculate the hyothetical daily removed rate
= (removal_rate * (pday_infected_noi))
daily_removed_noi
# Calculate our daily removed additional values
= max(daily_removed_noi - daily_removed, 0 )
noi_removed
# Keep track of our unintervened infected and removed values
-1] + (daily_infected_noi) - (noi_removed))
y_infected_noi.append(y_infected_noi[day-1] + (noi_removed))
y_removed_noi.append( y_removed_noi[day
else:
= 0
daily_infected_noi = 0
daily_removed_noi 0)
y_infected_noi.append(0)
y_removed_noi.append(
-1] - (daily_infected + (daily_infected_noi)))
y_susceptible.append( y_susceptible[day-1] + daily_infected - daily_removed )
y_infected.append( y_infected[day-1] + daily_removed )
y_removed.append( y_removed[day
return {'infected': np.rint(y_infected)
'removed' : np.rint(y_removed)
, 'susceptible': np.rint(y_susceptible)
, 'infected_noi': np.rint(y_infected_noi)
, 'removed_noi': np.rint(y_removed_noi) }
,
= calculate_values_SIR( removal_rate = rr.value
data = ir.value
, infection_rate = ip.value
, population = 150
, days = ii.value
, initial_infection = 15
, intervention_day
)
data
{'infected': array([ 1., 2., 3., 5., 9., 15., 26., 43., 71., 113., 169.,
231., 279., 285., 246., 161., 103., 64., 39., 23., 14., 8.,
5., 3., 2., 1., 1., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0.]),
'removed': array([ 0., 0., 1., 3., 5., 9., 16., 28., 48., 80., 132.,
210., 316., 444., 575., 688., 762., 810., 839., 857., 867., 874.,
877., 880., 881., 882., 882., 882., 882., 882., 882., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883., 883., 883., 883., 883., 883.,
883., 883., 883., 883., 883., 883.]),
'susceptible': array([1000., 999., 997., 993., 987., 977., 959., 930., 882.,
808., 700., 560., 406., 272., 180., 127., 99., 84.,
76., 71., 68., 66., 65., 65., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64., 64., 64., 64., 64.,
64., 64., 64., 64., 64.]),
'infected_noi': array([ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 24., 25., 21., 15., 11., 8., 5., 3., 2., 1., 1.,
1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 0., 0., 0.]),
'removed_noi': array([ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
0., 0., 0., 11., 23., 32., 39., 44., 48., 50., 52., 53., 53.,
54., 54., 54., 54., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55., 55.,
55., 55., 55., 55., 55., 55.])}
# This slider works with floating point numbers (hence being called Float Slider) and
# allows us to set a variable with this slider. This is going to be the way we set the
# infection rate.
= widgets.FloatSlider(
ir =1.187, # this is the initial value of our slider when it appears
valuemin=0.0, # the minimum value we'll allow
max=5.0, # the maximum value we'll allow
=.001, # by what increments the slider will change when we move it
step='Infection_rate:', # the name of the slider
description=False # Will this slider wait until it stops moving to
continuous_update# call it's update function or does it call the
# update function immediately?
)
= widgets.FloatSlider(
rr =.46,
valuemin=0.1, # this is set to greater than 0 because this is the denominator in the R0 calculation
max=2.0,
=.01,
step='Removal_Rate:',
description=False
continuous_update
)
= widgets.IntSlider(
ii =1,
valuemin=1,
max=50,
=1,
step='Initially Infected:',
description=False
continuous_update
)
= widgets.IntSlider(
ip =1000,
valuemin=500,
max=10_000_000,
=500,
step='Initial Population:',
description=False
continuous_update
)
= widgets.IntSlider(
iday =15,
valuemin=1,
max=150,
=1,
step='Day of intervention (reducing infection rate):',
description=False
continuous_update
)
= widgets.FloatSlider(
ie =.46,
valuemin=0.0,
max=1.0,
=.01,
step='Intervention effectiveness:',
description=False
continuous_update
)
= widgets.HBox(children=[ir, rr, ie])
first_slider_group = widgets.HBox(children=[ii, ip, iday])
second_slider_group
# First, we use the method created above to calculate a model using the initial
# values of the sliders we just created. Given that at this point we haven't
# displayed the sliders yet, their values will be the default values we set above.
= calculate_values_SIR( removal_rate = rr.value
data = ir.value
, infection_rate = ip.value
, population = 150
, days = ii.value
, initial_infection = iday.value
, intervention_day
)
# Next we add all the data traces to the chart
= go.Bar(x = list(range(1,len(data['infected'])))
infected_trace = data['infected']
,y ='Infected'
, name= dict(color='red')
, marker
)
= go.Bar(x = list(range(1, len(data['susceptible'])))
susceptible_trace = data['susceptible']
, y ='Susceptible'
, name= dict(color='rgba(0,0,255,0.5)')
, marker =0.5
, opacity
)
= go.Bar(x =list(range(1, len(data['removed'])))
removed_trace = data['removed']
,y ='Removed'
, name= dict(color='rgba(0,128,0,0.5)')
, marker =0.5)
, opacity
######
# ADDITION STARTS HERE
######
= go.Bar(x = list(range(1,len(data['infected_noi'])))
infected_trace_noi = data['infected_noi']
,y ='Infected No Intervention'
, name= dict(color='rgba(225,0,0,0.5)')
, marker
)
= go.Bar(x =list(range(1, len(data['removed_noi'])))
removed_trace_noi = data['removed_noi']
,y ='Removed No Intervention'
, name= dict(color='rgba(0,225,0,0.5)')
, marker
)
######
# ADDITION END HERE
######
# This trace is interesting as it's basically just drawing a straight line on the
# selected intervention day.
= go.Scatter(x = [iday.value, iday.value]
intervention_day = [0, ip.value]
,y ='Intervention day'
, name= dict(color='darkblue')
, marker = dict(width=5)
, line
)
# We create our figure adding all the traces we created to the data list, and setting some layout values in the layout parameter.
= go.FigureWidget(data=[ infected_trace
g
, infected_trace_noi
, removed_trace
, susceptible_trace
, intervention_day
, removed_trace_noi ],=go.Layout(
layout={
title'text': f'R0 = {ir.value / rr.value} <br /> Post-Intervention R0: {ir.value/ 2* rr.value} <br />Infection_rate={ir.value} Removal_rate={rr.value}',
'y':.95,
'x':0.5,
'xanchor': 'center',
'yanchor': 'top'}
='stack'
,barmode='x'
,hovermode=600
,height=dict(title='Number of Days')
,xaxis=dict(title='Number of People')
,yaxis
))
# This is to update the x-axis range to show only the days where we have cases
range=[0,np.where(data['infected']==0)[0][0]])
g.update_xaxes(
=1
x
# This method will be called any time one of the sliders is modified. It will re-run our model calculation
# with the new values and update the data for the 4 traces we added to the figure.
def response(change):
=150
num_days
# recalculate the model using the new values defined by the sliders
= calculate_values_SIR(removal_rate = rr.value
pop_values =ir.value
, infection_rate=ii.value
, initial_infection=ip.value
, population=num_days
, days= iday.value
, intervention_day = ie.value
, intervention_effectiveness
)
# Try to find the first day where we have no more infections,
# if that fails and we get an error, use the maximum number of
# days for which we've calculated the model. We use this later
# to update the x-axis range to keep our curve centered.
try:
= max( np.where(pop_values['infected_noi'][iday.value+1:]==0)[0][0] + iday.value +1,
end_infection 'infected']==0)[0][0] )
np.where(pop_values[except IndexError:
= num_days
end_infection
# plotly updates much faster if you send all the updates in one go,
# this construct allows us to open an update session on our chart
# and when we leave it, it will send all the updates to the figure
# at once.
with g.batch_update():
# update the y-axis values from the model
0].y = pop_values['infected']
g.data[1].y = pop_values['infected_noi']
g.data[2].y = pop_values['removed']
g.data[3].y = pop_values['susceptible']
g.data[5].y = pop_values['removed_noi']
g.data[
# Add the intervention day line
4].y = [0 , ip.value]
g.data[4].x = [iday.value, iday.value]
g.data[
# update the title to show the R values as well as the infection rate and the removal rate
= ir.value * max((1-ie.value),0)
ie_rate = {
g.layout.title 'text': f'R0 = {ir.value / rr.value} <br /> Post-Intervention R0: {ie_rate / rr.value} <br /> Infection_rate={ir.value} Removal_rate={rr.value}',
'x':0.5,
'xanchor': 'center',
'yanchor': 'top'}
# change the x-axis range to keep the infection curve in scope
= dict(range=[0,end_infection])
g.layout.xaxis
# Update each of the widgets and register our update method as the method to
# call when they change.
='value')
ir.observe(response, names='value')
rr.observe(response, names='value')
ii.observe(response, names='value')
ip.observe(response, names='value')
iday.observe(response, names='value')
ie.observe(response, names
# put the widgets, and our chart together into our layout
widgets.VBox([first_slider_group, second_slider_group, g])