レビンソン再帰
Recursive algorighm in linear algebra

レビンソン再帰法またはレビンソン・ダービン再帰法は、線形代数において、テプリッツ行列を含む方程式の解を再帰的に計算する手順である。このアルゴリズムはΘ ( n 2 )時間で実行され、Θ( n 3 )時間を要するガウス・ジョルダン消去法よりも大幅に改善されている。

レビンソン・ダービンアルゴリズムは、1947 年にノーマン・レビンソンによって最初に提案され、1960 年にジェームズ・ダービンによって改良され、その後、 WF トレンチと S. ゾハールによってそれぞれ4 n 2乗算と3 n 2乗算に改良されました。

データ処理の他の手法としては、シュアー分解やコレスキー分解などがあります。これらと比較すると、レビンソン再帰（特に分割レビンソン再帰）は計算速度が速い傾向がありますが、丸め誤差などの計算上の不正確さの影響を受けやすくなります。

テプリッツ行列に対するBareissアルゴリズム（一般的なBareissアルゴリズムと混同しないこと）は、レビンソン再帰とほぼ同程度の速度で実行されますが、O ( n2 )の空間を使用するのに対し、レビンソン再帰はO ( n )の空間しか使用しません。Bareissアルゴリズムは数値的に安定ですが[1] [2]、レビンソン再帰はせいぜい弱安定（つまり、条件付き線形システムに対して数値的に安定）です。[3]

漸近的に高速な、あるいは超高速なテプリッツアルゴリズムと呼ばれる新しいアルゴリズムは、様々なp（例えばp = 2、[4] [5] p = 3 [6] ）に対してΘ( n (log n ) p )で解くことができます。レビンソン再帰は、いくつかの理由から依然として人気があります。1つは、比較的理解しやすいこと、もう1つは、nが小さい場合（通常n  < 256）は超高速アルゴリズムよりも高速になる可能性があることです。[7]

導出

背景

行列方程式は次の形式に従う。

M x → = y → . {\displaystyle \mathbf {M} \,{\vec {x}}={\vec {y}}.}

レビンソン・ダービン法は、 M が主対角成分が非零である既知のテプリッツ行列である限り、このような任意の方程式に使用できます。ここでは既知のベクトルであり、はまだ決定されていない数x iからなる未知のベクトルです。 y → {\displaystyle {\vec {y}}} x → {\displaystyle {\vec {x}}}

この記事では、ê iはi番目の位置（値 1）を除いてすべてゼロのベクトルであるとします。その長さは周囲の文脈によって暗黙的に決定されます。Nという用語は上記の行列の幅を表します。MはN × N行列です。最後に、この記事では、上付き文字は帰納的インデックスを表し、下付き文字はインデックスを表します。例えば、この記事では、行列T nはMの左上のn × nブロックをコピーしたn × n行列です 。つまり、T n ij = M ijです。

T nもテプリッツ行列であり、次のように表される。

T n = [ t 0 t − 1 t − 2 … t − n + 1 t 1 t 0 t − 1 … t − n + 2 t 2 t 1 t 0 … t − n + 3 ⋮ ⋮ ⋮ ⋱ ⋮ t n − 1 t n − 2 t n − 3 … t 0 ] . {\displaystyle \mathbf {T} ^{n}={\begin{bmatrix}t_{0}&t_{-1}&t_{-2}&\dots &t_{-n+1}\\t_{1}&t_{0}&t_{-1}&\dots &t_{-n+2}\\t_{2}&t_{1}&t_{0}&\dots &t_{-n+3}\\\vdots &\vdots &\vdots &\ddots &\vdots \\t_{n-1}&t_{n-2}&t_{n-3}&\dots &t_{0}\end{bmatrix}}.}

導入手順

このアルゴリズムは2つのステップで進行します。最初のステップでは、順方向ベクトルと逆方向ベクトルと呼ばれる2組のベクトルが確立されます。順方向ベクトルは逆方向ベクトルの集合を取得するために使用され、その後すぐに破棄されます。逆方向ベクトルは2番目のステップで必要となり、そこでは目的の解を構築するために用いられます。

レビンソン・ダービン再帰法では、 n番目の「順方向ベクトル」を、次を満たす長さnのベクトルとして定義します。 f → n {\displaystyle {\vec {f}}^{n}}

T n f → n = e ^ 1 . {\displaystyle \mathbf {T} ^{n}{\vec {f}}^{n}={\hat {e}}_{1}.}

