2
$\begingroup$

I would appreciate it very much if someone would be able to provide some clarity, help or comment on this problem. I have been reading several papers on time series identification such as

and system identification

among several others.

The question is, given an unstable causal system, is it possible to turn it into an anticausal but stable system in order to simulate it? (assuming that I have the complete input sequence).

I think it might have to do with the time reversal property of the Z-Transform, i.e.

$\mathcal{Z} \left\{ x[n] \right \} = X(z)$ with it's ROC given by $r_2<|z|<r_1$, then

$\mathcal{Z} \left\{ x[-n] \right \} = X(z^{-1})$ with $\frac{1}{r_1}<|z|<\frac{1}{r_2}$

however i can't figure it out.

In the Forssell and Hjalmarsson paper says that backward filtering is common in off-line signal processing applications, like the zero-phase filtering (filtfilt in Matlab) and the implementation of non-causal Wiener filters. The relevant (at least for my question) equations are the (31), (32) and (33).

To illustrate the problem, please consider the following system given by

$y[n] = \frac{B(z)}{A(z)} u[n] = \frac{(z-0.1)}{(z-0.5)(z-5)} u[n]$;

assume that the complete input sequence $\left\{ u[n] \right\}^{N}_{n=0}$ is known (available).

Then, it is possible to separate the system as follows

$y[n] = \frac{B(z)}{A_s(z)} \frac{1}{A_a(z)} u[n],$

where $A_s(z)$ and $A_a(z)$ represents the stable and anti-stable monic parts of $A(z)$ respectively. Thus we can obtain the causal part of $y_{\text{C}}[n]$ by computing

$y_{\text{C}}[n] = \frac{B(z)}{A_s(z)} u[n] = \frac{(z-0.1)}{(z-0.5)} u[n]$

followed by

$y[n] = \frac{1}{A_a(z)} y_{\text{C}}[n] = \frac{1}{(z-5)} y_{\text{C}}[n]$

So here is my problem. The approach given in Forssell and Hjalmarsson is the following

in order to obtain $y[n]$ backward-filter the reversed sequence given by $\left\{ y_{\text{C}}^R[n] \right\}_{n=0}^{N}$ through the filter $\frac{1}{A_a(z)}$ as follows

$y^R[n] = \frac{1}{A_a(z)} y_{\text{C}}^R[n]$

the desired signal $y[n]$ is then finally obtained by reversing the sequence $\left\{ y^R[n] \right\}_{n=0}^{N}$

I apologize in advance for the length of the question and the orthographic mistakes. Any comments, suggestions and help will be greatly appreciated. Thanks

$\endgroup$
1
$\begingroup$

An unstable system can be stablized by all-pass decomposition for exampe, but simply constructing the anti-causal impulse respone extension will not provide the same output; so I think it's not possible, but I'm not rigorous at this point. Following just shows why.

Assume that your unstable and causal IIR filter has the impulse response $h_+[n]$ and by evaluating it into $ n < 0$ domain you construct a stable but anti-causal impulse response denoted by $h_{-}[n]$.

Your true output given by the unstable filtering is $$ y[n] = h_{+}[n] \star x[n]$$

Now you construct alternate output from a reversed input and the anti-causal filter as $$ y_{-}[n] = h_{-}[n] \star x[-n] $$

And by reversing $y_{-}[n]$ , do you get $y[n]$ ?

No as seen:

$$y_{-}[-n] = h_{-}[-n] \star x[n] \neq h_+[n] \star x[n] $$

that's to say unless $h_{-}[-n] = h_+[n]$ , they wont be the same. And indeed the reversed anti-causal impuse response will not be equal to the causal one, since one is stable and the other is not.

$\endgroup$
1
$\begingroup$

Here is one way to think about it: Let's say you have an unstable transfer function, which is unstable because of a single pole outside the unit circle. We can write this as $$H_0(z)=H_{stable}(z) \cdot \frac{1}{1-z^{-1} \cdot p}$$

That's a cascade of the stable and unstable part. Now we can multiply this with a unity filter that has a pole and zero at the inverse location of that pole. We get

$$H_0(z)=H_{stable}(z) \cdot \frac{1}{1-z^{-1} \cdot p} \cdot \frac{1-z^{-1}/ p}{1-z^{-1}/ p}$$

Reordering things a little we get: $$H_0(z)= \frac{H_{stable}(z)}{1-z^{-1} / p} \cdot \frac{1-z^{-1}/ p}{1-z^{-1} \cdot p} = H_1(z) \cdot H_2(z)$$

