統計数理（2006）
第 54 巻 第 1 号 177–190
c 2006 統計数理研究所
［原著論文］

von Mises 分布における経験 Bayes 推定
1

2

叶 雄 ・大西 俊郎（受付 2005 年 7 月 20 日；改訂 2005 年 10 月 20 日）

要旨

von Mises 分布の共役解析において，事後モードがエントロピー型損失関数の下で最適，すなわち，Bayes 推定量になっていることが知られている．本論文では，von Mises 事前分布の超パラメータのうち位置パラメータを既知とし拡散パラメータを周辺分布に対するモーメント法
で推定することにより経験 Bayes 推定量を構成した．様々な設定の下でシミュレーションを行い，経験 Bayes 推定量のパフォーマンスを最尤推定量と比較評価したところ，パラメータが高次元であるとき実際的な状況下で経験 Bayes 推定量が最尤推定量を圧倒的に優越することが分かった．

キーワード： 共役事前分布，事後モード，ピタゴラス関係．
1.

はじめに

von Mises 分布は，角度データを解析するとき標準的な分布として使用されている．ニーズの高さを考えると非常に重要な分布の一つである．円周上の正規分布と呼ばれることもあり，その密度関数はp(x; µ, τ ) ∝ exp{τ cos(x − µ)}で与えられ，位置パラメータ µ および拡散パラメータ τ をもつ．拡散パラメータ τ が既知定数のとき，von Mises 分布は曲指数型分布族に属する．
共役解析は，Bayes 統計学の一分野として重要な役割を果たしてきた．指数型分布族の共役解析が代表的であり，多くの研究者によって様々な研究がなされている．一方，いくつかの曲指数型分布族でも共役事前分布が存在する．実際，Mardia and El-Atoum（1976）が示したように，拡散パラメータ τ が既知の場合の von Mises 分布においても共役解析が可能である．vonMises 分布の共役解析を論じた論文としては，Guttorp and Lockhart（1988），Bagchi（1994）および Rodrigues et al.（2000）などがある．
本論文では，von Mises 分布において共役解析をベースにした経験 Bayes 推定量を構成する．
具体的には，超パラメータに対する推定方程式を明示的に与え，これを解くことにより経験
Bayes 推定量を導出する．また，そのパフォーマンスをシミュレーションにより評価する．

1
2

総合研究大学院大学 複合科学研究科統計科学専攻：〒106–8569 東京都港区南麻布 4–6–7
統計数理研究所：〒106–8569 東京都港区南麻布 4–6–7

178

2.

統計数理 第 54 巻

第 1 号 2006

共役解析

確率変数 x は次の密度関数をもつ von Mises 分布に従うとする．
1(2.1)
p(x; µ, τ0 ) =exp{τ0 cos(x − µ)}, 0 ≤ x < 2π.
2πI0 (τ0 )ここで，I0 (τ0 ) は第一種修正ベッセル関数であり，拡散パラメータ τ0 > 0 は既知とする．位置パラメータ µ ∈ [0, 2π) に対して，次の von Mises 事前分布を仮定する．
(2.2)π(µ; m, δ) =1
exp{δ cos(µ − m)}.
2πI0 (δ)

ここで，m および δ は超パラメータである．Mardia and El-Atoum（1976）は次の命題を証明した．
（i）事前分布（2.2）は標本分布（2.1）に対して共役で命題 2.1.（Mardia and El-Atoum（1976））
あり，事後密度は π(µ; µ̂MAP , δ∗ ) と書ける．ただし，µ̂MAP は事後モードであり，δ∗ = {τ02 + δ 2 +2τ0 δ cos(x − m)}1/2 である．
（ii）事後モード µ̂MAP は損失関数 L(µ̂ − µ) = 1 − cos(µ̂ − µ) の下で最適である．
事後モードは具体的には次式で与えられる：
(2.3)



µ̂MAP = µ̂MAP (x, m, δ)






τ0 cos x + δ cos m
−1


cos


=



 2π − cos−1δ∗

