Re: [問卦] 武漢肺炎也會數學？

看板Gossiping作者pulseaudio (pulseaudio)時間4年前 (2020/02/07 02:27)推噓17(17推 0噓 16→)

留言33則, 17人參與討論串6/12 (看更多)

1. 你算錯了，你把你的每一個vlist[n]印出來，加起來不等於1。程式裡有bug。我試着把每一項列出來： p(k) b(k) b(k) chi2(k) (0.210526-0.301030)**2/0.301030 = 0.027210 (0.315789-0.176091)**2/0.176091 = 0.110827 (0.052632-0.124939)**2/0.124939 = 0.041847 (0.105263-0.096910)**2/0.096910 = 0.000720 (0.105263-0.079181)**2/0.079181 = 0.008591 (0.105263-0.066947)**2/0.066947 = 0.021930 (0.052632-0.057992)**2/0.057992 = 0.000495 (0.000000-0.051153)**2/0.051153 = 0.051153 (0.052632-0.045757)**2/0.045757 = 0.001033 chi2(k)加起來乘於19 等於 5.012 p-value = 0.756 在合理範圍。 2. 那個公式是給N很大的時候用的。 chi2 = N * (p-b)^2 / b，只有N很大的時候分母才可以用 b。這個公式其實就是（省略sum符號）： chi2 = (N*p-N*b)^2 / 誤差^2。當N很大時，誤差=sqrt(N*b)。可是在本例裡面，9出現的機率的期望值是0.045757 樣本19，因此9出現的次數的期望值是0.87次，不到一次。這時候只能用Poisson誤差，不能用常態分佈的誤差。 Poisson誤差是不對稱的，不適合用chi-squared，要用likelihood。不過不必算，我們可以，估算，由於Poisson誤差比常態分佈大，真正的chi2值會比上面的 5.012要小，真正的p-value會比0.756要大。 3. 就跟hugoyu一樣，你使用的一系列數字不是獨立變量，而是下一個數字是前一個數字的累加。我懷疑不是獨立變量可以用這個定律。你應該用hugoyu提共的數字的第二列，那才是獨立變量。用那一組數字得到的 chi2是6.68，p-value是0.57，依然在合理的範圍內。 ※ 引述《Glamsight (安穩殘憶)》之銘言： : 生活中總是會遇到許許多多數據 : 也不知道究竟是真是假 : 所以現在檢定的時間到了 : 對於這種增量正比於存量的數據通常我們會使用一種神祕的定律 Benford's law : 簡單來說，該定律描述數據其第一個數據之分布情形 : 具體來說 (10 進位下)，認為數據第一位數字之機率有如下表 : 第一位數字 | 出現機率 : ------------------------ : 1 | 0.301029996 : 2 | 0.176091259 : 3 | 0.124938737 : 4 | 0.096910013 : 5 | 0.079181246 : 6 | 0.06694679 : 7 | 0.057991947 : 8 | 0.051152522 : 9 | 0.045757491 : 如，各國的地區人口數、各國 GDP、國土面積、放射性元素之半衰期等 : 其數據之第一位數字都成這樣的機率分布，十分之神奇 : 李永樂老師有視頻介紹，並附有該情況下的證明 : https://www.youtube.com/watch?v=CCo4k9Ax7cM