n番目の「後方ベクトル」も同様に定義され、次式を満たす長さnのベクトルです。 b → n {\displaystyle {\vec {b}}^{n}}

T n b → n = e ^ n . {\displaystyle \mathbf {T} ^{n}{\vec {b}}^{n}={\hat {e}}_{n}.}

Mが対称行列の場合には、重要な簡略化が起こります。その場合、2つのベクトルはb n i = f n n +1− iの関係を持ちます。つまり、それらは互いに行反転した関係にあります。この特殊なケースでは、これにより余分な計算を省くことができます。

後方ベクトルの取得

行列が対称でなくても、長さn − 1のベクトルからn番目の順方向ベクトルと逆方向ベクトル を次のように求めることができます。まず、順方向ベクトルをゼロで拡張して、次の式を得ます。

T n [ f → n − 1 0 ] = [       t − n + 1   T n − 1   t − n + 2       ⋮ t n − 1 t n − 2 … t 0 ] [   f → n − 1   0   ] = [ 1 0 ⋮ 0 ε f n ] . {\displaystyle \mathbf {T} ^{n}{\begin{bmatrix}{\vec {f}}^{n-1}\\0\\\end{bmatrix}}={\begin{bmatrix}\ &\ &\ &t_{-n+1}\\\ &\mathbf {T} ^{n-1}&\ &t_{-n+2}\\\ &\ &\ &\vdots \\t_{n-1}&t_{n-2}&\dots &t_{0}\\\end{bmatrix}}{\begin{bmatrix}\ \\{\vec {f}}^{n-1}\\\ \\0\\\ \\\end{bmatrix}}={\begin{bmatrix}1\\0\\\vdots \\0\\\varepsilon _{f}^{n}\end{bmatrix}}.}

T n −1からT nへ進む際、行列に追加された余分な列は、順方向ベクトルを拡張するために零点が用いられる場合、解を乱すことはありません。しかし、行列に追加された余分な行は解を乱し、最後に 不要な誤差項ε fを生み出します。上記の式は、その値を与えます。

ε f n   =   ∑ i = 1 n − 1   M n i   f i n − 1   =   ∑ i = 1 n − 1   t n − i   f i n − 1 . {\displaystyle \varepsilon _{f}^{n}\ =\ \sum _{i=1}^{n-1}\ M_{ni}\ f_{i}^{n-1}\ =\ \sum _{i=1}^{n-1}\ t_{n-i}\ f_{i}^{n-1}.}

この誤差はすぐに元に戻り、新しい順方向ベクトルから除去されます。しかし、まずは逆方向ベクトルを同様に（ただし逆方向に）拡張する必要があります。逆方向ベクトルについては、

T n [ 0 b → n − 1 ] = [ t 0 … t − n + 2 t − n + 1 ⋮       t n − 2   T n − 1   t n − 1     ] [   0   b → n − 1   ] = [ ε b n 0 ⋮ 0 1 ] . {\displaystyle \mathbf {T} ^{n}{\begin{bmatrix}0\\{\vec {b}}^{n-1}\\\end{bmatrix}}={\begin{bmatrix}t_{0}&\dots &t_{-n+2}&t_{-n+1}\\\vdots &\ &\ &\ \\t_{n-2}&\ &\mathbf {T} ^{n-1}&\ \\t_{n-1}&\ &\ &\end{bmatrix}}{\begin{bmatrix}\ \\0\\\ \\{\vec {b}}^{n-1}\\\ \\\end{bmatrix}}={\begin{bmatrix}\varepsilon _{b}^{n}\\0\\\vdots \\0\\1\end{bmatrix}}.}

前回と同様に、行列に追加された列はこの新しい逆ベクトルに影響を与えませんが、行が追加されると影響を受けます。ここで、別の望ましくない誤差ε b が次の値を持ちます。

ε b n   =   ∑ i = 2 n   M 1 i   b i − 1 n − 1   =   ∑ i = 1 n − 1   t − i   b i n − 1 .   {\displaystyle \varepsilon _{b}^{n}\ =\ \sum _{i=2}^{n}\ M_{1i}\ b_{i-1}^{n-1}\ =\ \sum _{i=1}^{n-1}\ t_{-i}\ b_{i}^{n-1}.\ }

これら2つの誤差項は、以下のように記述される高次の順方向ベクトルと逆方向ベクトルを形成するために用いられる。行列の線形性を用いると、すべてのに対して以下の恒等式が成立する。 ( α , β ) {\displaystyle (\alpha ,\beta )}

