Results and plotting¶

This tutorial covers how to inspect, export, and plot results from an STIsim simulation.

Setup¶

Let's run a simple gonorrhea simulation to work with:

In [1]:

import stisim as sti
import matplotlib.pyplot as plt

sim = sti.Sim(diseases='ng', n_agents=2000, start=2010, stop=2030)
sim.run(verbose=0)

import stisim as sti import matplotlib.pyplot as plt sim = sti.Sim(diseases='ng', n_agents=2000, start=2010, stop=2030) sim.run(verbose=0)

Initializing sim with 2000 agents

Out[1]:

Sim(n=2000; 2010—2030; networks=structuredsexual, maternalnet; diseases=ng)

Browsing available results¶

After running a sim, all results are stored in sim.results. Each disease has its own results object. You can list all available keys with sim.results.ng.keys(), or access individual results directly. STIsim stores results at multiple levels of detail: overall, by sex, and by age group.

In [2]:

# Access individual results
prev = sim.results.ng.prevalence
print(f'{prev.label}: min={prev.values.min():.3f}, mean={prev.values.mean():.3f}, max={prev.values.max():.3f}')

new_inf = sim.results.ng.new_infections
print(f'{new_inf.label}: min={new_inf.values.min():.0f}, mean={new_inf.values.mean():.1f}, max={new_inf.values.max():.0f}')

# There are many more -- including by sex and age group
print(f'\nTotal result keys: {len(sim.results.ng.keys())}')

# Access individual results prev = sim.results.ng.prevalence print(f'{prev.label}: min={prev.values.min():.3f}, mean={prev.values.mean():.3f}, max={prev.values.max():.3f}') new_inf = sim.results.ng.new_infections print(f'{new_inf.label}: min={new_inf.values.min():.0f}, mean={new_inf.values.mean():.1f}, max={new_inf.values.max():.0f}') # There are many more -- including by sex and age group print(f'\nTotal result keys: {len(sim.results.ng.keys())}')

Prevalence: min=0.007, mean=0.073, max=0.095
New infections: min=1, mean=40.3, max=61

Total result keys: 200

Each result is an ss.Result object. The raw data is in .values (a NumPy array) and the time axis is in .timevec:

In [3]:

print(f'Values (first 5): {prev.values[:5]}')
print(f'Timevec:          {prev.timevec[0]} to {prev.timevec[-1]} ({len(prev.values)} timesteps)')

print(f'Values (first 5): {prev.values[:5]}') print(f'Timevec: {prev.timevec[0]} to {prev.timevec[-1]} ({len(prev.values)} timesteps)')

Values (first 5): [0.00675219 0.00740741 0.0080429  0.0080429  0.00937082]
Timevec:          2010.01.01 to 2030.01.01 (241 timesteps)

Automatic plotting¶

The simplest way to visualize results is sim.plot(). By default, it shows only the high-level results (not every age/sex subgroup):

In [4]:

sim.plot()

sim.plot()

Figure(933.333x700)

Plotting specific results¶

To plot specific results, pass the key argument. You can use a single key or a list:

In [5]:

# Plot just prevalence and new infections
sim.plot(key=['ng.prevalence', 'ng.new_infections'])

# Plot just prevalence and new infections sim.plot(key=['ng.prevalence', 'ng.new_infections'])

Figure(800x600)

You can also plot results that are hidden from the default view. For example, sex-stratified prevalence:

In [6]:

sim.plot(key=['ng.prevalence_f', 'ng.prevalence_m'])

sim.plot(key=['ng.prevalence_f', 'ng.prevalence_m'])

Figure(800x600)

Exporting results¶

As a DataFrame¶

Export all results to a Pandas DataFrame with sim.to_df():

In [7]:

df = sim.to_df()
print(f'Shape: {df.shape}')
df.head()

df = sim.to_df() print(f'Shape: {df.shape}') df.head()

Shape: (241, 207)