τ0 cos x + δ cos mδ∗τ0 sin x + δ sin m ≥ 0 のとき,


その他．

後続の節で記述を簡単化するために，事後モードが x，m および δ に依存していることを明示的に表すときには，µ̂MAP (x, m, δ) と記す．
円周上のデータを扱っているので，損失関数 L(µ̂ − µ) = 1 − cos(µ̂ − µ) は自然な選択と思われる．この損失関数はいわゆるエントロピー型損失関数であることを示す．まず，2 つの標本密度 p(x; µ1 , τ0 ) および p(x; µ2 , τ0 ) の間の Kullback-Leibler 分離度を計算する．三角関数の加法定理によりcos(x − µ2 ) = cos(x − µ1 + µ1 − µ2 )= cos(x − µ1 ) cos(µ1 − µ2 ) − sin(x − µ1 ) sin(µ1 − µ2 )であることに注意する．これを用いると，
1
KL p(x; µ1 , τ0 ), p(x; µ2 , τ0 ) =
τ0

2π
0

{cos(x − µ1 ) − cos(x − µ2 )} p(x; µ1 , τ0 ) dx
2π

= {1 − cos(µ1 − µ2 )}0cos(x − µ1 ) p(x; µ1 , τ0 ) dx

2π

+ sin(µ1 − µ2 )である．等式0sin(x − µ1 ) p(x; µ1 , τ0 ) dx

2πp(x; µ, τ0 ) dx = 1 の両辺を µ および τ0 で偏微分すると，それぞれ
0
2π

(2.4)
0
2π

(2.5)
0

sin(x − µ)p(x; µ, τ0 ) dx = 0,cos(x − µ)p(x; µ, τ0 ) dx =I1 (τ0 )I0 (τ0 )von Mises 分布における経験 Bayes 推定

179を得る．ここで，I0 (τ0 ) = I1 (τ0 ) を用いた．したがって，


(2.6)



KL p(x; µ1 , τ0 ), p(x; µ2 , τ0 ) = τ0I1 (τ0 )L(µ1 − µ2 )I0 (τ0 )となる．
事後モードが他の推定量を優越する様子を明快に示す変形ピタゴラス関係を導く．
命題 2.2. 任意の推定量 µ̂ に対して，次の変形ピタゴラス関係が成り立つ．
2π

(2.7)
0

L(µ̂ − µ) − L(µ̂MAP − µ) π(µ; µ̂MAP , δ∗ ) dµ =I1 (δ∗ )L(µ̂ − µ̂MAP ).
I0 (δ∗ )

証明．三角関数の加法定理を用いると次の式変形ができる：
L(µ̂ − µ) − L(µ̂MAP − µ) = cos(µ̂MAP − µ) − cos(µ̂ − µ)
= L(µ̂ − µ̂MAP ) cos(µ̂MAP − µ) + sin(µ̂ − µ̂MAP ) sin(µ̂MAP − µ).

ここで，式（2.4）および式（2.5）において x，µ および τ0 をそれぞれ µ，µ̂MAP および δ∗ に置換したものを用いると，変形ピタゴラス関係（2.7）が得られる．
変形ピタゴラス関係（2.7）の右辺の意味を議論する．2 つの事前密度 π(µ; m1 , δ) および π(µ; m2 ,δ) の Kulback-Leibler 分離度は（2.6）と同様に
KL π(µ; m1 , δ), π(µ; m2 , δ) = δ

I1 (δ)
L(m1 − m2 )I0 (δ)となる．したがって，変形ピタゴラス関係（2.7）は次のように書き換えられる：
2π0L(µ̂, µ) − L(µ̂MAP , µ) π(µ; µ̂MAP , δ∗ ) dµ =1KL π(µ; µ̂MAP , δ∗ ), π(µ; µ̂, δ∗ ) .
δ∗

事後リスクの差が 2 つの密度の Kulback-Leibler 分離度を用いて表現できるというこの関係は，
von Mises 分布に限られたものではない．ある一定の曲指数型分布族の共役解析において成立することがわかっている（大西・柳本, 2005）．
3.

