improved comparison

src/covid19sim/log/track.py

@@ -1877,7 +1877,9 @@ def write_metrics(self):
     def write_for_training(self, humans, outfile, conf):
         """ Writes some data out for the ML predictor """
         data = dict()
-        data['hospitalization_per_day'] = self.hospitalization_per_day
+        data['hospitalization_per_day'] = self.hospitalization_per_day # how many new people get admitted to the hospital
+        import pdb; pdb.set_trace()
+        # TODO: write total number of _currently hospitalized people_, not how many people _go into the hospital_.
 
         # parse test results
         data['positive_test_results_per_day'] = []

src/covid19sim/plotting/compare_real_and_sim.py

@@ -8,13 +8,13 @@
 
 # Constants
 quebec_population = 8485000
-csv_path = "/Users/akpateln/Desktop/Personal/COVID19Tracker.ca Data - QC.csv"
-sim_dir_path = "/Users/akpateln/Desktop/Personal/covi-simulator/output/sim_v2_people-5000_days-30_init-0.002_uptake--1.0_seed-0_20200901-115548_654000"
-#sim_dir_path = "output/sim_v2_people-1000_days-150_init-0.004_uptake--1.0_seed-0_20200716-101640_592894/"
+csv_path = "COVID19Tracker.ca Data - QC.csv"
+# sim_dir_path = "output/sim_v2_people-20000_days-120_init-0.005_uptake--1.0_seed-0_20200909-143542_235853"
+sim_dir_path = "output/sim_v2_people-30000_days-120_init-0.005_uptake--1.0_seed-0_20200909-175946_676450"
 sim_priors_path = os.path.join(sim_dir_path, "train_priors.pkl")
 # sim_tracker_path = os.path.join(sim_dir_path, "tracker_data_n_100_seed_0_20200716-101538_.pkl")
 # ./output/sim_v2_people-5000_days-30_init-0.05_uptake--1.0_seed-0_20200831-105018_963000
-sim_tracker_path = os.path.join(sim_dir_path, "tracker_data_n_5000_seed_0_20200901-120248.pkl")
+sim_tracker_path = os.path.join(sim_dir_path, "tracker_data_n_30000_seed_0_20200909-232724.pkl")
 
 plot_real = False
 
@@ -28,15 +28,24 @@
 def parse_tracker(sim_tracker_data):
     sim_pop = sim_tracker_data['n_humans']
     dates = []
-    deaths = []
+    cumulative_deaths = []
     tests = Counter()
 
     # Get dates and deaths
     for k, v in sim_tracker_data['human_monitor'].items():
         dates.append(str(k))
         death = sum([x['dead'] for x in v])
-        deaths.append(float(death) * 100 / sim_pop)
-        # TODO: Add hospitalizations and case counts/positive tests
+        cumulative_deaths.append(death)
+
+    daily_deaths_prop = []
+    last_deaths = 0
+    for idx, deaths in enumerate(cumulative_deaths):
+        if idx == 0:
+            daily_deaths_prop.append(0)
+        else:
+            deaths = (float(deaths) * 100 / sim_pop)
+            daily_deaths_prop.append(deaths - last_deaths)
+            last_deaths = deaths
 
     # Get tests
     for test in sim_tracker_data['test_monitor']:
@@ -46,7 +55,7 @@ def parse_tracker(sim_tracker_data):
     # Get cases
     cases = [float(x) * 100 / sim_pop for x in sim_tracker_data['cases_per_day']]
 
-    return dates, deaths, tests, cases
+    return dates, daily_deaths_prop, tests, cases
 
 # Parse data
 sim_dates, sim_deaths, sim_tests, sim_cases = parse_tracker(sim_tracker_data)
@@ -55,32 +64,35 @@ def parse_tracker(sim_tracker_data):
 
 # Plot cases
 fig, ax = plt.subplots(figsize=(7.5, 7.5))
-real_dates = qc_data.loc[34:63, 'data » date'].to_numpy()
-real_cases = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:63, 'data » change_cases']]
-real_hospitalizations = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:63, 'data » total_hospitalizations']]
-real_deaths = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:63,'data » change_fatalities']]
-real_tests = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:63, 'data » change_tests']]
-
-ax.plot(real_dates, real_cases, label="QC Cases (per Day)")
-ax.plot(real_dates, real_hospitalizations, label="QC Hospital Utilization (by Day)")
-ax.plot(real_dates, real_deaths, label="QC Mortalities (per Day)")
+real_dates = qc_data.loc[34:133, 'dates'].to_numpy()
+real_cases = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:133, 'change_cases']]
+real_hospitalizations = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:133, 'total_hospitalizations']]
+real_deaths = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:133,'change_fatalities']]
+real_tests = [100 * float(x if str(x) != "nan" else 0) / quebec_population for x in qc_data.loc[34:133, 'change_tests']]
+
+ax.plot(real_dates, real_cases, label="Diagnoses per day")
+ax.plot(real_dates, real_hospitalizations, label="Hospital utilization")
+ax.plot(real_dates, real_deaths, label="Mortalities per day")
 # ax.plot(real_dates, real_tests, label="QC Tests (per Day)")
 
 ax.legend()
 plt.ylabel("Percentage of Population")
 plt.xlabel("Date")
-plt.yticks(plt.yticks()[0], [str(round(x, 2)) + "%" for x in plt.yticks()[0]])
+plt.yticks(plt.yticks()[0], [str(round(x, 3)) + "%" for x in plt.yticks()[0]])
 plt.xticks([x for i, x in enumerate(real_dates) if i % 10 == 0], rotation=45)
-plt.title("Quebec COVID Statistics")
-plt.savefig("qc_stats30.png")
+plt.title("Real Quebec COVID Statistics")
+plt.savefig("qc_stats.png")
 
 # Plot deaths and hospitalizations
 fig, ax = plt.subplots(figsize=(7.5, 7.5))
+ax.set_ylim([0.,0.03])
+sim_cases = np.array(sim_cases[1:])*.1 # adjust for unknown cases
 
-ax.plot(sim_dates, sim_deaths, label="Simulated Mortalities (per Day)")
-ax.plot(sim_dates, sim_hospitalizations[1:], label="Simulated Hospital Utilization (per Day)")
-ax.plot(sim_dates, sim_cases[1:], label="Simulated Cases (per Day)")
-ax.plot(list(sim_tests.keys()), list(sim_tests.values()), label="Simulated Tests (per Day)")
+ax.plot(sim_dates, sim_cases, label="Diagnoses per day")
+ax.plot(sim_dates, sim_hospitalizations[1:], label="Hospital admissions per day")
+ax.plot(sim_dates, sim_deaths, label="Mortalities per day")
+
+# ax.plot(list(sim_tests.keys()), list(sim_tests.values()), label="Simulated Tests (per Day)")
 
 # Goodness of Fit
 '''
@@ -94,9 +106,9 @@ def parse_tracker(sim_tracker_data):
 ax.legend()
 plt.ylabel("Percentage of Population")
 plt.xlabel("Date")
-plt.yticks(plt.yticks()[0], [str(round(x, 2)) + "%" for x in plt.yticks()[0]])
+plt.yticks(plt.yticks()[0], [str(round(x, 3)) + "%" for x in plt.yticks()[0]])
 plt.xticks([x for i, x in enumerate(real_dates) if i % 10 == 0], rotation=45)
-plt.title("Quebec COVID Statistics")
+plt.title("Simulated Quebec COVID Statistics")
 # plt.figtext(0.5, 0.01, fit_caption, wrap=True, horizontalalignment='center', fontsize=12)
 plt.savefig("sim_stats.png")