pypfilt.examples#

The pypfilt.examples module includes a chaotic dynamical system, a predator-prey system, a Gaussian random walk, and several epidemic SIR models.

Dynamical systems models#

class pypfilt.examples.lorenz.Lorenz63#

The Lorenz-63 system:

\[\begin{split}\frac{dx}{dt} &= \sigma (y - x) \\ \frac{dy}{dt} &= x (\rho - z) - y \\ \frac{dz}{dt} &= xy - \beta z \\ x_t &= [\sigma, \rho, \beta, x, y, z]^T\end{split}\]

This system exhibits chaotic behaviour in the neighbourhood of \(\sigma = 10\), \(\rho = 28\), \(\beta = \frac{8}{3}\).

field_types(ctx)#

Define the state vector \([\sigma, \rho, \beta, x, y, z]^T\).

d_dt(time, xt, ctx, is_forecast)#

The right-hand side of the ODE system.

Parameters:
  • time – The current time.

  • xt – The particle state vectors.

  • ctx – The simulation context.

  • is_forecast – True if this is a forecasting simulation.

can_smooth()#

Indicate which state vector fields can be smoothed.

class pypfilt.examples.lorenz.ObsLorenz63(obs_unit, settings)#

An observation model for the Lorenz-63 system:

\[y_t \sim N(\mu = x_t, \sigma = 1.5)\]

The observation unit must be the name of a field in the state vector. For example, the observation unit must be "y" to observe \(y(t)\):

[observations.y]
model = "pypfilt.examples.lorenz.ObsLorenz63"
distribution(ctx, snapshot)#

Return a frozen scipy.stats distribution that defines the observation model for each particle.

Parameters:
  • ctx – The simulation context.

  • snapshot (Snapshot) – The current particle states.

pypfilt.examples.lorenz.save_lorenz63_scenario_files()#

Save all of the example Lorenz-63 scenario files in the working directory.

This function creates (or overwrites) the following files:

  • 'lorenz63_simulate.toml': a scenario file for simulating observations.

  • 'lorenz63_forecast.toml': a scenario file for running forecasts, using the simulated observations.

  • 'lorenz63_forecast_regularised.toml': a scenario file for running forecasts with post-regularisation, using the simulated observations.

  • 'lorenz63_all.toml': a scenario file that defines all of the above scenarios.

pypfilt.examples.lorenz.lorenz63_simulate_toml()#

A scenario for the Lorenz63 model, which can be used for simulating observations.

Returns:

The scenario definition, represented as a TOML string.

Return type:

str

pypfilt.examples.lorenz.lorenz63_forecast_toml()#

A scenario for the Lorenz63 model, which can be used for forecasting.

Returns:

The scenario definition, represented as a TOML string.

Return type:

str

pypfilt.examples.lorenz.lorenz63_forecast_regularised_toml()#

A scenario for the Lorenz63 model, which can be used for forecasting, and enabled post-regularisation.

Returns:

The scenario definition, represented as a TOML string.

Return type:

str

pypfilt.examples.lorenz.lorenz63_all_scenarios_toml()#

All example scenarios for the Lorenz63 model.

Returns:

The scenario definitions, represented as a TOML string.

Return type:

str

Predator-prey system#

Models#

class pypfilt.examples.predation.LotkaVolterra#

An implementation of the (continuous) Lotka-Volterra equations.

field_types(ctx)#

Return a list of (field_name, field_dtype, field_shape) tuples that define the state vector.

The third element, field_shape, is optional and contains the shape of this field if it forms an array of type field_dtype.

These tuples must be in the same order as the state vector itself.

Parameters:

ctx – The simulation context.

d_dt(time, xt, ctx, is_forecast)#

Calculate the derivatives of x(t) and y(t).

can_smooth()#

Return the set of field names in the state vector that can be smoothed by the post-regularised particle filter (see post_regularise()).

Note

Models should only implement this method if the state vector contains parameters that can be smoothed.

class pypfilt.examples.predation.ObsModel(obs_unit, settings)#
distribution(ctx, snapshot)#

Return a frozen scipy.stats distribution that defines the observation model for each particle.

Parameters:
  • ctx – The simulation context.

  • snapshot (Snapshot) – The current particle states.

Example files#

pypfilt.examples.predation.write_example_files()#

Save the following example files to the working directory:

  • The forecast scenario file “predation.toml”;

  • The observations file “predation-counts-x.ssv”;

  • The observations file “predation-counts-y.ssv”;

  • The forecast scenario file “predation-datetime.toml”;

  • The observations file “predation-counts-x-datetime.ssv”; and

  • The observations file “predation-counts-y-datetime.ssv”;

pypfilt.examples.predation.example_toml_data()#

Return the contents of the example file “predation.toml”.

pypfilt.examples.predation.example_obs_x_data()#

Return the contents of the example file “predation-counts-x.ssv”.

pypfilt.examples.predation.example_obs_y_data()#

Return the contents of the example file “predation-counts-y.ssv”.

pypfilt.examples.predation.example_toml_datetime_data()#

Return the contents of the example file “predation-datetime.toml”.

pypfilt.examples.predation.example_obs_x_datetime_data()#

