高次元データの次元縮小について 一般化主成分分析から代数曲線当てはめ

高次元データの次元縮小について

一般化主成分分析から代数曲線当てはめ

水田正弘

1∗ 1

_{北海道大学情報基盤センター}

Abstract: 多次元データの次元縮小は、データ解析において基本的な手法である。線形な次元縮小 として、因子分析、主成分分析、射影追跡法、独立成分分析などがよく知られている。これらは、線 形射影が基本となっている。非線形な次元縮小としては、多くの可能性が考えられる。  本報告では、Gnanadesikan[3] による一般化主成分分析法の理論的考察、改良、さらには、その限 界から代数曲線当てはめに発展させていった過程を紹介する。

1

はじめに

昨今のデータ解析において、扱うべきデータの構造 は複雑化し、かつそのサイズも巨大化している。もっ とも基本的な構造である高次元データでは、その次元 縮小が大きなテーマとなっている。因子分析や主成分 分析は、線形射影に基づく次元縮小の手法であり、広 く知られている。線形な手法には、実用上多くの利点 がある。第一に、アルゴリズムの構築が比較的容易で あり、一般に計算時間が短い。第二に、適切な解が得 られた (と思われる) 場合には結果の解釈が分かりやす い。さらに、多くの手法では検定方法が開発されてい る。すなわち、統計理論に裏付けされ、かつ計算コス トの点からも実用的な手法群が存在する。しかし、線 形な手法で「曲がった構造」、すなわち非線形な構造を 見出すことは困難である。そこで、非線形な手法が注 目されている。非線形な次元縮小を、曲線当てはめの 観点からまとめてみる。 多次元のデータに曲線を当てはめることは、データ 解析において最も基本的な問題の一つである。変量が 説明変量と目的変量に分けることができる場合には、 回帰直線やその拡張である様々な回帰曲線のための手 法が利用できる。本報告では、各変量が同等な場合に おける曲線当てはめについて述べる。 表 1: データへの曲線当てはめの手法の分類 データ 外的基準あり 外的基準なし 曲線の種類 直線   回帰直線* 主成分分析*   解析曲線 回帰分析* 一般化主成分分析 代数曲線当てはめ* 自由曲線 平滑化* 主座標曲線  *データ点と曲線との「距離の 2 乗和」が最小 ∗_{連絡先：北海道大学情報基盤センター先端データ科学研究室}       〒 060-0811 札幌市北区北 11 条西 5 丁目        E-mail: [email protected] 当てはめる曲線を直線 (または、平面・超平面) に限 定すると、主成分分析で実行できる。通常、主成分分 析の説明では、分散が最大となる、変量の 1 次結合を 第 1 主成分と定義し、データの分散共分散行列または 相関行列の固有値問題に帰着できると述べている。す なわち、最大固有値に対応する固有ベクトルを使った 1 次結合によりデータの (線形な) 散らばりが最大とな る新たな変量が構成される。 図 1: 直線当てはめと曲線当てはめ [5] ここで、逆に最小固有値に対応する固有ベクトルに 注目すると、1 次結合によりデータの (線形な) 散らば りが最小となる新たな変量が構成される。これはデー タに当てはまる直線と解釈できる。このことを説明の 都合上、2 次元データに限定して検討する。図 1-b を参 照されたい。 直線の式を ax + by + c = 0 (a 6= 0 または、b 6= 0) と置く。これは、x と y が対等に表される式の表現で あり、y = d + ex などの表現とは異なる。データを (xi, yi) (i = 1, 2, · · · , n) としたとき、データ点 (xi, yi)

SIG-DMSM-A701-05 (7/25)

人工知能学会研究会資料