経験 Bayes 推定

確率変数 xi (1 ≤ i ≤ p) が互いに独立に von Mises 分布 p(xi ; µi , τ0 ) に従っているとし，パラメータ µi (1 ≤ i ≤ p) が互いに独立に von Mises 事前分布 π(µi ; mi , δ) に従っていると仮定する．
超パラメータ m = (m1 , . . . , mp )T は既知とし，δ を未知とする．周辺分布によるモーメント法を用いて超パラメータ δ を推定し，経験 Bayes 推定量を構成する．ここで，推定すべきパラメー
タ µ = (µ1 , . . . , µp )T の次元 p が大きい場合を想定している，すなわち，高次元パラメータの推定問題を考える．次節のシミュレーションでは，高次元の場合として p = 100 のケースを扱う．
また，比較のために低次元の例として p = 10 の場合を計算する．
von Mises 分布における経験 Bayes 推定において，Bagchi（1994）はモーメント法を用いて mおよび δ を推定している．本研究で m を既知としたのは，次の 2 つの理由による：（1）実際の推定では，事前データなどがありパラメータ µ のおおよその値が知れていることが多い．
（2）正規分布における James-Stein 推定量は，任意に与えた位置パラメータについて拡散パラメータを推定して構成される経験 Bayes 推定量として理解できるが，von Mises 分布においても Stein
現象と似たような現象が予想される．

180

統計数理 第 54 巻

第 1 号 2006

超パラメータ δ に対する推定方程式を導くとき，次の変形ピタゴラス関係が重要な役割を果たす．証明は命題 2.2 と本質的に同じなので省略する．
補題 3.1. 次の変形ピタゴラス関係が成り立つ：
2π

L(xi − mi ) − L(xi − µi ) −

(3.1)
0



I1 (τ0 )L(µi − mi ) p(xi ; µi , τ0 ) dxi = 0 (1 ≤ i ≤ p).
I0 (τ0 )超パラメータ δ をモーメント法で推定するときの基礎となる等式を命題の形で与える．
命題 3.1. 次の等式が成立する．


Emarg

(3.2)



p

1 
I1 (τ ) I1 (δ)
cos(xi − mi ) =.
p i=1
I0 (τ ) I0 (δ)

ただし，Emarg [ · ] は周辺分布による期待値を表す．
証明． 周辺分布を pmarg (x; τ0 , m, δ) =であり，

pmarg (xi ; τ0 , mi , δ) と書く．ただし，x = (x1 , . . . , xp )T2π

pmarg (xi ; τ0 , mi , δ) =p(xi ; µi , τ0 )π(µi ; mi , δ) dµi0である．積分の順序を交換することにより，
2π
0

L(xi − mi ) pmarg (xi ; τ0 , mi , δ) dxi
2π

=
0
2π

2π

=
0

0



2π

L(xi − mi )0p(xi ; µi , τ0 )π(µi ; mi , δ) dµidxiL(xi − mi ) p(xi; µi , τ0 ) dxi π(µi ; mi , δ) dµiを得る．補題 3.1 および等式
2π0

L(xi − µi ) p(xi; µi , τ0 ) dxi = 1 −

I1 (τ0 )I0 (τ0 )を用いると，
2π
0

L(xi − mi ) pmarg (xi ; τ0 , mi , δ) dxi
2π 

2π

=
0

0

2π 



L(xi − µi ) +

I1 (τ0 )L(µi − mi ) p(xi ; µi , τ0 ) dxi π(µi ; mi , δ) dµi
I0 (τ0 )

I1 (τ0 )I1 (τ0 )+L(µi − mi ) π(µi ; mi , δ) dµi
1−I0 (τ0 )I0 (τ0 )=02πとなる．ここで， 0 L(µi − mi ) π(µi ; mi , δ) dµi = 1 − I1 (δ)/I0 (δ) であるから，
2π
0

となる．
Emarg



p

L(xi − mi ) pmarg (xi ; τ0 , mi , δ) dxi = 1 −


p

