各種分布のパラメータの最尤推定（最尤推定大演習）

(1)

1

各種分布のパラメータの最尤推定（最尤推定大演習）

得られたデータx_ ( 1, ,N)が、特定の分布に従うかどうかを調べる際、分布のパラメータが既知であるとは限らない。そのため、多くの場合、与えられたデータを用いて各種分布のパラメータを推定し、その下で検定の問題を考えることになると思われる。そこで、メニュー［分析－基本統計－分布と検定］で表示される図 1 の分析実行メニューに、

パラメータを自動的に推定する機能を加えた。分布を選んで「推定」ボタンをクリックすると左のテキストボックスに推定値が表示される。

図1 分布と検定実行メニュー

ここでは分布毎にパラメータを推定するための方法を具体的に与えておく。

正規分布 (   x )

密度関数： ² ² ²

2

( | , ) 1 exp[ ( ) 2 ]

2

f x   x_  

   

尤度関数： ²

2 2 2

1

1 1

exp ( )

(2 ) 2

N

L N x_



  _ 

 

   







対数尤度： ² ²

2 1

log 1 ( ) log(2 )

2 2

N N

L x_



 

 

 



 

スコアベクトルUと情報行列

2





 

  

 

β ，

2

log log

L L





 

 

   

U ，

2 2 2 2

2 2 2 2 2

log log

log log ( )

L L

  

  

      

   

    

 

2 1

log 1 ( ) 0

N

L x_



 

 

  



 

1

1 ^N N  x^







2 2

4 2

1

log 1 ( ) 0

2 2

N N

L x_



 

  

  



   ² ²

1

1 ( )

N

N _ x^

 







以上で解析的に求めることが可能であるが、プログラムでは練習問題としてニュートン・

ラフソン法を用いて計算を試している。

(2)

2

2 2

log N2

L 

   

2 2

4 1

log 1 ( ) 0

N

L x_



  

 

    



 

2 2 2 2

6 4 4

1

log ( ) 1 ( )

2 2

N N N

L x_



 

 _  

   



   

初期値は₀ 0, ₀² 1を用いている。

注） ² ²

2

log 2N

L 

    ²logL ²4²²logL (^{2 2})

パラメータの推定にはどちらの値を用いるべきだろうか。

χ²分布 (0  x ) パラメータが離散的

密度関数： ^{2 1}

2

( ) 1 exp( 2)

2 ( 2)

n

f x n x x

n

  



尤度関数： ₂ ^{2 1}

1

1 exp( 2)

2 ( 2)

N n

Nn N

L x x

n _ ^ ^





 





対数尤度：

1 1

2 1

log log( ) log 2 log ( 2)

2 2 2

N N

n Nn

L x_ x_ N n

 

 







  

2 2

 n のとき、E(²)n の性質を用いて、

0.5

nx  ^注）  x ^はxを越えない最大の整数

これを元に(1 )  n 5 n_max  n 5の範囲で最大の対数尤度を与える自由度n_max^を求めて

いる。

Ｆ分布 (0  x )

密度関数：

1 1

1 2

2 2 1

1

( ) 2

1 2 2 1 2

( ) 1

( 2, 2) (1 )

n n

n x

f x B n n n xn n





 

   

尤度関数：

1 1

1 2

2 2 1

1

( ) 2

1 2 2 1 1 2

1

( 2, 2) (1 )

Nn N n

n n N

x L n

B n n n x n n



 



 

 

  





対数尤度：

1 1 2

1 2

1 1

1

1 2 1 2

log 1 log( ) log(1 )

2 2

log( ) log ( 2, 2)

2

N N

n n n

L x x n n

Nn n n N B n n

 

 

  

    

 

 

 

2 2

[ ] 2

E X n

n

 (n² 2)^， ²² ¹ ₂ ²

1 2 2

2 ( 2)

[ ] ( 2) ( 4)

n n n

V X n n n

  

  (n² 4) を利用して、

2

2 [ ] [X] 1 n E X

 E

 ^，

2

2 2

1 2 2

