diff --git a/README.md b/README.md index 0a1a894..cd91fe0 100644 --- a/README.md +++ b/README.md @@ -6,6 +6,10 @@ The official documentation can be found here https://memilio.readthedocs.io . ### Marimo notebook setup +We recommend using the editor VSCode and installing the Python and Marimo extensions. Only this this set up was +tested by us. If you encounter any issues with the installation process, please open an issue on our +[GitHub](https://github.com/SciCompMod/memilio/issues) page. + First, open a new terminal. Go to a directory for the new project, then create a virtual environment (here "venv") and install marimo: diff --git a/exercises/exercise07.py b/exercises/exercise07.py index d7b1dc6..5a7d0ed 100644 --- a/exercises/exercise07.py +++ b/exercises/exercise07.py @@ -1,6 +1,6 @@ import marimo -__generated_with = "0.20.4" +__generated_with = "0.24.0" app = marimo.App(width="medium") @@ -128,7 +128,8 @@ def _(AgeGroup, model, np, num_age_groups): model.parameters.CriticalPerSevere[AgeGroup(2)] = 0.25 * 2 # Set contact frequency - model.parameters.ContactPatterns.cont_freq_mat[0].baseline = np.ones((num_age_groups, num_age_groups)) * 10 + model.parameters.ContactPatterns.cont_freq_mat[0].baseline = np.ones( + (num_age_groups, num_age_groups)) * 10 return @@ -161,8 +162,10 @@ def _(AgeGroup, model, num_age_groups, osecir, total_population_per_region): # The population is equally distributed among the age groups for group in range(num_age_groups): # 1% of the population is initially infected, 0.5% Exposed and 0.5% in the pre- or asymptomatic state - model.populations[AgeGroup(group), osecir.InfectionState.Exposed] = 0.005 * total_population_per_region / num_age_groups - model.populations[AgeGroup(group), osecir.InfectionState.InfectedNoSymptoms] = 0.005 * total_population_per_region / num_age_groups + model.populations[AgeGroup(group), osecir.InfectionState.Exposed] = 0.005 * \ + total_population_per_region / num_age_groups + model.populations[AgeGroup(group), osecir.InfectionState.InfectedNoSymptoms] = 0.005 * \ + total_population_per_region / num_age_groups # The rest of the population is Susceptible model.populations.set_difference_from_group_total_AgeGroup( (AgeGroup(group), osecir.InfectionState.Susceptible), total_population_per_region / num_age_groups) @@ -180,7 +183,7 @@ def _(mo): @app.cell def _(graph, model, t0): # Add node with id 0 and copy beforehand initialized model to it - graph.add_node(id=0, model=model, t0=t0) + graph.add_node(id=0, model=model, t0=t0) return @@ -198,7 +201,8 @@ def _(AgeGroup, model, num_age_groups, osecir, total_population_per_region): for age in range(num_age_groups): # No infected individuals model.populations[AgeGroup(age), osecir.InfectionState.Exposed] = 0 - model.populations[AgeGroup(age), osecir.InfectionState.InfectedNoSymptoms] = 0 + model.populations[AgeGroup( + age), osecir.InfectionState.InfectedNoSymptoms] = 0 # The total population is Susceptible model.populations.set_difference_from_group_total_AgeGroup( (AgeGroup(age), osecir.InfectionState.Susceptible), total_population_per_region / num_age_groups) @@ -251,7 +255,8 @@ def _(mo): ### Exercise Set the mobility coefficients for `Dead` individuals to zero: - Hint: For the first age group this can be done via `mobility_coefficients[0 * int(osecir.InfectionState.Dead) + osecir.InfectionState.Dead] = 0`. + + Hint: For the first age group this can be done via `mobility_coefficients[0 * (int(osecir.InfectionState.Dead)+1) + int(osecir.InfectionState.Dead)] = 0`. """) return @@ -298,21 +303,26 @@ def _(mo): @app.cell -def _(dt_exchange, graph, osecir, t0): - # Create graph simulation and advance until tmax - sim = osecir.MobilitySimulation(graph, t0, dt=dt_exchange) +def _(dt_exchange, graph, osecir, t0, tmax): + import os + # Silence C++ log output; it can deadlock marimo's output pipe on Windows + devnull = os.open(os.devnull, os.O_WRONLY) + saved_stdout, saved_stderr = os.dup(1), os.dup(2) + os.dup2(devnull, 1) + os.dup2(devnull, 2) + try: + # Create graph simulation and advance until tmax + sim = osecir.MobilitySimulation(graph, t0, dt=dt_exchange) + sim.advance(tmax) + finally: + os.dup2(saved_stdout, 1) + os.dup2(saved_stderr, 2) + os.close(devnull) + os.close(saved_stdout) + os.close(saved_stderr) return (sim,) -@app.cell -def _(mo): - mo.md(r""" - ### Exercise - Please advance the simulation until `tmax`. - """) - return - - @app.cell def _(): # Insert code here @@ -344,8 +354,10 @@ def _(mo): @app.cell def _(osecir, result_region0, result_region1): - result_region0_interpolated = osecir.interpolate_simulation_result(result_region0) - result_region1_interpolated = osecir.interpolate_simulation_result(result_region1) + result_region0_interpolated = osecir.interpolate_simulation_result( + result_region0) + result_region1_interpolated = osecir.interpolate_simulation_result( + result_region1) return result_region0_interpolated, result_region1_interpolated @@ -368,8 +380,10 @@ def _(osecir, plt, result_region0_interpolated, result_region1_interpolated): # Plot the number of non-symptomatically infected for both regions fig, ax = plt.subplots() time = result_array_0[0, :] - InfectedNoSymptoms_0 = result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int(osecir.InfectionState.Dead) + 1, :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] - InfectedNoSymptoms_1 = result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int(osecir.InfectionState.Dead) + 1, :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] + InfectedNoSymptoms_0 = result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int( + osecir.InfectionState.Dead) + 1, :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] + InfectedNoSymptoms_1 = result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int( + osecir.InfectionState.Dead) + 1, :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] ax.plot(time, InfectedNoSymptoms_0, label='Infected No Symptoms Region 1') ax.plot(time, InfectedNoSymptoms_1, label='Infected No Symptoms Region 2') ax.set_xlabel('Time [days]') diff --git a/tutorial07.py b/tutorial07.py index da4dee3..bbdb1fe 100644 --- a/tutorial07.py +++ b/tutorial07.py @@ -1,6 +1,6 @@ import marimo -__generated_with = "0.19.11" +__generated_with = "0.24.0" app = marimo.App(width="medium") @@ -49,7 +49,7 @@ def _(): @app.cell def _(mo): mo.md(r""" - We set the simulation start time `t0`, the end time `tmax` and the initial step size `dt` as: + We set the simulation start time `t0` and the end time `tmax`.: """) return @@ -58,7 +58,6 @@ def _(mo): def _(): t0 = 0 tmax = 100 - dt = 0.1 return t0, tmax @@ -104,20 +103,20 @@ def _(mo): @app.cell def _(AgeGroup, model, np, num_age_groups): - for ag in range(num_age_groups): + for i in range(num_age_groups): # Set infection state stay times (in days) - model.parameters.TimeExposed[AgeGroup(ag)] = 3.2 - model.parameters.TimeInfectedNoSymptoms[AgeGroup(ag)] = 2. - model.parameters.TimeInfectedSymptoms[AgeGroup(ag)] = 6. - model.parameters.TimeInfectedSevere[AgeGroup(ag)] = 12. - model.parameters.TimeInfectedCritical[AgeGroup(ag)] = 8. + model.parameters.TimeExposed[AgeGroup(i)] = 3.2 + model.parameters.TimeInfectedNoSymptoms[AgeGroup(i)] = 2. + model.parameters.TimeInfectedSymptoms[AgeGroup(i)] = 6. + model.parameters.TimeInfectedSevere[AgeGroup(i)] = 12. + model.parameters.TimeInfectedCritical[AgeGroup(i)] = 8. # Set infection state transition probabilities - model.parameters.RelativeTransmissionNoSymptoms[AgeGroup(ag)] = 0.67 - model.parameters.TransmissionProbabilityOnContact[AgeGroup(ag)] = 0.1 - model.parameters.RecoveredPerInfectedNoSymptoms[AgeGroup(ag)] = 0.2 - model.parameters.RiskOfInfectionFromSymptomatic[AgeGroup(ag)] = 0.25 - model.parameters.DeathsPerCritical[AgeGroup(ag)] = 0.3 + model.parameters.RelativeTransmissionNoSymptoms[AgeGroup(i)] = 0.67 + model.parameters.TransmissionProbabilityOnContact[AgeGroup(i)] = 0.1 + model.parameters.RecoveredPerInfectedNoSymptoms[AgeGroup(i)] = 0.2 + model.parameters.RiskOfInfectionFromSymptomatic[AgeGroup(i)] = 0.25 + model.parameters.DeathsPerCritical[AgeGroup(i)] = 0.3 # The groups have an increasing risk of severe and critical infections model.parameters.SeverePerInfectedSymptoms[AgeGroup(0)] = 0.2 @@ -128,7 +127,8 @@ def _(AgeGroup, model, np, num_age_groups): model.parameters.CriticalPerSevere[AgeGroup(2)] = 0.25 * 2 # Set contact frequency - model.parameters.ContactPatterns.cont_freq_mat[0].baseline = np.ones((num_age_groups, num_age_groups)) * 10 + model.parameters.ContactPatterns.cont_freq_mat[0].baseline = np.ones( + (num_age_groups, num_age_groups)) * 10 return @@ -159,13 +159,15 @@ def _(mo): @app.cell def _(AgeGroup, model, num_age_groups, osecir, total_population_per_region): # The population is equally distributed among the age groups - for group in range(num_age_groups): + for j in range(num_age_groups): # 1% of the population is initially infected, 0.5% Exposed and 0.5% in the pre- or asymptomatic state - model.populations[AgeGroup(group), osecir.InfectionState.Exposed] = 0.005 * total_population_per_region / num_age_groups - model.populations[AgeGroup(group), osecir.InfectionState.InfectedNoSymptoms] = 0.005 * total_population_per_region / num_age_groups + model.populations[AgeGroup(j), osecir.InfectionState.Exposed] = 0.005 * \ + total_population_per_region / num_age_groups + model.populations[AgeGroup(j), osecir.InfectionState.InfectedNoSymptoms] = 0.005 * \ + total_population_per_region / num_age_groups # The rest of the population is Susceptible model.populations.set_difference_from_group_total_AgeGroup( - (AgeGroup(group), osecir.InfectionState.Susceptible), total_population_per_region / num_age_groups) + (AgeGroup(j), osecir.InfectionState.Susceptible), total_population_per_region / num_age_groups) return @@ -180,7 +182,7 @@ def _(mo): @app.cell def _(graph, model, t0): # Add node with id 0 and copy beforehand initialized model to it - graph.add_node(id=0, model=model, t0=t0) + graph.add_node(id=0, model=model, t0=t0) return @@ -195,13 +197,14 @@ def _(mo): @app.cell def _(AgeGroup, model, num_age_groups, osecir, total_population_per_region): # The population is equally distributed among the age groups - for age in range(num_age_groups): + for k in range(num_age_groups): # No infected individuals - model.populations[AgeGroup(age), osecir.InfectionState.Exposed] = 0 - model.populations[AgeGroup(age), osecir.InfectionState.InfectedNoSymptoms] = 0 + model.populations[AgeGroup(k), osecir.InfectionState.Exposed] = 0 + model.populations[AgeGroup( + k), osecir.InfectionState.InfectedNoSymptoms] = 0 # The total population is Susceptible model.populations.set_difference_from_group_total_AgeGroup( - (AgeGroup(age), osecir.InfectionState.Susceptible), total_population_per_region / num_age_groups) + (AgeGroup(k), osecir.InfectionState.Susceptible), total_population_per_region / num_age_groups) return @@ -216,7 +219,7 @@ def _(mo): @app.cell def _(graph, model, t0): # Add node with id 0 and copy beforehand initialized model to it - graph.add_node(id=1, model=model, t0=t0) + graph.add_node(id=1, model=model, t0=t0) return @@ -231,35 +234,50 @@ def _(mo): @app.cell -def _(graph, mio, model, np, osecir): +def _(graph, mio, model, np, num_age_groups, osecir): # One coefficient per (age group x compartment) mobility_coefficients = 0.1 * np.ones(model.populations.numel()) # Dead individuals do not commute - mobility_coefficients[osecir.InfectionState.Dead] = 0 + for l in range(num_age_groups): + mobility_coefficients[l * (int(osecir.InfectionState.Dead)+1) + + int(osecir.InfectionState.Dead)] = 0 mobility_params = mio.MobilityParameters(mobility_coefficients) # Add two edges to graph graph.add_edge(0, 1, mobility_params) graph.add_edge(1, 0, mobility_params) # Individuals are exchanged every half day dt_exchange = 0.5 - return + return (dt_exchange,) @app.cell def _(mo): mo.md(r""" ## Model simulation - + We now have finished initializing the metapopulation model. The graph-based simulation is created and advanced until `tmax` via: """) return @app.cell -def _(graph, osecir, t0, tmax): - # Create graph simulation and advance until tmax - sim = osecir.MobilitySimulation(graph, t0, dt=dt_exchange) - sim.advance(tmax) +def _(dt_exchange, graph, osecir, t0, tmax): + import os + # Silence C++ log output; it can deadlock marimo's output pipe on Windows + devnull = os.open(os.devnull, os.O_WRONLY) + saved_stdout, saved_stderr = os.dup(1), os.dup(2) + os.dup2(devnull, 1) + os.dup2(devnull, 2) + try: + # Create graph simulation and advance until tmax + sim = osecir.MobilitySimulation(graph, t0, dt=dt_exchange) + sim.advance(tmax) + finally: + os.dup2(saved_stdout, 1) + os.dup2(saved_stderr, 2) + os.close(devnull) + os.close(saved_stdout) + os.close(saved_stderr) return (sim,) @@ -288,8 +306,10 @@ def _(mo): @app.cell def _(osecir, result_region0, result_region1): - result_region0_interpolated = osecir.interpolate_simulation_result(result_region0) - result_region1_interpolated = osecir.interpolate_simulation_result(result_region1) + result_region0_interpolated = osecir.interpolate_simulation_result( + result_region0) + result_region1_interpolated = osecir.interpolate_simulation_result( + result_region1) return result_region0_interpolated, result_region1_interpolated @@ -297,7 +317,7 @@ def _(osecir, result_region0, result_region1): def _(mo): mo.md(r""" ## Visualization of model output - + Finally, we can compare the trajectories of all infection states for both regions. In the following, we plot the number of `InfectedNoSymptoms` aggregated over all age groups for both regions: """) return @@ -312,8 +332,10 @@ def _(osecir, plt, result_region0_interpolated, result_region1_interpolated): # Plot the number of non-symptomatically infected for both regions fig, ax = plt.subplots() time = result_array_0[0, :] - InfectedNoSymptoms_0 = result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int(osecir.InfectionState.Dead) + 1, :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] - InfectedNoSymptoms_1 = result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int(osecir.InfectionState.Dead) + 1, :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] + InfectedNoSymptoms_0 = result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int( + osecir.InfectionState.Dead) + 1, :] + result_array_0[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] + InfectedNoSymptoms_1 = result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms), :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + int( + osecir.InfectionState.Dead) + 1, :] + result_array_1[1 + int(osecir.InfectionState.InfectedNoSymptoms) + 2 * (int(osecir.InfectionState.Dead) + 1), :] ax.plot(time, InfectedNoSymptoms_0, label='Infected No Symptoms Region 1') ax.plot(time, InfectedNoSymptoms_1, label='Infected No Symptoms Region 2') ax.set_xlabel('Time [days]')