Werden wir Helden für einen Tag

Home | About | Archive

無聊統計學答案 #1

Posted on Jul 7, 2009 by Chung-hong Chan

這條問題幾近成了標準問題。由於市面上的統計學人都沒有貝氏分析的 training ,這條問題希望可起提綱挈領之效,讓大家去窺探貝氏分析之奧堂。

我們假定一個人能夠估中 25% 以上為有超能力。

傳統的頻率派分析的方法,是這樣的:

其實有關的,只是 6/20 這個實驗結果。我們要先設立兩個 Hypothesis

H0: p < = 0.25
HA: p > 0.25

故此只要推翻 H0 就有證據證明陳有超能力。
用最簡單的 Binomial test ,可知道沒有足夠證據推翻 H0 。用 R 的計法是這樣:


binom.test(6,20,(1/4),alternative="greater")$p.value # p-value = 0.3828

這個分析的結論會是:數據不能提供足夠證據證明陳振聰估中的比率是 25% 或以上,沒有足夠證據證明他的準繩度比亂估更高。故此,沒有足夠證明陳振聰有超能力。
但是,這問題是問「陳振聰沒有超能力的機會是幾多」,而不是有沒有證據證明陳振聰有超能力。這個是傳統的頻率派分析之限制。其中一個原因,是頻率派對機會率的定義,是經過反覆實驗後,再從結果獲得的頻率得出機會率。故此,機會率一定是實驗後才出現。之前本博解釋 Confidence Interval 時,就有講過 95% Confidence Interval 的 95% 不是指出 95% 機會包括母數。
要解決這個問題的另一個方法,是改變我們對機會率的看法。機會率,可以解讀為「知識的現狀」。故此,就算沒有進行實驗,我們都可以有一個機會率,稱為事前機會率( Prior probability );實驗過後,我們根據實驗結果修正我們的看法,又一個機會率,稱為事後機會率( Posterior probability )。這種思考方法,稱為貝(葉斯)氏思考法( Bayesian Thinking )。
貝氏和頻率派其實沒有優劣之分,只是頻率派較為常用,亦較為常教。但在解決「陳振聰沒有超能力的機會是幾多」這類問題,使用貝氏的結果會較好。但當然,正如一般的統計分析,都是 Garbage in Garbage out 。假如我們對「知識的現狀」了解太差,貝氏分析就不準確。
設 p 為陳氏估中的機會率。我們先要計算出事前機會率,我們可用實驗前的信心計算。其實總共有 40 個 X 。我們可以用 X / 40 得出事前機會率,我們用 g(p) 代表事前機會率的密度( Density )。以 R 計算事前機會率的方法為:


p < - seq(0.05,0.95,by=0.1)
prior <- c(2,4,16,8,4,2,1,1,1,1)
prior.p <- prior/sum(prior)
cbind(p,prior.p)

結果如下:


p prior.p
[1,] 0.05 0.050
[2,] 0.15 0.100
[3,] 0.25 0.400
[4,] 0.35 0.200
[5,] 0.45 0.100
[6,] 0.55 0.050
[7,] 0.65 0.025
[8,] 0.75 0.025
[9,] 0.85 0.025
[10,] 0.95 0.025

從 g(p) 我們看到,陳氏估中 25% 的事前機會率是四成。而陳氏估中 5%, 15% 或 25% 的機會率是 0.05+0.1+0.4 = 0.55 ,即五成五。

根據貝氏定律,事後機會率的密度 g(p|data) ((這個叫做 Density of p given the data 。是 Conditional Probability )) ,是等於事前機會率密度乘以 p 的似然值( Likelihood , L(p) )。即是

g(p|data) ∝ g(p) * L(p)

在今次 Discrete distribution ,似然值是這樣計算。

L(p) = ps(1-p)f

這個似然值是根據實驗數據計算出來的。只有兩個數要塞進去,就是 s 及 f 。就是 success 及 fail 的次數,即 6 及 14 。我們可以使用 R 計算出事後機會率:


library("LearnBayes") #pls install LearnBayes package
trial < - c(6,14)
post.p <- pdisc(p, prior.p, trial)

我們可列出 p 值, p 值的事前及事後機會率


cbind(p,prior.p,round(post.p,3)) #post.p have too many digits. round to 3

如下


p prior.p
[1,] 0.05 0.050 0.000
[2,] 0.15 0.100 0.040
[3,] 0.25 0.400 0.589
[4,] 0.35 0.200 0.299
[5,] 0.45 0.100 0.065
[6,] 0.55 0.050 0.007
[7,] 0.65 0.025 0.000
[8,] 0.75 0.025 0.000
[9,] 0.85 0.025 0.000
[10,] 0.95 0.025 0.000

從結果可見,事後機會率有所增降。例如陳氏估中 0.25 的機會率由四成修正至近六成。而陳氏估中 0.05, 0.15 或 0.25 的事後機會率,是 0 + 0.04 + 0.589 = 0.629 。根據貝氏分析法,陳氏估中 25% 或以下(亂估或更差)的可能性是 62.9% 。即是,陳有 62.9% 可能性沒有超能力。就只有從貝氏分析獲得的結果能夠這樣直接的解讀。當然,我們亦要「o岩 channel 」,這個 62.9% ,要用貝氏的方法去理解,只是反映「知識的現狀」。 ((亦因此,如果再有新的數據,我們可以再修改這個「現狀」。)) 而不是「 100 個陳振聰裡面有 62.9 個沒有超能力」那樣解釋。

參考書目: Introduction to Bayesian Thinking. In: Bayesian Computation with R. 2nd Eds.


Powered by Jekyll and profdr theme, a fork of true minimal theme