2-7-20
残差ベクトルの `1 ノルムを最小化する線形予測分析の解法∗◎亀岡弘和, 鎌本優, 原田登, 守谷健弘

(NTT コミュニケーション科学基礎研究所)
1序論

線形予測分析 (LPC) に基づく時系列信号の可逆圧
縮符号化方式では、LPC により生成される残差信号の符号化にエントロピー符号化 (Golomb-Rice 符号など) がしばしば用いられ、出現頻度の高い振幅値ほど短い符号が割り当てられる仕組みとなっている。
Golomb-Rice 符号 (以後、Rice 符号) の場合、残差振幅に割り当てられる符号量は、振幅の絶対値とほぼ比例関係にあるため、残差符号量は残差振幅の絶対値和である程度良く近似される。しかしながら、従来の LPC による可逆圧縮符号化方式では、残差振幅の二乗和を最小化するように予測係数を求めるため、残差符号量を直接的に最小化する規準での最適化となっていない。従って、予測係数の算出方法には符号量削減の余地が少なからず残されている。
LPC を含む線形回帰分析において、誤差ベクトルの要素の絶対値和を最小化する回帰係数の解法は、シンプレックス法に基づくもの [1] が古典的で有名である。しかし、この方法には、説明変数行列や回帰係数ベクトルのサイズが大きい場合に計算量が膨大となるという特徴があり、即時性が求められる圧縮符号化への応用には適していないという問題があった。
本稿では、誤差ベクトルの要素の絶対値和を最小化する回帰係数を効率的に求める手法の導出を出発点とし、Rice 符号量関数を絶対値和に比べてより良く表現するコスト関数をいくつか例示したのちに、同じ原理に基づいてこれらを最小化する回帰係数も同様に求められることを示す。

合を仮定した s の最尤解を求める問題である。すな
わち、`1 ノルム (絶対誤差)IX¯¯°°¯xi − (As)i ¯ (3)L(s) , kzk1 = °x − As°1 =i=1を最小化する s がこの場合の最尤解に対応する。この推定量のことは `1 推定量と呼び、x の外れ値に頑健なロバスト推定量の 1 つとして知られる。

2.2

線形予測分析 (LPC)

線形予測モデルは式 (1) の特殊例である。線形予測モデルは、あるディジタル信号の時刻 i におけるサンプル yi を、時刻 i − 1 から過去 J 個のサンプル
yi−1 , yi−2 , · · · , yi−J の線形結合で表現するモデル

yi ≈JXsj yi−jのことであり、自己回帰モデルとも呼ばれる。上記のモデルのもとで、I 個の連続するサンプル系列
y1 , · · · , yI を最も良く表す線形予測モデルの J 個の係数を推定する問題を線形予測分析と呼ぶ。式 (1) に対応させると、線形予測モデルは、
  
y1y0 y−1 . . . y1−J     s1  
y0 . . . y2−J   y2   y1..  =+ z (5).
..
..   
 ..   ..
.  .
.
.
   sJ
yIyI−1 yI−2 . . . yI−J2

絶対誤差を最小化する回帰係数の解法と置いたことに相当し、z が予測残差である。

2.1

線形回帰分析

2.3

x, s をそれぞれ I 次元および J 次元ベクトル、Aを I × J 行列とし、次のような線形回帰モデル
x = As + z(1)を仮定する。ただし、z は雑音を表す確率ベクトルである。x, A が与えられたもとで、s の最尤推定値を求める線形回帰分析の問題を考える。z にもし Gauss 性°°2
白色雑音を仮定すると、これは `2 ノルム °x − As°2を s に関して最小化する問題であり、I ≥ J ならば¡¢−1 Tbs = AT AA x(2)補助関数法の原理

提案法の特徴は、式 (3) の最小化に、補助関数法と呼ぶ原理を用いる点にある。補助関数法とは、適切に設計した補助関数を反復的に減少させていくことで目的関数を間接的に最小化する手法であり、EM アルゴリズムを一般化した原理である。EM アルゴリズム
同様に補助関数法はいかなる最適化問題に適用可能というわけではないが、ある要件を満たす補助関数が設計できれば大変効率の良い最適化アルゴリズムが導ける場合がある。以下に、補助関数法の原理を示し、補助関数が満たすべき要件を明らかにする。
目的関数 L(s) に対し、が最適解となる。これに対し、我々がここで考えるのは、z の各要素が互いに独立に Laplace 分布に従う場
∗(4)

j=1

L(s) = min L+ (s, B)
B(6)