1 
1 
cos(xi − mi ) =Emarg [ cos(xi − mi ) ]
p i=1
p i=1

I1 (τ0 ) I1 (δ)
I0 (τ0 ) I0 (δ)

von Mises 分布における経験 Bayes 推定

181

p

=

1  2π{1 − L(xi − mi )} pmarg (xi ; τ0 , mi , δ) dxi
p i=1 0に注意すると，命題の等式が得られる．
超パラメータ δ を推定するために次の推定方程式を考える．
p

I1 (δ) I0 (τ0 ) 1 =cos(xi − mi ).
I0 (δ) I1 (τ0 ) p i=1(3.3)

関数 I1 (δ)/I0 (δ) は δ ∈ Ê について狭義単調増加であり，値域は区間 (−1, 1) である．狭義単調増加性は次のように証明できる．
d I1 (δ) I1 (δ)I0 (δ) − I1 (δ)I0 (δ)=
dδ I0 (δ)
I02 (δ)であるから，右辺の分子が正であることを示せばよい．積分表現
In (δ) =

1
π

π
0

n は整数

cos(nθ) exp(δ cos θ) dθ ,（Abramowitz and Stegun, 1970, p. 376）を用いると，右辺の分子は


1
π2

π
0

2

cos θ eδ cos θ

π

dθ

e

0

δ cos θ



dθ −

π
0

cos θ eδ cos θ2 dθとなる．ここで微分と積分の順序交換を行っている．Schwarz の不等式からこれは正である．
積分表現において φ = π − θ と置換し積分を実行すると，関数 I1 (δ)/I0 (δ) は奇関数であることがわかる．奇関数なので I1 (0)/I0 (0) = 0 である．また，修正ベッセル関数の漸近展開（Abramowitz
and Stegun, 1970, p. 377）を用いると，δ が十分大きいときの近似式
1I1 (δ)≈1−
I0 (δ)
2δ

が得られる．
狭義単調性および値域が (−1, 1) であることから，



p


 I0 (τ0 ) 1 cos(xi − mi ) < 1

 I1 (τ0 ) p
i=1のときは，δ̂ は一意に定まる．パラメータ µi に対する経験 Bayes 推定量は，式（2.3）を用いて，µ̂EB, i = µ̂MAP (xi , mi , δ̂) (1 ≤ i ≤ p) となる．
推定方程式（3.3）が解を持たないときは次の手続きを行う．
p

I0 (τ0 ) 1 cos(xi − mi ) ≥ 1
I1 (τ0 ) p i=1のとき，µ̂EB, i = mi (1 ≤ i ≤ p) とする．この不等式が成立するのは，データ x と超パラメータm は十分近いときである．パラメータ µ の推定量として m を採用するのは自然である．一方，p

I0 (τ0 ) 1 cos(xi − mi ) ≤ −1
I1 (τ0 ) p i=1のとき，µ̂EB, i = mi + π とする．このケースでは，データ x と超パラメータ m は各成分が180 度近く離れていると考えられるので，m + π1 を推定量とすることは妥当である．ただし，

182

統計数理 第 54 巻

第 1 号 2006

1 = (1, . . . , 1) ∈ Êp である．ここで，lim µ̂MAP (xi , mi , δ̂) = mi ,δ̂→∞lim µ̂MAP (xi , mi , δ̂) = mi + πδ̂→−∞が成り立つ．したがって，推定方程式（3.3）が解を持たないときも含めて経験 Bayes 推定量は
形式的に µ̂MAP (xi , mi , δ̂) (1 ≤ i ≤ p) と書ける．
ここで，一つ注意すべきことは，δ̂ として負の値を許していることである．δ cos(µi − mi ) =−δ cos(µi − mi − π) であり，π(µi ; mi , δ) および π(µi ; mi + π, −δ) は分布としては同一であるので，これは自然な取扱いといえる．
4.

リスクの比較

経験 Bayes 推定量 µ̂EB のリスクを数値シミュレーションし，最尤推定量 µ̂ML = x のリスクと比較する．最尤推定量のリスクは次のように計算できる．