2 2 2

2 ( 2)

( 2) ( 4) [ ] 2 n n n

n n V X n

 

  

これを元に、ぶれが大きいので、

(1 ) n_i20n_imax  n_i 20の範囲で対数尤度を最大化す

(3)

3 るn_i_maxを求めている。

ｔ分布 (   x )

密度関数：

( 1) 2 1 2

2 2

( )

( ) 1

( )

n n

n

f x x

n n

   

     

尤度関数：

( 1) 2 1 2

2 2 1

( )

1 ( )

N n

n

L x n n



 

  



    

  



  

対数尤度：

2 1

2

1 2

( )

log 1 log 1 log

2 ( )

N n

n

x

L n N

n n



 





 

  

  



     

平均：E X[ ] 0

分散： [ ]

2 V X n

n

 ^{を利用して、}

2 [ ] [ ] 1 n V X

V X



これを元に

(1 )  n 5 nmax  n 5の範囲で最大の対数尤度を与える自由度

nmax^を求めて

いる。

ガンマ分布 (0  x )

密度関数： 1 ¹

( ) exp( )

( )

a

f x a x x b

b a

  



尤度関数： ¹

1

1 exp( )

[ ( )]

N a

a N

L x x b

b a _ ^ ^





 





対数尤度：

1 1

log ( 1) log 1 log log ( )

N N

L a x x Na b N a

 b 

 

 







  

1

log log log ( ) ( )

N

L a x_ N b N a a



  



   

2 1

log 1

N

L b b x_ N a b



  





2 2 2 2

logL a N ( )a ( )a ( )a ( )a 

         

2logL a b N b

    

2 2 3 2

1

log 2

N

L b b x_ N a b



   





初期値はa₀ 0.5, b₀ 0.5を用いている。

逆ガンマ分布 (0 x 1)

密度関数： ( ) ¹exp( ) ( )

a

b a

f x x b x

a

   



(4)

4

尤度関数： ¹

1

( ) exp( )

N N

a a

L b a x_ b x_



 



 

  





対数尤度：

1 1

log ( 1) log (1 ) log log ( )

N N

L a x_ b x_ Na b N a

 

  







  

1

log log log ( ) ( )

N

L a x_ N b N a a



   



   

1

log (1 )

N

L b x_ N a b



   





2 2 2 2

logL a N ( )a ( )a ( )a ( )a 

         

2logL a b N b

   

2 2 2

logL b N a b

   

初期値はa₀ 0.5, b₀ 0.5を用いている。

ベータ分布 (0 x 1)

密度関数：

1 1

(1 ) ( )

( ) (1 )

( , ) ( ) ( )

a b

x x a b

f x x x

B a b a b

 

  

  

 

尤度関数： ¹ ¹

1

1 (1 )

( , )

N

a b

L x x

B a b  ^ ^

 







対数尤度：

1 1

log ( 1) log ( 1) log(1 ) log ( , )

N N

L a x_ b x_ N B a b

 

 



 



 

1

( , )

log log

( , )

N

B a b a

L a x N

B a b

 

  





1

( , )

log log(1 )

( , )

N

B a b b

L b x N

B a b

 

  



 

2

2 2 ( , ) ( , )

log ( , )

aa a

B a b B a b

L a N

B a b B

 

      

 

2

( , ) ( , ) ( , )

log ( , ) ( , )

ab a b

B a b B a b B a b L a b N

B a b B a b

 

       

 

2

2 2 ( , ) ( , )

log ( , )

bb b

B a b B a b

L b N

B a b B

 

      

 

初期値の設定で、平均値が0に近い場合は1,5、1に近い場合は5,1、0.5に近い場合は0.5, 0.5などを使う。小さい方から大きい方へ近づけて行くことは問題ないが、大きい方から小さい方へ近づけて行く際にはエラーが出る。

ワイブル分布 (0  x )^{（失敗例）}

通常のa b, を使って最尤法を試みた。

(5)

5 密度関数： ^{f x}^{( )}^⁽^{a b x b}⁾

 

