【統計学】尤度関数（RとPython）

教科書有斐閣「統計学」第9章SECTION４の内容です。

尤度関数。正直苦手です。 ですので間違えている可能性も大なので炎上マーケ大歓迎です。

最尤推定法とは

自分で説明するよりも下のサイトのほうがはるかにわかりやすいです。

mathtrain.jp

下の実装もこのサイトのコインの例をそのままやってます。

表が出る確率が θ であるようなコインがある。このコインを 100 回投げたら 70 回表が出た。最尤法により θ を推定せよ。

Rで実装

#コインの成功確率の関数（＝尤度関数）
dice <- function(p,trial, success){
  (choose(trial, success))*(p^success)*((1-p)^(trial-success))
}

#試行100回、成功70回の時

#表になる確率が最大化するように尤度関数を最大化
optimize(dice, interval = c(0, 1),  trial = 100, success = 70 , maximum = TRUE)
#$maximum
#[1] 0.700014

#$objective
#[1] 0.08678386

Pythonで実装

import math

def combinations_count(n, r):
    return math.factorial(n) // (math.factorial(n - r) * math.factorial(r))

def coin(p, trial, success):
     cul = (combinations_count(trial, success))*(pow(p, success))*(pow((1-p),(trial-success)))
     return cul
 
m = 0
s = 0
for i in [0.001*x for x in range(1000)]:
    p = coin(i, 100, 70)
    if p>m:
        m = p
        s = i
    else:
        pass

print('コイントスの確率: {0:.4f}'.format(m))
#コイントスの確率: 0.0868

print('コインが表の確率: {0:.4f}'.format(s))
#コインが表の確率: 0.7000