Exercises

The LMMSE estimate is $\hat{d}[n] = \sum_\ell h[\ell] y[n-\ell]$. Orthogonality requires $\mathbb{E}\!\left[\left(d[n] - \sum_\ell h[\ell] y[n-\ell]\right) y^*[n-k]\right] = 0$ for every $k$.

Expanding the expectation: $\mathbb{E}[d[n] y^*[n-k]] = r_{dy}[k]$ and $\mathbb{E}[y[n-\ell] y^*[n-k]] = r_y[k-\ell]$. The equation becomes $\sum_\ell h[\ell] r_y[k-\ell] = r_{dy}[k]$ for all $k \in \mathbb{Z}$. $\blacksquare$

This is a discrete convolution equation on the full integer line. Taking the DTFT gives $\check{h}_{\text{nc}}(f) = S_{dy}(f)/S_y(f)$ directly.

$r_{dy}[k] = \mathbb{E}[d[n+k](d[n]+w[n])^*] = r_d[k]$, so $S_{dy}(f) = S_d(f)$. Similarly $S_y(f) = S_d(f)+S_w(f)$.

The non-causal Wiener filter is $\check{h}_{\text{nc}}(f) = \dfrac{S_d(f)}{S_d(f)+S_w(f)}$, and the MMSE is $\sigma^2_{nc} = \displaystyle\int_{-1/2}^{1/2} \frac{S_d(f) S_w(f)}{S_d(f)+S_w(f)}\,df$. $\blacksquare$

The filter is a frequency-selective attenuator: it passes frequencies where $\mathrm{SNR}(f) = S_d(f)/S_w(f) \gg 1$ and suppresses those with $\mathrm{SNR}(f) \ll 1$.

Taking the transfer function $H(z) = 1/(1-\alpha z^{-1})$ acting on white $u[n]$, $S_d(f) = \dfrac{\sigma_u^2}{|1-\alpha e^{-j2\pi f}|^2} = \dfrac{\sigma_u^2}{1-2\alpha\cos(2\pi f)+\alpha^2}$.

$S_y(f) = S_d(f)+\sigma_w^2$.

$\check{h}_{\text{nc}}(f) = \dfrac{S_d(f)}{S_d(f)+\sigma_w^2} = \dfrac{\sigma_u^2}{\sigma_u^2 + \sigma_w^2 |1-\alpha e^{-j2\pi f}|^2}$. $\blacksquare$

For $\alpha$ close to 1 the signal PSD concentrates near DC, and the Wiener filter becomes a low-pass filter. For $\alpha$ near 0 the signal is nearly white and the filter collapses to the scalar MMSE gain $\sigma_u^2/(\sigma_u^2+\sigma_w^2)$.

From orthogonality, $\sigma^2_{nc} = \mathbb{E}[e[n] d^*[n]] = r_d[0] - \sum_k h[k] r_{dy}^*[k]$.

$r_d[0] = \int S_d(f)\,df$ and $\sum_k h[k] r_{dy}^*[k] = \int H(f) S_{dy}^*(f)\,df$. Substituting $H(f) = S_{dy}(f)/S_y(f)$ yields $\sigma^2_{nc} = \int\!\left(S_d(f) - |S_{dy}(f)|^2/S_y(f)\right)df$. $\blacksquare$

The integrand is the residual PSD after projecting the desired signal onto the observation spectrum at each frequency. When $|S_{dy}(f)|^2 = S_d(f) S_y(f)$ (perfect coherence), the integrand vanishes and the estimate is exact.

$\displaystyle\int_{-1/2}^{1/2} |\log S_y(f)|\,df < \infty$ (equivalently, $\log S_y(f) \in L^1$).

Spectral factorization constructs $G(z) = \exp\!\big(\sum_{k=0}^\infty c_k z^{-k}\big)$ with Fourier coefficients of $\tfrac12 \log S_y$. These coefficients exist and decay only when $\log S_y \in L^1$. Without this, no causal $G$ with a causal inverse exists.

A PSD that vanishes on a positive-measure set (spectral null) violates Paley-Wiener; the process is then deterministic in the forward direction and cannot be represented as the output of a causal stable filter driven by white noise. $\blacksquare$

$S_y(z) = 5 - 2z - 2z^{-1} = -2 z^{-1}(z^2 - \tfrac{5}{2} z + 1) = -2 z^{-1}(z-2)(z-\tfrac12)$.

$(z-2)(z-\tfrac12) = \tfrac12 (2z-1)(z-\tfrac12)\cdot(\text{rearrange})$. Equivalently, $S_y(z) = 4 \cdot (1-\tfrac12 z^{-1})(1-\tfrac12 z)$.

