Abstract Model Classes

Introduction

A model in QInfer is a class that describes the probabilities of observing data, given a particular experiment and given a particular set of model parameters. The observation probabilities may be given implicitly or explicitly, in that the class may only allow for sampling observations, rather than finding the a distribution explicitly. In the former case, a model is represented by a subclass of Simulatable, while in the latter, the model is represented by a subclass of Model.

Simulatable - Base Class for Implicit (Simulatable) Models

Class Reference

class qinfer.Simulatable[source]

Bases: object

n_modelparams

Returns the number of real model parameters admitted by this model.

This property is assumed by inference engines to be constant for the lifetime of a Simulatable instance.

expparams_dtype

Returns the dtype of an experiment parameter array. For a model with single-parameter control, this will likely be a scalar dtype, such as "float64". More generally, this can be an example of a record type, such as [('time', 'float64'), ('axis', 'uint8')].

This property is assumed by inference engines to be constant for the lifetime of a Model instance.

is_n_outcomes_constant

Returns True if and only if the number of outcomes for each experiment is independent of the experiment being performed.

This property is assumed by inference engines to be constant for the lifetime of a Simulatable instance.

model_chain

Returns a tuple of models upon which this model is based, such that properties and methods of underlying models for models that decorate other models can be accessed. For a standalone model, this is always the empty tuple.

base_model

Returns the most basic model that this model depends on. For standalone models, this property satisfies model.base_model is model.

underlying_model

Returns the model that this model is based on (decorates) if such a model exists, or None if this model is independent.

sim_count
Q

Returns the diagonal of the scale matrix \(\matr{Q}\) that relates the scales of each of the model parameters. In particular, the quadratic loss for this Simulatable is defined as:

\[L_{\matr{Q}}(\vec{x}, \hat{\vec{x}}) = (\vec{x} - \hat{\vec{x}})^\T \matr{Q} (\vec{x} - \hat{\vec{x}})\]

If a subclass does not explicitly define the scale matrix, it is taken to be the identity matrix of appropriate dimension.

Returns:The diagonal elements of \(\matr{Q}\).
Return type:ndarray of shape (n_modelparams, ).
modelparam_names

Returns the names of the various model parameters admitted by this model, formatted as LaTeX strings.

n_outcomes(expparams)[source]

Returns an array of dtype uint describing the number of outcomes for each experiment specified by expparams.

Parameters:expparams (numpy.ndarray) – Array of experimental parameters. This array must be of dtype agreeing with the expparams_dtype property.
are_models_valid(modelparams)[source]

Given a shape (n_models, n_modelparams) array of model parameters, returns a boolean array of shape (n_models) specifying whether each set of model parameters represents is valid under this model.

simulate_experiment(modelparams, expparams, repeat=1)[source]
clear_cache()[source]

Tells the model to clear any internal caches used in computing likelihoods and drawing samples. Calling this method should not cause any different results, but should only affect performance.

experiment_cost(expparams)[source]

Given an array of experimental parameters, returns the cost associated with performing each experiment. By default, this cost is constant (one) for every experiment.

Parameters:expparams (ndarray of dtype given by expparams_dtype) – An array of experimental parameters for which the cost is to be evaluated.
Returns:An array of costs corresponding to the specified experiments.
Return type:ndarray of dtype float and of the same shape as expparams.
distance(a, b)[source]

Gives the distance between two model parameter vectors \(\vec{a}\) and \(\vec{b}\). By default, this is the vector 1-norm of the difference \(\mathbf{Q} (\vec{a} - \vec{b})\) rescaled by Q.

Parameters:
  • a (np.ndarray) – Array of model parameter vectors having shape (n_models, n_modelparams).
  • b (np.ndarray) – Array of model parameters to compare to, having the same shape as a.
Returns:

An array d of distances d[i] between a[i, :] and b[i, :].

update_timestep(modelparams, expparams)[source]

Returns a set of model parameter vectors that is the update of an input set of model parameter vectors, such that the new models are conditioned on a particular experiment having been performed. By default, this is the trivial function \(\vec{x}(t_{k+1}) = \vec{x}(t_k)\).

Parameters:
  • modelparams (np.ndarray) – Set of model parameter vectors to be updated.
  • expparams (np.ndarray) – An experiment parameter array describing the experiment that was just performed.
Return np.ndarray:
 

Array of shape (n_models, n_modelparams, n_experiments) describing the update of each model according to each experiment.

canonicalize(modelparams)[source]

Returns a canonical set of model parameters corresponding to a given possibly non-canonical set. This is used for models in which there exist model parameters \(\vec{x}_i\) and :math:vec{x}_j such that

\[\Pr(d | \vec{x}_i; \vec{e}) = \Pr(d | \vec{x}_j; \vec{e})\]

for all outcomes \(d\) and experiments \(\vec{e}\). For models admitting such an ambiguity, this method should then be overridden to return a consistent choice out of such vectors, hence avoiding supurious model degeneracies.

Note that, by default, SMCUpdater will not call this method.

Model - Base Class for Explicit (Likelihood) Models

If a model supports explicit calculation of the likelihood function, then this is represented by subclassing from Model. The likelihood function provided by a subclass is then used to implement Simulatable.simulate_experiment(), so that the primary method to be defined by a Model subclass is Model.likelihood().

Class Reference

class qinfer.Model[source]

Bases: qinfer.abstract_model.Simulatable

call_count
likelihood(outcomes, modelparams, expparams)[source]
is_model_valid(modelparams)[source]

Returns True if and only if the model parameters given are valid for this model.

simulate_experiment(modelparams, expparams, repeat=1)[source]
static pr0_to_likelihood_array(outcomes, pr0)[source]

Assuming a two-outcome measurement with probabilities given by the array pr0, returns an array of the form expected to be returned by likelihood method.

Parameters:
  • outcomes (numpy.ndarray) – Array of integers indexing outcomes.
  • pr0 (numpy.ndarray) – Array of shape (n_models, n_experiments) describing the probability of obtaining outcome 0 from each set of model parameters and experiment parameters.

DifferentiableModel - Base Class for Explicit Models with Differentiable Likelihoods

Class Reference

class qinfer.DifferentiableModel[source]

Bases: qinfer.abstract_model.Model

score(outcomes, modelparams, expparams, return_L=False)[source]

Returns the score of this likelihood function, defined as:

\[q(d, \vec{x}; \vec{e}) = \vec{\nabla}_{\vec{x}} \log \Pr(d | \vec{x}; \vec{e}).\]

Calls are represented as a four-index tensor score[idx_modelparam, idx_outcome, idx_model, idx_experiment]. The left-most index may be suppressed for single-parameter models.

If return_L is True, both q and the likelihood L are returned as q, L.

fisher_information(modelparams, expparams)[source]

Returns the covariance of the score taken over possible outcomes, known as the Fisher information.

The result is represented as the four-index tensor fisher[idx_modelparam_i, idx_modelparam_j, idx_model, idx_experiment], which gives the Fisher information matrix for each model vector and each experiment vector.

Note

The default implementation of this method calls score() for each possible outcome, which can be quite slow. If possible, overriding this method can give significant speed advantages.