T ( α [ f →   0 ] + β [ 0   b → ] ) = α [ 1 0 ⋮ 0 ε f ] + β [ ε b 0 ⋮ 0 1 ] . {\displaystyle \mathbf {T} \left(\alpha {\begin{bmatrix}{\vec {f}}\\\ \\0\\\end{bmatrix}}+\beta {\begin{bmatrix}0\\\ \\{\vec {b}}\end{bmatrix}}\right)=\alpha {\begin{bmatrix}1\\0\\\vdots \\0\\\varepsilon _{f}\\\end{bmatrix}}+\beta {\begin{bmatrix}\varepsilon _{b}\\0\\\vdots \\0\\1\end{bmatrix}}.}

αとβ を右辺がê 1またはê nとなるように選択すると、括弧内の量はそれぞれn番目の順方向ベクトルまたは逆方向ベクトルの定義を満たすことになります。このように α と β を選択すると、括弧内のベクトルの和は単純になり、目的の結果が得られます。

これらの係数を求めるには、次のようになります。 α f n {\displaystyle \alpha _{f}^{n}} β f n {\displaystyle \beta _{f}^{n}}

f → n = α f n [ f → n − 1 0 ] + β f n [ 0 b → n − 1 ] {\displaystyle {\vec {f}}^{n}=\alpha _{f}^{n}{\begin{bmatrix}{\vec {f}}^{n-1}\\0\end{bmatrix}}+\beta _{f}^{n}{\begin{bmatrix}0\\{\vec {b}}^{n-1}\end{bmatrix}}}

および は それぞれ次のようになります。 α b n {\displaystyle \alpha _{b}^{n}} β b n {\displaystyle \beta _{b}^{n}}

b → n = α b n [ f → n − 1 0 ] + β b n [ 0 b → n − 1 ] . {\displaystyle {\vec {b}}^{n}=\alpha _{b}^{n}{\begin{bmatrix}{\vec {f}}^{n-1}\\0\end{bmatrix}}+\beta _{b}^{n}{\begin{bmatrix}0\\{\vec {b}}^{n-1}\end{bmatrix}}.}

前の 2 つの式に 1 を掛けると、次の式が得られます。 T n {\displaystyle {\mathbf {T} }^{n}}

[ 1 ε b n 0 0 ⋮ ⋮ 0 0 ε f n 1 ] [ α f n α b n β f n β b n ] = [ 1 0 0 0 ⋮ ⋮ 0 0 0 1 ] . {\displaystyle {\begin{bmatrix}1&\varepsilon _{b}^{n}\\0&0\\\vdots &\vdots \\0&0\\\varepsilon _{f}^{n}&1\end{bmatrix}}{\begin{bmatrix}\alpha _{f}^{n}&\alpha _{b}^{n}\\\beta _{f}^{n}&\beta _{b}^{n}\end{bmatrix}}={\begin{bmatrix}1&0\\0&0\\\vdots &\vdots \\0&0\\0&1\end{bmatrix}}.}

ここで、上記の 2 つのベクトルの中央にあるすべてのゼロは無視され、縮小されるため、次の式だけが残ります。

[ 1 ε b n ε f n 1 ] [ α f n α b n β f n β b n ] = [ 1 0 0 1 ] . {\displaystyle {\begin{bmatrix}1&\varepsilon _{b}^{n}\\\varepsilon _{f}^{n}&1\end{bmatrix}}{\begin{bmatrix}\alpha _{f}^{n}&\alpha _{b}^{n}\\\beta _{f}^{n}&\beta _{b}^{n}\end{bmatrix}}={\begin{bmatrix}1&0\\0&1\end{bmatrix}}.}

これらを（クラメールの 2×2行列の逆公式を使用して）解くと、新しい前方ベクトルと後方ベクトルは次のようになります。

