X
トップページ⇒数学
1002コメント563KB
数学 統計に詳しい人が語るコロナウイルス
■ このスレッドは過去ログ倉庫に格納されています
0001１３２人目の素数さん
垢版 |
2020/02/29(土) 02:18:41.53ID:twdO677Q
東大数学科卒の元官僚はこう分析してるが、お前らはどうなると思う？

http://www.zakzak.co.jp/soc/news/200220/dom2002200003-n2.html

中国国外感染者の中国国内との比率をみると、
１月２０日の数字公表以降は、０・８〜２・６％で比較的安定している。
これは、新型肺炎の感染者のほとんどは中国国内、それも湖北省に集中しているからだ。
ちなみに中国国外での感染者数は、中国国内の１・１％だ（２月１６日現在）。
本コラムで紹介したが、現時点では、最終的な中国国内の感染者数は２０万人超と筆者は推計している。
となると、中国国外の感染者は数千人程度になるだろう。
中国国外のうち日本の比率は１割弱なので、日本の感染者数は数百人程度であろう。
その場合、死者も数人から１０人程度になるだろう。

こうした推計をすると、今の感染者は氷山の一角だと思われるが、今後の増加ペースはどうなるだろうか。
新型コロナウイルスの検査は簡単に行えるので、今後、日本での感染者数は増えていくだろう。
ある時点ではそれがネズミ算的に増えるかのように思える局面もあるだろうが、
筆者の推計が正しければ、現時点ではせいぜい数百人が一つのメドだ。
0503１３２人目の素数さん
垢版 |
2020/05/14(木) 14:10:15.04ID:kec+XbRE
upper_boundの制限を外して

