6.4. pypfilt.obs

While all observations models should derive from the Obs base class, the Univariate class provides a simple way to define observation models in terms of Scipy probability distributions.

class pypfilt.obs.Obs

The base class of observation models, which defines the minimal set of methods that are required.

Note

The observation model constructor (__init__) must accept two positional arguments: the observation unit (string) and the observation model settings (dictionary).

log_llhd(ctx, snapshot, obs, ixs=None)

Return the log-likelihood \(\mathcal{l}(y_t \mid x_t)\) for the observation \(y_t\) and every particle \(x_t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • obs – An observation for the current time-step, \(y_t\).
  • ixs – Sampling indices for the shape parameters (optional).
expect(ctx, snapshot, ixs=None)

Return the expected observation value \(\mathbb{E}[y_t]\) for every particle \(x_t\), at one or more times \(t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • ixs – Sampling indices for the shape parameters (optional).
quantiles(ctx, snapshot, probs, ixs=None)

Return the values \(y_i\) that satisfy:

\[y_i = \inf\left\{ y : p_i \le \sum_i w_i \cdot \mathcal{L}(y_t \le y \mid x_t^i)\right\}\]
Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • probs – The probabilities \(p_i\), which must be sorted in ascending order.
  • ixs – Sampling indices for the shape parameters (optional).
llhd_in(ctx, snapshot, y0, y1, ixs=None)

Return the weighted likelihood that \(y_t \in [y_0, y_1)\):

\[\sum_i w_i \cdot \mathcal{L}(y_0 \le y_t < y_1 \mid x_t^i)\]
Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • y0 – The (inclusive) minimum fraction of cases, \(y_0\).
  • y1 – The (exclusive) maximum fraction of cases, \(y_1\).
  • ixs – Sampling indices for the shape parameters (optional).
simulate(ctx, snapshot, rng, ixs=None)

Return a random sample of \(y_t\) for each particle \(x_t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • rng – The random number generator to use.
  • ixs – Sampling indices for the shape parameters (optional).
from_file(filename, time_scale)

Load observations from a space-delimited text file with column headers defined in the first line.

Parameters:
  • filename – The file to read.
  • time_scale – The simulation time scale.
Returns:

The data table of observations.

Return type:

numpy.ndarray

Note

Use read_fields() to implement this method. See the example implementation, below.

Example:
import numpy as np
import pypfilt.io

def from_file(self, filename, time_scale):
    fields = [pypfilt.io.time_field('date'), ('value', np.float_)]
    return pypfilt.io.read_fields(time_scale, filename, fields)
row_into_obs(row)

Convert a data table row into an observation dictionary.

Parameters:row – The data table row.
simulated_obs(ctx, snapshot, rng, ixs=None)

Return a simulated observation dictionary for each particle.

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • rng – The random number generator to use.
  • ixs – Sampling indices for the shape parameters (optional).
simulated_field_types(ctx)

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

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

Note

Use pypfilt.io.time_field() for columns that will contain time values. This ensures that the time values will be converted as necessary when loading and saving tables.

Parameters:ctx – The simulation context.
obs_into_row(obs, dtype)

Convert an observation dictionary into a data table row.

Parameters:
  • obs – The observation dictionary.
  • dtype – The NumPy dtype for the data table.
class pypfilt.obs.Univariate(obs_unit)

Define observation models in terms of a univariate scipy.stats distribution.

Implement the distribution() method and all of the Obs methods will be automatically implemented.

Each observation comprises a time ('date') and a floating-point value ('value').

Parameters:obs_unit – The observation unit, a unique identifier for the observations associated with this observation model.

Note

Override the log_llhd() method to handle, e.g., incomplete observations.

Note

The observation unit is stored in the unit attribute (see the example below).

Examples:
>>> from pypfilt.obs import Univariate
>>> # Define a Gaussian observation model with a known standard deviation.
>>> class MyObsModel(Univariate):
...     def __init__(self, obs_unit, settings):
...         super().__init__(obs_unit)
...         self.sdev = settings['parameters']['sdev']
...     def distribution(self, ctx, snapshot, ixs=None):
...         # Directly observe the state variable 'x'.
...         expect = snapshot.vec['state_vec']['x']
...         # Apply the sampling indices, if provided.
...         if ixs is not None:
...             expect = expect[ixs]
...         return scipy.stats.norm(loc=expect, scale=self.sdev)
...
>>> observation_unit = 'x'
>>> settings = {'parameters': {'sdev': 0.1}}
>>> obs_model = MyObsModel(observation_unit, settings)
>>> obs_model.unit
'x'

The observation model shown in the example above can then be used in a scenario definition:

[observations.x]
model = "my_module.MyObsModel"
file = "x-observations.ssv"
parameters.sdev = 0.2
descriptor.format.sdev = "0.1f"
descriptor.name.sdev = "sdev"
distribution(ctx, snapshot, ixs=None)

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.
  • ixs – Sampling indices for the shape parameters (optional).
quantiles_tolerance()

Return the minimum interval width when calculating quantiles for a continuous random variable.

Note

The default tolerance is 0.00001. Override this method to adjust the tolerance.

log_llhd(ctx, snapshot, obs, ixs=None)

Return the log-likelihood \(\mathcal{l}(y_t \mid x_t)\) for the observation \(y_t\) and every particle \(x_t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • obs – An observation for the current time-step, \(y_t\).
  • ixs – Sampling indices for the shape parameters (optional).
expect(ctx, snapshot, ixs=None)

Return the expected observation value \(\mathbb{E}[y_t]\) for every particle \(x_t\), at one or more times \(t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • ixs – Sampling indices for the shape parameters (optional).
quantiles(ctx, snapshot, probs, ixs=None)

Return the values \(y_i\) that satisfy:

\[y_i = \inf\left\{ y : p_i \le \sum_i w_i \cdot \mathcal{L}(y_t \le y \mid x_t^i)\right\}\]
Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • probs – The probabilities \(p_i\), which must be sorted in ascending order.
  • ixs – Sampling indices for the shape parameters (optional).
llhd_in(ctx, snapshot, y0, y1, ixs)

Return the weighted likelihood that \(y_t \in [y_0, y_1)\):

\[\sum_i w_i \cdot \mathcal{L}(y_0 \le y_t < y_1 \mid x_t^i)\]
Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • y0 – The (inclusive) minimum fraction of cases, \(y_0\).
  • y1 – The (exclusive) maximum fraction of cases, \(y_1\).
  • ixs – Sampling indices for the shape parameters (optional).
simulate(ctx, snapshot, rng, ixs=None)

Return a random sample of \(y_t\) for each particle \(x_t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • rng – The random number generator to use.
  • ixs – Sampling indices for the shape parameters (optional).
from_file(filename, time_scale, time_col='date', value_col='value')

Load count data from a space-delimited text file with column headers defined in the first line.

Parameters:
  • filename – The file to read.
  • time_scale – The simulation time scale.
  • time_col – The name of the observation time column; this will be renamed to 'date'.
  • value_col – The name of the observation value column; this will be renamed to 'value'.
row_into_obs(row)

Return an observation with fields 'date', 'value', and 'unit'.

obs_into_row(obs, dtype)

Convert an observation into a (date, value) tuple.

simulated_field_types(ctx)

Return the field types for simulated observations.

simulated_obs(ctx, snapshot, rng, ixs=None)

Return a simulated observation with fields 'date' and 'value'.

pypfilt.obs.expect(ctx, snapshot, unit)

Return the expected observation value \(\mathbb{E}[y_t]\) for every every particle \(x_t\) at time \(t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • unit (str) – The observation unit.
  • hist – The particle history matrix.
  • hist_ix – The index of time \(t\) in the history matrix.
pypfilt.obs.log_llhd(ctx, snapshot, obs)

Return the log-likelihood \(\mathcal{l}(y_t \mid x_t)\) for the observation \(y_t\) and every particle \(x_t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • obs – The list of observations for the current time-step.
pypfilt.obs.bisect_cdf(probs, cdf_fn, bisect_fn, y_lower, y_upper)

Use a bisection method to estimate the values \(y_i\) that satisfy:

\[y_i = \inf\left\{ y : p_i \le \sum_i w_i \cdot \mathcal{L}(y_t \le y \mid x_t^i)\right\}\]
Parameters:
  • probs – The probabilities \(p_i\), which must be sorted in ascending order.
  • cdf_fn – The CDF function \(f(y) = \sum_i w_i \cdot \mathcal{L}(y_t \le y \mid x_t^i)\).
  • bisect_fn – The bisection function f(a, b) that either returns the midpoint of the interval \([a, b]\) or None if the search should stop (e.g., because a tolerance has been reached).
  • y_lower – A lower bound for all \(y_i\).
  • y_upper – An upper bound for all \(y_i\).
pypfilt.obs.simulate(ctx, snapshot, unit, rng)

Return a random sample of \(y_t\) for each particle \(x_t\).

Parameters:
  • ctx – The simulation context.
  • snapshot (Snapshot) – The current particle states.
  • unit – The observation unit.
  • rng – The random number generator to use.