ODEPrior¶

class
probnum.diffeq.
ODEPrior
(driftmat: numpy.ndarray, dispmat: numpy.ndarray, ordint: int, spatialdim: int, precond_step: float = 1.0)[source]¶ Bases:
probnum.filtsmooth.statespace.continuous.linearsdemodel.LTISDEModel
Prior dynamic model for ODE filtering and smoothing.
 An ODE prior is an continuous LTI state space model with attributes:
 order of integration \(q\)
 spatial dimension of the underlying ODE
 projection to \(X_i(t)\) (the \((i1)\)th derivative estimate)
 A preconditioner \(P\) (see below)
Instead of the LTI SDE
\[d X(t) = [F X(t) + u] dt + L dB(t)\]the prior for the ODE Dynamics is given by
\[dX(t) = P F P^{1} X(t) dt + P L dB(t)\]where \(P\) is a preconditioner matrix ensuring stability of the iterations. Note that ODE priors do not have a drift term. By default, we choose \(P\) to be the matrix that maps to filtering iteration to the Nordsieck vector,
\[P = \text{diag }(h^{q}, h^{q+1}, ..., 1).\]Here, \(h\) is some expected average step size. Note that we ignored the factorials in this matrix. Our setting makes it easy to recover “no preconditioning” by choosing \(h=1\).
 If no expected step size is available we choose \(h=1.0\). This recovers \(P=I_{d(q+1)}\), hence no preconditioning.
 For fixed step size algorithms this quantity \(h\) is easy to choose
 For adaptive steps it is a bit more involved.
Since it doesn’t have to be exact, any more or less appropriate choice will do well. The main effect of this preconditioning is that the predictive covariances inside each filter iteration are wellconditioned: for IBM(\(q\)) priors, the condition number of the predictive covariances only depends on order of integration \(q\) and not on the step size anymore. Nb: this only holds if all required derivatives of the RHS vector field of the ODE are specified: None for IBM(1), Jacobian of \(f\) for IBM(2), Hessian of \(f\) for IBM(3). If this is not the case the preconditioner still helps but is not as powerful anymore.
Without preconditioning they can be numerically singular for small steps and higher order methods which especially makes smoothing algorithms unstable.
Another advantage of this preconditioning is that the smallest value that appears inside the algorithm is \(h^{q}\) (with preconditioning) instead of \(h^{2q+1}\) (without preconditioning).
The matrices \(F, u, L\) are the usual matrices for IBM(\(q\)), IOUP(\(q\)) or Matern(\(q+1/2\)) processes. As always, \(B(t)\) is sdimensional Brownian motion with unit diffusion matrix \(Q\).
Attributes Summary
diffusionmatrix
Evaluates Q. dimension
Spatial dimension (utility attribute). dispersionmatrix
driftmatrix
force
inverse_preconditioner
Convenience property to return the readilycomputed inverse preconditioner without having to remember abbreviations. preconditioner
Convenience property to return the readilycomputed preconditioner without having to remember abbreviations. Methods Summary
__call__
(arr_or_rv, RandomVariable], start, …)Transition a random variable or a realization of one. dispersion
(time, state, **kwargs)Evaluates l(t, x(t)) = L(t). drift
(time, state, **kwargs)Evaluates f(t, x(t)) = F(t) x(t) + u(t). jacobian
(time, state, **kwargs)maps t > F(t) precond2nordsieck
(step)Computes preconditioner inspired by Nordsieck. proj2coord
(coord)Projection matrix to \(i\)th coordinates. transition_realization
(real, start, stop, …)Transition a realization of a random variable from time \(t\) to time \(t+\Delta t\). transition_rv
(rv, start, stop, **kwargs)Transition a random variable from time \(t\) to time \(t+\Delta t\). Attributes Documentation

diffusionmatrix
¶ Evaluates Q.

dimension
¶ Spatial dimension (utility attribute).

dispersionmatrix
¶

driftmatrix
¶

force
¶

inverse_preconditioner
¶ Convenience property to return the readilycomputed inverse preconditioner without having to remember abbreviations.
Returns: Inverse preconditioner matrix \(P^{1}\) Return type: np.ndarray, shape=(d(q+1), d(q+1))