data {
// int<lower = K> upper_bound;
int upper_bound;

再生数の平均値を以下の出すブロックを加えて走らせてみた。
transformed parameters{
real mean_Rt;
real mean_Rt_adj;
mean_Rt = mean(Rt);
mean_Rt_adj = mean(Rt_adj);
}

その結果、
Inference for Stan model: fit_infection.
2 chains, each with iter=10000; warmup=2000; thin=5;
post-warmup draws per chain=1600, total post-warmup draws=3200.

mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff Rhat
mean_Rt 1.79 0 0.07 1.66 1.74 1.79 1.83 1.92 3122 1
mean_Rt_adj 1.79 0 0.07 1.66 1.74 1.79 1.84 1.93 3123 1
0504１３２人目の素数さん
垢版 |
2020/05/14(木) 16:35:35.75ID:kec+XbRE
再生算数を０〜１０人の一様分布にすると、収束しない。

> print(fit_u)
Inference for Stan model: fit_infection_u.
4 chains, each with iter=10000; warmup=5000; thin=5;
post-warmup draws per chain=1000, total post-warmup draws=4000.

mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff Rhat
mean_Rt 2.15 0.05 0.10 1.98 2.06 2.16 2.23 2.30 3 1.99
mean_Rt_adj 2.13 0.07 0.11 1.93 2.06 2.13 2.21 2.36 3 2.27
lp__ -789.19 7.18 14.43 -820.19 -797.69 -787.44 -779.42 -762.53 4 1.89

traceplotやchainの分布はこんな感じ、
https://i.imgur.com/z0RL1KW.png
https://i.imgur.com/VRrrsKw.png
0505１３２人目の素数さん
垢版 |
2020/05/14(木) 17:00:26.74ID:1Kd4DHQt
>>501
否定できないけど実証もできないことになるのでは？
単純なモデルだと現実の問題を解決できない事もあるし
0506１３２人目の素数さん
垢版 |
2020/05/14(木) 17:23:23.34ID:kec+XbRE
>>494
モデルで再生産数の事前分布は　平均2.5　標準偏差2.0の正規分布に設定されていたので
平均と標準偏差を変化させて、再生産数の事後分布を描出してみた。
かなり、事前分布の影響を受けるみたい。
https://i.imgur.com/OwqsFC1.png
0507１３２人目の素数さん
垢版 |
2020/05/14(木) 17:26:03.64ID:kec+XbRE
どうも、こういう境地だなぁ。

断定的な結論は出せないということを数字で確認しているだけｗ
0508１３２人目の素数さん
垢版 |
2020/05/14(木) 17:30:34.30ID:50GljeDM
基本再生算数はウイルスの特徴で決まるもので
実行再生算数は更に環境や人の行動の影響で変化するものだと理解している
実行再生算数を減らすようにするには
どう行動したらいいか
かつ経済活動を出来るだけ下げずに
0509１３２人目の素数さん
垢版 |
2020/05/14(木) 18:25:32.93ID:IN0uokww
西浦モデル、一度感染日ごとの人数を最尤法で確定して、そこからベイズを回してRtを求めてるんだよなあ。
本来は、日々の感染者数も確率的に揺れがあるはず。だからRtの誤差幅は発表よりも多く見積もるべきかもしれん。
0513１３２人目の素数さん
垢版 |
2020/05/15(金) 09:56:20.25ID:qjCTzgxb
【新型コロナ】 東京0.6％、東北6県0.4％陽性・・・抗体検査1000人実施 [影のたけし軍団★]
https://asahi.5ch.net/test/read.cgi/newsplus/1589502801/

複数の検査キットの性能評価と感染状況の確認が目的でしたが、東京都で献血した500人のうち3人、
東北6県で献血した500人のうち2人がいずれかの検査キットで陽性と判定されました。

満員電車など人との接触の多い東京とそうでない東北で陽性率に有意差はあるか？
0514１３２人目の素数さん
垢版 |
2020/05/15(金) 12:23:12.21ID:GeDvuMme
>>513
検査キットの特異度が99.5%以上あることが検証されました。
ってだけの話じゃね？

まぁ、東北６県と東京が同じってことは、特異度はそれより
高くはないってことだろうね。
0515１３２人目の素数さん
垢版 |
2020/05/15(金) 12:26:48.61ID:GeDvuMme
言い換えると、真の陽性率はわかりません、ってこと。
0517１３２人目の素数さん
垢版 |
2020/05/15(金) 18:43:40.23ID:qjCTzgxb
これは、味気がないな。

> Epi::twoby2(x)
2 by 2 table analysis:
------------------------------------------------------
Outcome : Col 1
Comparing : Row 1 vs. Row 2

Col 1 Col 2 P(Col 1) 95% conf. interval
Row 1 3 497 0.006 0.0019 0.0184
Row 2 2 498 0.004 0.0010 0.0158

95% conf. interval
Relative Risk: 1.5000 0.2517 8.9384
Sample Odds Ratio: 1.5030 0.2501 9.0339
Conditional MLE Odds Ratio: 1.5024 0.1713 18.0536
Probability difference: 0.0020 -0.0092 0.0139

Exact P-value: 1.0000
Asymptotic P-value: 0.6561
------------------------------------------------------
0522１３２人目の素数さん
垢版 |
2020/05/16(土) 12:36:46.54ID:cGFgpNwX
4/2の時点で感染者数を６０００くらいと見積もってるね。＞アヒル本の人
実数の三倍程度。いい線かもしれない。
0524１３２人目の素数さん
垢版 |
2020/05/16(土) 14:29:29.30ID:BxGcLzV+
int<lower = K> upper_bound;

↓

int　upper_bound;

にしてもエラーがでる。

stan_dataで
upper_bound = 147 にすると動くけど、何をやってんのか自分でもよくわからん。
0526１３２人目の素数さん
垢版 |
2020/05/16(土) 17:55:00.91ID:BxGcLzV+
>>525
名著として名高いStanの入門書

StanとRでベイズ統計モデリング (Wonderful R) (日本語) 単行本 ? 2016/10/25

表紙の色がアヒルのくちばしの色に似ているかららしい。
0529１３２人目の素数さん
垢版 |
2020/05/16(土) 19:28:52.51ID:BxGcLzV+
R に　EpiEstimというパッケージがあって、再生産数を算出する関数が搭載されている。
結局、infecterとinfecteeが発症するまでの期間serial intervalの分布をどう設定するかで結果が変わるみたいだなぁ。

Rのヘルプファイルを解読中。
Rのヘルプファイルは不親切設計で有名（理解できている人の備忘録みたいな性格だから）。
0530１３２人目の素数さん
垢版 |
2020/05/16(土) 20:49:50.96ID:5T7nYzW3
ここに居る人達って何者なの
学部の知識超えてるよね
統計でご飯食べてる人たち？
0531１３２人目の素数さん
垢版 |
2020/05/16(土) 21:03:16.59ID:YPW1s8+S
しかし、なんでも揃ってるなRって
0532１３２人目の素数さん
垢版 |
2020/05/16(土) 22:57:43.46ID:BxGcLzV+
>>529
Stanでの西浦モデルではinfecterとinfecteeが発症するまでの期間serial intervalの分布に

## Serial interval [Nishiura et al 2020 - only certain cases]
param1_SI = 2.305,
param2_SI = 5.452,

// serial interval
vector[K] gt = pweibull(param1_SI, param2_SI, K);

として使われているので、平均値などを出してみた。
乱数発生と理論値

> x=rweibull(1e5,param1_SI,param2_SI)
> mean(x) ; param2_SI*gamma(1+1/param1_SI)
[1] 4.829273
[1] 4.830129
> var(x) ; param2_SI^2*(gamma(1+2/param1_SI)-(gamma(1+1/param1_SI))^2)
[1] 4.907755
[1] 4.940682
> median(x) ; param2_SI*(log(2)^(1/param1_SI))
[1] 4.655777
[1] 4.6505
> density(x)$x[which.max(density(x)$y)] ; param2_SI*(1-1/param1_SI)^(1/param1_SI)
[1] 4.116837
[1] 4.259624
> optimise(function(x) dweibull(x,param1_SI,param2_SI),c(0,10),maximum = T)$max
[1] 4.259623

グラフにしてみた。
https://i.imgur.com/9vvCJuZ.png
正規分布で近似してもよさそうな感じだな。
0533１３２人目の素数さん
垢版 |
2020/05/17(日) 00:04:44.72ID:u5AEq3c8
>>532
ワイブル分布の平均 4.830129と標準偏差2.222765をそのまま正規分布のパラメータに使って、グラフを重ねてみる。
https://i.imgur.com/TnzGwWx.png

ワイブル分布で発生させた乱数をワイブルでフィットさせてAICを出してみた
Goodness-of-fit criteria
1-mle-weibull
Akaike's Information Criterion 438377.2
Bayesian Information Criterion 438396.2

ワイブル分布で発生させた乱数を正規分布でフィットさせてAICを出してみた。
Goodness-of-fit criteria
1-mle-norm
Akaike's Information Criterion 444280.9
Bayesian Information Criterion 444299.9

まぁ、許容範囲。

これで、
library(EpiEstim)の例にある、　mean_si std_siが求まった

## Estimate R with assumptions on serial interval
res <- estimate_R(incid, method = "parametric_si",
config = make_config(list(
mean_si = 4.83, std_si = 2.22)))

domestic ,　imported, unobserved の分類がよくわからんが、全部足してグラフを描いてみた。

https://i.imgur.com/rKBeWgq.png
0534１３２人目の素数さん
垢版 |
2020/05/17(日) 00:18:59.61ID:u5AEq3c8
別の論文だと対数正規分布がフィットすると西浦氏は記載している。

serila interval : infector と infecteeの発症間隔

https://www.ijidonline.com/article/S1201-9712(20)30119-3/pdf
その分布が平均4.7　標準偏差2.9の対数正規分布が最もフィットするのはいいんだが、
その分布を与えるパラメータの記述がほしい。
最小二乗法で求めてみた。　
$par
[1] 1.3862713 0.5679836

ワイブル分布にも似るとか書いてあるがパラメータ記載なし
この対数正規分布をワイブル分布で近似してみる。
Fitting of the distribution ' weibull ' by maximum likelihood
Parameters:
estimate Std. Error
shape 1.757488 0.00392072
scale 5.316986 0.01014664

https://i.imgur.com/Uzg6u84.png
点線が2項分布で実線がワイブル分布
0535１３２人目の素数さん
垢版 |
2020/05/17(日) 07:09:05.66ID:/ApxcPC4
“頭脳王”東大生・河野玄斗

基本的な数学でコロナウイルス検査を全員にしても意味がないことを証明してみた

https://www.youtube.com/watch?v=jMIScCb04qs
0536１３２人目の素数さん
垢版 |
2020/05/17(日) 07:28:25.69ID:mpH4uuPx
実効再生産数を感染者数の推移から推定する数理的原理をキチンと解説した本、または論文誰か知りませんか？
勉強してみたい。
0538１３２人目の素数さん
垢版 |
2020/05/17(日) 07:55:41.34ID:u5AEq3c8
Rで引数なしで関数を実行させようとすると、ソースが表示される。
その関数が呼び出した関数のソースを次々に辿っていけば内部計算がわかるので、そこから原理がつかめるかも（俺には無理なのでしないけど）

内部関数のときは:が3つとか、パッケージ名:::関数、パッケージ名:::関数.default（例、ｔ検定のソース表示は　stats:::t.test.default　）で表示される。

EpiEstim::estimate_R
EpiEstim:::process_si_data
EpiEstim:::process_config_si_from_data
coarseDataTools::dic.fit.mcmc
0539１３２人目の素数さん
垢版 |
2020/05/17(日) 08:10:41.41ID:u5AEq3c8
>>529
感染者に０が続くと再生産数の信頼区間幅がどんどん広くなってくる。
まあ、疫病用のソフトウェアと理解しておこう。

https://i.imgur.com/QbwNydN.png

infected=c(0,1,1,1,0,0,0,2,0,3,2,3,1,1,1,1,3,0,1,0,0,0,0,0,1,0,0,0,
0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0)
## Estimate R with assumptions on serial interval
res <- estimate_R(infected, method = "parametric_si",
config = make_config(list(
mean_si = 4.83, std_si = 2.22)))
0541１３２人目の素数さん
垢版 |
2020/05/17(日) 09:34:29.22ID:u5AEq3c8
基本コンセプトはこれだろうな。
Therefore, in practice, we apply our method to data consisting of daily counts of onset of symptoms where the infectivity profile ws is approximated by the distribution of the serial interval.

公表された感染者数で計算するために発症から診断/公表までの時間も考慮した点が西浦モデルの優れた点ではないかと感じている。
0542１３２人目の素数さん
垢版 |
2020/05/17(日) 12:07:27.73ID:S4iLXC97
>>537
ありがとう。
読んでみます。
でも"a new frame work"っていうならもっと大本のこの道の研究者なら「Rtの推定はコレ」みたいな標準的なものがあるんですかね？
どなたかご存知ないですか？
0545１３２人目の素数さん
垢版 |
2020/05/17(日) 15:31:13.14ID:u5AEq3c8
結局、これが核心部分

Rt =I t (number of new infections generated at time t) / Σ[s=1,t] I t-s * Ws ( = total infectiousness of infected individuals at time t)

Ws : an infectivity profile given by a probability distribution ws, dependent on time since infection of the case, s, but independent of calendar time, t.


E[I t] = Rt Σ[s=1,t] I t-s * Ws

Σ[s=1,t] I t-s * Wsの部分が畳み込み積分で

Ws ∝　serial interval　　ガンマ分布で近似するのが定石らしい。
0546１３２人目の素数さん
垢版 |
2020/05/17(日) 20:41:31.22ID:/XIaXEJI
>>545
それはそのモデルでのRtの定義ですね、
ではなく例えば4/1〜4/30までのデータx1,x2,‥,x30までが与えられた時、各日付のCIをどう定義してるのかがわからないんです。
統計量Xが一つあってその観測値xが一つある時そのレベルλのレベルλのCIがIであるとは

P(X<x|θ)>(1-λ)/2 & P(X>x|θ)>(1-λ)/2 ⇔ θ∈I

(θがIに入らないときはxが小さすぎてそんな観測値が得られる確率が(1-λ)/2以下になるか、xが大きすぎてそんな事が起こる確率が(1-λ)/2以下になってほとんど起こり得ない)

となりますが、統計量が複数になるとこの大きすぎて、小さすぎてと二つのハズレ領域を考えるだけでは済まなくなります。
ハズレ領域の設定の仕方は色々考えられるけど統計の問題なので自分で俺様流の領域を設定するわけにもいかないのでpublicにはどう処理してるのだろうと。
しかし疫学の教科書わざわざ買うほどには興味もないしw
知ってる人いたら教えてもらおうと。
まだ論文読んでないので書いてあるかもですけど。
0547１３２人目の素数さん
垢版 |
2020/05/18(月) 05:56:32.65ID:hVD/4gbI
>>546
>523氏が挙げてくれた
実効再生産数とその周辺　に　記述があったが、publicと呼べるのかどうかは門外漢なので知らない。

https://github.com/contactmodel/COVID19-Japan-Reff/blob/master/nishiura_Rt%E4%BC%9A%E8%AD%B0_12May2020.pdf

こんな記述があった
>>
Several non mathematical definitions.
A:
あるカレンダー時刻 t で起こっている 2 次感染数の 1 人あたり平均値
B:
あるカレンダー時刻 t で感染した 1 人がその後の経過で生み出す 2 次感染者数の平均値
C:
罹患率 有病割合比などから推定される 1 人あたりの 2 次感染者数（ actual reproduction number とか）
D:
予防接種など流行対策下での 1 人当たりが生み出す 2次感染者数
<<
と定義は一義的ではないみたいだが、

西浦モデルでのR1(t) R2(t)は　>540の論文では　Rt　Rct　（ｃはcohortの頭文字)として言及されているから、まあ、理論疫学者の間ではcommonな定義なんだろうと思う。
0548１３２人目の素数さん
垢版 |
2020/05/18(月) 06:16:59.93ID:hVD/4gbI
>>546
>4/1〜4/30までのデータx1,x2,‥,x30
x1,...,x30が感染者数なら非負整数で与えられるから、CIを考える必要があるかな？
集計ミスで19人が117人とかあったらしいから、信頼区間を考える必要があるのかもしれないとは思うけど。