f → n = 1 1 − ε b n ε f n [ f → n − 1 0 ] − ε f n 1 − ε b n ε f n [ 0 b → n − 1 ] {\displaystyle {\vec {f}}^{n}={1 \over {1-\varepsilon _{b}^{n}\varepsilon _{f}^{n}}}{\begin{bmatrix}{\vec {f}}^{n-1}\\0\end{bmatrix}}-{\varepsilon _{f}^{n} \over {1-\varepsilon _{b}^{n}\varepsilon _{f}^{n}}}{\begin{bmatrix}0\\{\vec {b}}^{n-1}\end{bmatrix}}}
b → n = 1 1 − ε b n ε f n [ 0 b → n − 1 ] − ε b n 1 − ε b n ε f n [ f → n − 1 0 ] . {\displaystyle {\vec {b}}^{n}={1 \over {1-\varepsilon _{b}^{n}\varepsilon _{f}^{n}}}{\begin{bmatrix}0\\{\vec {b}}^{n-1}\end{bmatrix}}-{\varepsilon _{b}^{n} \over {1-\varepsilon _{b}^{n}\varepsilon _{f}^{n}}}{\begin{bmatrix}{\vec {f}}^{n-1}\\0\end{bmatrix}}.}

これらのベクトルの和を求めると、前のベクトルからn番目の順方向ベクトルと逆方向ベクトルが得られます。あとはこれらのベクトルの最初のものを見つけ、簡単な和と乗算で残りのベクトルを求めるだけです。最初の順方向ベクトルと逆方向ベクトルは単純に次のようになります。

f → 1 = b → 1 = [ 1 M 11 ] = [ 1 t 0 ] . {\displaystyle {\vec {f}}^{1}={\vec {b}}^{1}=\left[{1 \over M_{11}}\right]=\left[{1 \over t_{0}}\right].}

後方ベクトルの使用

上記の手順により、 MのN個の逆ベクトルが得られます。そこから、より任意の式は次のようになります。

y → = M   x → . {\displaystyle {\vec {y}}=\mathbf {M} \ {\vec {x}}.}

解は、逆ベクトルの構築と同じ再帰的な方法で構築できます。したがって、 は中間値の列 へと一般化され、 となります。 x → {\displaystyle {\vec {x}}} x → n {\displaystyle {\vec {x}}^{n}} x → N = x → {\displaystyle {\vec {x}}^{N}={\vec {x}}}

次に、次のことに注意して再帰的に解を構築します。

T n − 1 [ x 1 n − 1 x 2 n − 1 ⋮ x n − 1 n − 1 ] = [ y 1 y 2 ⋮ y n − 1 ] . {\displaystyle \mathbf {T} ^{n-1}{\begin{bmatrix}x_{1}^{n-1}\\x_{2}^{n-1}\\\vdots \\x_{n-1}^{n-1}\\\end{bmatrix}}={\begin{bmatrix}y_{1}\\y_{2}\\\vdots \\y_{n-1}\end{bmatrix}}.}

次に、再度ゼロで拡張し、必要に応じてエラー定数を定義します。

T n [ x 1 n − 1 x 2 n − 1 ⋮ x n − 1 n − 1 0 ] = [ y 1 y 2 ⋮ y n − 1 ε x n − 1 ] . {\displaystyle \mathbf {T} ^{n}{\begin{bmatrix}x_{1}^{n-1}\\x_{2}^{n-1}\\\vdots \\x_{n-1}^{n-1}\\0\end{bmatrix}}={\begin{bmatrix}y_{1}\\y_{2}\\\vdots \\y_{n-1}\\\varepsilon _{x}^{n-1}\end{bmatrix}}.}

次に、 n番目の逆ベクトルを使用して誤差項を除去し、次のように目的の式に置き換えること ができます。

T n ( [ x 1 n − 1 x 2 n − 1 ⋮ x n − 1 n − 1 0 ] + ( y n − ε x n − 1 )   b → n ) = [ y 1 y 2 ⋮ y n − 1 y n ] . {\displaystyle \mathbf {T} ^{n}\left({\begin{bmatrix}x_{1}^{n-1}\\x_{2}^{n-1}\\\vdots \\x_{n-1}^{n-1}\\0\\\end{bmatrix}}+(y_{n}-\varepsilon _{x}^{n-1})\ {\vec {b}}^{n}\right)={\begin{bmatrix}y_{1}\\y_{2}\\\vdots \\y_{n-1}\\y_{n}\end{bmatrix}}.}

この方法をn = Nまで拡張すると、解が得られます。 x → {\displaystyle {\vec {x}}}

実際には、これらのステップは手順の残りの部分と同時に実行されることがよくありますが、一貫した単位を形成し、独自のステップとして扱う必要があります。

ブロックレビンソンアルゴリズム