preconditioner
¶ Convenience property to return the readilycomputed preconditioner without having to remember abbreviations.
Returns: Preconditioner matrix \(P\) Return type: np.ndarray, shape=(d(q+1), d(q+1))
Methods Documentation

__call__
(arr_or_rv: Union[numpy.ndarray, RandomVariable], start: float = None, stop: float = None, **kwargs) > ('RandomVariable', typing.Dict)¶ Transition a random variable or a realization of one.
The input is either interpreted as a random variable or as a realization. Accordingly, the respective methods are called:
transition_realization()
ortransition_rv()
.

dispersion
(time, state, **kwargs)¶ Evaluates l(t, x(t)) = L(t).

drift
(time, state, **kwargs)¶ Evaluates f(t, x(t)) = F(t) x(t) + u(t).

jacobian
(time, state, **kwargs)¶ maps t > F(t)

precond2nordsieck
(step: float) > (<class 'numpy.ndarray'>, <class 'numpy.ndarray'>)[source]¶ Computes preconditioner inspired by Nordsieck.
Computes the matrix \(P\) given by
\[P = I_d \otimes diag (1, h, h^2, ..., h^q)\]as well as its inverse \(P^{1}\).
Parameters: step (float) – Step size \(h\) used for preconditioning. If \(h\) is so small that \(h^q! < 10^{15}\), it is being set to \(h = (\cdot 10^{15})^{1/q}\). Returns:  precond (np.ndarray, shape=(d(q+1), d(q+1))) – Preconditioner matrix \(P\).
 invprecond (np.ndarray, shape=(d(q+1), d(q+1))) – Inverse preconditioner matrix \(P^{1}\).

proj2coord
(coord: int) → numpy.ndarray[source]¶ Projection matrix to \(i\)th coordinates.
Computes the matrix
\[H_i = \left[ I_d \otimes e_i \right] P^{1},\]where \(e_i\) is the \(i\)th unit vector, that projects to the \(i\)th coordinate of a vector. If the ODE is multidimensional, it projects to each of the \(i\)th coordinates of each ODE dimension.
Parameters: coord (int) – Coordinate index \(i\) which to project to. Expected to be in range \(0 \leq i \leq q + 1\). Returns: Projection matrix \(H_i\). Return type: np.ndarray, shape=(d, d*(q+1))

transition_realization
(real, start, stop, **kwargs)¶ Transition a realization of a random variable from time \(t\) to time \(t+\Delta t\).
For random variable \(x_t\), it returns the random variable defined by
\[x_{t + \Delta t} \sim p(x_{t + \Delta t}  x_t = r) .\]This is different to
transition_rv()
which computes the parametrization of \(x_{t + \Delta t}\) based on the parametrization of \(x_t\).Nb: Think of transition as a verb, i.e. this method “transitions” a realization of a random variable.
Parameters:  real – Realization of the random variable.
 start – Starting point \(t\).
 stop – End point \(t + \Delta t\).
Returns:  RandomVariable – Random variable, describing the state at time \(t + \Delta t\) based on realization at time \(t\).
 dict – Additional information in form of a dictionary, for instance the crosscovariance in the prediction step, access to which is useful in smoothing.
See also
transition_rv()
 Apply transition to a random variable.

transition_rv
(rv, start, stop, **kwargs)¶ Transition a random variable from time \(t\) to time \(t+\Delta t\).
For random variable \(x_t\), it returns the random variable defined by
\[x_{t + \Delta t} \sim p(x_{t + \Delta t}  x_t) .\]This returns a random variable where the parametrization depends on the paramtrization of \(x_t\). This is different to
transition_rv()
which computes the parametrization of \(x_{t + \Delta t}\) based on a realization of \(x_t\).Nb: Think of transition as a verb, i.e. this method “transitions” a random variable.
Parameters:  rv – Realization of the random variable.
 start – Starting point \(t\).
 stop – End point \(t + \Delta t\).
Returns:  RandomVariable – Random variable, describing the state at time \(t + \Delta t\) based on realization at time \(t\).
 dict – Additional information in form of a dictionary, for instance the crosscovariance in the prediction step, access to which is useful in smoothing.
See also
transition_realization()
 Apply transition to a realization of a random variable.