確かに、発症日に関してはファジーな部分があるとは思う。
いつから熱がありましたか？味がわからなくなりましたか？と問われても１〜３日位の幅はでるだろう。

RのEpiEstim::estimate_Rをヘルプファイルの実例を使って走らせてみた。
serial intervalの分布をデータから推定させるのに
　infecterの発症日の下限・上限
　infecteeの発症日の下限・上限
を設定する項目が(EL,ER,SL,SR)があった。
このデータに合致する分布関数(ガンマ、ワイブル、対数正規分布が指定可能）のパラメータを算出して計算させているみたい。

西浦モデルではワイブル分布で固定。
潜伏期間にも変動があるから、誰がinfecterで誰がinfecteeかを決定するのも難しいだろうとは思う。
後から発症した人間が感染源というのもありうるし。
0549１３２人目の素数さん
垢版 |
2020/05/18(月) 11:15:46.73ID:iMUOaQs/
>>548
統計データから真の罹患日を推定する方法もあるようですが
そこではないんです。
しりたいのはCIのハズレ領域の設定です。
1変数の場合、母数θに対して分布Fθが定まっている場合、レベルλに対して[0,1]の部分集合J(λ)が決まって、観測値xに対する信頼区間I(λ,x)は

θ∈I(λ,x)⇔Fθ(x)∈J(λ)

を満たす区間として定まります。
上の方でI(λ)が上下対称に取らないのはなぜという話題がありましたが、コレがその理由です。
J(λ)の方は上下の(1-λ)/2を削って((1-λ)/2,1-(1-λ)/2)をとり、上下対称に“ハズレ領域”をとりますが、それをもとに計算されるI(λ,x)は対称とはならないからです。
問題は観測値が2変数以上ある場合“ハズレ領域”をどう設定するものかわからないのです。
私が大学で勉強した時はそこまでやらなかったので。
普通に考えればI^3の中の立方体の体積がλになるように真ん中にとるんだろうなぁと思うんですけど。
0550１３２人目の素数さん
垢版 |
2020/05/18(月) 13:21:08.42ID:1P7V5xJn
職場でも最初に発症した人が感染源のように扱われるけど
潜伏期間の分布を考えたら断定はできない。

COVID19の潜伏期間の論文
https://www.nejm.org/doi/full/10.1056/NEJMoa2001316
から、

潜伏期間は対数正規分布で近似できてそのパラメータは
#--- incubation period ---
# from Li et al NEJM 2020
# lognormal mean = 5.2
ln.par1 = 1.434065
ln.par2 = 0.6612
という。


ある人物Xが新型コロナ肺炎に罹患したとする。
行動調査によって発症前にキャバクラに行っており接客したキャバ嬢がX発症の２日後に発症していたことがわかった。
Xがキャバ嬢から移された確率を求めよ。
0556１３２人目の素数さん
垢版 |
2020/05/18(月) 14:26:48.33ID:jxOnYm2h
>>554
ただ、それだとそもそもモードに近いとこをやってます。
信頼区間は密度関数を横に切るのではなく両裾を縦に切ってハズレ部分が1-λになるようにするので少しイメージが違うしきがします。
モードなのかメジアンなのかの違いです。
いずれにせよ、こうやればいいという拡張のための俺様ルールを設定するのはいくらでもできますが、統計の話なのでそんな俺様ルールについて語っても意味ありません。
0557１３２人目の素数さん
垢版 |
2020/05/18(月) 15:34:14.71ID:1P7V5xJn
>>556
単峰性の場合、信頼区間幅が最小になるのがHighest Density Interval

>550なら
HDI
> HDInterval::hdi(x)[1:2]
lower upper
0.5822687 12.5635525

分位数
> quantile(x,c(0.025,0.975))
2.5% 97.5%
1.148711 15.334698

HDIの方が幅が小さい。
0560１３２人目の素数さん
垢版 |
2020/05/18(月) 19:32:28.41ID:4BbElJo8
>>179
単にお勉強ができただけ。
頭が良くないのさ。
自分の頭で物事を考えるってことができない。
0561１３２人目の素数さん
垢版 |
2020/05/18(月) 19:34:11.31ID:4BbElJo8
>>181
その通り。
具体的には理学部の数学科と物理学科。
工学部にも時々もの凄いのがいる。
0562１３２人目の素数さん
垢版 |
2020/05/18(月) 19:49:54.23ID:1P7V5xJn
>>559
いや、Rのパッケージ　HDIntervalで計算しているから誤解していないと思う。
内部の処理コードをみると信頼区間幅が最小になるのを最尤法で出しているね。
0563１３２人目の素数さん
垢版 |
2020/05/18(月) 19:51:25.60ID:1P7V5xJn
>>561
ほんとうに頭のいい人は医学部でなく理学部か工学部に行く。
ほんとうに頭の悪い奴は底辺シリツ医大に行く。
シリツ医大には手先の器用な人もいるが、頭が器用な奴をみたことがない。
0564１３２人目の素数さん
垢版 |
2020/05/18(月) 20:05:42.34ID:1P7V5xJn
>>550(自答)
#
# 人物ｄが発症してdelay日後に濃厚接触したキャバ嬢ｃが発症
# ｃの感染がｄより先行していた確率は？
rm(list=ls())
stancode=
"
data{
real onset_delay;
real ln_par1;
real ln_par2;
}
parameters{
real <lower=0> d_incubation;
real <lower=0> c_incubation;
}
transformed parameters{
real infection_delay = onset_delay + d_incubation - c_incubation;
}
model{
d_incubation ~ lognormal(ln_par1,ln_par2);
c_incubation ~ lognormal(ln_par1,ln_par2);
}

"
model=stan_model(model_code = stancode)
#--- incubation period ---
# from Li et al NEJM 2020
# lognormal mean = 5.2
ln_par1 = 1.434065
ln_par2 = 0.6612
fn.stan <- function(delay, print=FALSE, ...){
dataList=list(onset_delay=delay,ln_par1=ln_par1,ln_par2=ln_par2)
fit=sampling(model,data=dataList, ...)
ms=rstan::extract(fit)
if(print) BEST::plotPost(ms$infection_delay,compVal=0,xlab='infection delay')
mean(ms$infection_delay < 0)
}
fn.stan(2,print=T,iter=5000,warmup=1000)
onset_delays=0:20
y=sapply(onset_delays,fn.stan)
plot(onset_delays,y, ylab='Pr[ Infected Later ])')

２日後の発症だと
> fn.stan(2,print=T,iter=5000,warmup=1000)
[1] 0.2945
３割くらいはあとから発症した方が先に感染していた可能性がある。
0565１３２人目の素数さん
垢版 |
2020/05/18(月) 20:48:18.33ID:1P7V5xJn
Temporal dynamics in viral shedding and transmissibility of COVID-19
https://www.nature.com/articles/s41591-020-0869-5
のRのコード
https://github.com/ehylau/COVID-19/blob/master/Fig1c_Rscript.R
と
西浦モデルのコード
https://nbviewer.jupyter.org/github/contactmodel/COVID19-Japan-Reff/blob/master/scripts/C.%20Calculating%20the%20Rt%20in%20Stan.ipynb

から発症間時間(serial interval)の分布を重ねてみた。

https://i.imgur.com/vrnra5F.png
0566１３２人目の素数さん
垢版 |
2020/05/18(月) 21:01:41.77ID:4BbElJo8
>>563
一方的で申し訳ないが私立大医学部は金持ちのバカ息子が行くイメージ。
0567１３２人目の素数さん
垢版 |
2020/05/18(月) 21:02:26.04ID:4BbElJo8
西浦さんはさんざん適当なことを言って世論を煽ってどう責任を取るのかな？
0569１３２人目の素数さん
垢版 |
2020/05/18(月) 21:20:19.41ID:AArRB0Ix
>>566
それは正しい認識。

凄いのはド底辺シリツ医の馬鹿さ加減だよ。
裏口バカと呼ばれるのがよくわかる。

http://imagizer.imageshack.com/img923/2715/RosCsf.jpg
http://i.imgur.com/XBFnEcU.jpg

馬鹿だという自覚がないので救いようがない。

ICU Bookの最終章の冒頭で著者がこう書いている。

In clinical matters, ignorance can be dangerous,
but ignorance of ignorance can be fatal.

「叱られないと勉強しない」の対偶を「勉強すると叱られる」
と答えるのはignorance can be dangerousの範疇だが、

ドヤ顔で
>対偶をとれば意味が逆になる例文。
というのは、まさに
ignorance of ignorance can be fatal.
0570１３２人目の素数さん
垢版 |
2020/05/21(木) 11:39:47.57ID:hAZkHjNF
西浦モデルの再生算数の事前分布を変化させてグラフにしてみた。

西浦モデルでのデフォルト
https://i.imgur.com/G1wVYgI.png

事前分布を一様分布(0,10)近似の正規分布で近似させた場合
https://i.imgur.com/doS5LEu.png

再生算数の平均０、標準偏差１の場合
https://i.imgur.com/doS5LEu.png


印象としては、西浦モデルは頑強性robustnessのあるモデルとは言えない。
事前分布に大きく影響されるモデルだと思う。
0571１３２人目の素数さん
垢版 |
2020/05/21(木) 11:42:59.41ID:hAZkHjNF
（url訂正)