p
p






E
L(xi − µi )
p(xi ; µi , τ0 ) =E L(xi − µi )  p(xi ; µi , τ0 )

i=1
i=1
i=1

 p


=p 1 −

I1 (τ0 )I0 (τ0 ).

ここで，E [ f (x)| p(x)] は密度関数 p(x) の下での f (x) の期待値を表す．パラメータ µ の値によらず一定であることに注意する．当然のことながら，m にも依存しない．
シミュレーションは以下の 4 つの場合に分けて行う．反復回数はすべて 50,000 である．固定したパラメータ µ に対して超パラメータ m を様々に変化させ，経験 Bayes 推定量のリスクを計算しグラフを描く．以下のグラフにおいて，経験 Bayes 推定量のリスクは ‘+’ でプロットし，
最尤推定量のリスクは実線で表示する．ケース 3 以外では τ0 = 1 とし，ケース 3 では拡散パラ
メータ τ0 を 2, 5, 10 と変化させて 2 つの推定量のリスクを比較している．
4.1 ケース 1
経験 Bayes 推定量のパフォーマンスの全体像をつかむために，まず最初に次のような最も簡単な設定の場合を扱う：µ = µ1 かつ m = m1.

図 1（a）は，p = 10, µ = µ1 に対して m = m1 の m を 0 から 2π まで変化させたときの経験 Bayes
推定量のリスクをグラフにしたものである．実際の計算では µ = 0 としているが，location familyがもつ対称性により結果は µ − m にのみ依存することが分かる．さらに，損失関数が偶関数であるからリスクは µ − m の符号によらず絶対値のみによることが分かる．この 2 点について略証を与える．確率変数 ei = xi − µ (1 ≤ i ≤ p) は互いに独立に von Mises 分布 vM(0, τ0 ) に従う．ここで，記号 vM(µ, τ0 ) は密度関数（2.1）をもつ von Mises 分布を表す．加法定理を用いると，経験 Bayes 推定量に対する損失の値は次のように表現できる．


−1
1 − cos µ̂EB, i − µ = 1 − δ̂∗i
τ0 cos ei + δ̂ cos(µ − m) .

ただし，



I1 (δ̂)1 I0 (τ0 )=p
I1 (τ0 )I0 (δ̂)cos(µ − m)

p

i=1

cos ei − sin(µ − m)

p

i=1



sin ei ,von Mises 分布における経験 Bayes 推定

183



1/2

δ̂∗i = τ02 + δ̂ 2 + 2τ0 δ̂{cos(µ − m) cos ei − sin(µ − m) sin ei }.

これらの等式から，適当な関数 fi を用いて



1 − cos µ̂EB, i − µ = fi (e1 , . . . , ep ; µ − m)








と表現できることが分かる．これはリスク pi=1 E 1 − cos µ̂EB, i − µ が µ − m にのみ依存することを意味している．
リスクが µ − m の絶対値のみに依存することは，µ − m < 0 のとき
cos(µ − m) cos ei = cos |µ − m| cos(−ei ),sin(µ − m) sin ei = sin |µ − m| sin(−ei ).

と書き換えられることおよび確率変数 −ei が von Mises 分布 vM(0, τ0 ) に従うことから証明できる．
グラフから読み取れることは，次の 2 点である：（1）|µ − m| ≈ 0 あるいは |µ − m| ≈ π のとき，圧倒的に経験 Bayes 推定量がよい．
（2）最尤推定量に対してリスクが 20%程度大きくなる|µ − m| ≈ π/2 のときを除けば，経験 Bayes 推定量のパフォーマンスは最尤推定量に勝っている．第 1 点については次のように説明できる：超パラメータ m はパラメータ µ の予想値という解釈ができる．予想と真の値が近い，すなわち |µ − m| ≈ 0 のとき，経験 Bayes 推定量はとてもよい振る舞いをする訳である．もっと詳しく言えば，δ の推定値が高い確率で大きな値になり，経験 Bayes 推定値も高い確率で m1 に近くなるのである．一方，予想と真の値が 180 度近く外れているとき，δ の推定値は高い確率で絶対値が大きな負の値になり，推定値が (m + π)1 に近くなるのである．経験 Bayes 推定量のパフォーマンスが必ずしも良くないのは |µ − m| ≈ π/2のときである．
図 1（b）は，図 1（a）において p = 100 の場合である．2 つのグラフを比較すると，高次元になれば経験 Bayes 推定量のパフォーマンスがよくなることが分かる．p = 10 のときでも最尤推定量を優越していた |µ − m| ≈ 0 および |µ − m| ≈ π の場合では，その優越度がさらに高まっている．また，p = 10 のときに最尤推定量に劣っていた |µ − m| ≈ π/2 の場合でも差がほとんど

