Skip to content

Commit e2b865b

Browse files
adding app.py
2 parents 4fb0ef0 + 6dac2d7 commit e2b865b

File tree

2 files changed

+291
-1
lines changed

2 files changed

+291
-1
lines changed

README.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -43,4 +43,4 @@ port, then set `PORT` in the `.env` file.
4343
## Deployment
4444
**Before you push your changes to master make sure that everything works in development mode.**
4545

46-
Changes merged to `master` will be automatically deployed to [pennchime.herokuapp.com](http://predictivehealthcare.pennmedicine.org/chime).
46+
Changes merged to `master` will be automatically deployed to [http://penn-chime.phl.io/](http://predictivehealthcare.pennmedicine.org/chime).

app.py

Lines changed: 290 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,290 @@
1+
import pandas as pd
2+
import streamlit as st
3+
import numpy as np
4+
import matplotlib
5+
matplotlib.use('Agg')
6+
import matplotlib.pyplot as plt
7+
8+
hide_menu_style = """
9+
<style>
10+
#MainMenu {visibility: hidden;}
11+
</style>
12+
"""
13+
st.markdown(hide_menu_style, unsafe_allow_html=True)
14+
15+
delaware = 564696
16+
chester = 519293
17+
montgomery = 826075
18+
bucks = 628341
19+
philly = 1581000
20+
S_default = delaware + chester + montgomery + bucks + philly
21+
known_infections = 31
22+
23+
# Widgets
24+
initial_infections = st.sidebar.number_input(
25+
"Currently Known Regional Infections", value=known_infections, step=10, format="%i"
26+
)
27+
current_hosp = st.sidebar.number_input(
28+
"Currently Hospitalized COVID-19 Patients", value=2, step=1, format="%i"
29+
)
30+
doubling_time = st.sidebar.number_input(
31+
"Doubling Time (days)", value=6, step=1, format="%i"
32+
)
33+
hosp_rate = (
34+
st.sidebar.number_input("Hospitalization %", 0, 100, value=5, step=1, format="%i")
35+
/ 100.0
36+
)
37+
icu_rate = (
38+
st.sidebar.number_input("ICU %", 0, 100, value=2, step=1, format="%i") / 100.0
39+
)
40+
vent_rate = (
41+
st.sidebar.number_input("Ventilated %", 0, 100, value=1, step=1, format="%i")
42+
/ 100.0
43+
)
44+
hosp_los = st.sidebar.number_input("Hospital LOS", value=7, step=1, format="%i")
45+
icu_los = st.sidebar.number_input("ICU LOS", value=9, step=1, format="%i")
46+
vent_los = st.sidebar.number_input("Vent LOS", value=10, step=1, format="%i")
47+
Penn_market_share = (
48+
st.sidebar.number_input(
49+
"Hospital Market Share (%)", 0.0, 100.0, value=15.0, step=1.0, format="%f"
50+
)
51+
/ 100.0
52+
)
53+
S = st.sidebar.number_input(
54+
"Regional Population", value=S_default, step=100000, format="%i"
55+
)
56+
57+
total_infections = current_hosp / Penn_market_share / hosp_rate
58+
detection_prob = initial_infections / total_infections
59+
60+
st.title("COVID-19 Hospital Impact Model for Epidemics")
61+
st.markdown(
62+
"""*This tool was developed by the [Predictive Healthcare team](http://predictivehealthcare.pennmedicine.org/) at Penn Medicine. For questions and comments please see our [contact page](http://predictivehealthcare.pennmedicine.org/contact/).*"""
63+
)
64+
65+
if st.checkbox("Show more info about this tool"):
66+
st.subheader(
67+
"[Discrete-time SIR modeling](https://mathworld.wolfram.com/SIRModel.html) of infections/recovery"
68+
)
69+
st.markdown(
70+
"""The model consists of individuals who are either _Susceptible_ ($S$), _Infected_ ($I$), or _Recovered_ ($R$).
71+
72+
The epidemic proceeds via a growth and decline process. This is the core model of infectious disease spread and has been in use in epidemiology for many years."""
73+
)
74+
st.markdown("""The dynamics are given by the following 3 equations.""")
75+
76+
st.latex("S_{t+1} = (-\\beta S_t I_t) + S_t")
77+
st.latex("I_{t+1} = (\\beta S_t I_t - \\gamma I_t) + I_t")
78+
st.latex("R_{t+1} = (\\gamma I_t) + R_t")
79+
80+
st.markdown(
81+
"""To project the expected impact to Penn Medicine, we estimate the terms of the model.
82+
83+
To do this, we use a combination of estimates from other locations, informed estimates based on logical reasoning, and best guesses from the American Hospital Association.
84+
85+
86+
### Parameters
87+
First, we need to express the two parameters $\\beta$ and $\\gamma$ in terms of quantities we can estimate.
88+
89+
- The $\\gamma$ parameter represents 1 over the mean recovery time in days. Since the CDC is recommending 14 days of self-quarantine, we'll use $\\gamma = 1/14$.
90+
- Next, the AHA says to expect a doubling time $T_d$ of 7-10 days. That means an early-phase rate of growth can be computed by using the doubling time formula:
91+
"""
92+
)
93+
st.latex("g = 2^{1/T_d} - 1")
94+
95+
st.markdown(
96+
"""
97+
- Since the rate of new infections in the SIR model is $g = \\beta S - \\gamma$, and we've already computed $\\gamma$, $\\beta$ becomes a function of the initial population size of susceptible individuals.
98+
$$\\beta = (g + \\gamma)/s$$
99+
100+
### Initial Conditions
101+
102+
- The total size of the susceptible population will be the entire catchment area for Penn Medicine entities (HUP, PAH, PMC, CCH)
103+
- Delaware = {delaware}
104+
- Chester = {chester}
105+
- Montgomery = {montgomery}
106+
- Bucks = {bucks}
107+
- Philly = {philly}
108+
- The initial number of infected will be the total number of confirmed cases in the area ({initial_infections}), divided by some detection probability to account for under testing {detection_prob:.2f}.""".format(
109+
delaware=delaware,
110+
chester=chester,
111+
montgomery=montgomery,
112+
bucks=bucks,
113+
philly=philly,
114+
initial_infections=initial_infections,
115+
detection_prob=detection_prob,
116+
)
117+
)
118+
119+
# The SIR model, one time step
120+
def sir(y, beta, gamma, N):
121+
S, I, R = y
122+
Sn = (-beta * S * I) + S
123+
In = (beta * S * I - gamma * I) + I
124+
Rn = gamma * I + R
125+
if Sn < 0:
126+
Sn = 0
127+
if In < 0:
128+
In = 0
129+
if Rn < 0:
130+
Rn = 0
131+
132+
scale = N / (Sn + In + Rn)
133+
return Sn * scale, In * scale, Rn * scale
134+
135+
136+
# Run the SIR model forward in time
137+
def sim_sir(S, I, R, beta, gamma, n_days, beta_decay=None):
138+
N = S + I + R
139+
s, i, r = [S], [I], [R]
140+
for day in range(n_days):
141+
y = S, I, R
142+
S, I, R = sir(y, beta, gamma, N)
143+
if beta_decay:
144+
beta = beta * (1 - beta_decay)
145+
s.append(S)
146+
i.append(I)
147+
r.append(R)
148+
149+
s, i, r = np.array(s), np.array(i), np.array(r)
150+
return s, i, r
151+
152+
153+
## RUN THE MODEL
154+
155+
S, I, R = S, initial_infections / detection_prob, 0
156+
157+
intrinsic_growth_rate = 2 ** (1 / doubling_time) - 1
158+
159+
recovery_days = 14.0
160+
# mean recovery rate, gamma, (in 1/days).
161+
gamma = 1 / recovery_days
162+
163+
# Contact rate, beta
164+
beta = (
165+
intrinsic_growth_rate + gamma
166+
) / S # {rate based on doubling time} / {initial S}
167+
168+
169+
n_days = st.slider("Number of days to project", 30, 200, 60, 1, "%i")
170+
171+
beta_decay = 0.0
172+
s, i, r = sim_sir(S, I, R, beta, gamma, n_days, beta_decay=beta_decay)
173+
174+
175+
hosp = i * hosp_rate * Penn_market_share
176+
icu = i * icu_rate * Penn_market_share
177+
vent = i * vent_rate * Penn_market_share
178+
179+
days = np.array(range(0, n_days + 1))
180+
data_list = [days, hosp, icu, vent]
181+
data_dict = dict(zip(["day", "hosp", "icu", "vent"], data_list))
182+
183+
projection = pd.DataFrame.from_dict(data_dict)
184+
185+
st.subheader("New Admissions")
186+
st.markdown("Projected number of **daily** COVID-19 admissions at Penn hospitals")
187+
188+
# New cases
189+
projection_admits = projection.iloc[:-1, :] - projection.shift(1)
190+
projection_admits[projection_admits < 0] = 0
191+
192+
plot_projection_days = n_days - 10
193+
projection_admits["day"] = range(projection_admits.shape[0])
194+
195+
fig, ax = plt.subplots(1, 1, figsize=(10, 4))
196+
ax.plot(
197+
projection_admits.head(plot_projection_days)["hosp"], ".-", label="Hospitalized"
198+
)
199+
ax.plot(projection_admits.head(plot_projection_days)["icu"], ".-", label="ICU")
200+
ax.plot(projection_admits.head(plot_projection_days)["vent"], ".-", label="Ventilated")
201+
ax.legend(loc=0)
202+
ax.set_xlabel("Days from today")
203+
ax.grid("on")
204+
ax.set_ylabel("Daily Admissions")
205+
st.pyplot()
206+
207+
admits_table = projection_admits[np.mod(projection_admits.index, 7) == 0].copy()
208+
admits_table["day"] = admits_table.index
209+
admits_table.index = range(admits_table.shape[0])
210+
admits_table = admits_table.fillna(0).astype(int)
211+
212+
if st.checkbox("Show Projected Admissions in tabular form"):
213+
st.dataframe(admits_table)
214+
215+
st.subheader("Admitted Patients (Census)")
216+
st.markdown(
217+
"Projected **census** of COVID-19 patients, accounting for arrivals and discharges at Penn hospitals"
218+
)
219+
220+
# ALOS for each category of COVID-19 case (total guesses)
221+
222+
los_dict = {
223+
"hosp": hosp_los,
224+
"icu": icu_los,
225+
"vent": vent_los,
226+
}
227+
228+
fig, ax = plt.subplots(1, 1, figsize=(10, 4))
229+
230+
census_dict = {}
231+
for k, los in los_dict.items():
232+
census = (
233+
projection_admits.cumsum().iloc[:-los, :]
234+
- projection_admits.cumsum().shift(los).fillna(0)
235+
).apply(np.ceil)
236+
census_dict[k] = census[k]
237+
ax.plot(census.head(plot_projection_days)[k], ".-", label=k + " census")
238+
ax.legend(loc=0)
239+
240+
ax.set_xlabel("Days from today")
241+
ax.grid("on")
242+
ax.set_ylabel("Census")
243+
st.pyplot()
244+
245+
census_df = pd.DataFrame(census_dict)
246+
census_df["day"] = census_df.index
247+
census_df = census_df[["day", "hosp", "icu", "vent"]]
248+
249+
census_table = census_df[np.mod(census_df.index, 7) == 0].copy()
250+
census_table.index = range(census_table.shape[0])
251+
census_table.loc[0, :] = 0
252+
census_table = census_table.dropna().astype(int)
253+
254+
if st.checkbox("Show Projected Census in tabular form"):
255+
st.dataframe(census_table)
256+
257+
st.markdown(
258+
"""**Click the checkbox below to view additional data generated by this simulation**"""
259+
)
260+
if st.checkbox("Show Additional Projections"):
261+
st.subheader(
262+
"The number of infected and recovered individuals in the hospital catchment region at any given moment"
263+
)
264+
fig, ax = plt.subplots(1, 1, figsize=(10, 4))
265+
ax.plot(i, label="Infected")
266+
ax.plot(r, label="Recovered")
267+
ax.legend(loc=0)
268+
ax.set_xlabel("days from today")
269+
ax.set_ylabel("Case Volume")
270+
ax.grid("on")
271+
st.pyplot()
272+
273+
# Show data
274+
days = np.array(range(0, n_days + 1))
275+
data_list = [days, s, i, r]
276+
data_dict = dict(zip(["day", "susceptible", "infections", "recovered"], data_list))
277+
projection_area = pd.DataFrame.from_dict(data_dict)
278+
infect_table = (projection_area.iloc[::7, :]).apply(np.floor)
279+
infect_table.index = range(infect_table.shape[0])
280+
281+
if st.checkbox("Show Raw SIR Similation Data"):
282+
st.dataframe(infect_table)
283+
284+
st.subheader("References & Acknowledgements")
285+
st.markdown(
286+
"""* AHA Webinar, Feb 26, James Lawler, MD, an associate professor University of Nebraska Medical Center, What Healthcare Leaders Need To Know: Preparing for the COVID-19
287+
* We would like to recognize the valuable assistance in consultation and review of model assumptions by Michael Z. Levy, PhD, Associate Professor of Epidemiology, Department of Biostatistics, Epidemiology and Informatics at the Perelman School of Medicine
288+
"""
289+
)
290+
st.markdown("© 2020, The Trustees of the University of Pennsylvania")

0 commit comments

Comments
 (0)