Out[7]:

	timevec	structuredsexual_n_edges	maternalnet_n_edges	ng_n_susceptible	ng_n_infected	ng_n_exposed	ng_n_asymptomatic	ng_n_symptomatic	ng_n_pid	ng_n_seeking_care	...	ng_new_symptomatic_m_50_65	ng_symp_prevalence_m_50_65	ng_new_symptomatic_m_65_100	ng_symp_prevalence_m_65_100	ng_rel_treat	n_alive	n_female	new_deaths	new_emigrants	cum_deaths
0	2010-01-01	353.0	0.0	1984.0	15.0	1.0	16.0	0.0	0.0	0.0	...	0.0	0.000000	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
1	2010-01-31	435.0	0.0	1983.0	16.0	1.0	9.0	8.0	0.0	4.0	...	1.0	0.005988	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
2	2010-03-03	497.0	0.0	1983.0	17.0	0.0	9.0	8.0	0.0	0.0	...	0.0	0.005917	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
3	2010-04-02	531.0	0.0	1981.0	17.0	2.0	11.0	8.0	0.0	1.0	...	0.0	0.005917	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
4	2010-05-03	568.0	0.0	1981.0	19.0	0.0	9.0	10.0	0.0	2.0	...	1.0	0.011561	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0

5 rows × 207 columns

You can also export a single result:

In [8]:

prev_df = sim.results.ng.prevalence.to_df()
prev_df.head()

prev_df = sim.results.ng.prevalence.to_df() prev_df.head()

Out[8]:

	timevec	value
0	2010-01-01	0.006752
1	2010-01-31	0.007407
2	2010-03-03	0.008043
3	2010-04-02	0.008043
4	2010-05-03	0.009371

Resampling to annual results¶

STIsim runs with monthly timesteps by default, but you'll often want annual results for reporting or comparison with data. Use the resample() method on any result:

In [9]:

# Monthly new infections (raw)
monthly = sim.results.ng.new_infections
print(f'Monthly: {len(monthly.values)} timesteps')

# Resample to annual (sums monthly counts)
annual = monthly.resample('year')
print(f'Annual:  {len(annual)} years')
print(annual)

# Monthly new infections (raw) monthly = sim.results.ng.new_infections print(f'Monthly: {len(monthly.values)} timesteps') # Resample to annual (sums monthly counts) annual = monthly.resample('year') print(f'Annual: {len(annual)} years') print(annual)

Monthly: 241 timesteps
Annual:  21 years
2010-12-31     80.0
2011-12-31    360.0
2012-12-31    453.0
2013-12-31    480.0
2014-12-31    510.0
2015-12-31    470.0
2016-12-31    482.0
2017-12-31    446.0
2018-12-31    490.0
2019-12-31    523.0
2020-12-31    531.0
2021-12-31    515.0
2022-12-31    498.0
2023-12-31    469.0
2024-12-31    572.0
2025-12-31    559.0
2026-12-31    564.0
2027-12-31    528.0
2028-12-31    574.0
2029-12-31    569.0
2030-12-31     47.0
Freq: YE-DEC, Name: value, dtype: float64

The resample() method automatically chooses the right aggregation: sum for counts (like new_infections) and mean for rates/proportions (like prevalence). You can also resample to a DataFrame:

In [10]:

annual_df = sim.results.ng.new_infections.resample('year', output_form='dataframe')
annual_df

annual_df = sim.results.ng.new_infections.resample('year', output_form='dataframe') annual_df

Out[10]:

	value
2010-12-31	80.0
2011-12-31	360.0
2012-12-31	453.0
2013-12-31	480.0
2014-12-31	510.0
2015-12-31	470.0
2016-12-31	482.0
2017-12-31	446.0
2018-12-31	490.0
2019-12-31	523.0
2020-12-31	531.0
2021-12-31	515.0
2022-12-31	498.0
2023-12-31	469.0
2024-12-31	572.0
2025-12-31	559.0
2026-12-31	564.0
2027-12-31	528.0
2028-12-31	574.0
2029-12-31	569.0
2030-12-31	47.0

You can also resample the entire sim's results at once via to_df():

In [11]:

annual_all = sim.to_df(resample='year', use_years=True)
annual_all.head()

annual_all = sim.to_df(resample='year', use_years=True) annual_all.head()