M が厳密にテプリッツ行列ではなくブロックテプリッツ行列である場合、ブロックテプリッツ行列を行列要素を持つテプリッツ行列とみなすことで、レビンソン再帰をほぼ同様に導くことができます (Musicus 1988)。ブロックテプリッツ行列は、信号処理アルゴリズムにおいて、複数の信号ストリーム（例えばMIMOシステム）や周期定常信号を扱う際に自然に出現します。

参照

注記

  1. ^ Bojanczykら（1995年）。
  2. ^ ブレント（1999年）。
  3. ^ クリシュナ＆ワン（1993）。
  4. ^ 「アーカイブコピー」（PDF） 。 2012年3月25日時点のオリジナル（PDF）からアーカイブ。 2013年4月1日閲覧。{{cite web}}: CS1 maint: archived copy as title (link)
  5. ^ 「アーカイブコピー」（PDF） 。 2009年11月15日時点のオリジナル（PDF）からアーカイブ。 2009年4月28日閲覧。{{cite web}}: CS1 maint: archived copy as title (link)
  6. ^ 「アーカイブコピー」（PDF）saaz.cs.gsu.edu . 2007年4月18日時点のオリジナル（PDF）からアーカイブ。2022年1月12日閲覧。{{cite web}}: CS1 maint: archived copy as title (link)
  7. ^ 「アーカイブコピー」（PDF） 。 2006年9月5日時点のオリジナル（PDF）からアーカイブ。2006年8月15日閲覧。{{cite web}}: CS1 maint: archived copy as title (link)

参考文献

ソースの定義

  • Levinson, N. (1947). 「フィルタ設計と予測におけるウィーナーRMS誤差基準」. J. Math. Phys. , v. 25, pp. 261–278.
  • ダービン, J. (1960). 「時系列モデルのフィッティング」Rev. Inst. Int. Stat. , v. 28, pp. 233–243.
  • Trench, WF (1964). 「有限テプリッツ行列の逆行列を求めるアルゴリズム」. J. Soc. Indust. Appl. Math. , v. 12, pp. 515–522.
  • Musicus, BR (1988). 「Toeplitz行列およびほぼToeplitz行列に対するLevinson法と高速Choleski法」RLE TR No. 538, MIT. [1]
  • Delsarte, P. および Genin, YV (1986). 「分割レビンソンアルゴリズム」. IEEE Transactions on Acoustics, Speech, and Signal Processing , v. ASSP-34(3), pp. 470–478.

さらなる研究

  • Bojanczyk, AW; Brent, RP; De Hoog, FR; Sweet, DR (1995). 「Bareiss分解アルゴリズムおよび関連するToeplitz分解アルゴリズムの安定性について」. SIAM Journal on Matrix Analysis and Applications . 16 : 40–57 . arXiv : 1004.5510 . doi :10.1137/S0895479891221563. S2CID  367586.
  • Brent RP (1999)、「構造化線形システムの高速アルゴリズムの安定性」、構造付き行列の高速信頼性アルゴリズム（編集者 T. Kailath、AH Sayed）、ch.4（SIAM）。
  • バンチ, JR (1985). 「テプリッツ方程式系を解く方法の安定性」SIAM J. Sci. Stat. Comput. , v. 6, pp. 349–364. [2]
  • Krishna, H.; Wang, Y. (1993). 「スプリット・レビンソン・アルゴリズムは弱安定である」 . SIAM Journal on Numerical Analysis . 30 (5): 1498– 1508. doi :10.1137/0730078.

要約

  • Bäckström, T. (2004). 「2.2. レビンソン・ダービン再帰」.音声の線形予測モデル化 ― 制約と線スペクトルペア分解. 博士論文. 報告書番号71 / ヘルシンキ工科大学音響・音声信号処理研究所. エスポー, フィンランド. [3]
  • クレアバウト、ジョン・F. (1976). 「第7章 最小二乗法の波形応用」.地球物理学的データ処理の基礎. パロアルト: ブラックウェル・サイエンティフィック・パブリケーションズ. [4]
  • Press, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007)「セクション2.8.2. Toeplitz行列」、Numerical Recipes: The Art of Scientific Computing (第3版)、ニューヨーク: Cambridge University Press、ISBN 978-0-521-88068-8
  • Golub, GH, Loan, CF Van (1996). 「セクション4.7: Toeplitzと関連システム」『行列計算』 、ジョンズ・ホプキンス大学出版局
Retrieved from "https://en.wikipedia.org/w/index.php?title=Levinson_recursion&oldid=1318699639"