A Solution to the Linear Predictive Coding Minimizing `1 Norm of the Residual Vector. by KAMEOKA
Hirokazu, KAMAMOTO Yutaka, HARADA Noboru and MORIYA Takehiro (NTT Communication Science Laboratories)

日本音響学会講演論文集- 707 -

2007年9月

が成り立つとき、L+ (s, B) を L(s) の補助関数、B を補助変数と定義すると、次の定理が成り立つ。
定理 1 (補助関数法). 補助関数 L+ (s, B) を、B に関して最小化するステップと、s に関して最小化するステップを繰り返すことで、目的関数 L(s) を単調減少させていくことができる。
証明： Bt+1 = argminB L+ (st , B) とし、st+1 =argmins L+ (s, Bt+1 ) とする。ただし、t は反復計
算 の ス テップ 数 を 表 す。s = st , B = Bt から s = st+1 , B = Bt+1 に更新されたときに、L(s) が非増加することを示す。明らかに L(st ) =L+ (st , Bt+1 ) であり、st+1 = argmins L+ (s, Bt+1 ) より、L+ (st , Bt+1 ) ≥ L+ (st+1 , Bt+1 ) である。さらに補助関数の定義より、L+ (st+1 , Bt+1 ) ≥ L(st+1 ) なのだから、結局、L(st ) ≥ L(st+1 ) である。
以上より、補助関数法を効果的に適用するには、次の 2 点を満たす補助関数を設計する必要がある。
+

1. argmin L (s, B) が解析的に求められる。
B

2. argmin L+ (s, B) が解析的に求められる。
s
これを踏まえ、式 (3) の最小化問題を解くための補助関数を以下で設計する。
2.4である。以上より、L+ (s, B) は補助関数の定義を満たす。次に、L+ (s, B) が補助関数の要件を満たすかどうかを確認する。

2.5

まず、L(s) = L+ (s, B) となる B は明らかに、

に b 6= 0 で接する原点に関して対称な 2 次関数はh0 (z) =z2|b|+2|b|2(8)で与えられる。∀z, f0 (z) ≤ h0 (z) であるので、
2¯ ¯¯z ¯ ≤ z + |b|2|b|2(9)が成り立つ。なお、言うまでもなく上不等式の等号成立は接点 z = b においてである。
この不等式を式 (3) に当てはめると、
!
Ã¡
¢2I
X
xi − (As)i|bi |+(10)L(s) ≤2|bi |2i=1¡¢T ¡¢= x − As B x − As + d ≡ L+ (s, B)のような不等式が立てられる。ただし、10
 2|b1 |
¢
1 ¡
..
B=
 , d = Tr B−1.
4
102|bI |日本音響学会講演論文集

(11)

1
 2|x1 −(As)1 |0


B=


..

.





(12)

12|xI −(As)I |0である。一方、B 固定のときに L+ (s, B) を最小化す
る s は、式 (12) より、大きい残差には小さい重みを、小さい残差には大きい重みを課す重みつき最小二乗解に相当し、正規方程式の解¡¢−1 T
s = AT BA
A Bx(13)で与えられる。AT BA は正定値対称行列であるので、¡ T¢−1A BAは Cholesky 分解で求めると良い。以上
より、L+ (s, B) は補助関数の要件を満たし、次の手順

1. s の初期値設定
2. s と B を収束するまで以下のように交互に更新
10
 2|x1 −(As)1 |


..
B←

.

102|xI −(As)I |¡ T¢−1 Ts ← A BAA Bx

2 次形式の補助関数の導出

式 (3) 最小化のための補助関数を設計するにあた
り、絶対値関数より小さくなることがなく (補助関数の定義より)、かつ、絶対値関数に接する (要件 1 より)2 次関数 (要件 2 より) を見つけたい。この中で 2次の係数が最小となるのは、原点に関して対称なものである。絶対値関数 f0 : R → R+¯ ¯
f0 (z) = ¯z ¯(7)補助関数法の更新則により、回帰係数の絶対誤差最小解を得ることができる。式 (3) は s に関して凸なので、収束解は大域最適解に一致する。ただし、大域最適解は一意とは限らず(大域最適解の集合は解空間上で凸集合をなし)、提案法はその中の 1 点を求めるものである。以上のアルゴリズムは、ロバスト推定の 1 つである Huber の M推定における反復重みつき最小二乗誤差推定 (IRLS)
法 [2] と同じ手続きをとる。`1 推定量を得るための上記手順においては Huber の IRLS 法との間に違いはないが、次章で両者の違いを明らかにする。

3