西浦モデルの再生算数の事前分布を変化させてグラフにしてみた。

西浦モデルでのデフォルト
https://i.imgur.com/G1wVYgI.png

事前分布を一様分布(0,10)近似の正規分布で近似させた場合
https://i.imgur.com/doS5LEu.png

再生算数の平均０、標準偏差１の場合
https://i.imgur.com/0J1RpDa.png


印象としては、西浦モデルは頑強性robustnessのあるモデルとは言えない。
事前分布に大きく影響されるモデルだと思う。
0572１３２人目の素数さん
垢版 |
2020/05/21(木) 11:59:00.68ID:hAZkHjNF
>>571
誤解されるのは不本意なの追記するけど、ソースを公開する西浦先生の姿勢は高く評価している。
隠蔽・改竄・破棄する安倍とは大違い。
0573１３２人目の素数さん
垢版 |
2020/05/21(木) 12:04:41.13ID:RwIzsMl5
だからベイズ推定が統計学の世界でメジャーにならんのだろ？
論理的根拠のない“事前分布”なるもので答えがひょいひょい変わるのでは社会的な影響が大きい防疫政策の決定には使えない。
普通の統計学の検定なら理論的に根拠のある数字、推論しか使わない。
計算は大変だけど。
0574１３２人目の素数さん
垢版 |
2020/05/21(木) 12:36:51.92ID:hAZkHjNF
>>573
成人の平均身長を１〜２ｍに事前分布にするのは納得できるし、
生まれる子供が男子である確率は０．４〜０．６というのも俺は納得できる。
PCRの感度が３０〜７０％として計算するのも納得できるからその設定で階層ベイズモデルを組むことには異論はないな。
0575１３２人目の素数さん
垢版 |
2020/05/21(木) 12:43:16.10ID:hAZkHjNF
こういう問題

