Skip to content

Commit

Permalink
Merge pull request #37 from twallema/master
Browse files Browse the repository at this point in the history
Modification to ideal network parameters in hpc scripts and demo
  • Loading branch information
twallema authored May 12, 2020
2 parents 603a07f + 78dffdc commit b902ed8
Show file tree
Hide file tree
Showing 3 changed files with 97 additions and 78 deletions.
44 changes: 25 additions & 19 deletions hpc/calibrate_stochastic.py
Original file line number Diff line number Diff line change
Expand Up @@ -24,48 +24,50 @@

# Construct the network G
# ~~~~~~~~~~~~~~~~~~~~~~~
numNodes = 30000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=7)
numNodes = 90000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=3)
# Baseline normal interactions:
G_norm = models.custom_exponential_graph(baseGraph, scale=200)
G_norm = models.custom_exponential_graph(baseGraph, scale=500)
models.plot_degree_distn(G_norm, max_degree=40)

# Construct the network G under social distancing
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
numNodes = 30000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=2)
numNodes = 90000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=1)
# Baseline normal interactions:
G_dist = models.custom_exponential_graph(baseGraph, scale=20000)
G_dist = models.custom_exponential_graph(baseGraph, scale=200000)
models.plot_degree_distn(G_dist, max_degree=40)

# Define model
# ~~~~~~~~~~~~
model = models.SEIRSNetworkModel(
# network connectivty
G = G_norm,
p = 0.6,
p = 0.51,
# clinical parameters
beta = 0.05,
sigma = 5.2,
beta = 0.20,
sigma = 4.0,
omega = 1.5,
zeta = 0,
a = 0.43, # probability of an asymptotic (supermild) infection
m = 1-0.43, # probability of a mild infection
h = 0.20, # probability of hospitalisation for a mild infection
c = 2/3, # probability of hospitalisation in cohort
mi = 1/6, # probability of hospitalisation in midcare
da = 14, # days of infection when asymptomatic (supermild)
dm = 14, # days of infection when mild
da = 6.5, # days of infection when asymptomatic (supermild)
dm = 6.5, # days of infection when mild
dc = 7,
dmi = 14,
dICU = 14,
dICUrec = 6,
dmirec = 6,
dhospital = 3, # days before reaching the hospital when heavy or critical
dhospital = 5, # days before reaching the hospital when heavy or critical
m0 = 0.49, # mortality in ICU
maxICU = 2000,
# testing
theta_S = 0,
theta_E = 0,
theta_I = 0,
theta_A = 0,
theta_M = 0,
theta_R = 0,
Expand All @@ -75,11 +77,13 @@
# back-tracking
phi_S = 0,
phi_E = 0,
phi_I = 0,
phi_A = 0,
phi_R = 0,
# initial condition
initN = 11.43e6, #results are extrapolated to entire population
initE = 100,
initE = 0,
initI = 3,
initA = 0,
initM = 0,
initC = 0,
Expand All @@ -89,13 +93,15 @@
initD = 0,
initSQ = 0,
initEQ = 0,
initIQ = 0,
initAQ = 0,
initMQ = 0,
initRQ = 0,
# monte-carlo sampling
monteCarlo = True,
monteCarlo = False,
repeats = 1
)

# Load data
# ~~~~~~~~~~~~
#[index,data] = model.obtainData()
Expand All @@ -111,19 +117,19 @@
# vector with dates
index=pd.date_range('2020-03-13', freq='D', periods=ICUvect.size)
# data series used to calibrate model must be given to function 'plotFit' as a list
idx = -26
idx = -23
index = index[0:idx]
data=[np.transpose(ICUvect[:,0:idx]),np.transpose(hospital[:,0:idx])]
# set optimisation settings
parNames = ['beta','p'] # must be a list!
parNames = ['beta'] # must be a list!
positions = [np.array([6]),np.array([4,5,6])] # must be a list!
bounds=[(10,100),(0.1,0.5),(0.4,0.8)] # must be a list!
bounds=[(10,100),(0.25,0.60)] # must be a list!
weights = np.array([0,1])
# run optimisation
theta = model.fit(data,parNames,positions,bounds,weights,setvar=True,maxiter=10,popsize=5)
theta = model.fit(data,parNames,positions,bounds,weights,setvar=True,maxiter=15,popsize=multiprocessing.cpu_count())
stop = timeit.default_timer()
print('Required time: ', stop - start,' seconds')