: 本文就考量下一次感染的人數亦可與現有已感染人數呈正比關係 : 及根據 hugoyo 提供之數據嘗試與 Benford's law 進行比較與檢定 : 首先，我們畫張圖 : https://s1.imgs.cc/img/aaaabbvbw.png?_w=750 : 看起來有點像，又有點不像 : 所以我們考慮具體一點的比較，即統計檢定 : 我們善意的假設數據沒有造假，故如下設計 : 虛無假設 H0: 數據沒有造假 : 即，該數據應與 Benford's law 分布相近 : 對立假設 H1: 數據存在造假 : 即，該數據與 Benford's law 分布不相似 : 檢驗方式邏輯大概如下 : 在善意考量沒有造假下，官方數據應該與 Benford's law 相似 : 那麼我們詢問在 Benford's law 的情況下，出現像是這次官方數據的機率大概有多高 : 具體的算法是使用 Chi-Square Test : 方法來自 Google 文獻 : https://evoq-eval.siam.org/Portals/0/Publications/SIURO/Vol1_Issue1/Testing : _for_the_Benford_Property.pdf : (有斷行，請注意) : 將本次事件之數據第一個數字 k 之機率定為 p(k) : 該 k 之數字於 Benford's law 下發生機率定為 b(k) : 則根據文獻可以知道一 Chi-Square Test 檢定量如下 : 檢定量 = 資料筆數 * sum[(p(k)-b(k))^2/b(k)] (k=1~9) : Chi-Square Test 自由度為 9-1=8 : 可以知道檢定量值約為 38.35 : 該檢定 p-value = 6.47 * 10^-6 約等於 0 : 也就是說不論設定多嚴苛的顯著水準 : 都會取得拒絕虛無假設的結論 (reject H0) : 也就是說，我們可以理性推論該數據存在疑慮。 : 以上如果有錯麻煩跟我說。 : 謝謝 : 下附 Python 程式碼 30 行 : def getHeader(value,base,onoff=True): : v=int(value/base) : if v!=0 and onoff: : v=getHeader(v,base,onoff=True) : else: : v=value : return v : from math import log : def benfordDistribution(n,base): : pr=log(n+1,base)-log(n,base) : return pr : from scipy.stats import chi2 : def getChi(vlist,base): : vsum=0 : for n in range(base-1): : pr=benfordDistribution(n+1,base) : vsum=vsum+((vlist[n]-pr)**2)/pr : v=len(vlist)*vsum : pvalue=1-chi2.cdf(v,base-2) : return {'Chi-Square 檢定量':v,'p-value':pvalue} : base=10 : vlist=[45,62,201,218,320,543,639,1356,2033,2744,4515,5974,7711,9692,11794, : 14384,17213,20448,24335] : v1list=[getHeader(v,base)/len(vlist) for v in vlist] : v=getChi(v1list,base) : for key in v: : print(key,v[key]) : ※ 引述《hugoyo (阿佑)》之銘言： : : 大概大概啦，他們公布的東西就是這樣 : : 天數日期確診累積確診死亡累積死亡死亡率 : : 1 1月18日 45 45 2 2 4.44% : : 2 1月19日 17 62 0 2 3.23% : : 3 1月20日 139 201 1 3 1.49% : : 4 1月21日 17 218 0 3 1.38% : : 5 1月22日 102 320 3 6 1.88% : : 6 1月23日 223 543 11 17 3.13% : : 7 1月24日 96 639 0 17 2.66% : : 8 1月25日 717 1356 24 41 3.02% : : 9 1月26日 677 2033 15 56 2.75% : : 10 1月27日 711 2744 24 80 2.92% : : 11 1月28日 1771 4515 26 106 2.35% : : 12 1月29日 1459 5974 26 132 2.21% : : 13 1月30日 1737 7711 38 170 2.20% : : 14 1月31日 1981 9692 43 213 2.20% : : 15 2月 1日 2102 11794 46 259 2.20% : : 16 2月 2日 2590 14384 45 304 2.11% : : 17 2月 3日 2829 17213 57 361 2.10% : : 18 2月 4日 3235 20448 64 425 2.08% : : 19 2月 5日 3887 24335 65 490 2.01% : : 猜猜看明天是不是 2.00% 左右?? : : 從一月底開始都穩穩地控制了，死亡率漸漸變低，這樣大家才會比較安心 : : 那麼，明天會多少確診呢? 既然大家覺得數據都是Garbage，再怎麼分析都是Garbage out : : 還是可以猜猜看明天的Garbage長什麼樣子呀~ : : 不想太認真用分佈的模型來積分成error function : : 用大家都能懂得多項式就好，抓個3次方。 : : https://i.imgur.com/pqua02g.png

還蠻平的~ : : https://i.imgur.com/ni4FgcQ.png

: : 所以明天中國大概確診 28297 人，增加 4046 人 : : 死亡 566 人，增加 76 人 : : 個人希望他們不要真的這樣玩數據。 -- ※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 203.77.33.103 (臺灣) ※ 文章網址: https://www.ptt.cc/bbs/Gossiping/M.1581013664.A.018.html

→

miku3920

02/07 02:28, 4年前 , 1^F

02/07 02:28, 1^F

推

Ponimp

02/07 02:28, 4年前 , 2^F

02/07 02:28, 2^F

推

tw15

02/07 02:30, 4年前 , 3^F

02/07 02:30, 3^F