図 1.

ケース 1 における経験 Bayes 推定量および最尤推定量のリスク．

184

統計数理 第 54 巻

第 1 号 2006

なくなっている．圧倒的に優越するか，劣るとしてもごく僅かということになって，提案した経験 Bayes 推定量のよい振る舞いが示されている．
高次元になったときのよい振る舞いを定性的に説明する．超パラメータ δ の推定方程式は，加法定理を用いると，I1 (δ) I0 (τ0 )=I0 (δ) I1 (τ0 )

p

p

1 
1 
cos(µ − m)cos(xi − µ) − sin(µ − m)sin(xi − µ)
p i=1
p i=1



と書ける．確率変数 ei = xi − µ (1 ≤ i ≤ p) は独立に同一に von Mises 分布 vM(0, τ0 ) に従う．大数の強法則を適用すると，確率 1 でpp

I1 (τ0 )1 
1 cos(xi − µ) =および limsin(xi − µ) = 0
p→∞ p
p→∞ p
I0 (τ0 )
i=1
i=1
lim

が成立することが分かる．したがって，確率 1 で
lim

I1 (δ̂)p→∞ I (δ̂)
0= cos(µ − m)が成り立つ．この等式は，高次元のとき高い確率で I1 (δ̂)/I0 (δ̂) ≈ cos(µ − m) が成り立つことを意味している．したがって，limδ→∞ {I1 (δ)/I0 (δ)} = 1 に注意すると |µ − m| ≈ 0 のとき高い確率で δ̂ は非常に大きな正の値になる．経験 Bayes 推定量の式 µ̂EB, i = µ̂MAP (xi , mi , δ̂) および（2.3）
から µ̂EB ≈ m1 となる．次に |µ − m| ≈ π の場合を考える．関数 I1 (δ)/I0 (δ) が奇関数であることに注意すると，上と同様の議論により δ̂ は絶対値が非常に大きな負の値になり，µ̂EB ≈ m + π1となることが分かる．最後に |µ − m| ≈ π/2 の場合を考える．等式 I1 (0)/I0 (0) = 0 に注意すると，高い確率で δ̂ ≈ 0 となることが分かる．したがって，µ̂EB ≈ µ̂ML となる．
4.2 ケース 2
実際的な設定で計算を行う．すなわち，
µi =

図 2.

π i
(1 ≤ i ≤ p) かつ m = m1
4 p

ケース 2 における経験 Bayes 推定量および最尤推定量のリスク．

von Mises 分布における経験 Bayes 推定