Return the contents of the example file “predation-counts-x-datetime.ssv”.

pypfilt.examples.predation.example_obs_y_datetime_data()#

Return the contents of the example file “predation-counts-y-datetime.ssv”.

Generating forecasts#

pypfilt.examples.predation.forecast(data_file)#

Run a suite of forecasts against generated observations, using a scalar time scale.

Parameters:

date_file – The name of the output HDF5 file.

pypfilt.examples.predation.plot(data_file, png=True, pdf=True)#

Save the plots produced by plot_params() and plot_forecasts().

This will save the plots to files whose names begin with “predation_params” and “predation_forecasts”.

Parameters:
  • png – Whether to save plots as PNG files.

  • pdf – Whether to save plots as PDF files.

pypfilt.examples.predation.plot_params(param_cints, pdf_file=None, png_file=None)#

Plot the parameter posteriors over the estimation run.

pypfilt.examples.predation.plot_forecasts(state_cints, x_obs, y_obs, pdf_file=None, png_file=None)#

Plot the population predictions at each forecasting date.

Other functions#

pypfilt.examples.predation.predation_instance(toml_file)#

Return an instance of the predation scenario from the specified TOML file.

pypfilt.examples.predation.predation_scalar_instance()#

Return an instance of the predation scenario using a scalar time scale.

pypfilt.examples.predation.predation_datetime_instance()#

Return an instance of the predation scenario using a datetime time scale.

pypfilt.examples.predation.apply_ground_truth_prior(instance)#

Define the predation model prior distribution for fixed ground truth.

pypfilt.examples.predation.save_scalar_observations(obs_tables, x_obs_file, y_obs_file)#

Save simulated observations to disk.

Gaussian random walk#

Models#

class pypfilt.examples.simple.GaussianWalk#

A Gaussian random walk.

\[\begin{split}x_t &= x_{t-1} + X_t \\ X_t &\sim N(\mu = 0, \sigma = 1)\end{split}\]

The initial values \(x_0\) are defined by the prior distribution for "x":

[prior]
x = { name = "uniform", args.loc = 10.0, args.scale = 10.0 }
field_types(ctx)#

Return a list of (field_name, field_dtype, field_shape) tuples that define the state vector.

The third element, field_shape, is optional and contains the shape of this field if it forms an array of type field_dtype.

These tuples must be in the same order as the state vector itself.

Parameters:

ctx – The simulation context.

update(ctx, time_step, is_fs, prev, curr)#

Perform a single time-step.

class pypfilt.examples.simple.GaussianObs(obs_unit, settings)#

A Gaussian observation model for the GaussianWalk model.

\[\mathcal{L}(y_t \mid x_t) \sim N(\mu = x_t, \sigma = s)\]

The observation model has one parameter: the standard deviation \(s\), whose value is defined by the "parameters.sdev" setting:

[observations.x]
model = "pypfilt.examples.simple.GaussianObs"
parameters.sdev = 0.2
distribution(ctx, snapshot)#

Return a frozen scipy.stats distribution that defines the observation model for each particle.

Parameters:
  • ctx – The simulation context.

  • snapshot (Snapshot) – The current particle states.

Support functions#

pypfilt.examples.simple.gaussian_walk_toml_data()#

Return the contents of the example file “gaussian_walk.toml”.

pypfilt.examples.simple.gaussian_walk_instance()#

Return an instance of the simple example scenario.

Epidemic SIR models#

Models#

class pypfilt.examples.sir.SirCtmc#

A continuous-time Markov chain implementation of the SIR model.

The model settings must include the following keys:

  • population_size: The number of individuals in the population.

field_types(ctx)#

Define the state vector structure.

can_smooth()#

The fields that can be smoothed by the post-regularisation filter.

init(ctx, vec)#

Initialise the state vectors.

update(ctx, time_step, is_forecast, prev, curr)#

Update the state vectors to account for all events that occur up to, and including, time.

active_particles(vec, stop_time)#

Return a Boolean array that identifies the particles whose most recent event occurred no later than stop_time.

select_next_event(ctx, vec, stop_time)#

Calculate the next event time and event type for each active particle.

class pypfilt.examples.sir.SirDtmc#

A discrete-time Markov chain implementation of the SIR model.

The model settings must include the following keys:

  • population_size: The number of individuals in the population.

field_types(ctx)#

Define the state vector structure.

can_smooth()#

The fields that can be smoothed by the post-regularisation filter.

init(ctx, vec)#

Initialise the state vectors.

update(ctx, time_step, is_forecast, prev, curr)#

Update the state vectors.

class pypfilt.examples.sir.SirOdeEuler#

An ordinary differential equation implementation of the SIR model, which uses the forward Euler method.

The model settings must include the following keys:

  • population_size: The number of individuals in the population.

field_types(ctx)#

Define the state vector structure.

can_smooth()#

The fields that can be smoothed by the post-regularisation filter.

init(ctx, vec)#

Initialise the state vectors.

update(ctx, time_step, is_forecast, prev, curr)#

Update the state vectors.

class pypfilt.examples.sir.SirOdeRk#

