snowdrop.src.numeric.filters package

Submodules

snowdrop.src.numeric.filters.diffuse module

Durbin-Koopman filter for non-stationary processes with missing and non-missing observations.

It is a modification of Standard Kalman filter which accounts for growth of a covariance matrix P as time progresses. The initial values of P matrix is assumed diffusive. During the first d diffusive periods matrix P is decomposed into Pstar and Pinf and separate equations for these matrices are used. After number d of diffusive periods standard Kalman filter is applied.

Translated from Dynare Matlab code to Python by A.Goumilevski.

snowdrop.src.numeric.filters.diffuse.diffuse_filter(T, Z, R, Q, H, Y, C, a0, pp, mm, Ts, Nd, n_shocks, data_index, Pinf1, Pstar1, decomp_flag=False, state_uncertainty_flag=False)[source]

Diffuse Kalman filter and smoother.

References:
  • S.J. Koopman and J. Durbin (2003, “Filtering and Smoothing of State Vector for Diffuse State Space Models”, Journal of Time Series Analysis, vol. 24(1), pp. 85-98).

  • Durbin/Koopman (2012): “Time Series Analysis by State Space Methods”, Oxford University Press, Second Edition, Ch. 5

Args:
T: mm*mm matrix

is the state transition matrix

Z: pp*mm matrix

is the selector matrix for observables in augmented state vector

R: mm*rr matrix

is the second matrix of the state equation relating the structural innovations to the state variables

Q: rr*rr matrix

is the covariance matrix of structural errors

H: pp*pp

is the matrix of variance of measurement errors

Y: Ts*pp array

is the vector of observations

C: mm array

is the vector of constants

a0: mm array

is the vector of initial values of endogenous variables

pp: int

is the number of observed variables

mm: int

is the number of state variables

Ts: int

is the sample size

Nd: int

is the shocks maximum lead number minus minimum lag number

n_shocks: int

is the number of shocks

data_index: 1*Tsarray

cell of column vectors of indices

nk: int

is the number of forecasting periods

Pinf1: mm*mm 2D array

is the diagonal matrix with with q ones and m-q zeros

Pstar1: mm*mm 2D array

is the variance-covariance matrix with stationary variables

kalman_tol: number, float

is the tolerance for reciprocal condition number

diffuse_kalman_tol: number, float

is the tolerance for reciprocal condition number (for Finf) and the rank of Pinf

decomp_flag: bool

if True, compute filter decomposition

state_uncertainty_flag: bool

if True, compute uncertainty about smoothed state estimate

Returns:
alphahat:

smoothed variables (a_{t|T})

epsilonhat:

smoothed measurement errors

etahat:

smoothed shocks

atilde:

matrix of updated variables (a_{t|t})

aK:

3D array of k step ahead filtered state variables (a_{t+k|t) (meaningless for periods 1:d)

P:

3D array of one-step ahead forecast error variance matrices

PK:

4D array of k-step ahead forecast error variance matrices (meaningless for periods 1:d)

decomp:

decomposition of the effect of shocks on filtered values

V:

3D array of state uncertainty matrices

N:

3D array of state uncertainty matrices

dLIK:

1D array of likelihood

log_likelihood:

1D array of cumulative likelihood

snowdrop.src.numeric.filters.diffuse.missing_obs_diffuse_multivariate_filter(T, Z, R, Q, H, Y, C, a0, pp, mm, Ts, Nd, n_shocks, data_index=None, nk=1, Pinf1=None, Pstar1=None, kalman_tol=1e-10, diffuse_kalman_tol=1e-06, decomp_flag=False, state_uncertainty_flag=False)[source]

Compute diffuse Kalman smoother without measurement error, in case of a non-singular variance-covariance matrix of observation errors.

Multivariate treatment of time series.

Note

Translated from Dynare Matlab code.

References:

See “Filtering and Smoothing of State Vector for Diffuse State Space Models”, S.J. Koopman and J. Durbin (2003, in Journal of Time Series Analysis, vol. 24(1), pp. 85-98). Durbin/Koopman (2012): “Time Series Analysis by State Space Methods”, Oxford University Press, Second Edition, Ch. 5

Args:
T: mm*mm matrix

is the state transition matrix

Z: pp*mm matrix

is the selector matrix for observables in augmented state vector

R: mm*rr matrix

is the second matrix of the state equation relating the structural innovations to the state variables

Q: rr*rr matrix

is the covariance matrix of structural errors

H: pp*pp

is the matrix of variance of measurement errors

Y: Ts*pp

is the vector of observations

C: mm

is the vector of constants

a0: mm

is the vector of initial values of endogenous variables

pp: int

is the number of observed variables

mm: int

is the number of state variables

Ts: int

is the sample size

Nd: int

is the shocks maximum lead number minus minimum lag number

n_shocks: int

is the number of shocks

data_index: 1*Ts

is the cell of column vectors of indices

nk: int

is the number of forecasting periods

Pinf1: mm*mm

is the diagonal matrix with with q ones and m-q zeros

Pstar1: mm*mm

is the variance-covariance matrix with stationary variables

kalman_tol:

is the tolerance for reciprocal condition number

diffuse_kalman_tol:

is the tolerance for reciprocal condition number (for Finf) and the rank of Pinf

decomp_flag:

if true, compute filter decomposition

state_uncertainty_flag:

if True, compute uncertainty about smoothed state estimate

Returns:
alphahat:

is the smoothed variables (a_{t|T})

epsilonhat:

is the smoothed measurement errors

etahat:

is the smoothed shocks

atilde:

is the matrix of updated variables (a_{t|t})

aK:

is the 3D array of k step ahead filtered state variables (a_{t+k|t) (meaningless for periods 1:d)

P:

is the 3D array of one-step ahead forecast error variance matrices

PK:

is the 4D array of k-step ahead forecast error variance matrices (meaningless for periods 1:d)

decomp:

is the decomposition of the effect of shocks on filtered values

V:

is the 3D array of state uncertainty matrices

N:

is the 3D array of state uncertainty matrices

dLIK:

is 1D array of likelihood

log_likelihood:

is 1D array of cumulative likelihood

snowdrop.src.numeric.filters.diffuse.missing_obs_diffuse_univariate_filter(T, Z, R, Q, H, Y, C, a0, pp, mm, Ts, Nd, n_shocks, data_index=None, nk=1, Pinf1=None, Pstar1=None, kalman_tol=1e-10, diffuse_kalman_tol=1e-06, decomp_flag=False, state_uncertainty_flag=False)[source]

Compute diffuse Kalman smoother in case of a singular var-cov matrix.

Univariate treatment of multivariate time series. It is applied for singular and non-singular matrix of observation errors.

References:
  • Durbin, Koopman (2012): “Time Series Analysis by State Space Methods”, Oxford University Press, Second Edition, Ch. 6.4 + 7.2.5

  • Koopman, Durbin (2000): “Fast Filtering and Smoothing for Multivariatze State Space Models”, Journal of Time Series Analysis, vol. 21(3), pp. 281-296.

  • S.J. Koopman and J. Durbin (2003): “Filtering and Smoothing of State Vector for Diffuse State Space Models”, Journal of Time Series Analysis, vol. 24(1), pp. 85-98.

Note

Translated from Dynare Matlab code.

Args:
T: mm*mm matrix

is the state transition matrix

Z: pp*mm matrix

is the selector matrix for observables in augmented state vector

R: mm*rr matrix

is the second matrix of the state equation relating the structural innovations to the state variables

Q: rr*rr matrix

is the covariance matrix of structural errors

H: pp

is the matrix of variance of measurement errors

Y: Ts*pp

is the vector of observations

C: mm

is the vector of constants

a0: mm

is the vector of initial values of endogenous variables

pp: int

is the number of observed variables

mm: int

is the number of state variables

Ts: int

is the sample size

Nd: int

is the shocks maximum lead number minus minimum lag number

n_shocks: int

is the number of shocks

data_index: array, dtype=float

is the Ts cell of column vectors of indices

nk: int

is the number of forecasting periods

Pinf1: mm*mm

is the diagonal matrix with with q ones and m-q zeros

Pstar1: mm*mm

is the variance-covariance matrix with stationary variables

kalman_tol:

is the tolerance for zero divider

diffuse_kalman_tol:

is the tolerance for zero divider

decomp_flag: if true,

is the compute filter decomposition

state_uncertainty_flag: i

is the f True, compute uncertainty about smoothed state estimate

Returns:
alphahat:

is the smoothed state variables (a_{t|T})

epsilonhat:

is the measurement errors

etahat:

is the smoothed shocks

a:

is the matrix of updated variables (a_{t|t})

aK:

is the 3D array of k step ahead filtered state variables a_{t+k|t} (meaningless for periods 1:d)

P:

is the 3D array of one-step ahead forecast error variance matrices

PK:

is the 4D array of k-step ahead forecast error variance matrices (meaningless for periods 1:d)

decomp:

is the decomposition of the effect of shocks on filtered values

V:

is the 3D array of state uncertainty matrices

N:

is the 3D array of state uncertainty matrices

dLIK:

is 1D array of likelihood

log_likelihood:

id 1D array of cumulative likelihood

snowdrop.src.numeric.filters.filter module

Created on Mon Feb 26 18:49:21 2018.

@author: agoumilevski

snowdrop.src.numeric.filters.filter.linear_filter(model, T, periods, y0, Qm=None, Hm=None, obs=None, MULT=1, skip=0, missing_obs=None, ind_non_missing=None)[source]

Apply Kalman filter to linear model.

Parameters:
param model:

Model object.

type model:

Model.

param T:

Time span.

type T:

int.

param periods:

Array of endogenous variables.

type periods:

numpy.array.

param y0:

Starting values of endogenous variables.

type y0:

numpy.array

param Qm:

Covariance matrix of errors of endogenous variables.

type Qm:

numpy.array.

param Hm:

Covariance matrix of errors of measurement variables.

type Hm:

numpy.array.

param obs:

Measurements.

type obs:

List.

param MULT:

Multiplier of the terminal time range. If set greater than one than the solution will be computed for this extended time range interval.

type MULT:

int.

param missing_obs:

Matrix of logical values with True or False value depending on observation missing or not.

type missing_obs:

Matrix of bools.

returns:

Filtered solution.

snowdrop.src.numeric.filters.filter.nonlinear_filter(model, T, periods, y0, Qm=None, Hm=None, obs=None, MULT=1, skip=0, missing_obs=None, ind_non_missing=None)[source]

Apply Kalman filter to non-linear model.

Parameters:
param model:

Model object.

type model:

Model.

param T:

Time span.

type T:

int.

param periods:

Array of endogenous variables.

type periods:

numpy.array.

param y0:

Starting values of endogenous variables.

type y0:

numpy.array

param Qm:

Covariance matrix of errors of endogenous variables.

type Qm:

numpy.array.

param Hm:

Covariance matrix of errors of measurement variables.

type Hm:

numpy.array.

param obs:

Measurements.

type obs:

list.

param MULT:

Multiplier of the terminal time range. If set greater than one than the solution will be computed for this extended time range interval.

type MULT:

int.

param missing_obs:

Matrix of logical values with True or False value depending on observation missing or not.

type missing_obs:

Matrix of bools.

returns:

Filtered solution.

snowdrop.src.numeric.filters.filters module

Created on Mon Feb 26 18:49:21 2018

@author: agoumilevski

snowdrop.src.numeric.filters.filters.BPASS(data, W1, W2, drift=True, bChristianoFitzgerald=True, cutoff=False)[source]

Bandpass filter to time series.

snowdrop.src.numeric.filters.filters.Bryson_Frazier_Smoother(x, u, F, H, K, Sinv, P, L, l)[source]

Implement modified Bryson-Frazier (BF) smoother for a state-space model.

BF algorithm is:

\[ \begin{align}\begin{aligned}L_{t} &= H'_{t}*S_{t}^{-1}*H_{t} + (I-K_{t}*H_{t})'*L_{t}*(I-K_{t}*H_{t})\\L_{t-1} &= F'_{t}*Lt_{t}*F_{t}\\L_{T} &= 0\\l_{t} &= -H'_{t}*S_{t}^{-1}*u_{t} + (I-K_{t}*H_{t})'*l_{t}\\l_{t-1} &= F'_{t}*lt_{t}\\l_{T} &= 0 \end{aligned}\end{align} \]

Smoothed variables and covariances

\[ \begin{align}\begin{aligned}xs_{t|T} &= x_{t} - P_{t}*l_{t}\\P_{t|T} &= P_{t} - P_{t}*L_{t}*P_{t}\end{aligned}\end{align} \]
Args:
x: array, dtype=float

is the state variables,

u: array, dtype=float

is the vector of one-step-ahead prediction errors

F: 2D array, dtype=float

is the state-transition matrix,

H: 2D array, dtype=float

is the observation matrix,

K: 2D array, dtype=float

is the Kalman gain,

Sinv: 2D array, dtype=float

is the inverse of matrix S,

P: 2D array, dtype=float

is the predicted error covariance matrix,

L: 2D array, dtype=float

is the covariance matrix,

l: 2D array, dtype=float

is the temporary matrix

For details, see https://en.wikipedia.org/wiki/Kalman_filter

snowdrop.src.numeric.filters.filters.DK_Filter(x, xtilde, y, v, T, Z, P, H, Q, K, C=0, ind_non_missing=None, t=-1)[source]

Implement non-diffusive Durbin & Koopman version of Kalman filter algorithm for state-space model.

Multivariate approach to Kalman filter.

State-space model

\[ \begin{align}\begin{aligned}x_{t} = T_{t} * x_{t-1} + C_{t} + w_{t}, w = Normal(0, Q)\\y_{t} = Z_{t} * x_{t} + v_{t}, v = Normal(0, H)\end{aligned}\end{align} \]

Durbin-Koopman algorithm

\[ \begin{align}\begin{aligned}v_{t} &= y_{t} - Z_{t} * x_{t}\\ x_{t} &= T_{t} * (x_{t-1} + K_{t} * v_{t}) + C_{t}\\ F_{t} &= H_{t} + Z_{t} * P_{t-1} * Z'_{t}\\ K_{t} &= P_{t} * Z'_{t} * F_{t}^{-1}\\L_{t} &= K_{t} * Z_{t} \\ P_{t+1} &= T_{t} * P_{t} *L'_{t-1} + Q_{t} \end{aligned}\end{align} \]
Args:
x: array, dtype=float

is the variables at t+1|t+1 step,

xtilde: array, dtype=float

is the variables at t|t+1 forecast step,

y: array, dtype=float

is the observation variables,

v: array, dtype=float

is the residual,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

P: 2D array, dtype=float

is the predicted error covariance matrix,

H: 2D array, dtype=float

is the covariance matrix of measurement errors,

Q: 2D array, dtype=float

is the covariance matrix of structural errors,

K: 2Darray, dtype=float

is the Kalman gain,

C: array, dtype=float

is the constant term matrix,

ind_non_missing: array, dtype=int

are the indices of non-missing observations,

tol1: number, float

Kalman tolerance for reciprocal condition number,

tol2: number, float

diffuse Kalman tolerance for reciprocal condition number (for Finf) and the rank of Pinf

For details, see Durbin & Koopman 2012, “Time series Analysis by State Space Methods” and

Koopman & Durbin 2000, “Fast Filtering and Smoothing For Multivariate State Space Models”

snowdrop.src.numeric.filters.filters.DK_Smoother(x, y, T, Z, F, iF, P, QRt, K, H, r, v, ind_non_missing=None, t=-1)[source]

Implement non-diffusive Durbin-Koopman smoother for state-space model.

Durbin-Koopman smoothing algorithm is

\[ \begin{align}\begin{aligned}r_{t-1} &= Z'_{t}*F_{t}^{-1}*v_{t} + L'_{t}**r_{t} \\r_{T} &= 0 \\L_{t} &= I - K_{t}*Z_{t} \end{aligned}\end{align} \]

Smoothed variables

\[s_{t} = x_{t} + P_{t}*r_{t-1}\]
Args:
x: array, dtype=float

is the updated state variables vector,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

F: 2D array, dtype=float

is the residual covariance matrix,

iF: 2D array, dtype=float

is the inverse of matrix F,

P: 2D array, dtype=float

is the updated error covariance matrix,

QRt: 2D array, dtype=float

is the product of the covariance matrix of structural errors and shock matrix

K: 2D array, dtype=float

is the Kalman gain matrix,

H: 2D array, dtype=float

is the covariance matrix of measurement errors,

r: array, dtype=float

is the vector of one-step-back errors,

v: array, dtype=float

is the vector of one-step-ahead prediction errors,

ind_non_missing: array, dtype=int

are the indices of non-missing observations.

For details, see Durbin & Koopman (2012), “Time series Analysis by State Space Methods” and Koopman & Durbin (2000), “Fast Filtering and Smoothing For Multivariate State Space Models” and Koopman & Durbin (2003), in Journal of Time Series Analysis, vol. 24(1), pp. 85-98

snowdrop.src.numeric.filters.filters.Durbin_Koopman_Diffuse_Filter(x, xtilde, y, T, Z, P, H, Q, C=0, bUnivariate=False, ind_non_missing=None, Pstar=None, Pinf=None, tol1=1e-06, tol2=1e-06, t=-1)[source]

Implement diffuse Durbin & Koopman version of Kalman filter algorithm for state-space model.

State-space model:

\[ \begin{align}\begin{aligned}x_{t} = T_{t} * x_{t-1} + C_{t} + w_{t}, w = Normal(0, Q)\\y_{t} = Z_{t} * x_{t} + v_{t}, v = Normal(0, H)\end{aligned}\end{align} \]

Durbin-Koopman algorithm

\[ \begin{align}\begin{aligned}v_{t} &= y_{t} - Z_{t} * x_{t}\\ x_{t} &= T_{t} * (x_{t-1} + K_{t} * v_{t}) + C_{t}\\ F_{t} &= H_{t} + Z_{t} * P_{t-1} * Z'_{t}\\ K_{t} &= P_{t} * Z'_{t} * F_{t}^{-1}\\L_{t} &= K_{t} * Z_{t}\\ P_{t+1} &= T_{t} * P_{t} *L'_{t-1} + Q_{t} \end{aligned}\end{align} \]

Note

Translated from Dynare Matlab code.

Args:
x: array, dtype=float

is the state variables,

xtilde: array, dtype=float

is the variables at t|t+1 forecast step,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

P: 2D array, dtype=float

is the predicted error covariance matrix,

H: 2D array, dtype=float

is the covariance matrix of measurement errors,

Q: 2D array, dtype=float

is the covariance matrix of structural errors,

C: array, dtype=float

is the constant term matrix,

ind_non_missing: array, dtype=int

are the indices of non-missing observations,

Pstar: 2D array, dtype=float

is the finite part of the predicted error covariance matrix,

Pinf: 2D array, dtype=float

is the infinite part of the predicted error covariance matrix,

tol1: number, float

Kalman tolerance for reciprocal condition number,

tol2: number, float

diffuse Kalman tolerance for reciprocal condition number (for Finf) and the rank of Pinf

For details, see Durbin & Koopman 2012, “Time series Analysis by State Space Methods” and

Koopman & Durbin 2000, “Fast Filtering and Smoothing For Multivariate State Space Models”

snowdrop.src.numeric.filters.filters.Durbin_Koopman_Diffuse_Smoother(x, y, T, Z, F, iF, QRt, K, Kstar, Kinf, H, P, P1, Pstar, Pinf, Pstar1, Pinf1, v, r, r0, r1, Linf=None, Lstar=None, ind_non_missing=None, diffuse=True, Finf_singular=False, bUnivariate=False, tol1=1e-06, tol2=1e-06, t=-1)[source]

Implement diffuse Durbin-Koopman smoother for state-space model.

Durbin-Koopman algorithm is

\[ \begin{align}\begin{aligned}r_{t-1} &= Z'_{t}*F_{t}^{-1}*v_{t} + L'_{t}**r_{t}\\r_{T} &= 0\\L_{t} &= I - K_{t}*Z_{t} \end{aligned}\end{align} \]

Smoothed variables

\[s_{t} = x_{t} + P_{t}*r_{t-1}\]

Note

Translated from Dynare Matlab code.

Args:
x: array, dtype=float

is the predicted state variables vector,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

F: 2D array, dtype=float

is the residual covariance matrix,

iF: 2D array, dtype=float

is the inverse of matrix F,

QRt: 2D array, dtype=float

is the product of the covariance matrix of structural errors and shock matrix

K: 2D array, dtype=float

is the Kalman gain matrix,

Kinf: 2D array, dtype=float

is the infinite part of the Kalman gain matrix,

Kstar: 2D array, dtype=float

is the finite part of the Kalman gain matrix,

H: 2D array, dtype=float

is the observation errors covariance matrix,

P: 2D array, dtype=float

is the updated error covariance matrix P,

P1:2D array, dtype=float

is the updated error covariance matrix P1,

Pstar: 2D array, dtype=float

is the finite part of P matrix,

Pinf: 2D array, dtype=float

is the infinite part of P matrix,

Pstar1: 2D array, dtype=float

is the finite part of P1 matrix,

Pinf1: 2D array, dtype=float

is the infinite part of P1 matrix,

ind_non_missing: array, dtype=int

are the indices of non-missing observations.

v: array, dtype=float

is the vector of one-step-ahead prediction errors,

r: array, dtype=float

is a smoothed state vector,

r0: array, dtype=float

is zero approximation in a Taylor expansion of a smoothed state vector,

r1: array, dtype=float

is the first approximation in a Taylor expansion of a smoothed state vector,

Pstar: 2D array, dtype=float

is the finite part of the predicted error covariance matrix,

Pinf: 2D array, dtype=float

is the infinite part of the predicted error covariance matrix,

Finf_singular: bool

True if Finf is a singular matrix and False otherwise,

bUnivariate: bool

True if univariate Kalman filter is used and False otherwise,

tol1: number, float

is Kalman tolerance for reciprocal condition number,

tol2: number, float

is diffuse Kalman tolerance for reciprocal condition number (for Finf) and the rank of Pinf

For details, see Durbin & Koopman (2012), “Time series Analysis by State Space Methods” and Koopman & Durbin (2000), “Fast Filtering and Smoothing For Multivariate State Space Models” and Koopman & Durbin (2003), in Journal of Time Series Analysis, vol. 24(1), pp. 85-98)

snowdrop.src.numeric.filters.filters.Durbin_Koopman_Non_Diffuse_Filter(x, xtilde, y, T, Z, P, H, Q, C=0, bUnivariate=False, ind_non_missing=None, tol1=1e-06, tol2=1e-06, t=-1)[source]

Implement non-diffusive Durbin & Koopman version of Kalman filter algorithm for state-space model.

State-space model:

\[ \begin{align}\begin{aligned}x_{t} = T_{t} * x_{t-1} + C_{t} + w_{t}, w = Normal(0, Q)\\y_{t} = Z_{t} * x_{t} + v_{t}, v = Normal(0, H)\end{aligned}\end{align} \]

Durbin-Koopman algorithm

\[ \begin{align}\begin{aligned}v_{t} &= y_{t} - Z_{t} * x_{t}\\ x_{t} &= T_{t} * (x_{t-1} + K_{t} * v_{t}) + C_{t}\\ F_{t} &= H_{t} + Z_{t} * P_{t-1} * Z'_{t}\\ K_{t} &= P_{t} * Z'_{t} * F_{t}^{-1}\\L_{t} &= K_{t} * Z_{t}\\ P_{t+1} &= T_{t} * P_{t} *L'_{t-1} + Q_{t} \end{aligned}\end{align} \]
Args:
x: array, dtype=float

is the variables at t+1|t+1 step,

xtilde: array, dtype=float

is the variables at t|t+1 forecast step,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

P: 2D array, dtype=float

is the predicted error covariance matrix,

H: 2D array, dtype=float

is the covariance matrix of measurement errors,

Q: 2D array, dtype=float

is the covariance matrix of structural errors,

C: array, dtype=float

is the constant term matrix,

ind_non_missing: array, dtype=int

are the indices of non-missing observations,

tol1: number, float

Kalman tolerance for reciprocal condition number,

tol2: number, float

diffuse Kalman tolerance for reciprocal condition number (for Finf) and the rank of Pinf

For details, see Durbin & Koopman 2012, “Time series Analysis by State Space Methods” and

Koopman & Durbin 2000, “Fast Filtering and Smoothing For Multivariate State Space Models”

snowdrop.src.numeric.filters.filters.Durbin_Koopman_Non_Diffuse_Smoother(x, y, T, Z, F, iF, P, P1, QRt, K, r, v, ind_non_missing=None, bUnivariate=False, tol=1e-06, t=-1)[source]

Implement non-diffusive Durbin-Koopman smoother for state-space model.

Durbin-Koopman algorithm is

\[ \begin{align}\begin{aligned}r_{t-1} &= Z'_{t}*F_{t}^{-1}*v_{t} + L'_{t}**r_{t}\\r_{T} &= 0\\N_{t-1} &= Z'_{t}*F_{t}^{-1}*Z_{t} + L'_{t}*N_{t}*L{t}\\N_{T} &= 0\\L_{t} &= I - K_{t}*Z_{t} \end{aligned}\end{align} \]

Smoothed variables

\[s_{t} = x_{t} + P_{t}*r_{t-1}\]
Args:
x: array, dtype=float

is the predicted variables vector,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

F: 2D array, dtype=float

is the residual covariance matrix,

iF: 2D array, dtype=float

is the inverse of matrix F,

P: 2D array, dtype=float

is the updated error covariance matrix,

QRt: 2D array, dtype=float

is the product of the covariance matrix of structural errors and shock matrix,

K:2D array, dtype=float

is the Kalman gain matrix,

r: array, dtype=float

is the vector of one-step-back errors,

v: array, dtype=float

is the vector of one-step-ahead prediction errors,

ind_non_missing: array, dtype=int

are the indices of non-missing observations.

bUnivariate: bool

True if univariate Kalman filter is used and False otherwise,

tol: number, float

is Kalman tolerance for reciprocal condition number.

For details, see Durbin & Koopman (2012), “Time series Analysis by State Space Methods” and Koopman & Durbin (2000), “Fast Filtering and Smoothing For Multivariate State Space Models” and Koopman & Durbin (2003), in Journal of Time Series Analysis, vol. 24(1), pp. 85-98

snowdrop.src.numeric.filters.filters.HPF(data, lmbd=1600)[source]

Hodrick-Prescott filter to a time series.

snowdrop.src.numeric.filters.filters.HPfilter(data, lmbd=1600)[source]

Apply a Hodrick-Prescott filter to a dataset.

The return value is the filtered data-set found according to:

\[\min_{T_{t}}[ \sum_{t=0}((X_{t} - T_{t})^2 + \lambda*((T_{t+1} - T_{t}) - (T_{t} - T_{t-1}))^2) ]\]
Args:
data: array, dtype=float

The data set for which you want to apply the HP filter. This mustbe a numpy array.

lmbd: array, dtype=float

This is the value for lambda as used in the equation.

Returns:
T: array, dtype=float

The solution to the minimization equation above (the trend).

Cycle: array, dtype=float

This is the ‘stationary data’ found by X - T.

Note

This function implements sparse methods to be efficient enough to handle very large data sets.

snowdrop.src.numeric.filters.filters.Kalman_Filter(x, y, T, Z, P, Q, H, C=0, ind_non_missing=None, bUnivariate=False, t=-1, tol=1e-06)[source]

Implement Kalman filter algorithm for state-space model.

Args:
x: array, dtype=float

is the state variables,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

P: 2D array, dtype=float

is the predicted error covariance matrix,

Q: 2D array, dtype=float

is the covariance matrix of state variables (endogenous variables),

H: 2D array, dtype=float

is the covariance matrix of space variables (measurement variables),

C: array, dtype=float

is the constant term matrix,

bUnivariate: bool

if True univariate Kalman filter is used and if False - multivariate,

ind_non_missing: array, dtype=int

are the indices of non-missing observations,

tol: number, dtype=float

is the tolerance parameter (iteration over the Riccati equation)

For details, see https://en.wikipedia.org/wiki/Kalman_filter

snowdrop.src.numeric.filters.filters.Kalman_Filter_SS(x, y, T, Z, K, F, iF, C=0, ind_non_missing=None)[source]

Implement Kalman filter algorithm for stationary state-space model.

Args:
x: array, dtype=float

is the state variables,

y: array, dtype=float

is the observation variables,

T: 2D array, dtype=float

is the state-transition matrix,

Z: 2D array, dtype=float

is the observation matrix,

K: 2D array, dtype=float

is the Kalman gain matrix,

F: 2D array, dtype=float

is the pre-fit residual covariance,

iF: 2D array, dtype=float

is the inverse of the pre-fit residual covariance,

C: array, dtype=float

is the constant term matrix,

bUnivariate: bool

if True univariate Kalman filter is used and if False - multivariate,

ind_non_missing: array, dtype=int

are the indices of non-missing observations,

tol: number, dtype=float

is the tolerance parameter

For details, see https://en.wikipedia.org/wiki/Kalman_filter

snowdrop.src.numeric.filters.filters.LRXfilter(y, lmbd=1600, w1=1, w2=0, w3=0, dprior=None, prior=None)[source]

Apply enchanced HP filter.

This is LRX filter named after Laxton-Rose-Xie.

A brief explanation of terms used in LRX filter

Syntax:

y_eq = LRXfilter(y,lmbd[,w1[,w2[,w3[,DPRIOR[,PRIOR]]]]])

Args:
lmbd: float

tightness constraint for HP filter,

y: array, dtype=float

series to be detrended,

w1: float

weight on (y - y_eq)^2,

w2: float

weight on [(y_eq - y_eq(-1)) - dprior]^2,

w3: float

weight on (y_eq - prior)^2,

dprior: array, dtype=float

growth rate of y in steady-state,

prior: array, dtype=float

previously computed y_eq.

Returns:

y_eq: The trend component,

gap: The cycle component

Note

the default value of lmbd is 1600,

the default value of w1 is a vector of ones,

the default value of w2 is a vector of zeros,

the default value of w3 is a vector of zeros,

the default value of dpr is a vector of zeros,

the default value of pr is a vector of zeros,

if any of w1, w2, w3, dpr, pr is of length 1, then it is extended to a vector of the appropriate length, in which all the entries are equal.

snowdrop.src.numeric.filters.filters.Rauch_Tung_Striebel_Smoother(Ytp, Xt, Xtp, Ft, Pt, Ptp, Ttp)[source]

Implement Rauch-Tung-Striebel (RTS) smoother for state-space model.

RTS algorithm is

\[ \begin{align}\begin{aligned}L_{t} &= P_{t|t} * F'_{t} * P_{t+1|t}^{-1}\\X_{t|T} &= X_{t|t} + L_{t} * (X_{t+1|T} - X_{t+1|t})\\P_{t|T} &= P_{t+1|t} + L_{t} * (P_{t+1|T} - P_{t+1|t}) * L'_{t}\end{aligned}\end{align} \]
Args:
Ytp: array, dtype=float

is the smoothed state variables at time t+1,

Xt: array, dtype=float

is the state variables at time t,

Xtp: array, dtype=float

is the state variables at time t+1,

Ft: 2D array, dtype=float

is the state-transition matrix at time t,

Pt: 2D array, dtype=float

is the predicted error covariance matrix at time t,

Ptp: 2D array, dtype=float

is the predicted error covariance matrix at time t+1,

Ttp: 2D array, dtype=float

is the smoothed covariance matrix of state variables.

For details, see https://en.wikipedia.org/wiki/Kalman_filter

snowdrop.src.numeric.filters.filters.Standard_Kalman_Filter(x, y, T, Z, P, Q, H, C=0, ind_non_missing=None, t=-1)[source]

Implement Kalman filter algorithm for state-space model.

This implementation assumes that measurement equation is specified at time t.

\[ \begin{align}\begin{aligned}x_{t+1} = T_{t} * x_{t} + C_{t} + w_{t}, w_{t} = Normal(0, Q),\\y_{t+1} = Z_{t} * x_{t+1} + v_{t}, v_{t} = Normal(0, H)\end{aligned}\end{align} \]

KF uses two phase algorithm

Predict

\[ \begin{align}\begin{aligned}x_{t|t} &= x_{t|t-1} + P_{t|t-1} * L_{t} * (y_{t} - Z_{t} * x_{t|t-1})\\P_{t|t} &= P_{t|t-1} - L * Z_{t} * P_{t|t-1}\\L{t} &= Z_{t} * [Z_{t}*P_{t|t-1}*Z'_{t} + H]^{-1}\end{aligned}\end{align} \]

Update

\[ \begin{align}\begin{aligned}v_{t} &= y_{t} - Z_{t} * x_{t|t-1}\\x_{t+1|t} &= T_{t} * x_{t|t} + C_{t}\\P_{t+1|t} &= T_{t} * P_{t|t} * T'_{t} + Q_{t}\end{aligned}\end{align} \]
Args:
xarray, dtype=float

is the state variables,

yarray, dtype=float

is the observation variables,

T2D array, dtype=float

is the state-transition matrix,

Z2D array, dtype=float

is the observation matrix,

P2D array, dtype=float

is the predicted error covariance matrix,

Q2D array, dtype=float

is the covariance matrix of state variables (endogenous variables),

H2D array, dtype=float

is the covariance matrix of space variables (measurement variables),

Carray, dtype=float

is the constant term matrix,

ind_non_missingarray, dtype=int

are the indices of non-missing observations.

For details, please see https://stanford.edu/class/ee363/lectures/kf.pdf

snowdrop.src.numeric.filters.filters.bandpass_filter(data, w1, w2, k=None)[source]

Apply a bandpass filter to data.

This band-pass filter is of kth order. It selects the band between w1 and w2.

Args:
data: array, dtype=float

The data you wish to filter

k: number, int

The order of approximation for the filter. A max value for this is data.size/2

w1: number, float

This is the lower bound for which frequencies will pass through.

w2: number, float

This is the upper bound for which frequencies will pass through.

Returns:
y: array, dtype=float

The filtered data.

snowdrop.src.numeric.filters.filters.fast_kalman_filter(Y, N, a, P, T, H, Z, kalman_tol=1e-10)[source]

Computes the likelihood of a stationary state space model, given initial condition for the states (mean and variance).

References:

Edward P. Herbst, (2012). ‘Using the ”Chandrasekhar Recursions” for Likelihood Evaluation of DSGE Models’’

Args:
Y: array, dtype=float

Measurement data,

N: scalar, int

Last period,

a: array, dtype=float

Initial mean of the state vector,

P: 2D array, dtype=float

Initial covariance matrix of the state vector,

T: 2D array, dtype=float

Transition matrix of the state equation,

H: 2D array, dtype=float

Covariance matrix of the measurement errors (if no measurement errors set H as a zero scalar),

Z: 2D array, dtype=float

Matrix relating the states to the observed variables or vector of indices,

kalman_tol: scalar, float

Tolerance parameter (rcond, inversibility of the covariance matrix of the prediction errors),

snowdrop.src.numeric.filters.kalman module

Basic implementation of the Kalman filter (and smoother).

Overview

The Kalman filter/smoother is a well-known algorithm for computing recursively the filtering/smoothing distributions of a linear Gaussian model, i.e. a model of the form:

\[\begin{split}X_0 & \sim N(\mu_0,C_0) \\ X_t & = F X_{t-1} + U_t, \quad U_t \sim N(0, C_X) \\ Y_t & = G X_t + V_t, \quad V_t \sim N(0, C_Y)\end{split}\]

Linear Gaussian models and the Kalman filter are covered in Chapter 7 of the book.

MVLinearGauss class and subclasses

To define a specific linear Gaussian model, we instantiate class MVLinearGauss (or one its subclass) as follows:

import numpy as np
from particles import kalman

ssm = kalman.MVLinearGauss(F=np.ones((1, 2)), G=np.eye(2), covX=np.eye(2),
                           covY=.3)

where the parameters have the same meaning as above. It is also possible to specify mu0`and `cov0 (the mean and covariance of the initial state X_0). (See the documentation of the class for more details.)

Class MVLinearGauss is a sub-class of StateSpaceModel in module state_space_models, so it inherits methods from its parent such as:

true_states, data = ssm.simulate(30)

Class MVLinearGauss implements methods proposal, proposal0 and logeta, which correspond respectively to the optimal proposal distributions and auxiliary function for a guided or auxiliary particle filter; see Chapter 11 and module state_space_models for more details. (That the optimal quantities are tractable is, of course, due to the fact that the model is linear and Gaussian.)

To define a univariate linear Gaussian model, you may want to use instead the more conveniently parametrised class LinearGauss (which is a sub-class of MVLinearGauss):

ssm = LinearGauss(rho=0.3, sigX=1., sigY=.2, sig0=1.)

which corresponds to model:

\[\begin{split}X_0 & \sim N(0, \sigma_0^2) \\ X_t|X_{t-1}=x_{t-1} & \sim N(\rho * X_{t-1},\sigma_X^2) \\ Y_t |X_t=x_t & \sim N(x_t, \sigma_Y^2)\end{split}\]

Another sub-class of MVLinearGauss defined in this module is MVLinearGauss_Guarniero_etal, which implements a particular class of linear Gaussian models often used as a benchmark (after Guarniero et al, 2016).

Kalman class

The Kalman filter is implemented as a class, Kalman, with methods filter and smoother. When instantiating the class, one passes as arguments the data, and an object that represents the considered model (i.e. an instance of MvLinearGauss, see above):

kf = kalman.Kalman(ssm=ssm, data=data)
kf.filter()

The second line implements the forward pass of a Kalman filter. The results are stored as lists of MeanAndCov objects, that is, named tuples with attributes ‘mean’ and ‘cov’ that represent a Gaussian distribution. For instance:

kf.filt[3].mean  # mean of the filtering distribution at time 3
kf.pred[7].cov  # cov matrix of the predictive distribution at time 7

The forward pass also computes the log-likelihood of the data:

kf.logpyt[5]  # log-density of Y_t | Y_{0:t-1} at time t=5

Smoothing works along the same lines:

kf.smoother()

then object kf contains a list called smooth, which represents the successive (marginal) smoothing distributions:

kf.smth[8].mean  # mean of the smoothing dist at time 8

It is possible to call method smoother directly (without calling filter first). In that case, the filtering step is automatically performed as a preliminary step.

Kalman objects as iterators

It is possible to perform the forward pass step by step; in fact a Kalman object is an iterator:

kf = kalman.Kalman(ssm=ssm, data=data)
next(kf)  # one step
next(kf)  # one step

If you run the smoother after k such steps, you will obtain the smoothing distribution based on the k first data-points. It is therefore possible to compute recursively the successive smoothing distributions, but (a) at a high CPU cost; and (b) at each time, you must save the results somewhere, as attribute kf.smth gets written over and over.

Functions to perform a single step

The module also defines low-level functions that perform a single step of the forward or backward step. Some of these function makes it possible to perform such steps in parallel (e.g. for N predictive means). The table below lists these functions. Some of the required inputs are MeanAndCov objects, which may be defined as follows:

my_predictive_dist = kalman.MeanAndCov(mean=np.ones(2), cov=np.eye(2))

Function (with signature)

predict_step(F, covX, filt)

filter_step(G, covY, pred, yt)

filter_step_asarray(G, covY, pred, yt)

smoother_step(F, filt, next_pred, next_smth)

class snowdrop.src.numeric.filters.kalman.Kalman(ssm=None, data=None)[source]

Bases: object

Kalman filter/smoother.

See the documentation of the module for more details.

filter()[source]

Forward recursion: compute mean/variance of filter and prediction.

next()[source]
smoother()[source]

Backward recursion: compute mean/variance of marginal smoother.

Performs the filter step in a preliminary step if needed.

property t
class snowdrop.src.numeric.filters.kalman.LinearGauss(**kwargs)[source]

Bases: MVLinearGauss

A basic (univariate) linear Gaussian model.

\[\begin{split}X_0 & \sim N(0, \sigma_0^2) \\ X_t|X_{t-1}=x_{t-1} & \sim N(\rho * X_{t-1},\sigma_X^2) \\ Y_t |X_t=x_t & \sim N(x_t, \sigma_Y^2)\end{split}\]

Note

If parameter sigma0 is set to None, it is replaced by the quantity that makes the state process invariant: \(\sigma_X^2 / (1 - \rho^2)\)

PX(t, xp)[source]

Law of X_t at time t, given X_{t-1} = xp

PX0()[source]

Law of X_0 at time 0

PY(t, xp, x)[source]

Conditional distribution of Y_t, given the states.

default_params = {'rho': 0.9, 'sigma0': None, 'sigmaX': 1.0, 'sigmaY': 0.2}
logeta(t, x, data)[source]
proposal(t, xp, data)[source]

Proposal kernel (to be used in a guided or auxiliary filter).

Parameter

t: int

time

x:

particles

data: list-like

data

proposal0(data)[source]
class snowdrop.src.numeric.filters.kalman.MVLinearGauss(F=None, G=None, covX=None, covY=None, mu0=None, cov0=None)[source]

Bases: StateSpaceModel

Multivariate linear Gaussian model.

\[X_0 & \sim N(\mu_0, cov_0) \ X_t & = F * X_{t-1} + U_t, \quad U_t\sim N(0, cov_X) \ Y_t & = G * X_t + V_t, \quad V_t \sim N(0, cov_Y)\]

The only mandatory parameters are covX and covY (from which the dimensions dx and dy of, respectively, X_t, and Y_t, are deduced). The default values for the other parameters are:

  • mu0:: an array of zeros (of size dx)

  • cov0: cov_X

  • F: Identity matrix of shape (dx, dx)

  • G: (dy, dx) matrix such that G[i, j] = 1[i=j]

Note

The Kalman filter takes as an input an instance of this class (or one of its subclasses).

PX(t, xp)[source]

Law of X_t at time t, given X_{t-1} = xp

PX0()[source]

Law of X_0 at time 0

PY(t, xp, x)[source]

Conditional distribution of Y_t, given the states.

check_shapes()[source]

Check all dimensions are correct.

logeta(t, x, data)[source]
proposal(t, xp, data)[source]

Proposal kernel (to be used in a guided or auxiliary filter).

Parameter

t: int

time

x:

particles

data: list-like

data

proposal0(data)[source]
class snowdrop.src.numeric.filters.kalman.MVLinearGauss_Guarniero_etal(alpha=0.4, dx=2)[source]

Bases: MVLinearGauss

Special case of a MV Linear Gaussian ssm from Guarnierio et al. (2016).

\[ \begin{align}\begin{aligned}G = cov_X = cov_Y = cov_0 = I_{d_x}\\F_{i, j} = lpha^ { 1 + |i-j|}\end{aligned}\end{align} \]

See MVLinearGauss for the definition of these quantities.

Parameters

alpha: float (default: 0.4)

value of alpha

dx: int (must be >1; default: 2)

dimension of state-space

Reference

Guarnierio et al (2016). The Iterated Auxiliary Particle Filter,

arxiv:1511.06286, JASA.

class snowdrop.src.numeric.filters.kalman.MeanAndCov(mean, cov)

Bases: tuple

cov

Alias for field number 1

mean

Alias for field number 0

snowdrop.src.numeric.filters.kalman.dotdot(a, b, c)[source]
snowdrop.src.numeric.filters.kalman.dotdotinv(a, b, c)[source]

a * b * c^{-1}, where c is symmetric positive

snowdrop.src.numeric.filters.kalman.filter_step(G, covY, pred, yt)[source]

Filtering step of Kalman filter.

Parameters

G: (dy, dx) numpy array

mean of Y_t | X_t is G * X_t

covX: (dx, dx) numpy array

covariance of Y_t | X_t

pred: MeanAndCov object

predictive distribution at time t

Returns

pred: MeanAndCov object

filtering distribution at time t

logpyt: float

log density of Y_t | Y_{0:t-1}

snowdrop.src.numeric.filters.kalman.filter_step_asarray(G, covY, pred, yt)[source]

Filtering step of Kalman filter: array version.

Parameters

G: (dy, dx) numpy array

mean of Y_t | X_t is G * X_t

covX: (dx, dx) numpy array

covariance of Y_t | X_t

pred: MeanAndCov object

predictive distribution at time t

Returns

pred: MeanAndCov object

filtering distribution at time t

logpyt: float

log density of Y_t | Y_{0:t-1}

Note

This performs the filtering step for N distinctive predictive means: filt.mean should be a (N, dx) or (N) array; pred.mean in the output will have the same shape.

snowdrop.src.numeric.filters.kalman.predict_step(F, covX, filt)[source]

Predictive step of Kalman filter.

Parameters

F: (dx, dx) numpy array

Mean of X_t | X_{t-1} is F * X_{t-1}

covX: (dx, dx) numpy array

covariance of X_t | X_{t-1}

filt: MeanAndCov object

filtering distribution at time t-1

Returns

pred: MeanAndCov object

predictive distribution at time t

Note

filt.mean may either be of shape (dx,) or (N, dx); in the latter case N predictive steps are performed in parallel.

snowdrop.src.numeric.filters.kalman.smoother_step(F, filt, next_pred, next_smth)[source]

Smoothing step of Kalman filter/smoother.

Parameters

F: (dx, dx) numpy array

Mean of X_t | X_{t-1} is F * X_{t-1}

filt: MeanAndCov object

filtering distribution at time t

next_pred: MeanAndCov object

predictive distribution at time t+1

next_smth: MeanAndCov object

smoothing distribution at time t+1

Returns

smth: MeanAndCov object

smoothing distribution at time t

snowdrop.src.numeric.filters.particle module

Created on Sat Jun 26 09:23:51 2021

@author: Alexei Goumilevski

class snowdrop.src.numeric.filters.particle.NonlinearStateSpaceModel(x=None, y=None, parameters=None, Z=None, sigmaX=None, sigmaY=None, model=None, ind_non_missing=None, n_shocks=None, *args, **kwargs)[source]

Bases: StateSpaceModel

State space model class.

PX(t, xp)[source]

Distribution of X_t given X_{t-1} = xp (p=past).

PX0()[source]

Distribution of X_0.

PY(t, xp, x)[source]

Distribution of Y_t given X_t=x, and X_{t-1}=xp.

class snowdrop.src.numeric.filters.particle.ParticleGibbs(niter=10, verbose=0, ssm_cls=None, prior=None, data=None, theta0=None, Nx=100, fk_cls=None, regenerate_data=False, backward_step=False, store_x=False)[source]

Bases: ParticleGibbs

update_theta(theta, x)[source]

Update model parameters.

class snowdrop.src.numeric.filters.particle.StateSpaceModel(x=None, y=None, theta0=None, data=None, T=None, R=None, Z=None, Qm=None, H=None, C=None, parameters=None, model=None, ind_non_missing=None, *args, **kwargs)[source]

Bases: StateSpaceModel

State space model class.

PX(t, xp)[source]

Distribution of X_t given X_{t-1} = xp (p=past).

PX0()[source]

Distribution of X_0.

PY(t, xp, x)[source]

Distribution of Y_t given X_t=x, and X_{t-1}=xp.

snowdrop.src.numeric.filters.particle.getMatrices(obj, par)[source]

Compute matrices.

snowdrop.src.numeric.filters.py_kalman module

class snowdrop.src.numeric.filters.py_kalman.MyBiermanKalmanFilter(transition_matrices=None, observation_matrices=None, transition_covariance=None, observation_covariance=None, transition_offsets=None, observation_offsets=None, initial_state_mean=None, initial_state_covariance=None, random_state=None, em_vars=['transition_covariance', 'observation_covariance', 'initial_state_mean', 'initial_state_covariance'], n_dim_state=None, n_dim_obs=None)[source]

Bases: BiermanKalmanFilter

This module implements Bierman’s version of the Kalman Filter. In particular, the UDU’ decomposition of the covariance matrix is used instead of the full matrix, where U is upper triangular and D is diagonal.

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_matrix=None, transition_offset=None, transition_covariance=None, observation_matrix=None, observation_offset=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_matrixoptional, [n_dim_state, n_dim_state] array

state transition matrix from time t to t+1. If unspecified, self.transition_matrices will be used.

transition_offsetoptional, [n_dim_state] array

state offset for transition from time t to t+1. If unspecified, self.transition_offset will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_matrixoptional, [n_dim_obs, n_dim_state] array

observation matrix at time t+1. If unspecified, self.observation_matrices will be used.

observation_offsetoptional, [n_dim_obs] array

observation offset at time t+1. If unspecified, self.observation_offset will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
class snowdrop.src.numeric.filters.py_kalman.MyCholeskyKalmanFilter(transition_matrices=None, observation_matrices=None, transition_covariance=None, observation_covariance=None, transition_offsets=None, observation_offsets=None, initial_state_mean=None, initial_state_covariance=None, random_state=None, em_vars=['transition_covariance', 'observation_covariance', 'initial_state_mean', 'initial_state_covariance'], n_dim_state=None, n_dim_obs=None)[source]

Bases: CholeskyKalmanFilter

This module implements the Kalman Filter in “Square Root” form (Cholesky factorization).

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_matrix=None, transition_offset=None, transition_covariance=None, observation_matrix=None, observation_offset=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate.

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_matrixoptional, [n_dim_state, n_dim_state] array

state transition matrix from time t to t+1. If unspecified, self.transition_matrices will be used.

transition_offsetoptional, [n_dim_state] array

state offset for transition from time t to t+1. If unspecified, self.transition_offset will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_matrixoptional, [n_dim_obs, n_dim_state] array

observation matrix at time t+1. If unspecified, self.observation_matrices will be used.

observation_offsetoptional, [n_dim_obs] array

observation offset at time t+1. If unspecified, self.observation_offset will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
class snowdrop.src.numeric.filters.py_kalman.MyKalmanFilter(transition_matrices=None, observation_matrices=None, transition_covariance=None, observation_covariance=None, transition_offsets=None, observation_offsets=None, initial_state_mean=None, initial_state_covariance=None, random_state=None, em_vars=['transition_covariance', 'observation_covariance', 'initial_state_mean', 'initial_state_covariance'], n_dim_state=None, n_dim_obs=None)[source]

Bases: KalmanFilter

This module implements Kalman Filter.

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_matrix=None, transition_offset=None, transition_covariance=None, observation_matrix=None, observation_offset=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate.

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_matrixoptional, [n_dim_state, n_dim_state] array

state transition matrix from time t to t+1. If unspecified, self.transition_matrices will be used.

transition_offsetoptional, [n_dim_state] array

state offset for transition from time t to t+1. If unspecified, self.transition_offset will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_matrixoptional, [n_dim_obs, n_dim_state] array

observation matrix at time t+1. If unspecified, self.observation_matrices will be used.

observation_offsetoptional, [n_dim_obs] array

observation offset at time t+1. If unspecified, self.observation_offset will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
class snowdrop.src.numeric.filters.py_kalman.MyUnscentedKalmanFilter(transition_functions=None, observation_functions=None, transition_covariance=None, observation_covariance=None, initial_state_mean=None, initial_state_covariance=None, n_dim_state=None, n_dim_obs=None, random_state=None)[source]

Bases: UnscentedKalmanFilter

This module contains “Square Root” implementations to the Unscented Kalman Filter. Square Root implementations typically propagate the mean and Cholesky factorization of the covariance matrix in order to prevent numerical error. When possible, Square Root implementations should be preferred to their standard counterparts.

References:

  • Terejanu, G.A. Towards a Decision-Centric Framework for Uncertainty Propagation and Data Assimilation. 2010. Page 108.

  • Van Der Merwe, R. and Wan, E.A. The Square-Root Unscented Kalman Filter for State and Parameter-Estimation. 2001.

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_function=None, transition_covariance=None, observation_function=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_functionoptional, function

state transition function from time t to t+1. If unspecified, self.transition_functions will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_functionoptional, function

observation function at time t+1. If unspecified, self.observation_functions will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
unscented_filter_correct(observation_function, moments_pred, points_pred, observation, points_observation=None, sigma2_observation=None)[source]

Integrate new observation to correct state estimates

Args:
observation_functionfunction

function characterizing how the observation at time t+1 is generated

moments_pred[n_dim_state] Moments

mean and covariance of state at time t+1 given observations from time steps 0…t

points_pred[2*n_dim_state+1, n_dim_state] SigmaPoints

sigma points corresponding to moments_pred

observation[n_dim_state] array

observation at time t+1. If masked, treated as missing.

points_observation[2*n_dim_state, n_dim_obs] SigmaPoints

sigma points corresponding to predicted observation at time t+1 given observations from times 0…t, if available. If not, noise is assumed to be additive.

sigma_observation[n_dim_obs, n_dim_obs] array

covariance matrix corresponding to additive noise in observation at time t+1, if available. If missing, noise is assumed to be non-linear.

Returns:

moments_filt : [n_dim_state] Moments mean and covariance of state at time t+1 given observations from time steps 0…t+1

snowdrop.src.numeric.filters.pykalman module

class snowdrop.src.numeric.filters.pykalman.MyBiermanKalmanFilter(transition_matrices=None, observation_matrices=None, transition_covariance=None, observation_covariance=None, transition_offsets=None, observation_offsets=None, initial_state_mean=None, initial_state_covariance=None, random_state=None, em_vars=['transition_covariance', 'observation_covariance', 'initial_state_mean', 'initial_state_covariance'], n_dim_state=None, n_dim_obs=None)[source]

Bases: BiermanKalmanFilter

This module implements Bierman’s version of the Kalman Filter. In particular, the UDU’ decomposition of the covariance matrix is used instead of the full matrix, where U is upper triangular and D is diagonal.

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_matrix=None, transition_offset=None, transition_covariance=None, observation_matrix=None, observation_offset=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_matrixoptional, [n_dim_state, n_dim_state] array

state transition matrix from time t to t+1. If unspecified, self.transition_matrices will be used.

transition_offsetoptional, [n_dim_state] array

state offset for transition from time t to t+1. If unspecified, self.transition_offset will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_matrixoptional, [n_dim_obs, n_dim_state] array

observation matrix at time t+1. If unspecified, self.observation_matrices will be used.

observation_offsetoptional, [n_dim_obs] array

observation offset at time t+1. If unspecified, self.observation_offset will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
class snowdrop.src.numeric.filters.pykalman.MyCholeskyKalmanFilter(transition_matrices=None, observation_matrices=None, transition_covariance=None, observation_covariance=None, transition_offsets=None, observation_offsets=None, initial_state_mean=None, initial_state_covariance=None, random_state=None, em_vars=['transition_covariance', 'observation_covariance', 'initial_state_mean', 'initial_state_covariance'], n_dim_state=None, n_dim_obs=None)[source]

Bases: CholeskyKalmanFilter

This module implements the Kalman Filter in “Square Root” form (Cholesky factorization).

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_matrix=None, transition_offset=None, transition_covariance=None, observation_matrix=None, observation_offset=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate.

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_matrixoptional, [n_dim_state, n_dim_state] array

state transition matrix from time t to t+1. If unspecified, self.transition_matrices will be used.

transition_offsetoptional, [n_dim_state] array

state offset for transition from time t to t+1. If unspecified, self.transition_offset will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_matrixoptional, [n_dim_obs, n_dim_state] array

observation matrix at time t+1. If unspecified, self.observation_matrices will be used.

observation_offsetoptional, [n_dim_obs] array

observation offset at time t+1. If unspecified, self.observation_offset will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
class snowdrop.src.numeric.filters.pykalman.MyKalmanFilter(transition_matrices=None, observation_matrices=None, transition_covariance=None, observation_covariance=None, transition_offsets=None, observation_offsets=None, initial_state_mean=None, initial_state_covariance=None, random_state=None, em_vars=['transition_covariance', 'observation_covariance', 'initial_state_mean', 'initial_state_covariance'], n_dim_state=None, n_dim_obs=None)[source]

Bases: KalmanFilter

This module implements Kalman Filter.

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_matrix=None, transition_offset=None, transition_covariance=None, observation_matrix=None, observation_offset=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate.

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_matrixoptional, [n_dim_state, n_dim_state] array

state transition matrix from time t to t+1. If unspecified, self.transition_matrices will be used.

transition_offsetoptional, [n_dim_state] array

state offset for transition from time t to t+1. If unspecified, self.transition_offset will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_matrixoptional, [n_dim_obs, n_dim_state] array

observation matrix at time t+1. If unspecified, self.observation_matrices will be used.

observation_offsetoptional, [n_dim_obs] array

observation offset at time t+1. If unspecified, self.observation_offset will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
class snowdrop.src.numeric.filters.pykalman.MyUnscentedKalmanFilter(transition_functions=None, observation_functions=None, transition_covariance=None, observation_covariance=None, initial_state_mean=None, initial_state_covariance=None, n_dim_state=None, n_dim_obs=None, random_state=None)[source]

Bases: UnscentedKalmanFilter

This module contains “Square Root” implementations to the Unscented Kalman Filter. Square Root implementations typically propagate the mean and Cholesky factorization of the covariance matrix in order to prevent numerical error. When possible, Square Root implementations should be preferred to their standard counterparts.

References:

  • Terejanu, G.A. Towards a Decision-Centric Framework for Uncertainty Propagation and Data Assimilation. 2010. Page 108.

  • Van Der Merwe, R. and Wan, E.A. The Square-Root Unscented Kalman Filter for State and Parameter-Estimation. 2001.

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_function=None, transition_covariance=None, observation_function=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_functionoptional, function

state transition function from time t to t+1. If unspecified, self.transition_functions will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_functionoptional, function

observation function at time t+1. If unspecified, self.observation_functions will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

index = None
unscented_filter_correct(observation_function, moments_pred, points_pred, observation, points_observation=None, sigma2_observation=None)[source]

Integrate new observation to correct state estimates

Args:
observation_functionfunction

function characterizing how the observation at time t+1 is generated

moments_pred[n_dim_state] Moments

mean and covariance of state at time t+1 given observations from time steps 0…t

points_pred[2*n_dim_state+1, n_dim_state] SigmaPoints

sigma points corresponding to moments_pred

observation[n_dim_state] array

observation at time t+1. If masked, treated as missing.

points_observation[2*n_dim_state, n_dim_obs] SigmaPoints

sigma points corresponding to predicted observation at time t+1 given observations from times 0…t, if available. If not, noise is assumed to be additive.

sigma_observation[n_dim_obs, n_dim_obs] array

covariance matrix corresponding to additive noise in observation at time t+1, if available. If missing, noise is assumed to be non-linear.

Returns:

moments_filt : [n_dim_state] Moments mean and covariance of state at time t+1 given observations from time steps 0…t+1

snowdrop.src.numeric.filters.unscented module

Inference for Non-Linear Gaussian Systems

This module contains the Unscented Kalman Filter (Wan, van der Merwe 2000) for state estimation in systems with non-Gaussian noise and non-linear dynamics.

class snowdrop.src.numeric.filters.unscented.AdditiveUnscentedKalmanFilter(transition_functions=None, observation_functions=None, transition_covariance=None, observation_covariance=None, initial_state_mean=None, initial_state_covariance=None, n_dim_state=None, n_dim_obs=None, random_state=None)[source]

Bases: UnscentedMixin

Implements the Unscented Kalman Filter with additive noise. Observations are assumed to be generated from the following process,

\[\begin{split}x_0 &\sim \text{Normal}(\mu_0, \Sigma_0) \\ x_{t+1} &= f_t(x_t) + \text{Normal}(0, Q) \\ z_{t} &= g_t(x_t) + \text{Normal}(0, R)\end{split}\]

While less general the general-noise Unscented Kalman Filter, the Additive version is more computationally efficient with complexity \(O(Tn^3)\) where \(T\) is the number of time steps and \(n\) is the size of the state space.

Args:
transition_functionsfunction or [n_timesteps-1] array of functions

transition_functions[t] is a function of the state at time t and produces the state at time t+1. Also known as \(f_t\).

observation_functionsfunction or [n_timesteps] array of functions

observation_functions[t] is a function of the state at time t and produces the observation at time t. Also known as \(g_t\).

transition_covariance[n_dim_state, n_dim_state] array

transition noise covariance matrix. Also known as \(Q\).

observation_covariance[n_dim_obs, n_dim_obs] array

observation noise covariance matrix. Also known as \(R\).

initial_state_mean[n_dim_state] array

mean of initial state distribution. Also known as \(\mu_0\).

initial_state_covariance[n_dim_state, n_dim_state] array

covariance of initial state distribution. Also known as \(\Sigma_0\).

n_dim_state: optional, integer

the dimensionality of the state space. Only meaningful when you do not specify initial values for transition_covariance, or initial_state_mean, initial_state_covariance.

n_dim_obs: optional, integer

the dimensionality of the observation space. Only meaningful when you do not specify initial values for observation_covariance.

random_stateoptional, int or RandomState

seed for random sample generation

filter(Z)[source]

Run Unscented Kalman Filter

Args:
Z[n_timesteps, n_dim_state] array

Z[t] = observation at time t. If Z is a masked array and any of Z[t]’s elements are masked, the observation is assumed missing and ignored.

Returns:
filtered_state_means[n_timesteps, n_dim_state] array

filtered_state_means[t] = mean of state distribution at time t given observations from times [0, t]

filtered_state_covariances[n_timesteps, n_dim_state, n_dim_state] array

filtered_state_covariances[t] = covariance of state distribution at time t given observations from times [0, t]

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_function=None, transition_covariance=None, observation_function=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_functionoptional, function

state transition function from time t to t+1. If unspecified, self.transition_functions will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_functionoptional, function

observation function at time t+1. If unspecified, self.observation_functions will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

sample(n_timesteps, initial_state=None, random_state=None)[source]

Sample from model defined by the Unscented Kalman Filter

Args:
n_timestepsint

number of time steps

initial_stateoptional, [n_dim_state] array

initial state. If unspecified, will be sampled from initial state distribution.

smooth(Z)[source]

Run Unscented Kalman Smoother

Args:
Z[n_timesteps, n_dim_state] array

Z[t] = observation at time t. If Z is a masked array and any of Z[t]’s elements are masked, the observation is assumed missing and ignored.

Returns:
smoothed_state_means[n_timesteps, n_dim_state] array

filtered_state_means[t] = mean of state distribution at time t given observations from times [0, n_timesteps-1]

smoothed_state_covariances[n_timesteps, n_dim_state, n_dim_state] array

filtered_state_covariances[t] = covariance of state distribution at time t given observations from times [0, n_timesteps-1]

class snowdrop.src.numeric.filters.unscented.Moments(mean, covariance)

Bases: tuple

covariance

Alias for field number 1

mean

Alias for field number 0

class snowdrop.src.numeric.filters.unscented.SigmaPoints(points, weights_mean, weights_covariance)

Bases: tuple

points

Alias for field number 0

weights_covariance

Alias for field number 2

weights_mean

Alias for field number 1

class snowdrop.src.numeric.filters.unscented.UnscentedKalmanFilter(transition_functions=None, observation_functions=None, transition_covariance=None, observation_covariance=None, initial_state_mean=None, initial_state_covariance=None, n_dim_state=None, n_dim_obs=None, random_state=None)[source]

Bases: UnscentedMixin

Implements the General (aka Augmented) Unscented Kalman Filter governed by the following equations,

\[x_0 = \sim \text{Normal}(\mu_0, \Sigma_0) x_{t+1} = f_t(x_t, \text{Normal}(0, Q)) z_{t} = g_t(x_t, \text{Normal}(0, R))\]

Notice that although the input noise to the state transition equation and the observation equation are both normally distributed, any non-linear transformation may be applied afterwards. This allows for greater generality, but at the expense of computational complexity. The complexity of UnscentedKalmanFilter.filter() is \(O(T(2n+m)^3)\) where \(T\) is the number of time steps, \(n\) is the size of the state space, and \(m\) is the size of the observation space.

If your noise is simply additive, consider using the AdditiveUnscentedKalmanFilter

Args:
transition_functionsfunction or [n_timesteps-1] array of functions

transition_functions[t] is a function of the state and the transition noise at time t and produces the state at time t+1. Also known as \(f_t\).

observation_functionsfunction or [n_timesteps] array of functions

observation_functions[t] is a function of the state and the observation noise at time t and produces the observation at time t. Also known as \(g_t\).

transition_covariance[n_dim_state, n_dim_state] array

transition noise covariance matrix. Also known as \(Q\).

observation_covariance[n_dim_obs, n_dim_obs] array

observation noise covariance matrix. Also known as \(R\).

initial_state_mean[n_dim_state] array

mean of initial state distribution. Also known as \(\mu_0\)

initial_state_covariance[n_dim_state, n_dim_state] array

covariance of initial state distribution. Also known as \(\Sigma_0\)

n_dim_state: optional, integer

the dimensionality of the state space. Only meaningful when you do not specify initial values for transition_covariance, or initial_state_mean, initial_state_covariance.

n_dim_obs: optional, integer

the dimensionality of the observation space. Only meaningful when you do not specify initial values for observation_covariance.

random_stateoptional, int or RandomState

seed for random sample generation

filter(Z)[source]

Run Unscented Kalman Filter

Args:
Z[n_timesteps, n_dim_state] array

Z[t] = observation at time t. If Z is a masked array and any of Z[t]’s elements are masked, the observation is assumed missing and ignored.

Returns:
filtered_state_means[n_timesteps, n_dim_state] array

filtered_state_means[t] = mean of state distribution at time t given observations from times [0, t]

filtered_state_covariances[n_timesteps, n_dim_state, n_dim_state] array

filtered_state_covariances[t] = covariance of state distribution at time t given observations from times [0, t]

filter_update(filtered_state_mean, filtered_state_covariance, observation=None, transition_function=None, transition_covariance=None, observation_function=None, observation_covariance=None)[source]

Update a Kalman Filter state estimate

Perform a one-step update to estimate the state at time \(t+1\) give an observation at time \(t+1\) and the previous estimate for time \(t\) given observations from times \([0...t]\). This method is useful if one wants to track an object with streaming observations.

Args:
filtered_state_mean[n_dim_state] array

mean estimate for state at time t given observations from times [1…t]

filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t given observations from times [1…t]

observation[n_dim_obs] array or None

observation from time t+1. If observation is a masked array and any of observation’s components are masked or if observation is None, then observation will be treated as a missing observation.

transition_functionoptional, function

state transition function from time t to t+1. If unspecified, self.transition_functions will be used.

transition_covarianceoptional, [n_dim_state, n_dim_state] array

state transition covariance from time t to t+1. If unspecified, self.transition_covariance will be used.

observation_functionoptional, function

observation function at time t+1. If unspecified, self.observation_functions will be used.

observation_covarianceoptional, [n_dim_obs, n_dim_obs] array

observation covariance at time t+1. If unspecified, self.observation_covariance will be used.

Returns:
next_filtered_state_mean[n_dim_state] array

mean estimate for state at time t+1 given observations from times [1…t+1]

next_filtered_state_covariance[n_dim_state, n_dim_state] array

covariance of estimate for state at time t+1 given observations from times [1…t+1]

sample(n_timesteps, initial_state=None, random_state=None)[source]

Sample from model defined by the Unscented Kalman Filter

Args:
n_timestepsint

number of time steps

initial_stateoptional, [n_dim_state] array

initial state. If unspecified, will be sampled from initial state distribution.

random_stateoptional, int or Random

random number generator

smooth(Z)[source]

Run Unscented Kalman Smoother

Args:
Z[n_timesteps, n_dim_state] array

Z[t] = observation at time t. If Z is a masked array and any of Z[t]’s elements are masked, the observation is assumed missing and ignored.

Returns:
smoothed_state_means[n_timesteps, n_dim_state] array

filtered_state_means[t] = mean of state distribution at time t given observations from times [0, n_timesteps-1]

smoothed_state_covariances[n_timesteps, n_dim_state, n_dim_state] array

filtered_state_covariances[t] = covariance of state distribution at time t given observations from times [0, n_timesteps-1]

class snowdrop.src.numeric.filters.unscented.UnscentedMixin(transition_functions=None, observation_functions=None, transition_covariance=None, observation_covariance=None, initial_state_mean=None, initial_state_covariance=None, n_dim_state=None, n_dim_obs=None, random_state=None)[source]

Bases: object

Methods shared by all Unscented Kalman Filter implementations.

snowdrop.src.numeric.filters.unscented.additive_unscented_filter(mu_0, sigma_0, f, g, Q, R, Z)[source]

Apply the Unscented Kalman Filter with additive noise

Args:
mu_0[n_dim_state] array

mean of initial state distribution

sigma_0[n_dim_state, n_dim_state] array

covariance of initial state distribution

ffunction or [T-1] array of functions

state transition function(s). Takes in an the current state and outputs the next.

gfunction or [T] array of functions

observation function(s). Takes in the current state and outputs the current observation.

Q[n_dim_state, n_dim_state] array

transition covariance matrix

R[n_dim_state, n_dim_state] array

observation covariance matrix

Returns:
mu_filt[T, n_dim_state] array

mu_filt[t] = mean of state at time t given observations from times [0, t]

sigma_filt[T, n_dim_state, n_dim_state] array

sigma_filt[t] = covariance of state at time t given observations from times [0, t]

snowdrop.src.numeric.filters.unscented.additive_unscented_smoother(mu_filt, sigma_filt, f, Q)[source]

Apply the Unscented Kalman Filter assuming additiven noise

Args:
mu_filt[T, n_dim_state] array

mu_filt[t] = mean of state at time t given observations from times [0, t]

sigma_filt[T, n_dim_state, n_dim_state] array

sigma_filt[t] = covariance of state at time t given observations from times [0, t]

ffunction or [T-1] array of functions

state transition function(s). Takes in an the current state and outputs the next.

Q[n_dim_state, n_dim_state] array

transition covariance matrix

Returns:
mu_smooth[T, n_dim_state] array

mu_smooth[t] = mean of state at time t given observations from times [0, T-1]

sigma_smooth[T, n_dim_state, n_dim_state] array

sigma_smooth[t] = covariance of state at time t given observations from times [0, T-1]

snowdrop.src.numeric.filters.unscented.augmented_points(momentses)[source]

Calculate sigma points for augmented UKF

Args:
momentseslist of Moments

means and covariances for multiple multivariate normals

Returns:
pointseslist of Points

sigma points for each element of momentses

snowdrop.src.numeric.filters.unscented.augmented_unscented_filter(mu_0, sigma_0, f, g, Q, R, Z)[source]

Apply the Unscented Kalman Filter with arbitrary noise

Args:
mu_0[n_dim_state] array

mean of initial state distribution

sigma_0[n_dim_state, n_dim_state] array

covariance of initial state distribution

ffunction or [T-1] array of functions

state transition function(s). Takes in an the current state and the process noise and outputs the next state.

gfunction or [T] array of functions

observation function(s). Takes in the current state and outputs the current observation.

Q[n_dim_state, n_dim_state] array

transition covariance matrix

R[n_dim_state, n_dim_state] array

observation covariance matrix

Returns:
mu_filt[T, n_dim_state] array

mu_filt[t] = mean of state at time t given observations from times [0, t]

sigma_filt[T, n_dim_state, n_dim_state] array

sigma_filt[t] = covariance of state at time t given observations from times [0, t]

snowdrop.src.numeric.filters.unscented.augmented_unscented_filter_points(mean_state, covariance_state, covariance_transition, covariance_observation)[source]

Extract sigma points using augmented state representation

Primarily used as a pre-processing step before predicting and updating in the Augmented UKF.

Args:
mean_state[n_dim_state] array

mean of state at time t given observations from time steps 0…t

covariance_state[n_dim_state, n_dim_state] array

covariance of state at time t given observations from time steps 0…t

covariance_transition[n_dim_state, n_dim_state] array

covariance of zero-mean noise resulting from transitioning from time step t to t+1

covariance_observation[n_dim_obs, n_dim_obs] array

covariance of zero-mean noise resulting from observation state at time t+1

Returns:
points_state[2 * n_dim_state + 1, n_dim_state] SigmaPoints

sigma points for state at time t

points_transition[2 * n_dim_state + 1, n_dim_state] SigmaPoints

sigma points for transition noise between time t and t+1

points_observation[2 * n_dim_state + 1, n_dim_obs] SigmaPoints

sigma points for observation noise at time step t+1

snowdrop.src.numeric.filters.unscented.augmented_unscented_smoother(mu_filt, sigma_filt, f, Q)[source]

Apply the Unscented Kalman Smoother with arbitrary noise

Args:
mu_filt[T, n_dim_state] array

mu_filt[t] = mean of state at time t given observations from times [0, t]

sigma_filt[T, n_dim_state, n_dim_state] array

sigma_filt[t] = covariance of state at time t given observations from times [0, t]

ffunction or [T-1] array of functions

state transition function(s). Takes in an the current state and the process noise and outputs the next state.

Q[n_dim_state, n_dim_state] array

transition covariance matrix

Returns:
mu_smooth[T, n_dim_state] array

mu_smooth[t] = mean of state at time t given observations from times [0, T-1]

sigma_smooth[T, n_dim_state, n_dim_state] array

sigma_smooth[t] = covariance of state at time t given observations from times [0, T-1]

snowdrop.src.numeric.filters.unscented.moments2points(moments, alpha=None, beta=None, kappa=None)[source]

Calculate “sigma points” used in Unscented Kalman Filter

Args:
moments[n_dim] Moments object

mean and covariance of a multivariate normal

alphafloat

Spread of the sigma points. Typically 1e-3.

betafloat

Used to “incorporate prior knowledge of the distribution of the state”. 2 is optimal is the state is normally distributed.

kappafloat

a parameter

Returns:
points[2*n_dim+1, n_dim] SigmaPoints

sigma points and associated weights

snowdrop.src.numeric.filters.unscented.points2moments(points, sigma_noise=None)[source]

Calculate estimated mean and covariance of sigma points

Args:
points[2 * n_dim_state + 1, n_dim_state] SigmaPoints

SigmaPoints object containing points and weights

sigma_noise[n_dim_state, n_dim_state] array

additive noise covariance matrix, if any

Returns:
momentsMoments object of size [n_dim_state]

Mean and covariance estimated using points

snowdrop.src.numeric.filters.unscented.unscented_correct(cross_sigma, moments_pred, obs_moments_pred, z)[source]

Correct predicted state estimates with an observation

Args:
cross_sigma[n_dim_state, n_dim_obs] array

cross-covariance between the state at time t given all observations from timesteps [0, t-1] and the observation at time t

moments_pred[n_dim_state] Moments

mean and covariance of state at time t given observations from timesteps [0, t-1]

obs_moments_pred[n_dim_obs] Moments

mean and covariance of observation at time t given observations from times [0, t-1]

z[n_dim_obs] array

observation at time t

Returns:
moments_filt[n_dim_state] Moments

mean and covariance of state at time t given observations from time steps [0, t]

snowdrop.src.numeric.filters.unscented.unscented_filter_correct(observation_function, moments_pred, points_pred, observation, points_observation=None, sigma_observation=None)[source]

Integrate new observation to correct state estimates

Args:
observation_functionfunction

function characterizing how the observation at time t+1 is generated

moments_pred[n_dim_state] Moments

mean and covariance of state at time t+1 given observations from time steps 0…t

points_pred[2*n_dim_state+1, n_dim_state] SigmaPoints

sigma points corresponding to moments_pred

observation[n_dim_state] array

observation at time t+1. If masked, treated as missing.

points_observation[2*n_dim_state, n_dim_obs] SigmaPoints

sigma points corresponding to predicted observation at time t+1 given observations from times 0…t, if available. If not, noise is assumed to be additive.

sigma_observation[n_dim_obs, n_dim_obs] array

covariance matrix corresponding to additive noise in observation at time t+1, if available. If missing, noise is assumed to be non-linear.

Returns:
moments_filt[n_dim_state] Moments

mean and covariance of state at time t+1 given observations from time steps 0…t+1

snowdrop.src.numeric.filters.unscented.unscented_filter_predict(transition_function, points_state, points_transition=None, sigma_transition=None)[source]

Predict next state distribution

Using the sigma points representing the state at time t given observations from time steps 0…t, calculate the predicted mean, covariance, and sigma points for the state at time t+1.

Args:
transition_functionfunction

function describing how the state changes between times t and t+1

points_state[2*n_dim_state+1, n_dim_state] SigmaPoints

sigma points corresponding to the state at time step t given observations from time steps 0…t

points_transition[2*n_dim_state+1, n_dim_state] SigmaPoints

sigma points corresponding to the noise in transitioning from time step t to t+1, if available. If not, assumes that noise is additive

sigma_transition[n_dim_state, n_dim_state] array

covariance corresponding to additive noise in transitioning from time step t to t+1, if available. If not, assumes noise is not additive.

Returns:
points_pred[2*n_dim_state+1, n_dim_state] SigmaPoints

sigma points corresponding to state at time step t+1 given observations from time steps 0…t. These points have not been “standardized” by the unscented transform yet.

moments_pred[n_dim_state] Moments

mean and covariance corresponding to time step t+1 given observations from time steps 0…t

snowdrop.src.numeric.filters.unscented.unscented_transform(points, f=None, points_noise=None, sigma_noise=None)[source]

Apply the Unscented Transform to a set of points

Apply f to points (with secondary argument points_noise, if available), then approximate the resulting mean and covariance. If sigma_noise is available, treat it as additional variance due to additive noise.

Args:
points[n_points, n_dim_state] SigmaPoints

points to pass into f’s first argument and associated weights if f is defined. If f is unavailable, then f is assumed to be the identity function.

f[n_dim_state, n_dim_state_noise] -> [n_dim_state] function

transition function from time t to time t+1, if available.

points_noise[n_points, n_dim_state_noise] array

points to pass into f’s second argument, if any

sigma_noise[n_dim_state, n_dim_state] array

covariance matrix for additive noise, if any

Returns:
points_pred[n_points, n_dim_state] SigmaPoints

points transformed by f with same weights

moments_pred[n_dim_state] Moments

moments associated with points_pred

snowdrop.src.numeric.filters.utils module

Compute initial error covariance matrices Pstar and Pinf.

The initial values of error covariance matrix P are assumed diffusive. During the first d diffusive periods matrix P is decomposed into Pstar and Pinf and separate equations for these matrices are applied. After the d diffusive periods the standard Kalman filter is used and only one error covariance matrix P is used. The initial condition for P(d) is Pstar(d).

Translated from Dynare Matlab code to Python by A.G.

snowdrop.src.numeric.filters.utils.compute_Pinf_Pstar(T, R, Q, N=0, n_shocks=-1, mf=None, qz_criterium=1.000001)[source]

Compute of Pinf and Pstar for Durbin and Koopman diffuse filter.

Kitagawa transformation of state space system with a quasi-triangular transition matrix with unit roots at the top, but excluding zero columns of the transition matrix.

For references please see G. Kitagawa, “An algorithm for solving the matrix equation”, International Journal of Control, Volume 25, 1977 - Issue 5

The transformed state space is y = [ss z x], were:

s: static variables (zero columns of T)

z: unit roots

x: stable roots

ss: s - z = stationarized static variables

Args:
T: double

matrix of transition equations

R: double

matrix of structural shock effects

Q: double

matrix of covariance of structural shocks

N: int

shocks maximum lead number minus minimum lag number

n_shocks: int

number pf shocks

mf: integer

vector of indices of observed variables in state vector

qz_criterium: double

numerical criterium for unit roots

Returns:
Pinf: double

matrix of covariance initialization for nonstationary part

Pstar: double

matrix of covariance of stationary part

Algorithm:

Real Schur transformation of transition equation

Note

Translated from Dynare version 4.5.1 to Python by AG

snowdrop.src.numeric.filters.utils.getCovarianceMatrix(Qm, shocks, equations, n)[source]

Return error covariance matrix of endogenous variables.

Parameters

Qmnumpy 2D ndarray.

Error covariance matrix.

shockslist.

Shocks names.

equationslist.

Measuremen equations.

nint.

Nimber of variables.

Returns

Qnumpy 2D array.

Full error covariance matrix of endogenous variables.

snowdrop.src.numeric.filters.utils.getEndogVarInMeasEqs(variables, measurement_equations)[source]

Retrun list of endogenous variables that are present in measurement equations.

Parameters

variableslist.

Names of variables.

measurement_equationslist.

Measurement equations.

Returns

ind_varlist.

Indices of endogenous variables.

var_measlist.

Names of endogenous variables.

snowdrop.src.numeric.filters.utils.getMissingInitialVariables(model, ind, x0)[source]

Return initial values of endogenous variables satisfying model equations. Parameters: :param model: Model object. :type model: `Model’. :param ind: Indices of variables with fixed values. :type ind: list. :param x0: Initial guess. :type x0: numpy ndarray.

snowdrop.src.numeric.filters.utils.getSteadyStateCovarianceMatrix(T, R, Qm, Hm, Z, n, Nd, n_shocks)[source]

Solve Lyapunov equation for stable part of error covariance matrix.

Predict:

\[P_{t|t-1} = T_{t} * P_{t-1|t-1} * F'_{t} + Q_{t}\]

For details, see https://en.wikipedia.org/wiki/Kalman_filter

Args:

T: is the state-transition matrix,

R: is the shock matrix,

Qm: is the covariance matrix of state variables (endogenous variables),

Hm: is the covariance matrix of space variables (measurement variables)

Z: is the observation matrix,

n: is the number of endogenous variables,

Nd: is the shocks maximum lead number minus minimum lag number,

n_shocks: is the number of shocks,

iterate: is the boolean variable. If True then iterative method is used, otherwise Lyapunov equation solver is applied.

For details, see https://en.wikipedia.org/wiki/Kalman_filter

snowdrop.src.numeric.filters.utils.getTimeShift(model)[source]

Get time shift between start of simulation and start of filtration.

snowdrop.src.numeric.filters.utils.sorter(x)[source]

Sort eigen values.

Module contents