Rice 符号量モデルとその最小化
残差振幅 z ∈ Z に対する Rice 符号量 R(z) は、(¥R(z) =¦z2r−1 + r + 1¥ −z−1¦+r+12r−1(z ≥ 0)(z < 0)(14)と表せる。ただし、b·c は切捨て整数化演算子を表し、r ∈ N は Rice パラメータである。2 章で最小化しようとした残差振幅の絶対誤差は、r = 1 と置いた場合の Rice 符号量のモデルに相当する。しかし Rice 符号化では実際のところ r > 1 と設定することが多く、式(14) は図 1 のように階段形の関数となる。ここで特- 708 -2007年9月

Fig. 1

Rice 符号量 R(z) (r = 4)に注目すべきは、z = 0 周辺の一定範囲内の値には等しい符号量が割り当てられる点であり、このことは、残差振幅の絶対値をできるだけ小さくすることだけに専念するよりも、残差振幅に適度な「遊び」を持たせ、その自由度の分だけ極端な外れ値を少なくする方が全体の符号量をより小さくできることを示している。残差振幅が通常 0 付近に集中する点を考慮
すると、コスト関数をモデル化する際、z = 0 周辺のR(z) を良く近似したものであった方が良い。
以上の問題意識に基づき、絶対値関数に代わる次のような 3 種類の Rice 符号量モデル f1 (z), f2 (z), f3 (z)( 2z+ βγ (|z| ≤ 2β)
f1 (z) , 4βγ|z|/γ(|z| > 2β)
pz2 + β2f2 (z) ,γ(|z| ≤ β − c)β/γ

f3 (z) ,(z−|β−c|)2+ βγ4cγ

|z|/γ(|z| > β + c)

IX¡¢
fn xi − (As)i (n = 1, 2, 3)(15)

i=1が立てられる。いずれの目的関数においても s の最適解は γ に依らないので、以後簡単の為 γ = 1 とする。
なお、これらはいずれも s に関してやはり凸であり、局所解と大域最適解は一致する。
各々のコスト関数 fn (z) (n = 1, 2, 3) に対し、2 章と同様の方針で、2 次関数 hn (z)4(|b|−β)2(z−sign(b)(β−c))+β4c (z−b)2 + β4(β−|b|)(|b| ≥ 2β)(β + c < |b| < 2β)(β − c < |b| < β + c)
(0 ≤ |b| < β − c)が導出できる (導出過程は略)。hn (z) (n = 1, 2, 3) はいずれも fn (z) と点 b で接し、fn (z) より決して小さくなることのない 2 次関数の中で 2 次の係数が最小のものである。これにより適切な補助関数が設計でき、

1. s の初期値設定
2. s と B, w を収束するまで次のように交互に更新(n)(n)
b1 (s)
0w1 (s)
 . ..
B← , w ←  .. .
(n)(n)
0bI (s)wI (s)¡ T¢−1 TA B(x − w)s ← A BA
により、目的関数 Ln (s) を最小化する回帰係数を得(n)ることができる。行列 B の対角要素 bi (s) およびベ(n)クトル w の要素 wi (s) の更新則は、符号量モデルの選択に応じてmax(|b|, 2β)z2+2 max(|b|, 2β)
222zb + 2β 2h2 (z) = p+ p2 b2 + β 22 b2 + β 2

1(1), wi = 02 max(|zi |, 2β)1(2)(2)bi = p 2, wi = 022 zi + β(1)h1 (z) =日本音響学会講演論文集

Rice 符号量モデル fn (z) (n = 1, 2, 3)

h3 (z) =
 2z

+ |b|
 2|b|2
 (z−(2sign(b)β−b))2 + β(β − c < |z| ≤ β + c)を導入する。β/γ = r + 1、γ = 2r−1 であり、β, γは Rice パラメータに依存して一意に決まる定数である。c は 0 < c ≤ β の範囲で適当に定める定数であり、c → 0 で f3 (z) は Vapnik’s ε-insensitive [3] に漸近する。また、c = β で f1 (z) = f3 (z) となる。各々の概形を図 2 に示す。以上の符号量モデルを用い、残差符号の全符号量を近似した 3 種類の目的関数

Ln (s) =Fig. 2bi- 709 -=2007年9月

4(|zi |−β)
1
4c

1

4(β−|zi |)(|zi | ≥ 2β)

71.4(β + c ≤ |zi | < 2β)71.3(β − c ≤ |zi | < β + c)

(0 ≤ |zi | < β − c)0(|zi | ≥ 2β)2sign(z )β − z (β + c ≤ |z | < 2β)iii(3)wi =sign(zi )(β − c) (β − c ≤ |zi | < β + c)zi
(0 ≤ |zi | < β − c)と与えられる。ただし、zi = xi − (As)i である。ここでは、B だけでなく w も補助変数であることに注意が必要である。
Huber の M 推定における IRLS 法では目的関数
L3 (s) の単調減少が保証される更新則を導くことができない。この点が、補助関数法により導出された提案法と IRLS 法の原理的な違いを明確に物語っている。

4

評価実験

圧縮後のファイルサイズ
× 100 (16)
元のファイルサイズ

図 3 より、予想通り、提案法 2 による圧縮率が最も高く、次いで提案法 1 による圧縮率が高かった。以上

日本音響学会講演論文集

Conventional
Devised 1

71.2

Devised 271.1
71.0
70.9
70.8
70.7715
Predictive order

31Fig. 3 `2 ノルム規準の線形予測符号化 (従来法) に対する、L2 (s) および L3 (s) を規準とした線形予測符
号化 (提案法 1&2) による圧縮率の比較の実験結果から、提案法の有効性と、提案法により実現された残差符号量削減の度合いを確認することができた。

5提案法の効果を確認するため、音響信号を対象として圧縮性能の評価実験を行った。本実験の目的は、提案法の原理を実証することだけでなく、LPC に基づく可逆圧縮符号化において残差の符号量削減の余地がどの程度残されているのかを見極めることである。
評価実験用に、RWC 研究用音楽データベース [4]に含まれる、サンプリング周波数 44.1kHz、16 ビットで収録されたステレオの音楽ファイル (WAV 形式)10個 (RWC-MDB-P-2001 No.33∼42) を実験データ (合計約 432MB) として用いた。LPC 分析フレーム長は2048 サンプルとした。
`2 ノルムを規準とした LPC を Levinson-Durbin アルゴリズムにより行う従来の方式を従来法、従来法で得られた予測係数を初期値として、L2 (s) およびL3 (s) を規準とした LPC を提案する反復最適化アルゴリズムにより行う方式をそれぞれ提案法 1、提案法2 と略称する。提案法 1 および提案法 2 の反復計算回数は 10 とした。なお、各方式ではいずれもチャネル間相関に基づく圧縮は行っていない。また、LPC の予測次数 (式 (5) における J) は固定とした。
L2 (s) は `2 ノルムに比べて Rice 符号量に則した目的関数となっているため、従来法に比べ、提案法 1 の圧縮率は向上するものと予想される。一方、L3 (s) はL2 (s) に比べてより Rice 符号量に則した目的関数となっているため、提案法 1 に比べ、提案法 2 の圧縮率は若干向上するものと予想される。
従来法、提案法 1、提案法 2 それぞれについての、予測次数が 7, 15, 31 の場合の圧縮率をを図 3 に示す。
圧縮率は、次のように算出した。
圧縮率 (%) =Compression ratio (%)(3)bi =

12|zi |1

ま と め

本稿では、従来の LPC による可逆圧縮符号化は残差符号量を直接的に小さくするものではなかった点に着目し、残差符号量をできるだけ小さくする線形予測分析法の確立を目的とし、(1) 残差の符号化に利用される Rice 符号における振幅と符号量との関係をモデル化し、(2)Rice 符号量モデルに基づいて立てられる目的関数の最小解を効率的に求める反復アルゴリ
ズム (補助関数法) を考案した。評価実験により、提案法の動作と提案法による圧縮性能を確認した。
提案法は、ある条件を満たす目的関数を規準とする場合には Huber の M 推定における IRLS 法と同一の反復アルゴリズムを与えることが明らかになった。
しかしながら、IRLS 法では収束性が保証されない目的関数に対しても、提案法の枠組では収束性が保証された反復アルゴリズムを導けることを例示した。このことで、提案法は IRLS 法より幅広いクラスの目的関数を扱える方法論であることが示される。

参考文献[1] A. Charnes, W. W. Cooper, and R. O. Ferguson, “Optimal Estimation of Executive Compensation by Linear Programming,” Management Science, Vol. 1, No. 2, pp. 138–151, 1955.
[2] P. J. Huber, Robust Statistics, New York: Wiley, 1981.
[3] V. Vapnik, Statistical Learning Theory, New
York: Wiley, 1998.
[4] 後藤, 橋口, 西村, 岡, “RWC 研究用音楽データ
ベース：研究目的で利用可能な著作権処理済み楽曲・楽器音データベース,” 情処論文, Vol. 45,
No. 3, pp. 728–738, 2004.

- 710 -2007年9月