\begin{equation*}
\newcommand\innerp[2]{\langle #1, #2 \rangle}
\newcommand\ve[1]{\boldsymbol{#1}}
\newcommand\parfrac[2]{\frac{\partial #1}{\partial #2}}
\newcommand\dfrac[2]{\frac{\mathrm{d} #1}{\mathrm{d} #2}}
\newcommand\mean[2]{\mathrm{E}_{#1} \left[ #2 \right]}
\end{equation*}
正則化が何故よいのか?というところで、解の表現力が落ちることで量子化誤差が減ったのがあり得るなと思い始めている。
とくに、最大絶対値に合わせたシフトが減るのは大きいかも。
Ill-Conditioning and Bandwidth Expansion in Linear Prediction of Speech を読んでいて、分散に定数加算することによる正則化で自己相関行列の条件数改善が理論的にできることが示されていたので、ここで確認する。
Heykin様のAdaptive Filter Theoryのp816に「自己相関行列 \(\ve{R}\) の固有値は、時系列のパワースペクトル密度 \(S(\omega)\) の最小値 \(S_{\min}\) と最大値 \(S_{\max}\) の間にある」という命題があるので示す。
(証明) \(\ve{R}\) の各固有値 \(\lambda_{i}\ (i = 1,...,M)\) に対する固有ベクトルを \(\ve{q}_{i}\) と書くと、レイリー商の性質から、
\begin{equation*}
\lambda_{i} = \frac{\ve{q}_{i}^{H} \ve{R} \ve{q}_{i}}{\ve{q}_{i}^{H} \ve{q}_{i}} \tag{1}
\end{equation*}
と書ける( \(\cdot^{H}\) :エルミート共役)。(1)の分子は二次形式で
\begin{equation*}
\ve{q}_{i}^{H} \ve{R} \ve{q}_{i} = \sum_{k=1}^{M} \sum_{l=1}^{M} \overline{q_{ik}} r(l-k) q_{il} \tag{2}
\end{equation*}
と書ける( \(\overline{\cdot}\) :共役、 \(r\) :自己相関関数)。ウィーナー・ヒンチンの定理より、 \(r\) は
\begin{equation*}
r(l-k) = \frac{1}{2\pi} \int_{-\pi}^{\pi} S(\omega) \exp[j\omega(l-k)] d\omega \tag{3}
\end{equation*}
と表せるから、(3)を(2)に代入すると、
\begin{align*}
\ve{q}_{i}^{H} \ve{R} \ve{q}_{i} &= \frac{1}{2\pi} \sum_{k=1}^{M} \sum_{l=1}^{M} \overline{q_{ik}} q_{il} \int_{-\pi}^{\pi} S(\omega) \exp[j\omega(l-k)] d\omega \\
&= \frac{1}{2\pi} \int_{-\pi}^{\pi} S(\omega) \left( \sum_{k=1}^{M} \overline{q_{ik}} \exp[-j\omega k] \right) \left( \sum_{l=1}^{M} q_{il} \exp[j\omega l] \right) d\omega \tag{4}
\end{align*}
ここで、 \(Q_{i}(\omega) = \sum_{k=1}^{M} \overline{q_{ik}} \exp[-j\omega k]\) ( \(\overline{q_{ik}}\) のDFT)とおくと、 \(\overline{Q_{i}(\omega)} = \overline{\sum_{k=1}^{M} \overline{q_{ik}} \exp[-j\omega k]} = \sum_{k=1}^{M} \overline{\overline{q_{ik}}} \overline{\exp[-j\omega k]} = \sum_{k=1}^{M} q_{ik} \exp[j\omega k]\) だから、
\begin{equation*}
\ve{q}_{i}^{H} \ve{R} \ve{q}_{i} = \frac{1}{2\pi} \int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} S(\omega) d\omega \tag{5}
\end{equation*}
また、インパルス関数 \(\delta(\cdot)\) のパワースペクトラムは平坦(全て1)で、再びウィーナー・ヒンチンの定理より
\begin{equation*}
\delta(l-k) = \frac{1}{2\pi} \int_{-\pi}^{\pi} \exp[j\omega(l-k)] d\omega
\end{equation*}
だから、これを使って \(\ve{q}_{i}^{H} \ve{q}_{i}\) を書き換えると、
\begin{align*}
\ve{q}_{i}^{H} \ve{q}_{i} &= \ve{q}_{i}^{H} \ve{I} \ve{q}_{i} = \sum_{k=1}^{M} \sum_{l=1}^{M} \overline{q_{ik}} q_{il} \delta(l-k) \\
&= \frac{1}{2\pi} \sum_{k=1}^{M} \sum_{l=1}^{M} \overline{q_{ik}} q_{il} \int_{-\pi}^{\pi} \exp[j\omega(l-k)] d\omega \\
&= \frac{1}{2\pi} \int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} d\omega \tag{6}
\end{align*}
(5), (6)を(1)に代入することで、
\begin{equation*}
\lambda_{i} = \displaystyle\frac{\int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} S(\omega) d\omega}{\int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} d\omega}
\end{equation*}
とまとめられる。ここで、
\begin{equation*}
S_{\min} \int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} d\omega \leq \int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} S(\omega) d\omega \leq S_{\max} \int_{-\pi}^{\pi} | Q_{i}(\omega)|^{2} d\omega
\end{equation*}
だから、
\begin{equation*}
S_{\min} \leq \lambda_{i} \leq S_{\max}
\end{equation*}
が成立し、命題が示された。(証明終)
命題からただちに、 \(\ve{R}\) の条件数 \(\kappa(\ve{R})\) は以下の上界を持つことがわかる。
\begin{equation*}
\kappa(\ve{R}) := \frac{\max_{i} \lambda_{i}}{\min_{i} \lambda_{i}} \leq \frac{S_{\max}}{S_{\min}}
\end{equation*}
未証明だが、 \(M \to \infty\) とすると、 \(\max_{i} \lambda_{i} \to S_{\max}\) かつ \(\min_{i} \lambda_{i} \to S_{\min}\) となる(要確認…)。このことから、 \(\lim_{M \to \infty} \frac{S_{\max}}{S_{\min}} = \kappa(\ve{R})\) が成立する。 \(M\) が十分大きければ、条件数の良い近似となっていることを示している。
正則化による影響を考える。行列 \(\ve{R}\) の対角要素に \(\varepsilon > 0\) を足すことは、自己相関関数としては、
\begin{equation*}
r^{\prime}(k) = r(k) + \varepsilon \delta(k)
\end{equation*}
という修正を行うことに対応する。 \(r^{\prime}\) に対応するパワースペクトル密度 \(S^{\prime}(\omega)\) は、ウィーナー・ヒンチンの定理より、
\begin{align*}
S^{\prime}(\omega) &= \int_{-\pi}^{\pi} r^{\prime}(k) \exp(-j\omega k) dk \\
&= \int_{-\pi}^{\pi} r(k) \exp(-j\omega k) dk + \varepsilon \int_{-\pi}^{\pi} \delta(k) \exp(-j\omega k) dk = S(\omega) + \varepsilon \exp(0) \\
&= S(\omega) + \varepsilon
\end{align*}
(全帯域が \(\varepsilon\) だけ持ち上がることになる)このことより、正則化を行った場合の条件数は以下の上界を持つ:
\begin{equation*}
\frac{S_{\max} + \varepsilon}{S_{\min} + \varepsilon}
\end{equation*}
なお、
\begin{align*}
\frac{S_{\max}}{S_{\min}} - \frac{S_{\max} + \varepsilon}{S_{\min} + \varepsilon} &= \frac{1}{S_{\min}(S_{\min} + \varepsilon)} \left\{ S_{\max}(S_{\min} + \varepsilon) - S_{\min}(S_{\max} + \varepsilon) \right\} \\
&= \frac{\varepsilon(S_{\max} - S_{\min})}{S_{\min}(S_{\max} + \epsilon)} \geq 0
\end{align*}
だから、
\begin{equation*}
\frac{S_{\max} + \varepsilon}{S_{\min} + \varepsilon} \leq \frac{S_{\max}}{S_{\min}}
\end{equation*}
となって、上界を抑えることができる。とくに \(\lim_{\varepsilon \to \infty} \frac{S_{\max} + \varepsilon}{S_{\min} + \varepsilon} = 1\) だから条件数も1に漸近する。これは、自己相関関数がほぼ単位行列に近づく(対角優位になる)ことに一致する。
条件数については理解が足りてない部分があった。端的に言って、Ax=bにおけるbの変化による解xの変化率の上限。以下の資料が分かりやすかった。