と直線との距離 diは、 di= | ax√i+ byi+ c | a2_{+ b}2 となる。従って、距離 diの 2 乗和 s(a, b, c) は、 s(a, b, c) = n X i=1 d2i = 1 a2_{+ b}2 n X i=1 (axi+ byi+ c)2 となる。この s(a, b, c) が最小となる a, b, c を求める。 ただし、a, b, c が同時に定数倍 (6= 0) されても同じ直線 となるので、a2+ b2_{= 1 (a ≥ 0) という制約をつけて} も問題はない。従って、 s(a, b, c) = n X i=1 (axi+ byi+ c)2 を最小とする a, b, c を求める。上式を展開したあと、c について整理すると、 s(a, b, c) = n X i=1 {(axi+ byi) + c}2 = nc2+ 2c à a n X i=1 xi+ b n X i=1 yi ! + n X i=1 (axi+ byi)2 = n{c + (a¯x + b¯y)}2 + n X i=1 (axi+ byi)2− n(a¯x + b¯y)2 となる。ただし、¯x =Pn_i=1xi/n, ¯y = P_n i=1yi/n とす る。従って、任意の固定した a, b に対して、s(a, b, c) は c = −(a¯x + b¯y) のとき最小値をとる。そこで、 s(a, b, −(a¯x + b¯y)) = (n − 1)(a2σ2x+ 2abσxy+ b2σ2y) となる。ただし、 σ2x = 1 n − 1 n X i=1 (xi− ¯x)2, σxy = 1 n − 1 n X i=1 (xi− ¯x)(yi− ¯y), σ2 y = 1 n − 1 n X i=1 (yi− ¯y)2 とする。ここで、a2+ b2_{= 1 (a ≥ 0) より、} ( a = cos θ b = sin θ (−π/2 ≤ θ < π/2) と置き換えると、

s(a, b, −(a¯x + b¯y))

= (n − 1)(σx2cos2θ + 2σxysin θ cos θ + σ2ysin2θ) = (n − 1)(σ 2 x+ σy2) 2 +(n − 1) r (σ2 x− σ2y)2 4 + σxy2 sin (2θ + α) ただし、    sin α = _√ (σ2x−σy2)/2 (σ2 x−σy2)/4+σ2xy cos α = _√ σxy (σ2 x−σ2y)/4+σxy2 (−π ≤ α < π) とする。従って、sin (2θ + α) = −1 のとき、s(a, b, c) は最小となる。−π/2 ≤ θ < π/2, −π ≤ α < π より、 −2π ≤ 2θ+α < 2π である。ゆえに、sin (2θ + α) = −1 となる θ は θ = −π 4 − α 2, 3 4π − α 2 のうち、−π/2 ≤ θ < π/2 を満たすものである。そこ で、θ が得られる。従って、a, b, c が得られる。以上に より、データ点に対して当てはまる直線 ax + by + c = 0 が求まる。この直線と直交する方向、すなわちベクト ルで表現すると (a, b) の方向は第 2 主成分である。以 上の計算は、データの分散共分散行列の固有値を求め ることと同じである。 線形の場合には、「分散が最大である」ことと、「直 線と点との距離の総和が最小」であることは同値であ るが、非線形 (曲線) の場合には「分散が最大である」 ことと、「曲線と点との距離の総和が最小」であること は同値ではない (図 2)。従って、以上の議論を、図 1-d のように、直線ではなく曲線に適用しようとするとか なり難しい問題となる。

2

一般化主成分分析法

主成分分析法を非線形な手法に拡張する研究として、 Gnanadesikan and Wilks [4] は、一般化主成分分析法 を提案した。これは、p 次元データを k 個の関数によ り、k 次元空間に写像し、その空間で通常の主成分分 析を実行するものである (図 2)。 2 5 1 2 PCA Z Z R R y x F(x) 図 2: 2 次の一般化主成分分析の概念図

以下、p 変量 x = (x1, x2, · · · , xp) で表現された n 個 の個体からなるデータを扱う。fi(x)(i = 1, 2, · · · , k) を p 次元で定義された k 個の実数値関数とする。 一般化主成分分析法では、k 個の互いに独立な新た な変数 (x の関数)z1, z2, · · · , zk を見つけることからは じめる。ここで、各 zj(j = 1, 2, · · · , k) は、fi(x)(i = 1, 2, · · · , k) の線形結合であり、 zj= k X i=1 lijfi(x) = l>jf (x), とする。ここで、lj= (l1j, l2j, · · · , lkj)>は、k 個のベク トルで、l>jlj= 1 および f (x) = (f1(x), f2(x), · · · , fk(x))> を満たすものとする。このとき、通常の主成分分析法と 同様な考察により、ベクトル l1, l2, · · · , lkは、(f1(x), f2(x), · · · , fk(x)) の分散共分散行列の固有ベクトルと なる。特に、最小固有値に対応した関数 zkはデータに 対する当てはめの関数と見なすことができる。 通常の主成分分析法は、実数値関数 fi(x) を xi(i = 1, 2, · · · , p) としたときの一般化主成分分析法である。 p 次元データに対する 2 次の一般化主成分分析法と は、以下の実数値関数で定義される。 ( fi(x) = xi (i = 1, 2, · · · , p) fi(x) = xjxm (i = p + 1, · · · , (p2+ 3p)/2), ここで、j と m は、i(i = p + 1, · · · , (p2+ 3p)/2) に対 して i = {(2p − j + 3)j/2} + m − 1, 1 ≤ j ≤ m ≤ p, によって一意的に決定される。 例えば、2 次元データに対する 2 次の一般化主成分 分析では f1(x, y) = x f2(x, y) = y f3(x, y) = x2 f4(x, y) = xy f5(x, y) = y2. を利用して、2 次元データ (x, y) を (x, y, x2, xy, y2_{) の} 5 次元データに拡張した後に、通常の主成分分析を適 用する。すなわち、 t1x + t2y + t3x2+ t4xy + t5y2, の分散を制約条件、 t2 1+ t22+ t23+ t24+ t25= 1, のもとで最小にする t1, t2, t3, t4, t5を求める。これは、 (x, y, x2_{, xy, y}2_{) の分散共分散行列の固有値問題を解く} ことで得られる。そのとき 2 次曲線 t1x + t2y + t3x2+ t4xy + t5y2= t, を当てはめの曲線と見なすことができる。ただし、t は、 t1x + t2y + t3x2+ t4xy + t5y2 の平均値である。

2.1

座標系の直交変換に対して不変な手法

の決定

一般化主成分分析では、座標系の回転 (直交変換) や 平行移動に対して「結果が不変」とは限らない。ここ で、「結果が不変」であるとは、新たな座標系で手法を 適用した後、元の座標系に戻した場合、元の座標系で手 法を適用した結果と一致することを言う。これは、デー タの次元縮小や曲線当てはめの観点から自然な性質で ある。先の 2 次元データに対する 2 序の一般化主成分 分析法でも、座標系の回転や平行移動に対して結果が 不変ではない。そこで、座標系の直交変換に対して結 果が不変である一般化主成分分析法を検討する。また、 平行移動に関して結果が不変な手法に関する結果につ いては、以下では省略する。 以下では、以下の条件を仮定する。 A1 f1(x), f2(x), · · · , fk(x) は、x の関数として線形 独立である。 A2 任意の直交行列 T に対して、f (T x) ≡ W f (x) を満たす行列 W が存在する。 A3 fi(x) は連続関数である。 条件 A1 と A3 は、一般化主成分分析においては、自 然な条件である。また、条件 A2 は、座標変換を考察 する上で、必要最小限の条件である。(一般化主成分分 析法の特殊な場合としての) 通常の主成分分析法、2 次 の一般化主成分分析法では、これらの 3 つの条件を満 たしていることは明らかである。 「結果が不変」であることを数学的に記述する。任 意の直交変換を実行したデータ x∗= T x について、一 般化主成分分析を提要した結果を z∗ j = l∗>j f (x∗) = l∗>j f (T x)(j = 1, 2, · · · , k) とする。ここで、l∗jを、Cov(f (x∗)) の固有ベクトルと する。f で規定される一般化主成分分析法の「結果が 不変」であるとは、任意の直交行列 T に関して、 l> jf (x) ≡ ±l∗>j f (T x)(j = 1, 2, · · · , k)

が、x の関数として成立することである。ただし、正 負の符号は、固有ベクトルの向きの自由性を表わして いるだけで本質的ではない。 f (x) による一般化主成分分析法が直交変換に関して 不変である必要十分条件は、任煮の直交行列 T に対し て W が直交行列となることである。すなわち、もし手 法が不変であれば、行列 W は、 (l∗₁, l∗₂, · · · , l∗_k)(l1, l2, · · · , lk)>, となり、これは直交行列である。逆に、もし、W が直 交行列ならば、W>l∗jは、Cov(f (x)) の固有ベクトル である。従って、 l> j = ±l∗>j W が得られる。 Mizuta[8] は、一般化主成分分析法が直交変換に関し て不変である必要十分条件を示した。 定理 条件 A1,A2,A3 のもとで、2 次元データ (x, y) に関す る一般化主成分分析法で直交変換に関して不変である ものは、以下の関数により規定される手法に限られ、か つこれ以外のものは存在しない。 (1) s 組の関数: f2i−1(x, y) = gi( p

x2 + y2)¡xNi −¡ Ni₂ ¢y2xNi−2 +¡ Ni₄ ¢y4xNi−4 − · · ·¢ −hi(

p

x2 + y2)¡NiyxNi−1 −¡ Ni3 ¢ y3xNi−3 +¡ Ni₅ ¢y5xNi−5 − · · ·¢ f2i(x, y) = gi( p x2 + y2) ¡ NiyxNi−1 − ¡ Ni 3 ¢ y3xNi−3 + ¡ Ni 5 ¢ y5xNi−5 − · · · ¢ +hi( p x2 + y2) ¡ xNi − ¡ Ni 2 ¢ y2xNi−2 + ¡ Ni 4 ¢ y4xNi−4 − · · · ¢ (i = 1, 2, · · · , s), ただし、gi, hiは、 p x2_{+ y}2_{に関する任意の連続関数} であり、Niは任意の正の整数とする。 (2)px2_{+ y}2_{に関する連続関数。} p 次元空間の座標系の直交変換は任意の 2 つの座標 軸に関する回転の合成で表現できるので、上記の定理 は、一般の p 次元データに対して拡張できる。 直交変換に関して不変な一般化主成分分析法を規定 するいくつかの代表的な関数をまとめておく。 (1) 3 次元における 1 次関数 (通常の主成分分析法): x, y, z. (2) 3 次元における 2 次関数: x2_{, y}2_{, z}2_,√_2xy,√_2yz,√_2zx. (3) 3 次元における 3 次関数:

x3_{, y}3_{, z}3_,√_3x2_y,√_3y2_z,√_3z2_x,√_3xy2_,√_3yz2_,√_3zx2_,√_6xyz.

(4) 3 次元における q 次関数: q q! i!j!k!xiyjzk (i + j + k = q; 0 ≤ i, j, k). (5) p 次元における q 次関数: r q! Qp i=1kt! Qp t=1(xt)kt Pp t=1kt= q; 0 ≤ kt.

3

代数曲線当てはめについて

一般化主成分分析や後述する Principal Curves (主 座標曲線) はデータ点と曲線との（集合としての）距離 の 2 乗和が最小になっているわけではない。すなわち、 これらの手法によって、データに曲線を当てはめるこ とは可能であるが、得られた曲線は必ずしも垂直な距 離の 2 乗和を最小にしていない。Gnanadesikan[3] も、 「“ 垂直な ”距離の 2 乗和を最小にすることにより、特 定の型の非線形関数（たとえば 2 次多項式）でもよい から、データに当てはめるためのアルゴリズムを開発 しておくことは非常に有益であろう」と指摘している。 そこで、2 次元データに対して垂直な距離の 2 乗和 が最小となる代数曲線を求めるアルゴリズムを検討す る。ここで、代数曲線とは多項式の値が 0 となる集合 として定義される。すなわち、多項式 f (x) = f (x, y) で定義される代数曲線は、 Z(f ) = {x ∈ R2: f (x) = 0} = {(x, y) : f (x, y) = 0} と表される。例えば、f (x, y) = x2+ y2− 1 とすると、 半径 1 の円となり、f (x, y) = x2− y とすると放物線と なる。 Z(x2_{+ 2y}2_{− 1)} Ellipse -10 -5 5 10 -10 -5 5 10 Z(y2_{− x}2_{+ 1)} Hyperbola

-1 -0.5 0.5 1 -0.4 -0.2 0.2 0.4 Z(y2_{− x}2_{+ x}4₎ Eight Curve 0.2 0.4 0.6 -0.4 -0.2 0.2 0.4 Z((x2_{+ y}2_{) − (4x}2_{− 2x + 4y}2₎2₎ Limacon of Pascal 図 3: 2 次元空間における代数曲線の例 また、3 次元空間においては、2 つの多項式が同時に 0 となる集合として代数曲線が定義される (図 4)。 図 4: 3 次元空間における代数曲線

3.1

点と曲線との距離および近似距離

n 個のデータ点を (αi, βi)(i = 1, 2, · · · , n) とする。多 項式 f (x, y) によって定義される代数曲線を Z(f ) = {(x, y) : f (x, y) = 0} = {x : f (x) = 0} とする。このとき、ある点 a = (α, β) と、曲線 Z(f ) との集合としての距離は、 dist(a, Z(f )) = inf{k a − y k: y ∈ Z(f )} (1) によって定義される。しかし、一般の代数曲線に対し て、この距離を求めることは困難であるとされてきた。 Z(f)={(x,y) : f(x,y)=0} (x,y) (α,β) 図 5: 曲線と点との距離 そこで、Taubin(1991) は、点 a に最も近い Z(f ) の 点 y の近似として ˆ y = a − (∇f (a)t)+f (a) (2) = (α, β) − 1 k ∇f (a) k2( ∂f ∂x, ∂f ∂y) |(α,β)f (a) を利用し、 dist(a, Z(f ))2≈ f (a) 2 k ∇f (a) k2 を提案した。ここで、(∇f (a)t)+_{は ∇f (a)}t_の一般化 逆行列である。さらに Taubin は、この近似距離を基準 として、データ点からの近似距離の 2 乗の総和が最小 となる代数曲線の導出法を示した。 Z(f) f(x) Tangent Space Exact Distance Data Point Approximate Distance 図 6: 近似距離 しかし、この近似距離と (1) による距離には差があ る。これを簡単な例で示す。x2+ 2y2_{= 1 で定義され} る楕円と、x 軸上の点 (0, e), e > 1 の (1) による距離は、 明らかに e − 1 となる。これを近似距離で求めると、表 2 のようになる。e が 1 に近い場合には比較的よい近似 となっているが、e の値が大きくなるにつれ近似の精度 が悪くなっていく。 そこで、(1) による距離を求める方法および、その距 離に基づいた曲線当てはめについて考察する。

図 7: 楕円と点との距離 表 2: 正確な距離と近似距離 e 正確な 近似 距離 距離 1.1 0.100000 0.095455 1.2 0.200000 0.183333 1.3 0.300000 0.265385 1.4 0.400000 0.342857 1.5 0.500000 0.416667 1.6 0.600000 0.487500 1.7 0.700000 0.555882 1.8 0.800000 0.622222 1.9 0.900000 0.686842

3.2

点と曲線との距離の導出法

まず、点 a = (α, β) と代数曲線 Z(f ) との (集合とし ての) 距離を求める方法を提案する。多項式 f (x, y) の 次数が２次の場合には、通常のラグランジュの乗数法 により、 L(x, y, λ) = (x − α)2+ (y − β)2+ λf (x, y) の x, y, λ による偏微分から得られる連立方程式、      2(x − α) + λ∂f_∂x = 0, 2(y − β) + λ∂f_∂y = 0, f (x, y) = 0 (3) を変形することにより、λ の 4 次方程式が得られ、正 確な解を求めることができる。 多項式 f (x, y) の次数が 3 次以上の場合には、上記の 連立方程式を直接解くことは困難である。そこで、拡 張ラグランジュ乗数法 ([17],p.68 など) を利用して求め ることを考察した。 Q(x, µ, r) = (x−α)2_{+(y −β)}2_{+µf (x, y)+}1 2rf (x, y) 2 と置き、以下のアルゴリズムを実行する。 Step 0   µ0, r0, x0を初期値として与え、k = 0 とす る。 Step 1   Q(x, µk, rk) を x について最小化し、その 点を xkとする。 Step 2   f(x, y)2_{が十分小さくなったら停止する。} Step 3   rk+1= crk (cは事前に与えられた1以上の定 数) Step 4   µk+1= µk+ rkf (x, y) Step 5   k = k + 1 として、Step 1 へいく。 一般に、点との距離が極小になる曲線上の点は複数 存在する。最小の距離を得るためには、初期値 x0の与 え方が重要になる。ここでは (2) 式の値を採用するこ とにより、多くの場合に最小な距離の有する点を得る ことができる。 また、Step 1 は無制約最小化問題であり、Q の x に 関する微分は、 ∂ ∂xQ(x, µ, r) = 2(x − α) + (µ + rf (x, y)) ∂f ∂x ∂

∂yQ(x, µ, r) = 2(y − β) + (µ + rf (x, y)) ∂f ∂y から容易に得られるので、準ニュートン法 [18] を用い ることができる。

3.3

2 乗距離の総和が最小になる曲線の導出

前節により、データ (αi, βi)(i = 1, 2, · · · , n) と多項 式 f に対して、データ点と代数曲線 Z(f ) との (集合と しての) 距離を求めることができる。そこで、2 乗距離 の総和 R が最小になる代数曲線を求めることが目標と なる。R の多項式の係数に関する微分を用いない直接 法によっても解くことができる。しかし、計算の効率 化のために微分を求めておく。 p 次多項式 f を f (x, y) = AtX とおく。ただし、 X = (1, x, y, x2_{, xy, · · · , y}p₎t_, A = (a1, a2, · · · , aq)t, q = (p + 2)(p + 1) 2 とする。また、各データ点 (αi, βi) から最も近い Z(f) 上の点を (xi, yi) とおく。 R = n X i=1 ¡ (xi− αi)2+ (yi− βi)2 ¢ , (4) より、 ∂R ∂a = 2 n X i=1 µ (xi− αi)∂xi ∂a + (yi− βi) ∂yi ∂a ¶ , (5)

となる。 ここで、f (xi, yi) = AtXi = 0 の両辺を aj (j = 1, 2, · · · , q) で微分すると、 At∂Xi ∂aj + ( ∂A ∂aj) t_X i= 0 (6) が得られる。 (αi, βi) は (xi, yi) からの法線上にあるので、 (yi− βi) ∂f ∂x ¯ ¯ ¯ ¯ (xi,yi) = (xi− αi) ∂f ∂y ¯ ¯ ¯ ¯ (xi,yi) (7) となる。以下、添字 i を省略する。 (y − β)At∂X ∂x = (x − α)A t∂X ∂y (8) の両辺を ajで微分すると、 ∂y ∂ajA t∂X ∂x + (y − β) µ (∂A ∂aj) t∂X ∂x + A t ∂2X ∂x∂aj ¶ = ∂x ∂ajA t∂X ∂y + (x − α) µ (∂A ∂aj) t∂X ∂y + A t ∂2X ∂y∂aj ¶ (9) となる。ここで、 ∂X ∂aj = ∂X ∂x ∂x ∂aj +∂X ∂y ∂y ∂aj ∂2_X ∂x∂aj = ∂2X ∂x2 ∂x ∂aj + ∂2X ∂x∂y ∂y ∂aj ∂2_X ∂y∂aj = ∂2_X ∂x∂y ∂x ∂aj + ∂2_X ∂y2 ∂y ∂aj ∂X ∂x = (0, 1, 0, 2x, y, 0, 3x 2_{, 2xy, y}2_{, 0, · · ·)}t ∂X ∂y = (0, 0, 1, 0, x, 2y, 0, x 2_{, 2xy, 3y}2_{, · · ·)}t ∂2_X ∂x2 = (0, 0, 0, 2, 0, 0, 6x, 2y, 0, 0, · · ·) t ∂2_X ∂x∂y = (0, 0, 0, 0, 1, 0, 0, 2x, 2y, 0, · · ·) t ∂2_X ∂y2 = (0, 0, 0, 0, 0, 2, 0, 0, 2x, 6y, · · ·) t

より、(6) 式と (9) 式は、_∂a∂x_j と _∂a∂y_j に関する連立 1 次 方程式                  At ∂X ∂x∂a∂xj + A t ∂X ∂y ∂y ∂aj = −( ∂A ∂aj) t_X At₍₋∂X ∂y − (x − α)∂ 2_X ∂x∂y + (y − β)∂ 2_X ∂x2)_∂a∂x_j +At³∂X ∂x − (x − α)∂ 2_X ∂y2 + (y − β)∂ 2_X ∂x∂y ´ ∂y ∂aj = (∂A ∂aj) t³_{(x − α)}∂X ∂y − (y − β)∂X∂x ´ (10) となる。従って、(xi, yi) における_∂a∂x_j と_∂a∂y_j の値が求

まり、(5) 式より R の ajに関する偏微分が得られる。 この偏微分を利用して、準ニュートン法により R が 最小になる代数曲線を求めることができる。準ニュー トン法の初期値 (初期曲線) としては、一般化主成分分 析の解などを利用できる。 例えば、代数曲線の次数 p = 2 の代数曲線 f (x, y) = g + ax + by + cx2_{+ dxy + ey}2 _{= 0 に対して距離の 2} 乗の総和 R の a による偏微分は、 (xi − αi)∂xi

∂a + (yi − βi) ∂yi

∂a =

¡

−((ax − by + bβ + 2cx2 + dxy − 2x2e + 2xαe − 2y2e + 2yβe)(x − α)

+ (ay − aβ + bx + 2dx2 − dxα + dy2 − dyβ + 2xye)(y − β)) ¢

/

¡a2 + 4acx + 3ady − adβ − 2axe + 2aαe + b2 − 2bcy + 2bcβ

+ 3bdx − bdα + 4bye + 4c2x2 + 4cdxy − 4cx2e + 4cxαe − 4cy2e + 4cyβe + 2d2x2 − d2xα + 2d2y2 − d2yβ + 4dxye + 4y2e2¢

を (5) 式に代入することで得られる。

3.4

数値例

56 個の 2 次元人工データに対して本手法を適用した。 初期曲線として、一般化主成分分析による 2 次曲線を 与え、1 次∼4 次の代数曲線を求めた。R の値を表 3 に まとめる。また、2 次の曲線の形状を図 9 に、3 次の曲 線を図 10 に示す。 表 3: 人工データに対する代数曲線の当てはめ 次数 R の値 1 0.13891422 2 0.00484526 3 0.00397380 4 0.00396130 0.2 0.4 0.6 0.8 0.6 0.8 1.0 1.2 1.4 • • •• • • •• •• •• • • ••• • •• •• ••• • •• • •••• •• ••• • ••• •••• •_{• •••}• • • • • 図 8: 2 次の代数曲線 (人工データ)

0.2 0.4 0.6 0.8 0.6 0.8 1.0 1.2 1.4 • • • • • • •• •• • •• • ••• • •• ••• •• • •• • •••• •• ••• • ••• •••_••_{• •••}• • •_{• •} 図 9: 3 次の代数曲線 (人工データ) 以上の議論は、3 次元空間におけるデータに対して も適用できる。3 次元空間における代数曲線は 2 つの 多項式によって定義できる。また、代数曲面は 1 つの 多項式によって定まる。さらに、代数曲面が有界であ るとの制約条件のもとで当てはめを行なう方法も研究 されている [11][15]。以下に例を示す。 (a) データ点 -1 0 1 -1 0 1 -1 0 1 -1 0 1 -1 0 1 (b) 代数曲面による当てはめ -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 (c) 有界な代数曲面による当てはめ 図 10: 3 次元データに対する代数曲面当てはめ

4

おわりに

当初、「非線形」という言葉に引かれ、一般化主成分 分析を研究を推進してきた。これは、最近の用語でい うと、カーネル法である。Gnanadisikan[3] に書かれて いる一般化主成分分析法では、座標系の変換に対して 不変性を持たないことに気が付き、その解決を模索し た。座標系の直交変換に関しては、本報告で述べたと おり、ほぼ完全な定理を導出することができた。平行 移動に関しても、目的関数を少し変形することで、自 分なりに納得できる結果を得ることができた。しかし、 一般化主成分分析法により得られた結果の解釈方法に ついては、疑問が残ったままであった。 いろいろ検討した結果、データ点との「距離」の 2 乗和が最小となる曲線の導出が自然であるとの結論を 得た。そこで、一般化主成分分析法の結果を代数曲線 とみなすことにより、解釈可能な定式化をすることが できた。すなわち代数曲線当てはめに課題を移した。2 次の代数曲線ですら、この定式化のもとで問題を解く のは時間がかかった。その後、数値最適化法などによ り、一般的な代数曲線当てはめがほぼ計算可能となっ た。しかし、代数曲線 (代数曲面) の不安定さなど多く の問題が残ったままである。今回は、現時点までの研 究過程のみを報告させていただいた。

参考文献

[1] Berry, M. J. A.(2000). Mastering Data Mining. John Wiley & Sons.

[2] Cleveland, W. S.(1979). Robust Locally Weighted Regression & Smoothing Scatterplots, J. Amer. Stat. Assoc., Vol.74, No. 368, 829-836.

[3] Gnanadesikan, R. (1977). Methods for Statistical Data Analysis of Multivariate Observations. John Wiley & Sons (丘本 正・磯貝恭史(1979).統計的多

変量データ解析.日科技連).

[4] Gnanadesikan, R. and Wilk, M. B. (1969), Data Ana-lytic Methods in Multivariate Statistical Analysis, In Multivariate Analysis II, Krishnaiah, P. R. ed. Aca-demic Press, 593–638.

[5] Hastie, T. and Stuetzle, W. (1989). Principal Curves. J.Amer.Statist.Assoc.,84, 502–516.

[6] Keren, D., Cooper, D. and Subrahmonia, J. (1994). Describing Complicated Objects by Implicit Polyno-mials, IEEE trans. Patt. Anal. Machine Intell., 16, 1, 38–53.

[7] Kriegman, D. J. and Ponce, J.(1990). On Recogniz-ing and PositionRecogniz-ing Curved 3D Objects from Image Contours, IEEE trans. Patt. Anal. Machine Intell., 12, 12, 1127–1137.

[8] Mizuta, M. (1984). Generalized Principal Compo-nents Analysis Invariant under Rotations of a Co-ordinate System, J. Japan Statist. Soc.,14, 1–9. [9] Mizuta, M. (1995). A Derivation of the Algebraic

Curve for Two-Dimensional Data using the Least-Squares Distance. In Data Science and Its Applica-tion, Hayashi, C. et al. (eds.), 167–176, Academic Press, Tokyo.

[10] Mizuta, M.(1996). Algebraic Curve Fitting for Multi-dimensional Data with Exact Squares Distance. Pro-ceedings of 1996 IEEE International Conference on Systems, Man and Cybernetics, 516–521.

[11] Mizuta, M.(1998). Bounded Algebraic Curve Fitting for Multidimensional Data Using the Least-Squares Distance. Data Science, Classification, and Related Methods, Springer, 610–616.

[12] Mizuta, M.(1999). Sliced Inverse Regression with Projection Pursuit. Applied Stochastic Models and Data Analysis – Quantitative Methos in Business and Industry Society – edited by: H.Bacelar-Nicolau, F.Costa Nicolau and J.Janssen INSTITUTO NA-CIONAL DE ESTAT´ISTICA), 51–56.

[13] Taubin, G.(1991). Estimation of Planar Curves, Sur-faces, and Nonplanar Space Curves Defined by Im-plicit Equations with Applications to Edge and Range Image Segmentation, IEEE Trans. Patt. Anal. Machine Intell., Vol.13, No.11, 1115–1138.

[14] Taubin, G.(1994). Distance Approximations for Ras-tering Implicit Curves. ACM Trans. on Graphics, 13, 1, 3–42.

[15] Taubin, G., Cukierman, F., Sullivan, S., Ponce, J. and Kriegman, D. J. (1994). Parameterized Fami-lies of Polynomials for Bounded Algebraic Curve and Surface Fitting, IEEE trans. Patt. Anal. Machine In-tell., 16, 3, 287–303.

[16] Zhou, Y. Z.(1990) Fitting Smooth Curves, Proc. 10th Int. Conf. on Pattern Recognition, 455–459.

[17] ASNOP研究会(1991).非線形最適化プログラミング, 日刊工業新聞社. [18] 茨木俊秀・福島雅夫(1991).最適化プログラミング,岩 波書店. [19] 西川 一・三宮信夫・茨木俊秀(1982).最適化,岩波書店. [20] 水田正弘(1992). 最小2乗距離を有する2次曲線の導 出法,計算機統計学,第5巻第2号, 107-115. [21] 水田正弘(1996).データ点と代数曲線との正確な距離の 効率的な導出法,日本行動計量学会第24回大会発表論 文抄録集, 210–211. [22] 水田正弘・馬場康維(1990). Principal Curvesについ て,第58回日本統計学会講演報告集, 244–246. [23] 水田正弘・馬場康維(1993). Principal Curvesと数量化 3類を用いた質的データの1次元構造の抽出,統計数理, 第41巻第1号, 1–11. [24] 山下信之・南 弘征・水田正弘・佐藤義治(1992).プリ ンシパル曲線のアルゴリズムの改良とその計算量の評 価,計算機統計学,第5巻第1号, 33–43.