|
18 | 18 | class SimSirModel:
|
19 | 19 |
|
20 | 20 | def __init__(self, p: Parameters) -> SimSirModel:
|
| 21 | + # TODO missing initial recovered value |
| 22 | + susceptible = p.susceptible |
| 23 | + recovered = 0.0 |
| 24 | + recovery_days = p.recovery_days |
| 25 | + |
| 26 | + rates = { |
| 27 | + key: d.rate |
| 28 | + for key, d in p.dispositions.items() |
| 29 | + } |
| 30 | + |
| 31 | + lengths_of_stay = { |
| 32 | + key: d.length_of_stay |
| 33 | + for key, d in p.dispositions.items() |
| 34 | + } |
21 | 35 |
|
22 | 36 | # Note: this should not be an integer.
|
23 | 37 | # We're appoximating infected from what we do know.
|
24 | 38 | # TODO market_share > 0, hosp_rate > 0
|
25 |
| - self.infected = infected = ( |
| 39 | + infected = ( |
26 | 40 | p.current_hospitalized / p.market_share / p.hospitalized.rate
|
27 | 41 | )
|
28 | 42 |
|
29 |
| - self.detection_probability = ( |
| 43 | + detection_probability = ( |
30 | 44 | p.known_infected / infected if infected > 1.0e-7 else None
|
31 | 45 | )
|
32 | 46 |
|
33 |
| - # TODO missing initial recovered value |
34 |
| - self.recovered = recovered = 0.0 |
35 |
| - |
36 |
| - self.intrinsic_growth_rate = intrinsic_growth_rate = \ |
| 47 | + intrinsic_growth_rate = \ |
37 | 48 | (2.0 ** (1.0 / p.doubling_time) - 1.0) if p.doubling_time > 0.0 else 0.0
|
38 | 49 |
|
39 |
| - self.gamma = gamma = 1.0 / p.recovery_days |
| 50 | + gamma = 1.0 / recovery_days |
40 | 51 |
|
41 | 52 | # Contact rate, beta
|
42 |
| - self.beta = beta = ( |
| 53 | + beta = ( |
43 | 54 | (intrinsic_growth_rate + gamma)
|
44 |
| - / p.susceptible |
| 55 | + / susceptible |
45 | 56 | * (1.0 - p.relative_contact_rate)
|
46 | 57 | ) # {rate based on doubling time} / {initial susceptible}
|
47 | 58 |
|
48 | 59 | # r_t is r_0 after distancing
|
49 |
| - self.r_t = beta / gamma * p.susceptible |
| 60 | + r_t = beta / gamma * susceptible |
50 | 61 |
|
51 | 62 | # Simplify equation to avoid division by zero:
|
52 | 63 | # self.r_naught = r_t / (1.0 - relative_contact_rate)
|
53 |
| - self.r_naught = (intrinsic_growth_rate + gamma) / gamma |
54 |
| - |
55 |
| - # doubling time after distancing |
56 |
| - # TODO constrain values np.log2(...) > 0.0 |
57 |
| - self.doubling_time_t = 1.0 / np.log2( |
| 64 | + r_naught = (intrinsic_growth_rate + gamma) / gamma |
| 65 | + doubling_time_t = 1.0 / np.log2( |
58 | 66 | beta * p.susceptible - gamma + 1)
|
59 | 67 |
|
60 |
| - self.raw_df = raw_df = sim_sir_df( |
| 68 | + raw_df = sim_sir_df( |
61 | 69 | p.susceptible,
|
62 | 70 | infected,
|
63 | 71 | recovered,
|
64 | 72 | beta,
|
65 | 73 | gamma,
|
66 | 74 | p.n_days,
|
67 | 75 | )
|
68 |
| - |
69 |
| - rates = { |
70 |
| - key: d.rate |
71 |
| - for key, d in p.dispositions.items() |
72 |
| - } |
73 |
| - |
74 |
| - lengths_of_stay = { |
75 |
| - key: d.length_of_stay |
76 |
| - for key, d in p.dispositions.items() |
77 |
| - } |
78 |
| - |
79 |
| - i_dict_v = get_dispositions(raw_df.infected, rates, p.market_share) |
80 |
| - r_dict_v = get_dispositions(raw_df.recovered, rates, p.market_share) |
81 |
| - |
82 |
| - self.dispositions = { |
83 |
| - key: value + r_dict_v[key] |
84 |
| - for key, value in i_dict_v.items() |
85 |
| - } |
86 |
| - |
87 |
| - self.dispositions_df = pd.DataFrame(self.dispositions) |
88 |
| - self.admits_df = admits_df = build_admits_df(p.n_days, self.dispositions) |
89 |
| - self.census_df = build_census_df(admits_df, lengths_of_stay) |
| 76 | + dispositions_df = build_dispositions_df(raw_df, rates, p.market_share) |
| 77 | + admits_df = build_admits_df(dispositions_df) |
| 78 | + census_df = build_census_df(admits_df, lengths_of_stay) |
| 79 | + |
| 80 | + self.susceptible = susceptible |
| 81 | + self.infected = infected |
| 82 | + self.recovered = recovered |
| 83 | + |
| 84 | + self.detection_probability = detection_probability |
| 85 | + self.recovered = recovered |
| 86 | + self.intrinsic_growth_rate = intrinsic_growth_rate |
| 87 | + self.gamma = gamma |
| 88 | + self.beta = beta |
| 89 | + self.r_t = r_t |
| 90 | + self.r_naught = r_naught |
| 91 | + self.doubling_time_t = doubling_time_t |
| 92 | + self.raw_df = raw_df |
| 93 | + self.dispositions_df = dispositions_df |
| 94 | + self.admits_df = admits_df |
| 95 | + self.census_df = census_df |
90 | 96 |
|
91 | 97 |
|
92 | 98 | def sir(
|
@@ -119,55 +125,49 @@ def gen_sir(
|
119 | 125 |
|
120 | 126 |
|
121 | 127 | def sim_sir_df(
|
122 |
| - s: float, i: float, r: float, beta: float, gamma: float, n_days |
| 128 | + s: float, i: float, r: float, beta: float, gamma: float, n_days: int |
123 | 129 | ) -> pd.DataFrame:
|
124 | 130 | """Simulate the SIR model forward in time."""
|
125 | 131 | return pd.DataFrame(
|
126 | 132 | data=gen_sir(s, i, r, beta, gamma, n_days),
|
127 | 133 | columns=("day", "susceptible", "infected", "recovered"),
|
128 | 134 | )
|
129 | 135 |
|
130 |
| - |
131 |
| -def get_dispositions( |
132 |
| - patients: np.ndarray, |
| 136 | +def build_dispositions_df( |
| 137 | + sim_sir_df: pd.DataFrame, |
133 | 138 | rates: Dict[str, float],
|
134 | 139 | market_share: float,
|
135 |
| -) -> Dict[str, np.ndarray]: |
| 140 | +) -> pd.DataFrame: |
136 | 141 | """Get dispositions of patients adjusted by rate and market_share."""
|
137 |
| - return { |
138 |
| - key: patients * rate * market_share |
139 |
| - for key, rate in rates.items() |
140 |
| - } |
141 |
| - |
142 |
| - |
143 |
| -def build_admits_df(n_days, dispositions) -> pd.DataFrame: |
144 |
| - """Build admits dataframe from Parameters and Model.""" |
145 |
| - days = np.arange(0, n_days + 1) |
146 |
| - projection = pd.DataFrame({ |
147 |
| - "day": days, |
148 |
| - **dispositions, |
| 142 | + patients = sim_sir_df.infected + sim_sir_df.recovered |
| 143 | + return pd.DataFrame({ |
| 144 | + "day": sim_sir_df.day, |
| 145 | + **{ |
| 146 | + key: patients * rate * market_share |
| 147 | + for key, rate in rates.items() |
| 148 | + } |
149 | 149 | })
|
150 |
| - # New cases |
151 |
| - admits_df = projection.iloc[:-1, :] - projection.shift(1) |
152 |
| - admits_df["day"] = range(admits_df.shape[0]) |
| 150 | + |
| 151 | + |
| 152 | +def build_admits_df(dispositions_df: pd.DataFrame) -> pd.DataFrame: |
| 153 | + """Build admits dataframe from dispositions.""" |
| 154 | + admits_df = dispositions_df.iloc[:-1, :] - dispositions_df.shift(1) |
| 155 | + admits_df.day = dispositions_df.day |
153 | 156 | return admits_df
|
154 | 157 |
|
155 | 158 |
|
156 | 159 | def build_census_df(
|
157 |
| - admits_df: pd.DataFrame, lengths_of_stay |
| 160 | + admits_df: pd.DataFrame, |
| 161 | + lengths_of_stay: Dict[str, int], |
158 | 162 | ) -> pd.DataFrame:
|
159 |
| - """ALOS for each category of COVID-19 case (total guesses)""" |
160 |
| - n_days = np.shape(admits_df)[0] |
161 |
| - census_dict = {} |
162 |
| - for key, los in lengths_of_stay.items(): |
163 |
| - census = ( |
164 |
| - admits_df.cumsum().iloc[:-los, :] |
165 |
| - - admits_df.cumsum().shift(los).fillna(0) |
166 |
| - ).apply(np.ceil) |
167 |
| - census_dict[key] = census[key] |
168 |
| - |
169 |
| - census_df = pd.DataFrame(census_dict) |
170 |
| - census_df["day"] = census_df.index |
171 |
| - census_df = census_df[["day", *lengths_of_stay.keys()]] |
172 |
| - census_df = census_df.head(n_days) |
173 |
| - return census_df |
| 163 | + """ALOS for each disposition of COVID-19 case (total guesses)""" |
| 164 | + return pd.DataFrame({ |
| 165 | + 'day': admits_df.day, |
| 166 | + **{ |
| 167 | + key: ( |
| 168 | + admits_df[key].cumsum().iloc[:-los] |
| 169 | + - admits_df[key].cumsum().shift(los).fillna(0) |
| 170 | + ).apply(np.ceil) |
| 171 | + for key, los in lengths_of_stay.items() |
| 172 | + } |
| 173 | + }) |
0 commit comments