185の場合を考える．このケースでは，第 3 節で述べたように事前データなどから m ≈ π/8 と知れていることが多い．図 2 は，µi = (π/4) · (i/p) (1 ≤ i ≤ p) に対して m = m1 の m を 0 から 2π まで変化させたときの経験 Bayes 推定量のリスクをグラフにしたものである．大雑把に言えば，µ = (π/8)1 の場合と似たような結果になることが予想される．図 1 と比較するとそのことが確認できる．p = 10 のとき，m = 5π/8 あるいは m = 13π/8 あたり以外では経験 Bayes 推定量のパフォーマンスは最尤推定量に勝っている．高次元になると，これまでと同様にほとんどの mの値に対して経験 Bayes 推定量は最尤推定量を優越する．最尤推定量に劣る m = 5π/8 あるいは m = 13π/8 においてもその差は極めて小さい．
4.3 ケース 3
拡散パラメータ τ0 以外はケース 1 と同じ設定で，τ0 を 2, 5, 10 と変化させる．図 3 はそれぞれの τ0 について経験 Bayes 推定量および最尤推定量のリスクをグラフにしたものである．これまでと同様に，高次元のとき経験 Bayes 推定量のパフォーマンスはよい．一方，10 次元の場
合でも τ0 が大きくなるにつれて経験 Bayes 推定量の最尤推定量に対する相対的なパフォーマ
ンスはよくなっていくことが分かる．この事実の理由として考えられるのは，τ0 が大きくなると von Mises 分布が正規分布に近づくことである．計算は省略するが，第 3 節で導出した経験
Bayes 推定量は τ0 が十分大きく，かつ，µ − m ≈ 0 のとき，James-Stein 推定量と近い形になる．第 3 節で少し述べたように，Stein 現象に似た状況が現れたと言える．
4.4 ケース 4
経験 Bayes 推定量にとってかなり不利な設定で計算してみる．具体的には次のような設定に
する：


µi =

µ (1 ≤ i ≤ p/2)
かつ m = m1.
µ (p/2 + 1 ≤ i ≤ p)

図 4（a）および図 4（b）は，µ = 0，µ = π に対して m を 0 から 2π まで変化させたときの経験 Bayes 推定量のリスクをグラフにしたものである．この設定ではすべての m に対して経験Bayes 推定量は最尤推定量に劣っている．パラメータの配置がこのような場合，経験 Bayes 推
定量のパフォーマンスはよくない．しかし，そのような不利な状況でも次元が高くなるとリスクの差はほとんどないことが分かる．
第 4.1 節と同様に，高次元のときの振る舞いを調べてみる．推定方程式は


p/2

p/2

1 
1 
I1 (δ) I0 (τ0 )=cos(µ − m)cos(xi − µ) − sin(µ − m)sin(xi − µ)I0 (δ) I1 (τ0 )
p j=1
p j=1
p
p
1 
1 
+ cos(µ − m)cos(xi − µ ) − sin(µ − m)sin(xi − µ )
p
p


j=p/2

j=p/2

となる．大数の強法則から，確率 1 でp/2

I1 (τ0 )2 cos(xi − µ) =p→∞ p
I0 (τ0 )
j=1
lim

p/2

2 
sin(xi − µ) = 0，
p→∞ p
j=1

および lim

p
p
I（
2 
2 
1 τ0 ）
および limcos(xi − µ ) =sin(xi − µ ) = 0
p→∞ p
p→∞ p
I（
0 τ0 ）

lim

j=p/2

j=p/2



186

統計数理 第 54 巻

図 3.

第 1 号 2006

ケース 3 における経験 Bayes 推定量および最尤推定量のリスク．

von Mises 分布における経験 Bayes 推定

図 4.

ケース 4 における経験 Bayes 推定量および最尤推定量のリスク．

187

188

統計数理 第 54 巻

第 1 号 2006が成立する．したがって，確率 1 で
(4.1)

lim

I1 (δ̂)p→∞ I (δ̂)
0=1{cos(µ − m) + cos(µ − m)}2が成り立つ．|µ − µ | = π のとき，右辺は m の値によらず 0 になる．つまり，p が十分大きい
とき，経験 Bayes 推定量は確率 1 で最尤推定量に一致する．これは図 4（b）のグラフはほぼ一定であることに合致する．
図 4（c）および図 4（d）は，µ = 0.2π ，µ = 0.8π に対して m = m1 の m を 0 から 2π まで変化させたときの経験 Bayes 推定量のリスクをグラフにしたものである．図 4（a）および図 4（b）に比べると少し経験 Bayes 推定量にとって有利な設定である．10 次元の場合，経験 Bayes 推定量のパフォーマンスは最尤推定量に優越するときもあれば逆のときもある．しかし，100 次元になるとたいていの場合で経験 Bayes 推定量のパフォーマンスが最尤推定量を優越する．劣る場
合でも差はごく僅かである．
式（4.1）を用いて高次元のときの振る舞いを調べてみる．加法定理から，（4.1）の右辺において cos(µ − m) + cos(µ − m) = (cos µ + cos µ ) cos m + (sin µ + sin µ ) sin m と変形できる．今の場
合，µ + µ = π であるから，確率 1 で
lim