^a^¹^exp^_^

 

^{x b} ^a^_

尤度関数：

 

¹

 

¹

1 1

( ) exp exp[ ]

N N

a a

N N Na a a a

L a b x b_ x b_ a b x_ x b_

 

   

 

 





 





対数尤度：

1 1

log ( 1) log log log

N N

a a

L a x_ b x_ N a Na b

 



 

 







 

1 1 1

log log log log log

N N N

a a a a

L a x_ b b x_ b x_ x_ N a N b

  

 

  

  



 







  

1 1

log

N

a a

L b ab x_ Na b



 



  





2 2 2

1 1

2 2

1

log (log ) 2 log log

(log )

N N

a a a a

N

a a

L a b b x b b x x

b x x N a

  

 



 

 





      

  

 



2 1 1

1 1

log (1 log ) log

N N

a a a a

L a b a b b x_ ab x_ x_ N b

 

   

 

    







 

2 2 2 2

1

log ( 1)

N

a a

L b a a b x_ Na b



 



    





この方法は、収束が思うように行かず、エラーとなった。

ワイブル分布 (0  x ) 再度

上記の失敗を踏まえ、生存時間分析で用いたパラメータの推定法を利用する。

密度関数： ^{f x}^{( )}^⁽^{a b x b}⁾

 

^a^¹^exp^_^

 

^{x b} ^a^_

尤度関数：

 

¹

 

¹

1 1

( ) exp exp[ ]

exp[ ]

N N

a a

N N Na a a a

N

N N a a

L a b x b x b a b x x b

a e x x e

   

 

 

 



   

 





 

   

 

 



 

1

1 1

log ( , ) log ( 1) log

( 1) log log

N

a

N N

a

L a b a a x x e

a x e x N a N

  



  

 





 

 

      

    



 

1 1

log log log

N N

L x e x xa N a

a



  

 

    



 

1

log

N

L e^ x_a N

 

   





2

2 2

2

1

log (log )

N

L e x xa N a

a

  



   





(6)

6

1

log log

N

L e x xa

a

  

 

  

 



2 2

1

log

N

L e^ x_a

 _

  





初期値はa₀ 2,  2を用いている。

指数分布 (0  x )

密度関数： f t( )aexp(ax) ^（x0）尤度関数：

1 1

exp( ) exp

N N

i

L a ax_ a a x_





 

    









対数尤度：

1

log log

N

i

L N a a x_



 



1

log 0

N N

L x

a a  ^

   





1 N

a N x_







2

2 log N2

a L a

  



この逆数は、推定値の分散を与える。（今回は使わない）

推定値は解析的に求まるが、練習問題として最尤法を用いてみる。

初期値はa0.1を用いている。

ポアソン分布 (0  x )^，整数

確率関数：P x( )e a x^^a ^x ! 尤度関数：

1 1

! !

N N

x x

a Na

L e a ^ x_ e a ^ x_

 

 

 









対数尤度：

1 1

log log log !

N N

L Na a x_ x_

 

  







1

log 1

N

L a N x

a ^

    



2 2

2 1

log 1

N

L a x

a  ^

   



初期値はa0.1を用いている。

2 項分布 (0  x )^，整数

まず以下の関係を使って、度数n^{を求める。}

[ ]

E X np^，V X[ ]npq [ ]2

[ ] [ ] n E X

E X V X

 

(7)

7 次に最尤法を使って、確率pを求める。

確率関数：P x( ) _nC p_x ^x(1p)^{n x}^ 尤度関数：

1

(1 )

N

x n x

n x

L C p ^ p ^

 









対数尤度：

1 1 1

log log log log(1 ) ( )

N N N

n x

L C p x p n x

  

  









 





1 1

log ( )

1

N N

L p x n x

p  ^ p  ^

    







2 2

1 1

log ( )

(1 )

N N

L p x n x

p __ ^ p __ ^

     







これも解析的に解を求めることができるが、最尤法の演習とする。

初期値はp0.5を用いている。

各種分布のパラメータの最尤推定（最尤推定大演習）