Take $G(z) = 2(1-\tfrac12 z^{-1})$ (zero inside unit disk at $z=\tfrac12$, causal and stable), and $\sigma_\nu^2 = 1$. Check: $|G(e^{j2\pi f})|^2 \cdot \sigma_\nu^2 = 4|1-\tfrac12 e^{-j2\pi f}|^2 = 4(1 - \cos 2\pi f + \tfrac14) = 5 - 4\cos 2\pi f$. $\blacksquare$

Let $\nu[n] = (1/G) * y$. In frequency, $S_\nu(f) = S_y(f)/|G(f)|^2 = \sigma_\nu^2$, constant. So $\nu[n]$ is white with variance $\sigma_\nu^2$.

Because $\nu[n]$ is white, projecting $d[n]$ onto $\{\nu[k]: k\leq n\}$ is a sequence of independent scalar projections. The causal Wiener filter is then the causal part of the cross-spectrum, normalized by $\sigma_\nu^2$: $\check{h}_c(f) = \dfrac{1}{G(f)}\left[\dfrac{S_{dy}(f)}{G^*(f)}\right]_+$.

The bracketed $[\cdot]_+$ operation keeps only causal (non-negative index) Fourier coefficients; it is the formal projection onto the past.

$S_y(z) = \dfrac{\sigma_u^2 + \sigma_w^2 (1-\alpha z^{-1})(1-\alpha z)}{(1-\alpha z^{-1})(1-\alpha z)}$. The numerator is a second-order symmetric polynomial; write it as $c^2 (1-\beta z^{-1})(1-\beta z)$ with $|\beta|<1$ and $c>0$. Matching coefficients gives $c^2 (1+\beta^2) = \sigma_u^2 + \sigma_w^2 (1+\alpha^2)$ and $c^2 \beta = \sigma_w^2 \alpha$.

After cancelling the common causal pole, one obtains $\check{h}_c(z) = \dfrac{K(1-\alpha z^{-1})}{1 - \beta z^{-1}}$, where $K = (\alpha-\beta)/(c(1-\alpha\beta))$ arises from the $[\cdot]_+$ operation. This is a first-order IIR smoother.

The causal MMSE is strictly larger than the non-causal one unless $\beta=0$ (i.e. the signal is white). The gap is quantified by the Kolmogorov-Szego formula applied to the prediction-error spectrum. $\blacksquare$

For a purely non-deterministic WSS process with PSD $S_y(f)>0$ a.e. satisfying Paley-Wiener, the minimum one-step prediction-error variance is $\sigma_\nu^2 = \exp\!\Big(\displaystyle\int_{-1/2}^{1/2} \log S_y(f)\,df\Big)$.

The arithmetic mean of $S_y(f)$ is $r_y[0]$ (signal power). The geometric mean is $\sigma_\nu^2$. The ratio is the achievable prediction-gain: $r_y[0]/\sigma_\nu^2 \geq 1$, with equality iff the process is white.

Concentrating PSD energy (colored spectra) lowers $\sigma_\nu^2$: predictable structure reduces the innovations floor. $\blacksquare$

$\begin{pmatrix} 4 & 2 \\ 2 & 4 \end{pmatrix}\begin{pmatrix}a_1\\a_2\end{pmatrix} = \begin{pmatrix}2\\1\end{pmatrix}$.

Determinant $= 12$. $a_1 = (4\cdot 2 - 2\cdot 1)/12 = 1/2$, $a_2 = (4\cdot 1 - 2\cdot 2)/12 = 0$.

$\sigma_e^2 = r_y[0] - a_1 r_y[1] - a_2 r_y[2] = 4 - 1 - 0 = 3$. $\blacksquare$

One can show (Kailath-Sayed-Hassibi, Ch. 7) that $\sigma^2_{c} - \sigma^2_{nc} = \int_{-1/2}^{1/2} \Big|\Big[\frac{S_{dy}}{G^*}\Big]_-\Big|^2 / \sigma_\nu^2\,df$, i.e. the energy of the anti-causal part of the whitened cross-spectrum.

The right-hand side is non-negative and vanishes iff $[S_{dy}/G^*]_-\equiv 0$, that is, iff the whitened cross-correlation is purely causal. For the additive-noise model this fails unless $S_d(f) S_w(f) = 0$ almost everywhere. $\blacksquare$

Access to the future strictly improves estimation whenever signal and noise overlap in frequency; the smoothing gap is a direct measure of that overlap.

$y[n] = \sum_{k\geq 0} g[k]\nu[n-k]$. The MMSE $m$-step predictor of $y[n+m]$ from $\{y[k]:k\leq n\}$ equals the projection onto $\{\nu[k]:k\leq n\}$, namely $\hat y[n+m\mid n] = \sum_{k\geq m} g[k] \nu[n+m-k]$.