(4.2)

I1 (δ̂)p→∞ I (δ̂)
0= sin µ sin mとなる．したがって，高次元であり，かつ，m = 0 あるいは m = π のとき，経験 Bayes 推定量と最尤推定量の振る舞いは似ている．
図 4（e）および図 4（f）は，µ = 0.4π ，µ = 0.6π に対して m = m1 の m を 0 から 2π まで変化させたときの経験 Bayes 推定量のリスクをグラフにしたものである．図 4（c）および図 4（d）に比べてさらに経験 Bayes 推定に有利な設定である．また，近似的に µ = 0.5π1 という設定と理解することもできる．10 次元のときでもすでに経験 Bayes 推定のパフォーマンスはよい．高次元になると経験 Bayes 推定は最尤推定量よりよい推定量だと考えてほぼ差し支えない．
先と同様に，高次元のとき（4.2）を用いて定性的な振る舞いを理解できる．経験 Bayes 推定
量の振る舞いが最もよくなるのは，m = π/2 あるいは m = 3π/2 のときである．このとき，δ̂ の絶対値が大きくなっている．図 4（d）と図 4（f）の場合を比較すると，前者では（4.2）の右辺の値
が sin 0.2π であるのに対し，後者では sin 0.4π となる．これは図 4（d）および図 4（f）のグラフ
における m = π/2 あるいは m = 3π/2 近辺での谷の深さに合致している．
謝辞

本論文をまとめるにあたり，統計数理研究所の柳本武美教授に有益なコメントをいただいた．

参 考 文 献Abramowitz, M. and Stegun, I. A.（1970）. Handbook of Mathematical Functions, Dover, New York.
Bagchi, P.（1994）. Empirical Bayes estimation in directional data, Journal of Applied Statistics, 21,
317–326.
Guttorp, P. and Lockhart, R. A.（1988）. Finding the location of a signal: A Bayesian analysis, Journalof the American Statistical Association, 83, 322–330.
Mardia, K. V. and EL-Atoum, S. A. M.（1976）. Bayesian inference for the von Mises-Fisher distribution, Biometrika, 63, 203–206.

von Mises 分布における経験 Bayes 推定

189

大西俊郎，柳本武美（2005）．Dual structure in the conjugate analysis of curved exponential families,
京都大学数理解析研究所講究録，1439, 202–220.
Rodrigues, J., Leite, J. G. and Milan, L. A.（2000）. An empirical Bayes inference for the von Mises
distribution, Australian & New Zealand Journal of Statistics, 42, 433–440.

190

Proceedings of the Institute of Statistical Mathematics Vol. 54, No. 1, 177–190(2006)

Empirical Bayes Estimation in von Mises Distribution
Xiong Ye1 and Toshio Ohnishi2
1 Department of Statistical Science, The Graduate University for Advanced Studies
2 The Institute of Statistical Mathematics

In the conjugate analysis of the von Mises distribution, it is known that the posteriormode is the optimal estimator, that is, the Bayes estimator under the Kullback-Leibler
loss function. This paper derives an empirical Bayes estimator by estimating the dispersion hyper-parameter of the von Mises prior distribution under the assumption that thelocation hyper-parameter is known. The estimation of the dispersion hyper-parameter isbased on the method of moment for the marginal density. We compare the risk of theempirical Bayes estimator with that of the maximum likelihood estimator (MLE) by anumerical-simulation method. In the practical situation, we find that the empirical Bayes
estimator is far superior to the MLE in high-dimension cases.

Key words: Conjugate prior, posterior mode, Pythagorean relationship.