<sup font-size: smaller> a Mathematics and computational analysis, University of Illinois Urbana-Champaign
b Environmental health and epidemiology, National Environmental Public Health Tracking Program, Centers for Disease Control and Prevention
c Spatial analysis and geovisualization, University of Maine
d Geography, Louisiana State University
e Human geography, Hunter College
f Spatial analysis and modeling, geovisualization, University of Southern California
g Spatially explicit agent-based modeling, University of Illinois Urbana-Champaign
</sup>
Influenza is a contagious respiratory illness caused by influenza viruses. The flu is spread from person to person, often by droplets when those who are infected cough, sneeze, or talk. Sometimes people can be infected by touching a surface that has the virus on it, although this is not as common. Symptoms include fever, cough, sore throat, muscle aches, headaches, and fatigue. Influenza has a high burden of disease despite vaccination and public health policies that have been implemented. The CDC estimates that more than 900,000 people were hospitalized and over 80,000 people died from the flu in the 2017-2018 season (CDC, 2018).
Spatial simulation modeling (also known as agent-based modeling or ABM) serves as a tool to simulate and capture various spatiotemporal phenomena. ABM can be used to explore the epidemiology and spread of certain diseases such as influenza. The phenomena being modeled are the results from dynamic interactions and behaviors of individual agents and the environments. The modeling process, however, is not always simple—it requires the conceptualization and understanding of the spatiotemporal process at an individual level. An ABM has several key components, including its agents, properties, environment, and rules, which specify how the agents interact and what happens when they interact. This work will apply a spatially explicit ABM to simulate influenza spread.
To what degree can a spatially explicit agent-based model simulate influenza spread in the Queen Anne neighborhood in Seattle, Washington?
1.1 Geographic data
The geographic data of households/workplaces/schools were obtained from the City of Seattle Open Data Hub. Then, the data were extracted based on the zoning types that overlay the building
footprints. The building footprint and zoning data can be downloaded
by clicking the respective links.
1.2 Neighborhood selection
For this workshop, we have selected Queen Anne neighborhood in Seattle as a use case example, and datasets were subsetted by choosing x and y limits (in UTM 10N coordinate system) of the Queen Anne neighborhood.
xmin = 546302; xmax = 549610; ymin = 5274394; ymax = 5277978;
nhouses_old = len(houses); houses = houses[(houses.X >= xmin) & (houses.X <= xmax) & (houses.Y >= ymin) & (houses.Y <= ymax)]; nhouses_new = len(houses); print(nhouses_new," of ",nhouses_old," houses were selected.")
1.3 Human agents
Interactions between the human agents and their environments (attributes and methods) are included in the Human.py and Env.py programs. The datasets needed for this include population distribution and household size. All of these datasets are at the block group level (specific to Queen Anne neighborhood) and were obtained from the American Community Survey for 2013-2017 5-year estimates.
Population distribution data - Table B01001 - Sex by Age
Household size data - Table B11016 - Household Type by Household Size
These data were joined to the Queen Anne neighborhood block groups in order to calculate the weights for population distribution and household size across the entire neighborhood. Age groups were applied using standard 5-year age groups across 18 categories (e.g., starting with 0-4 and ending with 85+). Household size used only family households (i.e., excluded non-family households) and used six categories (size 2, size 3, size 4, size 5, and size 6 or more).
1.4 Influenza data
Influenza data were obtained from the 2018-2019 Weekly Flu Surveillance Reports from the Seattle & King County Department of Public Health. The starting week used was October 7 - 13 2018 and the ending week used was April 28 - May 4 2019. All weekly data were used in this time period to sum the total influenza cases (e.g., summing influenza A cases, influenza B cases, etc.) over time. These data are used to generate a plot of reference data against which to compare the results of the simulation.
1.5 Additional data
Race and median household income data were also obtained from the American Community Survey (2013-2017 5-year estimates) for the Queen Anne neighborhood block groups. This information was used for visual analysis to examine the distribution of influenza cases in relation to different sociodemographic characteristics.
Race data - Table B02001 - Race
Median household income data - Table B19013 - Median household income in the past 12 months
2.1 The locations of households, schools, and workplaces are obtained by overlaying building footprints data with zoning data (e.g., residential, commercial, educational).
2.2 Age distributions and the composition of households
are synthesized as similar to American Community Survey (ACS) 2013-17 5-year estimates data (see Find the Data Section for more information).
2.3 Individual humans have age and age-dependent daily
behaviors (e.g., commuting to schools (aged 6 to 19) and workplaces (aged 20 to 64), and staying at home all day (the rest of people)).
3.1 A random sample of people get infected at the beginning of the simulation.
3.2 The disease process is explained through a susceptible-exposed-infectious-recovered (SEIR) model.
3.3 An individual human can spread influenza
to his/her classmates/co-workers and household members.
4.1 The model output will be saved as 'results.csv' file.
4.2 The file consists of three columns, including 'time', 'x-coordinate', and 'y-coordinate' of influenza outbreaks.
Import necessary libraries to run this model.
!pip install --user folium
import numpy as np
import folium, itertools, math, pickle, random, time
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.pylab as pylab
import matplotlib.image as mpimg
import matplotlib.animation as animation
from os import path
from pyproj import Proj, transform
import os
os.environ["PROJ_LIB"] = '/opt/conda/envs/python2/share/proj'
import Env
import Human
Import python codes containing the classes and functions for environment and human behaviors.
Files of households, workplaces, and schools for the Queen Anne neighborhood are imported.
Choose x and y limits (in UTM) and filter houses to select those in the Queen Anne neighborhood.
# Files of households, workplaces, and schools are imported.
CITY = "QueenAnne" # name of the city you're loading data for
houses = pd.read_csv('./data/{}/houses_points.csv'.format(CITY))
schools = pd.read_csv('./data/{}/schools_points.csv'.format(CITY))
works = pd.read_csv('./data/{}/works_points.csv'.format(CITY))
# Choose x and y limits (in UTM) and filter houses
xmin = 546302;
xmax = 549610;
ymin = 5274394;
ymax = 5277978;
nhouses_old = len(houses);
houses = houses[(houses.X >= xmin) & (houses.X <= xmax) & (houses.Y >= ymin) & (houses.Y <= ymax)];
nhouses_new = len(houses);
print(nhouses_new," of ",nhouses_old," houses filtered.")
Inspect first five rows of the house, school, and workplace datasets, then create figures to plot the house, school, and workplace datasets.
house_cat = houses.DESCR.unique().tolist() # Obtain the type of households
house_dict = {}
rd_color =["#"+''.join([random.choice('0123456789ABCDEF') for j in range(6)])
for i in range(len(house_cat))] # generated random color based on the types of households
for i in range(len(house_cat)):
houses.loc[houses['DESCR'] == house_cat[i], 'color'] = rd_color[i] # set value for color column
house_dict[house_cat[i]] = rd_color[i]
housemyProj = Proj("+proj=utm +zone=10N, +north +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
houselonlat = housemyProj(houses['X'].values.tolist(), houses['Y'].values.tolist(), inverse=True)
houses['lat'] = houselonlat[1]
houses['lon'] = houselonlat[0]
houses.head()
school_cat = schools.CODE.unique().tolist() # Obtain the type of schools
school_dict = {}
rd_color =["#"+''.join([random.choice('0123456789ABCDEF') for j in range(6)])
for i in range(len(school_cat))] # generate random color based on the types of schools
for i in range(len(school_cat)):
schools.loc[schools['CODE'] == school_cat[i], 'color'] = rd_color[i] # set value for color column
school_dict[school_cat[i]] = rd_color[i]
schoolmyProj = Proj("+proj=utm +zone=10N, +north +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
schoollonlat = schoolmyProj(schools['X'].values.tolist(), schools['Y'].values.tolist(), inverse=True)
schools['lat'] = schoollonlat[1]
schools['lon'] = schoollonlat[0]
schools.head()
work_cat = works.DESCR.unique().tolist() # obtain the type of workplaces
work_dict = {}
rd_color =["#"+''.join([random.choice('0123456789ABCDEF') for j in range(6)])
for i in range(len(work_cat))] # generate random color based on the types of workplaces
for i in range(len(work_cat)):
works.loc[works['DESCR'] == work_cat[i], 'color'] = rd_color[i] # set value for color column
work_dict[work_cat[i]] = rd_color[i]
workmyProj = Proj("+proj=utm +zone=10N, +north +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
worklonlat = workmyProj(works['X'].values.tolist(), works['Y'].values.tolist(), inverse=True)
works['lat'] = worklonlat[1]
works['lon'] = worklonlat[0]
works.head()
# count the number by house types
for key,value in house_dict.items():
sub_houses = houses[houses.DESCR == key]
print("The number of", key, "households is ",len(sub_houses)) # counting the number of houeholds
print("The number of households in total:", len(houses)) # counting the number of households
m = folium.Map(location=[47.638972, -122.361579], tiles='cartodbpositron', zoom_start=14) for i in range(0, len(houses)): folium.CircleMarker( location=[houses.iloc[i]['lat'], houses.iloc[i]['lon']], popup=houses.iloc[i]['DESCR'], radius=2, color='grey', weight=0.1, fill=True, fill_color=houses.iloc[i]['color'], fill_opacity=0.6, legend_name = 'Number of incidents per district' ).add_to(m) legend_html = '''
Legend''' m.get_root().html.add_child(folium.Element(legend_html)) m.save('maps/houses_interactive_map.html')
''' color_leg = "" for key, value in house_dict.items(): color_leg += ' ' +' ' color_leg += key + '
' legend_html += color_leg + '
# display the interactive map for houses
from IPython.display import IFrame
IFrame(src='maps/houses_interactive_map.html', width=900, height=700)
# count the number by school types
for key,value in school_dict.items():
sub_schools = schools[schools.CODE == key]
print("The number of", key, "schools is ",len(sub_schools)) # counting the number of schools
print("The number of schools in total:", len(schools)) # counting the number of schools
m = folium.Map(location=[47.638972, -122.361579], tiles='cartodbpositron', zoom_start=14) for i in range(0, len(schools)): folium.CircleMarker( location=[schools.iloc[i]['lat'], schools.iloc[i]['lon']], popup=schools.iloc[i]['NAME'], radius=5, color='grey', weight=0.1, fill=True, fill_color=schools.iloc[i]['color'], fill_opacity=0.6, legend_name = 'Number of incidents per district' ).add_to(m)
legend_html = '''Legend''' m.get_root().html.add_child(folium.Element(legend_html)) m.save('maps/schools_interactive_map.html')
''' color_leg = "" for key, value in school_dict.items(): color_leg += ' ' +' ' color_leg += key + '
' legend_html += color_leg + '
# display the interactive map
from IPython.display import IFrame
IFrame(src='maps/schools_interactive_map.html', width=900, height=700)
# count the number by workplace types
for key, value in work_dict.items():
sub_works = works[works.DESCR == key]
print("The number of", key, "workplace is ",len(sub_works)) # counting the number of schools
print("The number of workplace in total:", len(works)) # counting the number of schools
m = folium.Map(location=[47.638972, -122.361579], tiles='cartodbpositron', zoom_start=14) for i in range(0, len(works)): folium.CircleMarker( location=[works.iloc[i]['lat'], works.iloc[i]['lon']], popup=works.iloc[i]['DESCR'], radius=2, color='grey', weight=0.1, fill=True, fill_color=works.iloc[i]['color'], fill_opacity=0.6 ).add_to(m) legend_html = '''
Legend''' m.get_root().html.add_child(folium.Element(legend_html)) m.save('maps/works_interactive_map.html')
''' color_leg = "" for key, value in work_dict.items(): color_leg += ' ' +' ' color_leg += key + '
' legend_html += color_leg + '
# display the interactive map
from IPython.display import IFrame
IFrame(src='maps/works_interactive_map.html', width=900, height=700)
Before running the simulation, the environments (schools, workplaces, households) and individual human agents needs to be set up.
The functions of settingSchools (), settingWorks(), and settingHouseholds() are included in Env.py.
The function of settingHumanAgent() is included in Human.py.
try: # try to open the preprocessed data
with open("data/{}/schoolList.pickle".format(CITY), "rb") as f:
schoolList = pickle.load(f)
print("...loaded school data from pickle...")
except: # if data does not exist, write it
schoolList = Env.settingSchools('./data/{}'.format(CITY))
with open("data/{}/schoolList.pickle".format(CITY), "wb") as f:
pickle.dump(schoolList, f)
try:
with open("data/{}/workList.pickle".format(CITY), "rb") as f:
workList = pickle.load(f)
print("...loaded work list from pickle...")
except:
workList = Env.settingWorks('./data/{}'.format(CITY))
with open("data/{}/workList.pickle".format(CITY), "wb") as f:
pickle.dump(workList, f)
try:
with open("data/{}/houseList.pickle".format(CITY), "rb") as f:
houseList = pickle.load(f)
print("...loaded house list from pickle...")
except:
houseList = Env.settingHouseholds('./data/{}'.format(CITY), schoolList, workList, houses = houses)
with open("data/{}/houseList.pickle".format(CITY), "wb") as f:
pickle.dump(houseList, f)
try:
with open("data/{}/peopleList.pickle".format(CITY), "rb") as f:
peopleList = pickle.load(f)
print("...loaded people list from pickle...")
except:
peopleList = Human.settingHumanAgent(houseList)
with open("data/{}/peopleList.pickle".format(CITY), "wb") as f:
pickle.dump(peopleList, f)
Inspect the age distribution of the human agents.
# make a histogram for age distributions of human agents
age = [peopleList[i].age for i in range(0,len(peopleList))] # make a list for people's age
plt.figure(figsize=(10,10))
plt.hist(age, bins=20, color='orange', edgecolor='black', linewidth=1.2)
plt.xlabel("Age", size=15)
plt.ylabel("Number of people", size=15)
print("Total population:", len(peopleList))
We are working with a very large simulated data set that networkx is incapable of visualizing, but it is able to visualize subsets (or rather subgraphs)! We have provided two ways to subset your graphs, randomly subsampling by person ID and by taking all people in a geographic subset of our study area.
To help ease the computational strain on networkx, we are using something called a line graph rather than a graph. A line graph of an undirected graph G, is a graph that represents the adjacencies between the edges of G. In our case, this means instead of the nodes being people and the edges being places they visit, the nodes are places and the edges are people that visit both places. Because of our model assumptions that people are tied to their houses, schools, and places of work, a line graph reduces the graph while preserving all connectivity information.
This is the line graph for Queen Anne by randomly selecting every 10th person:
This is the line graph for Queen Anne by taking a geographic subset (people within the bounding box x in (547,956 - 548,783) and y in (5,277,186 - 5,277,082))
introRate: the number of people who are randomized as exposed at the beginning of the simulation.
reproduction: the expected average number of secondary infection cases taking place per a single infection case.
infRate: the probability that an individual human becomes infectious when he/she contacts an infectious person.
days: the duration for the simulation to be run.
Attention!!! Now, this model is assumed to be run for 30 days, but to be fully completed, the days parameter should be 155.
import ipywidgets as widgets
introRateW = widgets.IntSlider(min=0,max=10,value=5, description="introRate")
reproW = widgets.FloatSlider(min=0, max=5, step=0.1,value=2.4, description="reproduction")
infW = widgets.FloatSlider(min=0, max=1, step=.01, value=0.18, description="infRate")
daysW = widgets.IntSlider(min=0, max=155, value = 30, description="days")
display(introRateW, reproW, infW, daysW)
introRate = float(introRateW.value)/10000.0
reproduction = reproW.value
infRate = infW.value
days = daysW.value # the list for counting infectious people over simulation days.
#initial infections based on the introRate
peopleList = Env.initialInfection(peopleList, introRate)
from tqdm.notebook import tqdm
total_simulation_time=0
infectiousCount=[]
exposedX = []
exposedY = []
exposedT = []
for t in tqdm(range(1,days), desc="Time Steps"):
start = time.time()
susceptibleList=[]
exposedList=[]
infectiousList=[]
recoveredList=[]
for p in tqdm(range(0,len(peopleList)), desc="Update States", leave=False):
person = peopleList[p]
person.incubating()
person.recovering()
if (person.S==True):
susceptibleList.append(person)
if (person.E==True):
exposedList.append(person)
if (person.I==True):
infectiousList.append(person)
if (person.R==True):
recoveredList.append(person)
if (person.recT is None):
person.recT = t;
if (person.infT==t):
exposedX.append(person.infX)
exposedY.append(person.infY)
exposedT.append(person.infT)
infectiousCount.append(len(infectiousList))
for p in tqdm(peopleList, desc="Infection", leave=False):
p.infecting(peopleList, infRate, reproduction, t) #infecting function is included in Human.py
iteration_time = time.time()-start
total_simulation_time+=iteration_time
print("Simulation Day: {:3d}, Remaining Days: {:3d}, # of infectious people: {:5d} took {:4.2f} seconds".format(t, days-t,len(exposedList),iteration_time))
total_simulation_time/=60.0 # get minutes
hours = False
daysBool = False
if total_simulation_time >= 60:
total_simulation_time/=60.0 # get hours
hours = True
if total_simulation_time >= 24:
total_simulation_time/=24.0
daysBool = True
if not (hours and days):
print("\n\n Total Simulation Time: {:.1f} minutes.".format(total_simulation_time))
elif not daysBool:
print("\n\n Total Simulation Time: {:.1f} hours.".format(total_simulation_time))
else:
print("\n\n Total Simulation Time: {:.1f} days.".format(total_simulation_time))
This saves the simulation results to a .csv file to test for stochasticity later. Each time the simulation is run, every function starting with "Synthesize a Population" must be rerun.
## Visualization of the simulated flu cases
realdata = pd.read_csv('./data/{}/fludata.csv'.format(CITY))
plt.figure(figsize=(14,7))
plt.plot(realdata['Date'], realdata['Positive'])
infCases=pd.DataFrame({'time':range(1,days),'infectious':infectiousCount})
plt.figure(figsize=(14,7))
plt.plot(infCases['time'],infCases['infectious'])
plt.title('The number of infectious people over day of year', size=20)
plt.xlabel("Days", size=15)
plt.ylabel("Number of people", size=15)
## Visualization of actual flu data from the Queen Anne neighbhorhood in Seattle
realdata = pd.read_csv('./data/{}/fludata.csv'.format(CITY))
plt.figure(figsize=(14,7))
plt.plot(realdata['Date'], realdata['Positive'])
plt.title('The number of infectious people over day of year', size=20)
plt.xlabel("Date", size=15)
plt.ylabel("Number of infectious people", size=15)
import matplotlib.ticker as ticker
ax = plt.axes()
ax.xaxis.set_major_locator(ticker.MultipleLocator(5))
# Make the list of information of influenza outbreaks.
houseX=[peopleList[i].houseX for i in range(len(peopleList)) if peopleList[i].infX != None]
houseY=[peopleList[i].houseY for i in range(len(peopleList)) if peopleList[i].infY != None]
exposedX=[peopleList[i].infX for i in range(len(peopleList)) if peopleList[i].infX != None]
exposedY=[peopleList[i].infY for i in range(len(peopleList)) if peopleList[i].infY != None]
exposedT=[peopleList[i].infT for i in range(len(peopleList)) if peopleList[i].infX != None]
recoveredX=[peopleList[i].houseX for i in range(len(peopleList)) if peopleList[i].recT != None]
recoveredY=[peopleList[i].houseY for i in range(len(peopleList)) if peopleList[i].recT != None]
recoveredT=[peopleList[i].recT for i in range(len(peopleList)) if peopleList[i].recT != None]
Output file will be saved as results.csv
results = pd.DataFrame({'time':exposedT, 'xCoord':exposedX, 'yCoord':exposedY})
results.to_csv('data/{}/results.csv'.format(CITY))
results = pd.read_csv('data/{}/results.csv'.format(CITY))
xmin = 546302;
xmax = 549610;
ymin = 5274394;
ymax = 5277978;
for t in range(days):
fig = plt.figure() # size of figure
ilistX = [houseX[i] for i in range(len(exposedX)) if exposedT[i] <= t]
ilistY = [houseY[i] for i in range(len(exposedY)) if exposedT[i] <= t]
rlistX = [recoveredX[i] for i in range(len(recoveredX)) if recoveredT[i] <= t]
rlistY = [recoveredY[i] for i in range(len(recoveredY)) if recoveredT[i] <= t]
plt.scatter(houses['X'],houses['Y'], s = 0.2, c = '0.8') # drawing scatter plots of households' locations
plt.scatter(ilistX,ilistY, s = 10, c = 'r', label ='Infectious cases')
plt.scatter(rlistX,rlistY, s = 10, c = 'b', label ='Recovered cases')
plt.xlabel("UTM X-Coordinates", size=15)
plt.ylabel("UTM Y-Coordinates", size=15)
plt.xlim(xmin,xmax)
plt.ylim(ymin,ymax)
plt.title('Time = '+str(t)+' days, I(t) = '+str(len(ilistX)-len(rlistX)), size = 20)
plt.axes().set_aspect('equal', 'datalim')
plt.legend(markerscale = 2.0)
plt.tight_layout()
fig.savefig('img/Figures/figure%d.png' % t, dpi=300)
pylab.close(fig);
frames = []
fig = plt.figure(figsize=(15,15))
plt.xticks([])
plt.yticks([])
plt.tight_layout()
for t in range(days):
img = mpimg.imread('img/Figures/figure%d.png' % t)
frames.append([plt.imshow(img, animated=True)])
pylab.close(fig);
ani = animation.ArtistAnimation(fig, frames, interval=500, blit=False)
ani.save('img/movie.mp4')
import io
import base64
from IPython.display import HTML
video = io.open("img/movie.mp4", 'r+b').read()
encoded = base64.b64encode(video)
HTML(data='''<video width="960p" height="720p" alt="test" controls>
<source src="data:video/mp4;base64,{0}" type="video/mp4" />
</video>'''.format(encoded.decode('ascii')))
import geopandas as gpd
county_geo = gpd.read_file('./data/shapefiles/qanne.geojson')
county_geo_data = county_geo.to_crs(epsg='4326').to_json()
result_data = pd.read_csv("./data/QueenAnne/results.csv")
myProj = Proj("+proj=utm +zone=10N, +north +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
lonlat = myProj(result_data['xCoord'].values.tolist(), result_data['yCoord'].values.tolist(), inverse=True)
result_data['lat'] = lonlat[1]
result_data['lon'] = lonlat[0]
from folium.plugins import HeatMap
heat_data = [[row['lat'],row['lon']] for index, row in result_data.iterrows()]
heat_data_2 = [[[row['lat'],row['lon']] for index,
row in result_data[result_data['time'] == i].iterrows()] for i in range(0,155)]
from folium import plugins
m = folium.Map(location=[47.648972, -122.361579], tiles = 'cartodbpositron', zoom_start=13)
heatmap = plugins.HeatMapWithTime(heat_data_2,auto_play=True,max_opacity=1,
radius = 20,
name='Influenza Cases Time Series',
gradient={.6: 'PINK', .95: 'RED', 1: 'CRIMSON'})
heatmap.add_to(m)
m.save('./maps/HeatMap.html')
from IPython.display import IFrame
IFrame(src='./maps/HeatMap.html', width=900, height=1000)
Replicability
Reproducibility
Explore future publications
Florida Health. (2019). Florida Influenza Surveillance Reports. Accessed July 2019. http://www.floridahealth.gov/diseases-and-conditions/influenza/florida-influenza-surveillance-report-archive/index.html.
Guo, D., Li, K. C., Peters, T. R., Snively, B. M., Poehling, K. A., & Zhou, X. (2015). Multi-scale modeling for the transmission of influenza and the evaluation of interventions toward it. Scientific reports, 5, 8980. Accessed July 2019. https://www.nature.com/articles/srep08980.
Helmholtz. (2016). The ODD Protocol for describing individual-based and agent-based models. Accessed July 2019. https://www.ufz.de/export/data/2/100067_ODD_Update_template.doc.
Kang, J. Y., & Aldstadt, J. (2019). Using multiple scale space-time patterns in variance-based global sensitivity analysis for spatially explicit agent-based models. Computers, Environment and Urban Systems, 75, 170-183.
Lowen, A. C., Mubareka, S., Steel, J., & Palese, P. (2007). Influenza virus transmission is dependent on relative humidity and temperature. PLoS pathogens, 3(10), e151.
Mao, L., & Bian, L. (2010). Spatial–temporal transmission of influenza and its health risks in an urbanized area. Computers, Environment and Urban Systems, 34(3), 204-215.
QUEBS. (2013). The ODD protocol for documenting agent-based models. Accessed July 2019. https://qubeshub.org/community/groups/odd.
Seattle & King County Public Health. (2019). King County Influenza Update: May 4, 2019. Accessed July 2019. https://www.kingcounty.gov/depts/health/communicable-diseases/disease-control/~/media/depts/health/communicable-diseases/documents/influenza/2019/week-18.ashx.
Tracy, Melissa, Magdalena Cerda, and Katherine M. Keyes. (2018). Agent-Based Modeling in Public Health: Current Applications and Future Directions. Annual Review of Public Health 39:77–94. Accessed July 2019. https://doi.org/10.1146/annurev-publhealth040617-014317.
Willem, Lander. (2015). Agent-Based Models For Infectious Disease Transmission: Exploration, Estimation & Computational Efficiency. PhD dissertation. University of Antwerp. Accessed July 2019. https://ibiostat.be/publications/phd/landerwillem.pdf.