はじめに
https://yuhi-sa.github.io/posts/20260310_adaptive_filter/1/ と https://yuhi-sa.github.io/posts/20260501_lms_nlms/1/ では、LMS・NLMSアルゴリズムを中心に適応フィルタを解説しました。両記事とも最後に「RLSはLMSより高速に収束するが計算コストが高い」「RLSはカルマンフィルタと数学的に等価」と触れるにとどめ、導出には踏み込んでいませんでした。本記事ではその続きとして、**RLS（Recursive Least Squares、逐次最小二乗法）**を行列反転補題から導出し、Pythonでフルスクラッチ実装します。さらに、RLSとカルマンフィルタが特定の条件下で数値的に完全一致することを実験で確かめ、忘却係数 \(\lambda\) が定常誤差と追従性のトレードオフをどう支配するかを検証します。
問題設定：指数重み付き最小二乗
RLSは、LMSのような瞬時勾配の確率的近似ではなく、過去の全観測に対する重み付き二乗誤差を厳密に最小化します。長さ \(M\) のFIRフィルタ係数 \(\mathbf{w}\) に対し、時刻 \(n\) でのコスト関数を
\[ J(\mathbf{w}, n) = \sum_{i=1}^{n} \lambda^{n-i} \bigl[d(i) - \mathbf{x}(i)^\top \mathbf{w}\bigr]^2 \tag{1} \]と定義します。\(\lambda \in (0, 1]\) は忘却係数（forgetting factor）で、過去のサンプルの寄与を指数的に減衰させます。\(\lambda = 1\) ならすべてのサンプルを等しく扱う通常の最小二乗、\(\lambda < 1\) なら直近のサンプルを重視する適応的な推定になります。
\(J(\mathbf{w}, n)\) を \(\mathbf{w}\) で微分してゼロと置くと、正規方程式
\[ \mathbf{R}(n) \mathbf{w}(n) = \mathbf{r}(n), \qquad \mathbf{R}(n) = \sum_{i=1}^n \lambda^{n-i} \mathbf{x}(i)\mathbf{x}(i)^\top, \quad \mathbf{r}(n) = \sum_{i=1}^n \lambda^{n-i} \mathbf{x}(i) d(i) \tag{2} \]が得られます。\(\mathbf{R}(n)\) を毎回逆行列計算すると \(O(M^3)\) かかりますが、次の行列反転補題を使えば \(O(M^2)\) の再帰更新に落とし込めます。
行列反転補題によるO(M²)再帰更新
\(\mathbf{R}(n)\) は前時刻から次のように更新できます。
\[ \mathbf{R}(n) = \lambda \mathbf{R}(n-1) + \mathbf{x}(n)\mathbf{x}(n)^\top \tag{3} \]\(P(n) = \mathbf{R}(n)^{-1}\) とおき、Sherman-Morrison公式（行列反転補題）
\[ (A + \mathbf{u}\mathbf{v}^\top)^{-1} = A^{-1} - \frac{A^{-1}\mathbf{u}\mathbf{v}^\top A^{-1}}{1 + \mathbf{v}^\top A^{-1}\mathbf{u}} \tag{4} \]を \(A = \lambda \mathbf{R}(n-1)\) 、\(\mathbf{u} = \mathbf{v} = \mathbf{x}(n)\) に適用すると、
\[ P(n) = \frac{1}{\lambda}\left[P(n-1) - \frac{P(n-1)\mathbf{x}(n)\mathbf{x}(n)^\top P(n-1)}{\lambda + \mathbf{x}(n)^\top P(n-1)\mathbf{x}(n)}\right] \tag{5} \]が得られます。逆行列を直接計算せず、既知の \(P(n-1)\) から行列積とスカラー除算だけで \(P(n)\) を求められるのがポイントです。ここでカルマンゲインにあたる量
\[ \mathbf{k}(n) = \frac{P(n-1)\mathbf{x}(n)}{\lambda + \mathbf{x}(n)^\top P(n-1)\mathbf{x}(n)} \tag{6} \]を定義すると、式(5)は \(P(n) = \frac{1}{\lambda}\bigl[P(n-1) - \mathbf{k}(n)\mathbf{x}(n)^\top P(n-1)\bigr]\) と簡潔に書けます。重み更新は
\[ e(n) = d(n) - \mathbf{x}(n)^\top \mathbf{w}(n-1), \qquad \mathbf{w}(n) = \mathbf{w}(n-1) + \mathbf{k}(n) e(n) \tag{7} \]です。1ステップあたりの演算量は \(\mathbf{k}(n)\) の計算が支配的で \(O(M^2)\) 、LMS/NLMSの \(O(M)\) より重くなります。
Python実装
import numpy as np
def rls_filter(xs, d, M, lam=1.0, delta=1.0):
"""RLSアルゴリズムによる適応フィルタ
Parameters
----------
xs : ndarray, shape (N, M)
各時刻の入力ベクトル
d : ndarray, shape (N,)
所望信号
M : int
フィルタ長
lam : float
忘却係数（0 < lam <= 1）
delta : float
初期共分散の正則化パラメータ（P(0) = I / delta）
Returns
-------
W : ndarray, shape (N, M)
各時刻のフィルタ係数の推移
"""
w = np.zeros(M)
P = np.eye(M) / delta
W = np.zeros((len(d), M))
for n in range(len(d)):
x = xs[n]
Px = P @ x
k = Px / (lam + x @ Px)
e = d[n] - x @ w
w = w + k * e
P = (P - np.outer(k, x @ P)) / lam
W[n] = w
return W
np.outer(k, x @ P) が式(5)のランク1更新項、Px / (lam + x @ Px) が式(6)のゲインに対応します。
数値実験1：LMSとの収束速度比較
\(M=8\) タップのFIRシステム同定問題（SNR 20dB）で、LMS（ステップサイズ \(\mu=0.05\) ）とRLS（\(\lambda=1\) ）を200試行モンテカルロ平均した学習曲線を比較しました。
| サンプル数 \(n\) | LMSのMSE | RLSのMSE |
|---|---|---|
| 20 | 0.2916 | 0.0279 |
| 50 | 0.0310 | 0.0112 |
| 100 | 0.0145 | 0.0120 |
雑音床の2倍（MSE ≈ 0.0207）に到達するまでのサンプル数は、LMSが57サンプル、RLSが27サンプルでした。RLSはLMSの約半分のサンプル数で収束しており、理論通り「LMSは線形収束、RLSは二次収束（入力の統計に依存せず最短時間で正規方程式の解に到達）」という特性が数値的に確認できます。一方でRLSは1ステップ \(O(M^2)\) のためLMSの \(O(M)\) より計算コストが高く、この収束速度と計算量のトレードオフが用途選択の基準になります。
RLSとカルマンフィルタの等価性
https://yuhi-sa.github.io/posts/20260224_kalman_filter/1/ で導出したカルマンフィルタの更新式（予測・イノベーション・ゲイン・更新）を思い出してください。RLSは、以下のようにカルマンフィルタの特殊ケースとして定式化できます。
状態空間モデルとして見たRLS：
| カルマンフィルタの記号 | RLSでの対応 | 意味 |
|---|---|---|
| 状態 \(\mathbf{x}_k\) | フィルタ係数 \(\mathbf{w}\) | 推定したいパラメータ |
| 状態遷移 \(A\) | \(I\) （単位行列） | 係数は時間的に一定と仮定 |
| プロセスノイズ \(Q\) | \(0\) | 係数の真値は変化しない |
| 観測行列 \(H\) | \(\mathbf{x}(n)^\top\) | 時刻ごとに変わる入力ベクトル |
| 観測ノイズ \(R\) | \(\lambda\) （定数） | 忘却係数がそのまま観測分散に対応 |
| 初期共分散 \(P_0\) | \(I / \delta\) | RLSの正則化パラメータと一致 |
このマッピングでカルマンフィルタの更新式（\(K_k = P_{k|k-1}H^\top S_k^{-1}\) 等）を書き下すと、式(6)(7)のRLS更新式と完全に一致します。実際に両者を実装して同じデータに適用し、係数の推移を比較しました。
def kalman_param_filter(xs, d, M, R=1.0, delta=1.0):
w = np.zeros(M)
P = np.eye(M) / delta
W = np.zeros((len(d), M))
for n in range(len(d)):
x = xs[n]
Px = P @ x
S = x @ Px + R # イノベーション共分散
K = Px / S # カルマンゲイン
e = d[n] - x @ w
w = w + K * e
P = P - np.outer(K, x @ P)
W[n] = w
return W
\(M=4\) タップ、500サンプルのシステム同定問題で \(\lambda = 1\) 、\(R = 1\) 、\(\delta\) を揃えて両者を実行した結果：
max |RLS - KF| over all n, all taps (lambda=1, R=1): 0.0
final RLS weights: [ 0.99811908 -0.49835497 0.29976788 0.19937108]
final KF weights: [ 0.99811908 -0.49835497 0.29976788 0.19937108]
true weights : [ 1. -0.5 0.3 0.2]
**全時刻・全タップで誤差が完全にゼロ（機械精度で一致）**となり、RLSがカルマンフィルタの特殊ケースであることが数値的に裏付けられました。この等価性は実務上も有用で、RLSの忘却係数を調整する代わりに、より一般的なカルマンフィルタの枠組み（プロセスノイズ \(Q\) を明示的に導入する）に切り替えれば、非定常性のモデル化を柔軟に行えます。
忘却係数λの効果：定常性 vs 追従性
\(\lambda < 1\) の近似的な直感は「プロセスノイズ \(Q > 0\) を許容するカルマンフィルタに近づく」というものです。これを検証するため、系の途中（1000サンプル目）でフィルタ係数の真値が急変する非定常システムを用意し、\(\lambda\) を変えてRLSを適用しました。
| \(\lambda\) | ジャンプ直前の定常誤差 | ジャンプ後100サンプルの誤差 |
|---|---|---|
| 1.000 | 0.0013 | 2.0197 |
| 0.995 | 0.0006 | 1.3605 |
| 0.980 | 0.0010 | 0.2947 |
| 0.900 | 0.0030 | 0.0053 |
\(\lambda = 1\) は定常状態での誤差が最小（過去の全データを均等に使うため統計的に最も効率的）ですが、システムが変化した後の追従が最も遅く、100サンプル経過してもほとんど新しい値に収束していません。逆に \(\lambda = 0.9\) は定常誤差がやや大きい（直近10サンプル程度しか実質的に使っていないのと同等）反面、変化への追従は非常に速く、100サンプル後にはほぼ収束しています。\(\lambda\) は「定常時の精度」と「非定常環境への追従性」を直接トレードオフする設計パラメータであることが、この実験から定量的に確認できます。
計算量とLMS/NLMSとの比較
https://yuhi-sa.github.io/posts/20260310_adaptive_filter/1/ で示した比較表を再掲します。
| 特性 | LMS | NLMS | RLS |
|---|---|---|---|
| 計算量（1ステップ） | \(O(M)\) | \(O(M)\) | \(O(M^2)\) |
| メモリ | \(O(M)\) | \(O(M)\) | \(O(M^2)\) |
| 収束速度 | 遅い（線形） | 中程度 | 速い（二次収束） |
| 定常誤差 | ステップサイズ依存 | ステップサイズ依存 | 小さい（\(\lambda\) 依存） |
| 数値安定性 | 高い | 高い | \(P\) 行列の発散に注意 |
RLSは行列 \(P\) の更新を反復するため、有限精度演算では丸め誤差が蓄積して \(P\) が非対称・非正定値になり発散することがあります（実務では \(P \leftarrow (P + P^\top)/2\) で対称化する、あるいはQR分解ベースのRLS変種を使うなどの対策が取られます）。本記事のPython実装は教育目的の素朴な実装であり、長時間運用する組み込み用途ではこの数値安定性対策が必要な点に注意してください。
用途別ガイド
| 用途 | 推奨 | 理由 |
|---|---|---|
| タップ数が小さく（\(M \lesssim 20\) ）高速収束が必須 | RLS | 二次収束で最短時間で正規方程式の解に到達 |
| タップ数が大きい・組み込み低リソース | LMS / NLMS | \(O(M)\) で計算・メモリコストを抑制 |
| 環境が緩やかに非定常 | RLS（\(\lambda < 1\) ） | 忘却係数で定常誤差と追従性を調整可能 |
| 環境の変化が急激・大きい | カルマンフィルタ | プロセスノイズ \(Q\) を明示的にモデル化できる |
| 長時間連続運用（数値安定性が重要） | NLMS | \(P\) 行列の発散リスクがなく実装も単純 |
関連記事
- 適応フィルタ（LMS/RLS）の理論とPython実装 - RLSの位置づけとノイズキャンセレーション応用を含む包括解説です。本記事はRLSの導出と等価性に焦点を絞った発展編です。
- LMS／NLMSアルゴリズムの理論とPython実装 - 姉妹記事。LMS/NLMSの収束条件・定常誤差の詳細な導出はこちらを参照してください。
- カルマンフィルタの理論とPython実装 - 本記事で数値的に一致を示した状態空間フィルタです。
- 非線形カルマンスムーザ：Extended RTS SmootherとUnscented RTS Smootherの理論とPython実装 - カルマン系の発展編。RLSとの等価性を踏まえると、非線形拡張への理解も深まります。
- ウィーナーフィルタの理論とPython実装 - RLS・LMSが逐次近似する最適解そのものです。
- 自己相関関数の理論とPython実装 - 正規方程式 \(\mathbf{R}(n)\mathbf{w}(n) = \mathbf{r}(n)\) の \(\mathbf{R}\) が入力の自己相関行列であることの理解に直結します。
- SGD・Adamの理論と比較 - LMSと同じ確率的勾配降下法の枠組みの深層学習版です。RLSは二次法（自然勾配法に近い発想）である点が対照的です。
- デジタルフィルタ設計指針ハブ：選定の3軸と特性比較 - 固定設計フィルタとの比較・選定指針として活用できます。
おすすめ書籍
適応フィルタを含むディジタルフィルタの基礎を体系的に学べる定番書です。本記事のRLS導出の理論的背景を補完します。
※ 上記は Amazon アソシエイトのリンクです。
参考文献
- Haykin, S. (2014). Adaptive Filter Theory (5th ed.). Pearson. Chapters 9-13.
- Sayed, A. H. (2008). Adaptive Filters. Wiley-IEEE Press. Chapter 5 (RLS).
- Kalman, R. E. (1960). A new approach to linear filtering and prediction problems. Journal of Basic Engineering, 82(1), 35-45.
- Haykin, S. (2001). Kalman Filtering and Neural Networks. Wiley. Chapter 1 (RLS-Kalman equivalence).