In [1], the author shows an efficient way of implementing the forward and backward filter using matrices. One can also implement this using filtfilt command in MATLAB. However, I am interested in implementing this using matrices.


Any digital filter can be represented in its state-space format as follows: \begin{align} x(t+1) &= F x(t) + G u(t) \notag \\ y(t) &= H x(t) + D u(t) \notag \end{align}

This filter can be expressed as \begin{align} \mathcal{Y} = \mathcal{H} U + \mathcal{O}x_0 \end{align} where $\mathcal{Y} = [y_0, y_1,\ldots,y_{N-1}]^T$ is the vector of outputs, $U = [u_0,u_1,\ldots,u_{N-1}]^T$ is the vector of inputs, $x_0$ is the initial state, \begin{align} \mathcal{H} = \begin{bmatrix} D & 0 & \ldots & 0 \\ HG & D & 0 & \vdots \\ \vdots & & \ddots & \\ HF^{N-2}G & \ldots & HG & D \end{bmatrix}, \end{align} and \begin{align} \mathcal{O} = \begin{bmatrix} H \\ HF \\ \vdots \\ HF^{N-1} \end{bmatrix}. \end{align}

Problem Statement

In forward-backward filtering proposed in [1], the author mentions that the forward and backward filters are different (generally speaking). However, if we look at the implementation of the forward-backward filter derivation, the same filter is used, i.e., only the time-reversal of the input and outputs are performed but the filter transfer function does not change. What should be an efficient representation of the state-space filter while applying the forward filter and backward filter, i.e., determine $(F_f, G_f, H_f, D_f)$ and $(F_b, G_b, H_b, D_b)$.


[1] F. Gustafsson, "Determining the initial states in forward-backward filtering," in IEEE Transactions on Signal Processing, vol. 44, no. 4, pp. 988-992, Apr 1996. http://www.diva-portal.org/smash/get/diva2:315708/FULLTEXT02


1 Answer 1


The backwards filter, as the name suggests filters the input signal with the same filter $(F,G,H,D)$ backwards in time

$$ \begin{align} x(t-1) &= F\,x(t) + G\,u(t) \\ y(t) &= H\,x(t) + D\,u(t) \end{align} \tag{1} $$

But in order to define a filter using a state space model, the update law does have to go forwards in time, which can be done by solving the first equation for $x(t)$

$$ x(t) = F^{-1}\,\left(x(t-1) - G\,u(t)\right) \tag{2} $$

however $x(t)$ can only a function of $x(t-1)$ and $u(t-1)$. This can be solved by using a new definition for the state $z(t) = x(t) + F^{-1}\,G\,u(t)$, so

$$ z(t) = x(t) + F^{-1}\,G\,u(t) = F^{-1}\,x(t-1) = F^{-1}\,\left(z(t-1) - F^{-1}\,G\,u(t-1)\right) \tag{3} $$

which can be rewritten in the following more standard form

$$ z(t+1) = F^{-1}\,z(t) - F^{-2}\,G\,u(t). \tag{4} $$

Substituting the definition of $z(t)$ in to the expression for $y(t)$ gives

$$ y(t) = H\,x(t) + D\,u(t) = H\,\left(z(t) - F^{-1}\,G\,u(t) \right) + D\,u(t) \tag{5} $$

which can be simplified down to

$$ y(t) = H\,z(t) + \left(D - H\,F^{-1}\,G\right)\,u(t). \tag{6} $$

Combining equation $(4)$ and $(6)$ gives the following (forwards) state space model, which is equivalent to the (backwards) state space model in equation $(1)$

$$ \begin{align} z(t+1) &= F^{-1}\,z(t) - F^{-2}\,G\,u(t) \\ y(t) &= H\,z(t) + \left(D - H\,F^{-1}\,G\right)\,u(t) \end{align} \tag{7} $$

So the backwards filter based on the forwards filter $(F,G,H,D)$ can be written as $(F^{-1},-F^{-2}\,G,H,D - H\,F^{-1}\,G)$.

  • $\begingroup$ The problem with backward filter obtained in (7) is that it is not stable? For instance, if we obtain a Hankel matrix using the state-space representation in (7) and obtain the outputs, i.e, $\mathcal{Y} = \mathcal{H}U$, then the outputs seem to be diverging. $\endgroup$
    – Maxtron
    May 16, 2018 at 15:55
  • $\begingroup$ I am not sure if $z(t) = x(t) + F^{-1}Gu(t)$ is the right way to represent the state. We are starting with the backward filter $(F,G,H,D)$ and determining the forward filter $(F^{-1}, -F^{-2}G, H, D-HF^{-1}G)$. The new states of the original forward filter were only dependent on previous inputs and state, but in case of (7), the new states are dependent on current input $u(t+1)$ (appears in $z(t+1)$) and the previous states and inputs. $\endgroup$
    – Maxtron
    May 16, 2018 at 16:05
  • 1
    $\begingroup$ @Maxtron1 If the forward filter is stable, then indeed the backwards filter will be unstable. In order to avoid this it would be better to actually simulate backwards in time, which can be done applying the input backwards, by starting at the last sample and end at the first sample. The output of this result should then be reversed in time as well. $\endgroup$
    – fibonatic
    May 16, 2018 at 16:43
  • $\begingroup$ There seems to be some problem with the reconstruction of the original signal using forward-backward or backward-forward method when the filters are designed as matrices. In particular, I am trying to implement the two equations that appear at the end of Section 2 in reference [1]. I compute the Hankel matrix ($\mathcal{H}$) and observability matrix ($\mathcal{O}$) for the forward and backward filters, respectively, based on your derivation. $\endgroup$
    – Maxtron
    May 16, 2018 at 17:29

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.