Reindexing, the predictor has transfer function, in $z$-domain, $\hat Y(z)/Y(z) = \dfrac{1}{G(z)}\Big[z^m G(z)\Big]_{+, \text{shifted}}$, which reduces to $\hat{h}_m(z) = \dfrac{[z^m G(z)]_+}{G(z)}$, where $[\cdot]_+$ keeps non-negative powers of $z^{-1}$.

$\sigma_m^2 = \sigma_\nu^2 \sum_{k=0}^{m-1} |g[k]|^2$, which grows with $m$ and saturates at $r_y[0]$. $\blacksquare$

The true recursion $y[n]=\sum_k a_k y[n-k]+\nu[n]$ shows that the AR($p$) predictor $\hat y[n\mid n-1] = \sum_k a_k y[n-k]$ uses only the past $p$ samples and is optimal.

The prediction error is $\nu[n]$, hence $\sigma_1^2 = \sigma_\nu^2$. This is also what the Kolmogorov-Szego formula yields: $\sigma_\nu^2 = \exp \int \log S_y(f)\,df$. $\blacksquare$

For $|f|<B$, $S_d(f)\gg \sigma_w^2$ so $H(f)\approx 1$ (signal passes).

For $|f|>B$, $S_d(f)\approx 0$ so $H(f)\approx 0$ (noise suppressed).

Near the band edge the filter smoothly transitions; its sharpness depends on the roll-off of $S_d(f)$. The Wiener filter behaves as a signal-matched low-pass filter. $\blacksquare$

For $y = d+w$ with $d,w$ independent, $H(f) = S_d/(S_d+S_w)\in[0,1]$.

$\mathbb{E}|\hat d[n]|^2 = \int |H(f)|^2 S_y(f) df \leq \int S_y(f) df = r_y[0]$. Since $r_y[0] \geq r_d[0]$ (adding noise adds variance), in general this bounds $\hat d$ by the observation power. When additionally $S_d \leq S_y$ one also gets $\int|H|^2 S_d\leq r_d[0]$.

The Wiener filter is never a gain; it attenuates uniformly and preserves no noise power at frequencies where the signal is absent. $\blacksquare$

Taking expectations yields $\bar{\mathbf{h}}_{n+1} = \bar{\mathbf{h}}_n + \mu(\mathbf{r}_{yd} - \mathbf{R}\bar{\mathbf{h}}_n)$, which is steepest descent on $J(\mathbf{h}) = \mathbb{E}|d-\mathbf{h}^H\mathbf{y}|^2$. Its fixed point solves $\mathbf{R}\bar{\mathbf{h}}_\star = \mathbf{r}_{yd}$, the Wiener-Hopf normal equations.

Convergence in the mean requires $0 < \mu < 2/\lambda_{\max}(\mathbf{R})$. Slow-converging modes are dictated by $\lambda_{\min}(\mathbf{R})$, so the condition number $\lambda_{\max}/\lambda_{\min}$ controls the misadjustment-speed tradeoff. $\blacksquare$

Estimating $d[n]$ from data up to $n+L$ is equivalent (up to a delay) to causally estimating $d'[n]=d[n-L]$ from $y[n]$. Apply causal Wiener with cross-spectrum $e^{-j2\pi f L} S_{dy}(f)$.

$H_L(f) = \dfrac{e^{j2\pi f L}}{G(f)}\Big[\dfrac{S_{dy}(f) e^{-j2\pi f L}}{G^*(f)}\Big]_+$.

As $L\to\infty$ the bracketed causal projection absorbs the entire spectrum $S_{dy}/G^*$, yielding $H_\infty(f) = \dfrac{S_{dy}(f)}{|G(f)|^2} = \dfrac{S_{dy}(f)}{S_y(f)}$, which is the non-causal Wiener filter. $\blacksquare$

Every extra unit of allowable delay lets the smoother exploit more of the anti-causal cross-correlation. The curve of MMSE vs. $L$ is the canonical smoothing tradeoff.

Imposing $h[k]=0$ for $k<0$ and $k\geq M$ and restricting orthogonality to $k=0,\ldots,M-1$: $\sum_{\ell=0}^{M-1} h[\ell] r_y[k-\ell] = r_{dy}[k]$, i.e. $\mathbf{R}\mathbf{h}=\mathbf{r}$ with $\mathbf{R}_{k\ell}=r_y[k-\ell]$ Hermitian Toeplitz.

$\mathbf{R}$ is singular iff $y[n]$ lies in a proper subspace (e.g. a deterministic sum of sinusoids with fewer than $M$ components). Then the estimator is underdetermined on the null space.

Use diagonal loading $\mathbf{R}+\delta\mathbf{I}$ (ridge), or drop toward a shorter filter. Diagonal loading is the MAP solution under a Gaussian prior on $\mathbf{h}$. $\blacksquare$