# Make a graphical representation of results
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
model.plotFit(index,data,positions,modelClr=['red','orange'],legendText=('ICU (model)','ICU (data)','Hospital (model)','Hospital (data)'),titleText='Belgium',filename-'calibration.svg')
model.plotFit(index,data,positions,modelClr=['red','orange'],legendText=('ICU (model)','ICU (data)','Hospital (model)','Hospital (data)'),titleText='Belgium',filename-'calibration90K.svg')
34 changes: 19 additions & 15 deletions hpc/sim_stochastic.py
Original file line number Diff line number Diff line change
Expand Up @@ -24,43 +24,44 @@

# Construct the network G
# ~~~~~~~~~~~~~~~~~~~~~~~
numNodes = 60000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=7)
numNodes = 90000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=3)
# Baseline normal interactions:
G_norm = models.custom_exponential_graph(baseGraph, scale=200)
G_norm = models.custom_exponential_graph(baseGraph, scale=500)
models.plot_degree_distn(G_norm, max_degree=40)

# Construct the network G under social distancing
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
numNodes = 60000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=2)
numNodes = 90000
baseGraph = networkx.barabasi_albert_graph(n=numNodes, m=1)
# Baseline normal interactions:
G_dist = models.custom_exponential_graph(baseGraph, scale=20000)
G_dist = models.custom_exponential_graph(baseGraph, scale=200000)
models.plot_degree_distn(G_dist, max_degree=40)

# Define model
# ~~~~~~~~~~~~
model = models.SEIRSNetworkModel(
# network connectivty
G = G_norm,
p = 0.6,
p = 0.51,
# clinical parameters
beta = 0.05,
sigma = 5.2,
beta = 0.03,
sigma = 4.0,
omega = 1.5,
zeta = 0,
a = 0.43, # probability of an asymptotic (supermild) infection
m = 1-0.43, # probability of a mild infection
h = 0.20, # probability of hospitalisation for a mild infection
c = 2/3, # probability of hospitalisation in cohort
mi = 1/6, # probability of hospitalisation in midcare
da = 14, # days of infection when asymptomatic (supermild)
dm = 14, # days of infection when mild
da = 6.5, # days of infection when asymptomatic (supermild)
dm = 6.5, # days of infection when mild
dc = 7,
dmi = 14,
dICU = 14,
dICUrec = 6,
dmirec = 6,
dhospital = 3, # days before reaching the hospital when heavy or critical
dhospital = 5, # days before reaching the hospital when heavy or critical
m0 = 0.49, # mortality in ICU
maxICU = 2000,
# testing
Expand All @@ -75,11 +76,13 @@
# back-tracking
phi_S = 0,
phi_E = 0,
phi_I = 0,
phi_A = 0,
phi_R = 0,
# initial condition
initN = 11.43e6, #results are extrapolated to entire population
initE = 100,
initE = 0,
initI = 1,
initA = 0,
initM = 0,
initC = 0,
Expand All @@ -89,12 +92,13 @@
initD = 0,
initSQ = 0,
initEQ = 0,
initIQ = 0,
initAQ = 0,
initMQ = 0,
initRQ = 0,
# monte-carlo sampling
monteCarlo = True,
repeats = 10
monteCarlo = False,
repeats = 3
)

# Create checkpoints dictionary
Expand Down
Loading

0 comments on commit b902ed8

Please sign in to comment.