Why Wiener filter is not optimal?

Question

I am investigating Wiener filter in the frequency domain, given by $$ W(\omega) = \frac{S_{ud}(\omega)}{S_{uu}(\omega)}$$ where $S_{xy}(\omega)$ is the cross-power spectral density (CPSD), see scipy.signal.csd or MATLAB cpsd. I used the Python's scipy.signal.csd for estimating $S_{xy}$.

Here $d[n]$ is the desired signal ("mic") and $u[n]$ ("ref") is the reference signal.

Now, here I describe how I generate the signal (in Python). I generate a 10 seconds white noise, and also generate a FIR transfer function true_h (there are 2 ways to generate it: one is designing a FIR bandpass filter and using it, the second is measured from a real room experiment).

from scipy import signal

# Generate filter
true_L = 150 # filter size for simulation
true_h = signal.firwin(true_L, cutoff=[1000, 3000], pass_zero="bandpass", fs=FS)
fft_h = np.fft.rfft(true_h)

# Generate 10 seconds of white noise
np.random.seed(0)
ref = np.random.randn(10*FS)/100
mic = signal.lfilter(true_h, 1, ref)

The goal is to estimate $w = \mathcal{F}^{-1}(W)$ where $W$ is the estimated Wiener filter in frequency and $\mathcal{F}^{-1}$ is the inverse discrete Fourier transform.

The Wiener in the frequency domain is calculated as:

from scipy import signal
# ref = reference signal
# mic = desired signal
freqs, Pxy = signal.csd(ref, mic, fs=FS, nperseg=nfft)
freqs, Pxx = signal.csd(ref, ref, fs=FS, nperseg=nfft)
   
wiener = Pxy / Pxx

The application of the estimated filter $w$ is done as follows:

# Inverse FFT for Wiener filter
wiener_taps = np.fft.irfft(wiener)

# Apply Wiener
error_wiener = mic - signal.lfilter(wiener_taps, 1, ref)

The desired outcome is that $e[n] = d[n] - w[n] \ast u[n]$ will be almost $0$. But this is not the case. I get a residual which is only around -10 to -20 dB.

When I use frequency/fast block LMS (FBLMS or FLMS), I obtain a residual $e[n]$ which is about -40 dB, i.e. the estimation of true_h is almost perfect.

I also plot both Wiener and FLMS estimates of true_h compared to the true transfer function (both in the taps and frequency domain) and get a very good match, but the match is better in FLMS than in Wiener).

This is odd, since Wiener is considered as the optimal filter in the MSE (mean square error) sense. But in this example, Wiener is performing worse than FLMS. See the results below:

The results:

Fig 1: Echo Return Loss Enhancement (ERLE). Black = original, Green = after Wiener, Blue = after LMS.

Fig 2: residue (e = d - w * u). Green: Wiener filter, Red: FLMS filter. Blue: adaptive FLMS (irrelevant to this question).

Added on edit: Here I further investigate Wiener. I generated a 150-taps FIR, and apply it on a white noise signal (generated by normal distrubtion):

# Generate filter
true_L = 150 # filter size for simulation
true_h = signal.firwin(true_L, cutoff=[1000, 3000], pass_zero="bandpass", fs=FS)
# Generate 10 seconds of white noise
np.random.seed(0)
ref = np.random.randn(10*FS)/100
mic = signal.lfilter(true_h, 1, ref)

Then calculated Wiener by:

1. The default CSD parameters of scipy.

f1, Pxy = signal.csd(ref, mic, fs=FS)
f1, Pxx = signal.csd(ref, ref, fs=FS)
wiener_freq = Pxy/Pxx

Fig 3: Misalignment and residual after Wiener (default framing)

2. Taking one long FFT on the entire signal (one "frame"):

f1, Pxy = signal.csd(ref, mic, fs=FS, nperseg=len(ref))
f1, Pxx = signal.csd(ref, ref, fs=FS, nperseg=len(ref))
wiener_freq = Pxy/Pxx

Fig 4: Misalignment and residual after Wiener (one big frame)

As we can see, the second method yielded far better result, but still not a perfect 0 residual (i.e. residual < -100 dB).

Added in Edit: It seem that this is the case when the signals are synthetic: i.e. the reference is a white noise generated by numpy.random.randn and the transfer function is generated by scipy.signal.firwin with about 50-250 taps. For real data, it seems that Wiener performs better.

Added in another edit: Figure 5: Wiener (orange) vs. FBLMS (green) for synthetic white noise signal and FIR transfer function.

Figure 6: Wiener (orange) vs. FBLMS (green) for real signal (recorded using an audio device)

As it can be seen, for the synthetic case FBLMS is better than Wiener while for the real case it is vice versa, Wiener is better and FBLMS converges to Wiener.

References:

• https://webee.technion.ac.il/people/shimkin/Estimation09/ch3_Wiener.pdf

Answer 1

The Wiener filter minimizes the MSE under a series of assumptions, which usually don’t match reality. That means that, in practice, it will not give the solution with the smallest MSE possible.

Using the Wiener filter for deconvolution is the simplest thing you can do, and there are many other, more complex deconvolution algorithms that produce better results.