$H_1(z)$ is now stable since we "flipped" the pole into the unit circle. You can actually show that $H_1(z)$ has the same magnitude as the original transfer function $H_0(z)$, just a different phase and it's also causal.

We have isolated the "unstable" pole in $H_2(z)$ and paired it with a zero at the inverse location. That makes it an allpass filter. The inverse of $H_2(z)$ is actually a perfectly well behaved allpass filter.

So what happens if you inverse and allpass filter? Since the magnitude is unity, it remains unaffected and you only flip the phase, , i.e. $$H_a(z) = e^{j \cdot \phi(\omega)}, \frac{1}{H_a(z)}= e^{-j \cdot \phi(\omega)}$$

Flipping the phase (without changing the magnitude) in the frequency domain corresponds to time reversal in the time domain. Hence: the impulse response of an inverse allpass filter is the time reverse of the original allpass filter.

Putting this all together we get: A unstable system can be expressed as the cascade of a stable system and some time reversed allpass filters (one for each pole outside the unit circle).

The allpass filter has an infinite length impulse response, so it is not possible to exactly time reverse it. In practice, you would need to find a delay that captures "enough" of the reversed impulse response to meet your application's requirements.

$\endgroup$
1
$\begingroup$

First of all, I want to thank Fat32 and Hilmar for answering the question. I'm sorry I couldn't answer before, because I had some external problems to solve.

Finally I was able to solve the problem by using the reverse time property of the Zeta-Transform (despite having some difficulties with the initial conditions).

Consider a non-minimum phase system given by (with all zeros outside the unit circle)

$y[n] = H(z) u[n]$

We want to obtain the unknown input sequence $u[n]$ from the output sequence $y[n]$ and the filter $H(z)$. Clearly, since $H(z)$ is a non-minimum phase system, it is not possible to use its inverse $H^{-1}(z)$ (because it is unstable)

$u[n] = H^{-1}(z) y[n]$

However, by applying the reverse time property we have

$u[-n] = H^{-1}(\frac{1}{z}) y[-n],$

which is now stable, since we mapped all the poles into the unit circle (remember that $H(z)$ has all of its zeros outside the unit circle)

We can recover the sequence $u[n]$ by filtering the reversed sequence $y[-n]$ through the stable filter $H^{-1}(\frac{1}{z})$ and finally reversing its result.

It has some issues with the initial conditions of the filter, however for system (or time series) identification it is possible to ignore (or to cut) the last samples.

Note that it is possible to extend this result to any non-minimum phase system $H(z)$ by separating the stable an unstable parts of $H(z)$

Here you have a simple example in Matlab:

$\texttt{clear}$

$\texttt{close all}$

$\texttt{clc}$

$\texttt{a = [1 -2];}$

$\texttt{b = [1 -0.3];}$

$\texttt{G=tf(a,b,-1,'Variable','z')}$

$\texttt{pzmap(G)}$

$\texttt{u=randn(1,500);}$

$\texttt{y=lsim(G,u);}$

$\texttt{G1=tf(flip(a),flip(b),-1,'Variable','z') %equal to evaluate G(1/z)}$

$\texttt{u2=flipud(lsim(G1^-1,flipud(y)));}$

$\texttt{figure}$

$\texttt{plot(u,'b')}$

$\texttt{hold on}$

$\texttt{plot(u2,'r')}$

Regards and thanks a lot

$\endgroup$
  • $\begingroup$ This caused me a great headache to figure out what seems wrong in your explanation. Indeed the solution works perfectly but your explanation is misleading: stabilization has nothing to do with time-reversal property of Z-transform. It's about anti-causal & stable impulse response associated the inverse system $1/H(z)$ by treating its ROC as $|z| < p_{min}$. This stable (but acausal) impulse response can be obtained from evaluating the causal inverse-Z of $1/H(1/z)$. Hence your time-reversal argument. Code: ui = filter( a(end:-1:1), b(end:-1:1), y(end:-1:1)) ur = ui(end:-1:1) $\endgroup$ – Fat32 Mar 22 at 2:15
  • $\begingroup$ Nevertheless, since all the signals are reversed somewhow, it seems pretty ok to address this issue into time-reversal property. But keep in mind that, you can obtain the exact same result without using the time-reversal concept: just evaluate the acausal impulse respone for $n<0$, window it to a proper length, and shift it to make it causal and use filtering or convolution... of course above code is a much clever means of getting that result. $\endgroup$ – Fat32 Mar 22 at 2:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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