Proof of Lemma 1
3.4 アルゴリズムのまとめ
以上説明してきたサンプリングアルゴリズムをまとめ ると以下のようになる．
1. 初期化：z1, . . . , znを生成する((6)式)．また，(13) 式のグラム行列 K のうち，K1, K3, K5 を計算 する．
2. i= 1,2, . . . , nについて以下を実行する (a) zi を生成する((6)式)
(b) グラム行列K のk2, k3,k4 を計算する．
(c) µ, V を計算する((14),(16)式)．
(d) N[µ, V]に従ってf を生成する．
(e) (7)式の確率でzi を採択し，採択されればグ ラム行列を更新しzi をzi で置き換える．
3. ステップ2のサイクルを反復し，グリッド点にお けるf を生成する．
4 数値実験
人工的に生成したデータを用いて数値実験を行った．
データはまず(−6,6)区間に一様分布でz をn= 20点 生成し，各ziから分散σx2のノイズを加えてxi を生成 した．一方yi は関数
f(zi) = exp(−z2i/5) sin(zi) (17) に分散σ2y のノイズを加えて生成した．カーネル関数に はいずれもβ = 1のガウスカーネルを用いた．関数を 評価するグリッド点は(−6,6)の間に50個等間隔に置 いた．
σ2x=σ2y = 0.12 のときの実行例(図 1)を示す．提案 手法と比較する手法として，入力へのノイズを仮定しな い通常のカーネル最小二乗回帰によるあてはめを行った．
提案手法の結果は繰り返しのサイクルを 100 サイク ル行った単純期待値であり，burn-in等の処理は行って いない．ノイズレベルが比較的小さいのでカーネル回帰 との違いがわかりにくい(提案手法の実線とカーネル回 帰の破線が重なっている部分が多い)が，この図では左 端のデータがない部分でカーネル回帰のあてはまりが 悪い．そこで，フィッティングを定量的に評価するため に，データを生成した曲線との区間(−6,6)での二乗距
離(グリッド点における関数値の二乗平均)の変化をプ
ロットしたのが図2 である．この例では二乗距離で25
%程度の改善が見られる．
これらのデータはあくまでワンショットのあてはめな ので，乱数の種を変えて50回行ったときの提案手法と
−6 −4 −2 0 2 4 6
−2−1012
x
y
Target Regression PostMean PostStd
図 1: あてはめの例．σx2=σ2y= 0.12．提案手法はモン テカルロ法によるグリッド点 ({ui}mi=1)におけるf(ui) の期待値と標準偏差． 点線 (Target): データを生成 した f, 破線 (Regression): 通常のカーネル回帰, 実線 (PostMean): 提案手法（期待値），細い実線(PostStd):
エラーバー（期待値±標準偏差）
通常のカーネル回帰との比較を行った．結果は図3に示 したとおり，ほとんどの場合でカーネル回帰の結果を上 回っている．．
次に，ノイズレベルをσx2=σ2y= 0.32 とかなり大き くした場合にどうなるかについて同様の実験をやった結
果が図4,5,6である．図3と図6を比べてみるとノイズ
レベルが高い方が提案手法による改善が大きい様子がわ かる．
5 おわりに
測定誤差モデルにおいて，再生核ヒルベルト空間上の 関数の事後分布を生成するようなマルコフ連鎖モンテカ ルロ法を提案した．サンプリングの順序を入れ替えるト リックによって無限次元空間のサンプリングを有限次元 のサンプリングに変換することができた．
測定誤差モデルは非常に推定が困難なモデルとされて
いる[4]．今後は理論的あるいは実験的に収束性や近似
精度の分析が必要である．また，位相応答曲線の推定な ど実データを用いた有効性の検証も今後の課題である．
謝辞 本研究は科研費18500184の助成を受けたもので ある．
0 20 40 60 80 100
0.000.010.020.03
iteration
MSE
図2: 二乗誤差の推移(σx2=σy2= 0.32)．破線は通常の カーネル回帰による結果
0.02 0.04 0.06 0.08
0.020.040.060.08
regression
proposed
図 3: カーネル回帰と提案手法との比較 (50 回, σx2 = σy2= 0.12)
−6 −4 −2 0 2 4 6
−2−1012
x
y
Target Regression PostMean PostStd
図4: あてはめの例(σx2=σy2= 0.32)
0 20 40 60 80 100
0.000.020.040.060.080.100.12
iteration
MSE
図5: 二乗誤差の推移(σx2=σy2= 0.32)
0.05 0.10 0.15
0.050.100.15
regression
proposed
図 6: カーネル回帰と提案手法との比較 (50 回, σx2 = σy2= 0.32)
参考文献
[1] K. Nakae, Y. Iba, and T. Aoyagi. Statistical esti-mation of phase response curves and application in neural science. In Int. Soc. for Bayesian Anal-ysis 9th World Meeting, 2008.
[2] 中江, 伊庭, 青柳. レプリカ交換モンテカルロ法に よる位相応答曲線の統計的推定. 第11回情報論的 学習理論ワークショップ(IBIS), 2008.
[3] P. J. Bickel, C. A. J. Klaassen, Y. Ritov, and J. A.
Wellner. Efficient and Adaptive Estimation for Semiparametric Models. Springer-Verlag, 1993.
[4] 甘利,川鍋. 線形関係の推定−最小2乗法は最良で あるのか? 応用数理, Vol. 6, No. 2, pp. 96–109, 1996.
[5] S. M. Berry, R. J. Carroll, and D. Ruppert.
Bayesian Smoothing and Regression Splines for Measurement Error Problems. J. of the Ameri-can Statistical Association, Vol. 97, No. 457, pp.
160–169, 2002.
[6] 赤穂. カーネル多変量解析—非線形データ解析の 新しい展開. 確率の情報と科学.岩波書店, 2008.
[7] C. E. Rasmussen and C. K. I. Williams. Gaussian Processes for Machine Learning. MIT Press, 2006.
[8] A. J. Smola, S. V. N. Vishwanathan, and T. Hof-mann. Kernel methods for missing variables. In Proc. of 10th Int. Workshop on AI& Statistics 2005, pp. 325–334, 2005.
[9] 藤木,赤穂. 入力空間での計量に基づいた核主成分 分析.信学技報PRMU2008-121, Vol. 108, No. 327, pp. 69–74, 2008.
[10] 伊庭,種村,大森,和合,佐藤,高橋. 計算統計II. 統 計科学のフロンティア6.岩波書店, 2005.
情報論的学習理論テクニカルレポート2009
Technical Report on Information-Based Induc-tion Sciences 2009 (IBIS2009)
一次元正規分布のなす空間への曲線あてはめ
Curve fitting in the space of one-dimensional normal distribution
藤木 淳
∗Jun Fujiki
赤穂 昭太郎
†Shotaro Akaho
Abstract: We propose a method for extracting a one-dimensional structure from a set of parameters of one-dimensional normal distributions. In this paper, the one-dimensional structure is represented as a curve, which contains a linear parameter, on the manifold of one-dimensional normal distributions. And the fitting error is measured by metric tensor, which is the second order approximation of e-divergence and/or m-divergence. In this formulation, the estimation of curve fitting is represented by the framework of Jacobian kernel principal component analysis, which is the extension of kernel principal component.
Keywords: information geometry, manifold fitting, Euclideanization, kernel principal component analysis, Jacobian kernel
1 概要
1次元正規分布の平均と分散の組がが複数個与えられ たとき，それらは2次元空間をなすが，この2次元空間 の中から構造を抽出するという問題を考える．2次元空 間の部分空間であるから，構造としては0次元空間であ る点及び1次元空間である曲線が考えられる．
0次元構造は重心である．複数の1次元正規分布の重 心として，複数の1次元正規分布から測ったe-ダイバー ジェンスの和またはm-ダイバージェンスの和を最小に
するm-重心及びe-重心が提案されている[3]．1次元構
造は曲線であり，この構造を抽出する手法として，単 純には主成分分析（principle component analysis；
PCA）などが適用可能と考えられるが，PCAには以下
のような問題点がある．PCAは平均と分散の組が従う 分布がユークリッド（欧幾里得；欧氏） 空間中の正規 分布，つまり1次元正規分布の平均と分散の組が与える 計量が一定であることを仮定しているが，1次元正規分 布の平均と分散の組が与える自然な計量は，情報幾何学 的には一定ではない[4, 5]．また，PCAは線型構造は抽 出できるが，非線型構造は抽出に不向きである．
そこで本稿では，1次元正規分布の平均と分散をパラ メータとする2次元空間の点列への曲線あてはめによ
∗産業技術総合研究所 脳神経情報研究部門,〒305-8568茨城県つ くば市梅園1-1-1中央第2, e-mail [email protected],
The National Institute of Advanced Industrial Science and Tech-nology (AIST) Neuroscience Research Institute, Central 2, 1-1-1 Umezono, Tsukuba, Ibaraki 305-8568
†e-mail [email protected]
り1次元構造を探すことを考える．その際，確率分布 とあてはめ曲線の距離をダイバージェンスそのものでは なく，ダイバージェンスを2次形式で近似したもので測 る．このあてはめ問題は核主成分分析（kernel PCA； KPCA）[9]におけるあてはめ誤差を入力空間の計量で 近似したヤコビ核主成分分析（Jacobian kernel PCA； JKPCA）[8]の枠組みで捉えることができる．核関数 であてはめ問題を表現することにより多項式核に対応す る直線や2次曲線等の多項式曲線だけでなく，ガウス核 などの様々な核関数に対応する曲線があてはめ可能とな る．また，この結果はダイバージェンスを最小化するこ とによってe-平坦空間やm-平坦空間をあてはめる手法 [3]を含めたダイバージェンスを最小化する曲線あては めの初期値として利用することが可能である．