Out[11]:

	timevec	structuredsexual_n_edges	maternalnet_n_edges	ng_n_susceptible	ng_n_infected	ng_n_exposed	ng_n_asymptomatic	ng_n_symptomatic	ng_n_pid	ng_n_seeking_care	...	ng_new_symptomatic_m_50_65	ng_symp_prevalence_m_50_65	ng_new_symptomatic_m_65_100	ng_symp_prevalence_m_65_100	ng_rel_treat	n_alive	n_female	new_deaths	new_emigrants	cum_deaths
0	2010	557.250000	0.0	1980.416667	16.416667	3.166667	11.250000	8.333333	0.0	1.833333	...	3.0	0.004332	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
1	2011	658.250000	0.0	1913.250000	71.000000	15.750000	47.666667	39.083333	0.0	6.166667	...	12.0	0.021894	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
2	2012	685.500000	0.0	1863.750000	117.750000	18.500000	73.666667	62.583333	0.0	7.833333	...	22.0	0.042883	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
3	2013	714.833333	0.0	1846.583333	133.666667	19.750000	89.166667	64.250000	0.0	7.250000	...	18.0	0.038823	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0
4	2014	737.833333	0.0	1840.916667	137.833333	21.250000	83.666667	75.416667	0.0	9.333333	...	19.0	0.046988	0.0	0.0	0.0	2000.0	0.0	0.0	0.0	0.0

5 rows × 207 columns

Custom plots with Matplotlib¶

For publication-quality figures or custom layouts, use the result values directly with Matplotlib:

In [12]:

fig, axes = plt.subplots(1, 2, figsize=(10, 4))

# Left: prevalence by sex
axes[0].plot(sim.timevec, sim.results.ng.prevalence_f, label='Female')
axes[0].plot(sim.timevec, sim.results.ng.prevalence_m, label='Male')
axes[0].set_xlabel('Year')
axes[0].set_ylabel('Prevalence')
axes[0].set_title('NG prevalence by sex')
axes[0].legend()

# Right: annual new infections (resampled)
annual_inf = sim.results.ng.new_infections.resample('year', output_form='dataframe', use_years=True)
axes[1].bar(annual_inf.index, annual_inf['value'])
axes[1].set_xlabel('Year')
axes[1].set_ylabel('New infections')
axes[1].set_title('Annual NG infections')

plt.tight_layout()
fig

fig, axes = plt.subplots(1, 2, figsize=(10, 4)) # Left: prevalence by sex axes[0].plot(sim.timevec, sim.results.ng.prevalence_f, label='Female') axes[0].plot(sim.timevec, sim.results.ng.prevalence_m, label='Male') axes[0].set_xlabel('Year') axes[0].set_ylabel('Prevalence') axes[0].set_title('NG prevalence by sex') axes[0].legend() # Right: annual new infections (resampled) annual_inf = sim.results.ng.new_infections.resample('year', output_form='dataframe', use_years=True) axes[1].bar(annual_inf.index, annual_inf['value']) axes[1].set_xlabel('Year') axes[1].set_ylabel('New infections') axes[1].set_title('Annual NG infections') plt.tight_layout() fig

Out[12]:

Summary¶

Task	Method
List all result keys	sim.results.ng.keys()
Access a result's values	sim.results.ng.prevalence.values
Auto-plot high-level results	sim.plot()
Plot specific results	sim.plot(key=['ng.prevalence', 'ng.new_infections'])
Plot hidden results	sim.plot(key='ng.prevalence_f')
Export all to DataFrame	sim.to_df()
Export one result	sim.results.ng.prevalence.to_df()
Resample to annual	result.resample('year')
Annual DataFrame	sim.to_df(resample='year', use_years=True)

Exercises¶

1. Run a sim with both 'ng' and 'ct'. Use sim.plot(key=...) to compare the prevalence of the two diseases on the same figure.
2. Export the annual new infections for gonorrhea to a DataFrame and save it to a CSV file.
3. Make a Matplotlib figure that plots monthly symptomatic prevalence (ng.symp_prevalence) alongside total prevalence (ng.prevalence). What fraction of infections are symptomatic?