An ordinary differential equation implementation of the SIR model, which uses the explicit Runge-Kutta method of order 5(4).

The model settings must include the following keys:

  • population_size: The number of individuals in the population.

field_types(ctx)#

Define the state vector structure.

can_smooth()#

The fields that can be smoothed by the post-regularisation filter.

init(ctx, vec)#

Initialise the state vectors.

d_dt(time, xt, ctx, is_forecast)#

The right-hand side of the system.

class pypfilt.examples.sir.SirSde#

A stochastic differential equation implementation of the SIR model.

The model settings must include the following keys:

  • population_size: The number of individuals in the population.

field_types(ctx)#

Define the state vector structure.

can_smooth()#

The fields that can be smoothed by the post-regularisation filter.

init(ctx, vec)#

Initialise the state vectors.

update(ctx, time_step, is_forecast, prev, curr)#

Update the state vectors.

class pypfilt.examples.sir.SirObs(obs_unit, settings)#

A binomial observation model for the example SIR models.

\[ \begin{align}\begin{aligned}\mathcal{L}(y_t \mid x_t) &\sim B(n, p)\\n &= S(t-\Delta) - S(t)\end{aligned}\end{align} \]
Parameters:
  • obs_unit – A descriptive name for the data.

  • settings – The observation model settings dictionary.

The settings dictionary should contain the following keys:

  • observation_period: The observation period \(\Delta\).

For example, for daily observations that capture 80% of new infections:

[observations.cases]
model = "pypfilt.examples.sir.SirObs"
observation_period = 1
parameters.p = 0.8
new_infections(ctx, snapshot)#

Return the number of new infections \(S(t-\Delta) - S(t)\) that occurred during the observation period \(\Delta\) for each particle.

distribution(ctx, snapshot)#

Return the observation distribution for each particle.

Example files#

pypfilt.examples.sir.sir_toml_data()#

Return the contents of the example file “sir.toml”.

Epidemic SIRV models#

Models#

class pypfilt.examples.sirv.SirvOde#

An ordinary differential equation implementation of the SIRV model, using the forward Euler method, which supports two types of interventions:

  • Reducing R0 by some fraction (R0_reduction) once daily incidence reaches some threshold (R0_inc_threshold); and

  • Moving individuals from S to V at some rate (Vaccine_rate), beginning at some time (Vaccine_start).

The model settings must include the following keys:

  • population_size: The number of individuals in the population.

field_types(ctx)#

Define the state vector structure.

can_smooth()#

The fields that can be smoothed by the post-regularisation filter.

init(ctx, vec)#

Initialise the state vectors.

update(ctx, time_step, is_forecast, prev, curr)#

Update the state vectors.

class pypfilt.examples.sirv.Incidence(obs_unit, settings)#

A perfect observation model for incidence in the SIRV model.

Parameters:
  • obs_unit – A descriptive name for the data.

  • settings – The observation model settings dictionary.

The settings dictionary should contain the following keys:

  • observation_period: The observation period \(\Delta\).

For example:

[observations.cases]
model = "pypfilt.examples.sir.SirvCumIncidence"
observation_period = 1
new_infections(ctx, snapshot)#

Return the number of new infections \(S(t-\Delta) - S(t)\) that occurred during the observation period \(\Delta\) for each particle.

distribution(ctx, snapshot)#

Return the observation distribution for each particle.

class pypfilt.examples.sirv.FinalSize#

Calculate the final epidemic size for each particle.

field_types(ctx, obs_list, name)#

Return the column names and data types, represented as a list of (name, data type) tuples. See the NumPy documentation for details.

Note

To ensure that time and string values are handled appropriately when loading and saving tables, use pypfilt.io.time_field() to define time columns and pypfilt.io.string_field() to define string columns. For example:

fields = [time_field('fs_time'), time_field('time'),
          string_field('name), ('value', float)]
Parameters:
  • ctx – The simulation context.

  • obs_list – A list of all observations.

  • name – The table’s name.

n_rows(ctx, forecasting)#

Return the number of rows required for a single simulation.

Parameters:
  • ctx – The simulation context.

  • forecastingTrue if this is a forecasting simulation, otherwise False.

add_rows(ctx, fs_time, window, insert_fn)#

Record rows of summary statistics for some portion of a simulation.

Parameters:
  • ctx – The simulation context.

  • fs_time – The forecasting time; if this is not a forecasting simulation, this is the time at which the simulation ends.

  • window – A list of Snapshot instances that capture the particle states at each summary time in the simulation window.

  • insert_fn – A function that inserts one or more rows into the underlying data table; see the examples below.

The row insertion function can be used as follows:

# Insert a single row, represented as a tuple.
insert_fn((x, y, z))
# Insert multiple rows, represented as a list of tuples.
insert_fn([(x0, y0, z0), (x1, y1, z1)], n=2)
finished(ctx, fs_time, window, insert_fn)#

Record rows of summary statistics at the end of a simulation.

The parameters are as per add_rows().

Derived classes should only implement this method if rows must be recorded by this method; the provided method does nothing.

Example files#

pypfilt.examples.sirv.sirv_toml_data()#

Return the contents of the example file “sirv.toml”.