あるタクシー会社のタクシーには１から通し番号がふられている。
タクシー会社の規模から保有タクシー台数は１００台以下とわかっている（弱情報事前分布）。
この会社のタクシーを５台みかけた。最大の番号が６０であった。
この会社の保有するタクシー台数の期待値と９５％信用区間(信頼区間)を求めよ。

をベイズで解くときは、
６０台〜１００台である確率を一様分布として処理しているから
これに異論があるのは理解できるけど

日本人成人の平均身長を推定に１〜２ｍを事前分布に想定するのには俺は異論はないね。
一様分布ではなくてガンマ分布にすべきだというのは議論になるとは思うけど。
0576１３２人目の素数さん
垢版 |
2020/05/21(木) 13:00:03.35ID:hAZkHjNF
>>573
ベイズ信奉者から、ベイズ論者を採用したGAFA(Google等）が成長したと教わった。
迷惑メールのフィルタリングとか雑音の除去とか日常的に役立っているというぞ。
0577１３２人目の素数さん
垢版 |
2020/05/21(木) 13:04:49.04ID:hAZkHjNF
>普通の統計学の検定なら理論的に根拠のある数字、推論しか使わない。
>計算は大変だけど。

普通の統計学こそ、無理やり既知の分布に当てはめようとするんだよね。
MCMCの方の方が応用が広い。
0578１３２人目の素数さん
垢版 |
2020/05/21(木) 13:12:59.01ID:hAZkHjNF
プロ野球選手の打率は？と問われたら選手次第で異なる、と誰でもわかるのに
確定不能の平均値が存在していると妄想して計算を始めるのが古典主義統計。
つまり、値は存在するけど確定できないという信仰の世界。

昨今の新コロナでいえば、PCRの検査の感度・特異度が一定と考えるのが古典主義。
プロ野球選手の打率と同じでそんなのは場面場面で変化するよ、と考えて計算するのがベイズ。
0579１３２人目の素数さん
垢版 |
2020/05/21(木) 13:16:04.83ID:EdpRqUGW
>>576
もちろん推定の有力な方法であるにせよ、元の仮定に何の根拠もないわけだからそれから得られる結論には論理的根拠はない、ないが、数学的に伝統的な手法で与えるた結論と大差ない事がなんらかの保証があるなら、有力になる。
それが論理的に“大差ない結論が得られる”事が示されてるなら単なる計算手法に過ぎないし、示されていなくても経験的に“よい結果ぎ得られる事が多い”ならそのジャンルではそこそこ信頼するに値するんだろう。
しかしなんらかのモデルでは答えが一意に定まらず、事前分布の選び方により大きく答えが違ってしまう場合があっても不思議はないし、そのような場合ではやはり、“ではどう計算するのが正しいかのか”の論証を待たなければ信頼するのは危険になる。
0580１３２人目の素数さん
垢版 |
2020/05/21(木) 13:19:12.39ID:hAZkHjNF
古典的（頻度主義）統計信者って、この計算はどうやるんだ？
俺は乱数発生させて計算できるけど、
そればベイズなのかどうかは知らんが、条件付き確率なのでベイズなんだろうな。
（開業医スレに投稿したけど、回答できるやつは０）

COVID19の潜伏期間の論文

https://www.nejm.org/doi/full/10.1056/NEJMoa2001316

結論は
#--- incubation period ---
# from Li et al NEJM 2020
# lognormal mean = 5.2
ln.par1 = 1.434065
ln.par2 = 0.6612


ある開業医が新型コロナ肺炎に罹患したとする。
行動調査によって発症前にキャバクラに行っており接客したキャバ嬢が開業医発症の２日後に発症していたことがわかった。
キャバ嬢は開業医から移されたと主張して１億円の賠償を求めている。
潜伏期間には幅がありキャバ嬢から移された可能性もあると主張してその確率を計算して賠償金を値切りたい。
いくら値切れるか計算せよ。
0581１３２人目の素数さん
垢版 |
2020/05/21(木) 13:31:57.01ID:3Y6HuNza
>>580
潜伏期間なんて数学的に決定できるハズないやん？
数学がやるのは例えば感染日と実際発症した日が確定できるような症例がある程度以上あって、それが従うと病理学的に信頼できる分布の族があって、その中からデータに最も“沿う”分布を選び出すことしかできない。
それは何件も実際の患者のウイルス量の統計をとったり、ウイルスが体内でどのように増えていくのかの病理学的研究データがあって初めて可能になる。
0582１３２人目の素数さん
垢版 |
2020/05/21(木) 17:20:41.64ID:Du+rgnHD
>>571
よく理解できていないので質問ですけど
事前分布とは具体的に何の分布ですか？
基本再生算数の推定値の分布？
実行再生算数の推定値の分布？

実行再生算数の事前分布は基本再生算数の分布としたらいいのかなと思いますけど
0585１３２人目の素数さん
垢版 |
2020/05/21(木) 18:22:49.17ID:hAZkHjNF
潜伏期は病原体と罹患者のパワーバランスで決まるだろうから
定数でなくてばらつきはあると思うね。
0586１３２人目の素数さん
垢版 |
2020/05/21(木) 18:25:39.76ID:8+K/mYXQ
>>572
>隠蔽・改竄・破棄する安倍とは大違い。
安倍ちゃんがやってるわけじゃないだろ。
指示もされてねーのに、官僚がやってんだよ。

そんなこともわからんようでは、あんたの解析もあてにならんな。
0588１３２人目の素数さん
垢版 |
2020/05/21(木) 18:28:21.45ID:hAZkHjNF
>>586
安倍じゃなければ官僚もまともだったんじゃねえの？
安倍らの集団に組み込まれると東大卒も馬鹿になるようだぞ。
0589１３２人目の素数さん
垢版 |
2020/05/21(木) 18:31:24.22ID:8+K/mYXQ
>>581
その通り。
結局、そういう本質的なデータや理論抜きでは、ベイズ推定やったって
限界があるし、結果の説得力もない。
まあ、適用限界ってものがあるのは何やっても同じだけど。
0590１３２人目の素数さん
垢版 |
2020/05/21(木) 18:32:55.80ID:8+K/mYXQ
>>588
>安倍らの集団に組み込まれると東大卒も馬鹿になるようだぞ。

統計的に確認してみれば？ベイズ推定でもっともらしい結果が出るんじゃね？ｗ
0596１３２人目の素数さん
垢版 |
2020/05/21(木) 19:01:50.90ID:hAZkHjNF
>>581
いやウイルスの振る舞いも素粒子の振る舞いもばらつきがあるんじゃないの？
存在も位置も確率的にしかわからない、という議論になると思う。
0597１３２人目の素数さん
垢版 |
2020/05/21(木) 21:14:17.00ID:hAZkHjNF
>>571
最初の２つでは頑強と言えそう。
平均身長のたとえだと
事前分布に1〜２ｍを選んでも１〜１０ｍを選んでも同様の結果だが、０〜１ｍを選ぶと現実離れした結果が返ってくるという感じだな。
0599１３２人目の素数さん
垢版 |
2020/05/22(金) 01:26:24.95ID:MSJJjK3u
>>591
馬鹿の定義は？単なるお前の主観じゃねーかｗ

ベイズ統計にのめり込むと馬鹿になる事後確率分布でも求めてろ。
■ このスレッドは過去ログ倉庫に格納されています
ニューススポーツなんでも実況