20/03/11 20:21:10.84 hVKkfTiV.net
>>51
ご指摘通り。
それどころか陰性陽性の順番も考慮してないのも数字遊びではある。
55:132人目の素数さん
20/03/11 20:26:08.38 pvvBKa5t.net
有病率が異なる集団がたとえば自覚症状やスクリーニング後だったりするなら、ppvを求めるときに同一感度特異度を使うのは強い仮定かもね
興味がある集団に対しての真の感度と特異度はわからないだろうけれども
参照した感度はどういう集団に対してのものなのかを把握するのは大切かと思います
56:132人目の素数さん
20/03/12 17:06:50.22 tGwY7wTH.net
体調が悪くなった時コロナウイルスに感染している確率を計算してみた(確率的思考入門)
URLリンク(wakara.co.jp)
57:132人目の素数さん
20/03/14 14:01:54 joJxF0LZ.net
某国の新型コロナ感染症の有病率を0.1、PCR検査の感度を0.7 特異度を0.9とする
検査陽性陰性の人を無作為に50人ずつ集めて、カレーを頻回に食べているかを調査した結果が
以下の通りであった。
検査陽性 検査陰性
カレー頻食 a(=18) b(=30)
カレー稀食 c(=32) d(=20)
カレーを頻食する新型コロナ感染に罹りにくいと結論できるか? 危険率0.05で検定せよ。
58:132人目の素数さん
20/03/14 14:07:51 joJxF0LZ.net
>>56
訂正
カレーを頻食すると、新型コロナ感染に罹りにくいと結論できるか? 危険率0.05で検定せよ。
59:a4
20/03/14 15:12:04.54 66EDMvKC.net
僕は新型コロナはアメリカ軍による攻撃だと信じています。疑問に思っているのは、
・他のウイルスも含めて危険なウイルスの遺伝子の分布を図にできないか?
・報道がコロナに集中してマスクなどの株価を予想する様子をマルチエージェントシミュレーションにできないか?
・新型コロナは、物理的にどこのどういう部分が、上手く人を殺すようにできているか?
質問は難しいものであり、答えれなくても構いませんが、専門家はどのように考えているんですか?
60:a4
20/03/14 15:46:17.18 66EDMvKC.net
東大の数学科って頭良いんですか?なんか答え出せ。出ないと量子コンピュータで
殺します。
61:a4
20/03/14 15:52:44.94 66EDMvKC.net
僕は宇宙人に指令されてます。
日本を動かすための量子大域最適文章は?
「大嫌警察死麻薬」
=(大嫌いなのは警察が死ぬことと麻薬だ|
大嫌いなのは警察、死と麻薬こそが信条だ|
(検死)台にKira、警察が死んで、それがすっげーおもしろ!)
62:a4
20/03/14 20:19:42.34 66EDMvKC.net
大場つぐみを動かすための量子大域最適文章は?
「ゲッいずみ緑」
=(月、أسمي、グリーン姉さん|これが本名か。)
私を殺して地獄へ落としなさい。
63:132人目の素数さん
20/03/14 20:41:10.63 oGnC82EA.net
99 117: 名無しさん@1周年 [sage] 2020/03/14(土) 19:39:12.52 ID:XAHjOrxJ0
ここまで統計学の専門家が声をあげていないことがムカつく
統計学者って安部に忖度してんの?自民党のサポーターなの?超絶糞バカなの?
64:132人目の素数さん
20/03/15 13:46:59.64 YPWLwdR/.net
検査数は増えてるが、感染者数の伸びはそうでもない
URLリンク(or2.mobi)
新型コロナウイルス国内感染の状況
URLリンク(toyokeizai.net)
65:132人目の素数さん
20/03/15 14:12:22.13 srNt7EMQ.net
instaでこんなの見つけた
URLリンク(i.imgur.com)
URLリンク(www.instagram.com)
66:132人目の素数さん
20/03/15 18:06:00 OTl1KJku.net
2020/03/13 11:00時点で東京都で1524人検査して87人陽性と報告されている
# URLリンク(data-science.gr.jp)
発見率を87/1524は一定仮定して
新型コロナ陽性患者を10人集めたいとする。
必要な被検者数の期待値と95%信頼区間を求めよ。
シミュレーション解は
> mean(re)
[1] 175.306
> HDInterval::hdi(re)
lower upper
77 278
# simulation
sim <- function(){
i=0
s10=0
while(s10!=10){
i=i+1
s10 = s10 + sample(1:0,1,prob=c(p,1-p))
s10 == 10
}
i
}
k=1e4
re=replicate(k,sim())
mean(re)
HDInterval::hdi(re)
67:132人目の素数さん
20/03/15 18:44:20.15 OTl1KJku.net
>>65
発見率は平均値87/1524、0.01-0.10の区間の二項分布に従うと仮定する。
新型コロナ陽性患者を10人集めたい。
必要な被検者数の期待値と95%信頼区間を求めよ、とすると。
期待値175人 95%CIは142-210人
68:132人目の素数さん
20/03/15 20:00:49.03 1+9jqr1h.net
URLリンク(blogos.com)
■英国政府の分析と方針 名無しさんが作成したまとめ
・ピークは10~14週間後にやってくると分析
・海外への修学旅行の中止や基礎疾患のある人のクルーズ船旅行の自粛を呼びかけ
・イギリスの実際の感染者は他国の検査数と陽性率を見ると5000人から1万人
・もはや感染を止めることは不可能なのでゆっくり感染して集団免疫を獲得すること目指す
・来シーズンにはワクチンが開発されている可能性があるため
感染して免疫を獲得した人とワクチン接種を受けた人の割合をコントロールしながら60%に持っていく戦略
・今封じ込めすぎて感染のピークがNHS(国民医療サービス)の繁忙期である冬にやって来ると大変なので
一番暇な夏に来るようにコントロールする
・最悪シナリオの罹患率はドイツの70%を上回る80%に設定して実行計画を立てている
・一斉休校はしない
休校は効果はあるものの最小限であり、効果を上げるのは13~16週間以上の休校が必要になる
効果より害の方が大きい
・フライト制限についても中国便を95%削減してもエピデミックを1~2日遅らすぐらいの効果しかなく
50%ぐらい減らすのがちょうど良い
・イギリスの感染はイタリアより4週間遅れて進んでいる
・基本戦略
(1)ハッピーバースデーを2回歌いながら石鹸と温水で手を洗う
(2)熱や咳の症状がある人は1週間、自宅で自己隔離→ピークを20~25%削減
(3)家族全員を自宅で隔離→ピークを25%削減(未実施)
(4)新型コロナウイルスに脆弱なお年寄りをケア→死亡率を20~30%削減
・集団免疫を獲得するまで先は長いので今からあまり頑張りすぎないように
69:132人目の素数さん
20/03/16 06:09:25.25 CVVw1pKV.net
国内で患者数が大幅に増えたときに備えた医療提供体制の確保について
今後、国内で新型コロナウイルス感染症患者数が大幅に増えたときに備え、各都道府県、保健所設置市及び特別区(以下「都道府県等」という。)における外来を受診する患者数等について、以下の数式を用いて計算いただき、ピーク時の医療需要の目安としてご活用の上、必要な医療提供体制を確保していただくようお願いいたします。
(1)(ピーク時において1 日あたり新たに新型コロナウイルス感染症を疑って外来を受診する患者数)=(0-14 歳人口)×0.18/100+(15-64 歳人口) ×0.29/100+(65 歳以上人口) ×0.51/100
(2)(ピーク時において1 日あたり新型コロナウイルス感染症で入院治療が必要な患者数)=(0-14 歳人口)×0.05/100+(15-64 歳人口)×0.02/ 100+(65 歳以上人口) ×0.56/100
(3)(ピーク時において1 日あたり新型コロナウイルス感染症で重症者として治療が必要な患者数)=(0-14 歳人口)×0.002/100+(15-64 歳人口) ×0.001/100+(65 歳以上人口) ×0.018/100
注1)ピーク時は、各都道府県等において疫学的関連性が把握できない程度に感染が拡大した時点から概ね3か月後に到来すると推計されている。ただし、公衆衛生上の対策を行うことにより、ピークが下がるとともに後ろ倒しされる。
注2)重症者とは、集中治療や人工呼吸器を要する管理が必要な患者を指す。
注3)当該計算式は、都道府県等の単位以下における医療提供体制を確保するためのものであるとともに、各都道府県等によってピークを迎える時期が異なるため、全国の人口を用いて計算することや単純に各自治体が算出するピークの数値を足し合
70:わせることは、不適切な取扱いとなることに留意いただきたい。なお、当該計算式については、今後新たな知見等により変更される可能性がある。 注4)実際には、ピーク時に至るまでの日々の患者数の増加はばらつきがあり、増加曲線は推計通りの形にならない可能性が高いため、現実の患者の発生動向も踏まえて適切に体制を確保することが必要。 注5)当該計算式については、今後新たな知見等により変更される可能性がある。
71:132人目の素数さん
20/03/16 06:50:44.79 CVVw1pKV.net
この行列を使って何か引き出せるだろうか?
> v=c(0.18,0.29,0.51,
+ 0.05,0.02,0.56,
+ 0.002,0.001,0.018)
> (mat=matrix(v,3,byrow=T))
[,1] [,2] [,3]
[1,] 0.180 0.290 0.510
[2,] 0.050 0.020 0.560
[3,] 0.002 0.001 0.018
> mat*100
[,1] [,2] [,3]
[1,] 18.0 29.0 51.0
[2,] 5.0 2.0 56.0
[3,] 0.2 0.1 1.8
>
72:132人目の素数さん
20/03/16 08:14:56.44 CVVw1pKV.net
単なる連立方程式を解くだけだから、つまらんね。
問題 : ある都市でピーク時に外来1000人、入院600人、重症20人であったとすると、この都市の14歳以下の人口は何人と推測されるか?
v=c(0.18,0.29,0.51,
0.05,0.02,0.56,
0.002,0.001,0.018)
(mat=matrix(v*100,3,3,byrow=T))
# URLリンク(www.toukei.metro.tokyo.lg.jp)
(x=matrix(c(1601348,9035668,3103714)/1e4,ncol=1))
(y=round(mat%*%x))
solve(mat,y)
y=matrix(c(1000,600,20))
solve(mat,y)
73:a4 ◆L1L.Ef50zuAv
20/03/16 09:58:31 yKeJ8TdS.net
「安倍晋三」を動かすための量子大域最適文章は?
「アベシトミナミヘイケ」
=(安倍氏と低学歴な大学である南カリフォルニア大学へ行け。|
安倍氏、富、並、へ行け。
あべし、と、南アフリカ共和国へ行け。)
「新型肺炎の強制入院」を動かすための量子大域最適文章は?
「哈佛死吧」
=(ハーバード大学は死ねばいいんじゃね?|統合失調症。)
74:a4 ◆L1L.Ef50zuAv
20/03/16 10:10:58 yKeJ8TdS.net
僕は星籍を地球からアルファ星に変えました。地球を量子で軍事威嚇します。
75:a4
20/03/16 11:20:55.58 yKeJ8TdS.net
申し訳ございません。僕の統合失調症の陽性症状の急性期でした。数学板的には
ジョン・ナッシュとして定義があるでしょうか。でも僕の中国の友達も困ってるみたい
なので書き込みさせていただきました。情報としては、流行ってる地域と報道され
てるのに、その友達の友達に感染した人はいるか?と聞いても、いない、と返って
来ました。僕は一
76:旦落ちますね。
77:132人目の素数さん
20/03/16 18:07:14.25 KM2jIN+Z.net
【y=X^2】イタリアさんの死者数と感染者数、指数関数的に増えていた
スレリンク(news板)
78:132人目の素数さん
20/03/18 15:35:44 shuQMR+Y.net
前川喜平氏が説く「数学必修廃止論」に疑問 出会い系バーでの「貧困調査」の具体的な成果なぜ示さないのか
前川喜平・前文部科学事務次官が週刊東洋経済4月14日号で、貧困対策の一つとして、「高校中退をなくすには数学の必修を廃止するのがいい」と発言している。
高校中退を防ぐという方向性はいいだろう。低学歴者は低所得になりがちで、犯罪率も高いことは各種の調査で示されている。
せっかく「貧困調査」で出会い系バーに通ったのだから、前川氏には延べ何人の女の子を調査し、その中で高校中退の数は何人だったのかを示してもらいたかった。
さらに、高校中退の理由はどうだったのか。これも貧困調査を行ったのならば当然示せるだろう。そうした調査の成果が具体的に示されていないので、実のところ、前川氏が説く数学必修廃止論の理由はよく分からない。
数学をある程度知らないと、自然科学のみならず多くの社会科学を習得することはできない。数学の必修廃止は日本の国力を低下させることにつながるのではないか。特に、社会に必要なエンジニアの育成にも支障が出るだろう。
文部科学省による調査をみても、高校の中退理由は、「学業不振」が1割弱、「学校生活への不適応」が4割、「進路変更」が3割強である。
つまり、数学必修を廃止しても、中退理由の1割も除くことができない。数学の必修化をやめれば中退が少なくなるとの結論を導き出すことはできないだろう。
これらの理由の推移をみると、かつては学業不振が多かったが、最近は低下しており、学校生活への不適応が徐々に増えている。
他校への転校などの進路変更はいいとして、学校生活の不適応をいかに減らすかが、中退を防ぐためには重要だろう。前川氏の出会い系バーにおけるフィールドワークに基づく具体的な対策を聞いてみたい。
URLリンク(www.zakzak.co.jp)
79:132人目の素数さん
20/03/19 21:03:10.68 2Eqy074l.net
ソースTBS
上昌広先生と森永先生
「ランダムに1000人にPCR検査すれば統計学的に市中感染率はわかるんですよ」
Twitterの東大生
「PCR検査が6割程度しか正確じゃないのに、なんでそんなことが言えるんだ」
↑
どっちがただしい?
東大生が正しい場合、6割だとして、
何人に検査すれば信頼区間99%に収まる?
80:132人目の素数さん
20/03/19 22:28:38.86 mW+UfAsE.net
>>68
URLリンク(www.mhlw.go.jp)
これか。厚労省もけっこう悲惨な事態を想定してんだな
81:132人目の素数さん
20/03/20 06:19:35 p5Mf5Wxl.net
>>76
PCRで診断確定するのに偽陰性ってどういうこと?
偽陰性と言われた人はどうやって新型コロナ肺炎だと確定されたんだ?
82:132人目の素数さん
20/03/21 00:26:49.38 Ivwmz468.net
都道府県ごとのシミュレーションによる検討
URLリンク(www.fttsus.jp)
83:132人目の素数さん
20/03/21 13:59:14 bagTkMOY.net
>>76
感染率をθとして検体をn回とった時の陽性者の数Xは
感度をp、(感染者が陽性と判定される確率)
特異率をq (非感染者が陽性と判定される確率)
とするとXの分布は平均が
pθ+q(1-θ)
となるのでこの値を推定すればθの値も推定できます。
感度-特異率
が正の値なら回数を増やせばいくらでも小さい99%信頼区間に入れることができます。
千回のとき実際どれくらいの市中感染者数をどれくらいの信頼度で測れるかは計算機ないとわかんね
84:132人目の素数さん
20/03/21 18:03:20.37 MXwMv06E.net
確率微分方程式ですね
85:132人目の素数さん
20/03/21 21:43:31 RyI2Q/uv.net
>>80
特異度qはTN/(TN+FP)だから
偽陽性率は1-qなので
X ~ pΘ+(1-q)(1-Θ)
では?
86:132人目の素数さん
20/03/21 23:12:32.18 hCC4s83x.net
>>80
特異度は「陰性のものを正しく陰性と判定する確率」です。(wikiより引用)
認識が間違っていたため、式も間違っているパターンですね。
82さんの式が正しい。
87:132人目の素数さん
20/03/21 23:25:44.12 bagTkMOY.net
ありゃ?
記憶と反対だったか。
88:132人目の素数さん
20/03/21 23:29:01.08 bagTkMOY.net
まぁそれなら
感度+特異度>1
のときは検査数を上げていいなら望むだけ市中感染率を正しく推定できる
ですな。
なので上先生の勝ち。
東大生の負け。
89:132人目の素数さん
20/03/21 23:54:36.88 RyI2Q/uv.net
カレーで免疫ができるかの検定
某国の新型コロナ感染症の有病率を0.1、PCR検査の感度を0.7 特異度を0.9とする
検査陽性陰性の人を無作為に各々100人ずつ集めて、カレーを頻回に食べているかを調査した結果が
以下の通りであった。
カレーを頻食する新型コロナ感染に罹りにくいと結論できるか?
検査陽性 検査陰性
カレー頻食 a(=36) b(=60)
カレー稀食 c(=64) d(=40)
PCR(+) PCR(-)
Exposed 36 60 96
Nonexposed 64 40 104
Total 100 100 200
にPPV(陽性的中率),NPV(陰性的中率)を使って計算すると
Disease Nondisease Total
Exposed 17.89286 78.10714 96
Nonexposed 29.42857 74.57143 104
Total 47.32143 152.67857 200
検査陽性は有意にカレー暴露が少ないといえるが、
感染に関しては有意にカレー暴露が少ないとはいえない、
という結論になった。
p値は 各々 0.0007032002 と0.1092375975
達人の検算を希望。
90:132人目の素数さん
20/03/22 01:10:50.10 liILqu/N.net
r=x/nとして
Θ = (q + r - 1)/(p + q - 1)
91:132人目の素数さん
20/03/22 08:39:00.57 liILqu/N.net
>>85
プログラムを組んで計算させてみた(有病率の事前分布は一様分布を仮定)
感度0.7 特異度0.9として
100人に1人陽性だった場合の有病率の期待値、中央値、最頻値、95%信頼区間
> sim(100,1)
mean median mode lower upper
0.0196 0.0166 0.0101 0.0004 0.0464
1000人に10人陽性だった場合の有病率の期待値、中央値、最頻値、95%信頼区間
> sim(1000,10)
mean median mode lower upper
0.0110 0.0107 0.0099 0.0050 0.0175
10000人に100人陽性だった場合の有病率の期待値、中央値、最頻値、95%信頼区間
> sim(10000,100)
mean median mode lower upper
0.0101 0.0101 0.0100 0.0082 0.0121
100000人に1000人陽性だった場合の有病率の期待値、中央値、最頻値、95%信頼区間
> sim(100000,1000)
mean median mode lower upper
0.0100 0.0100 0.0100 0.0094 0.0106
92:132人目の素数さん
20/03/22 08:53:48.86 liILqu/N.net
>>85
感度+特異度>1 の条件は必要?
感度0.4 特異度0.5でシミュレーションしてみたけど
nを増やせば信頼区間が狭くなっていく
> sim(100,1,sen=0.4,spc=0.5)
mean median mode lower upper
0.0196 0.0166 0.0102 0.0004 0.0463
> sim(1000,10,sen=0.4,spc=0.5)
mean median mode lower upper
0.0110 0.0107 0.0099 0.0050 0.0176
> sim(10000,100,sen=0.4,spc=0.5)
mean median mode lower upper
0.0101 0.0101 0.0101 0.0082 0.0121
> sim(100000,1000,sen=0.4,spc=0.5)
mean median mode lower upper
0.0100 0.0100 0.0100 0.0094 0.0106
93:132人目の素数さん
20/03/22 09:13:41.93 t2LOPpzl.net
>>89
ああ、p+q<1だと逆にあてにならなすぎて逆張りしてθが推定できるんだな。
94:132人目の素数さん
20/03/22 14:31:30.77 liILqu/N.net
プログラムにバグがあったので修正というより、stanでのMCMCに変更。
> PCR(100,1)
mean lower upper
0.01497496587 0.00004299248 0.04480331292
> PCR(1000,10)
mean lower upper
0.001642898836 0.000001299175 0.004982795701
> PCR(10000,100)
mean lower upper
0.000174692810 0.000000052897 0.000560041026
> PCR(10000,1000)
mean lower upper
0.003876080828 0.000002574352
95: 0.010176384444 まあ、nを増やすほど信頼区間が狭くなるという結論には変わりない。
96:132人目の素数さん
20/03/22 20:15:28 liILqu/N.net
【富山県最強伝説】新型コロナウイルスPCR検査件数 54人 陽性0人
スレリンク(newsplus板)
ある集団から54人を無作為に選んでPCR検査したら陽性0であった。PCR検査の感度0.7 特異度0.9としてこの集団の有病率の期待値と95%信頼区間を求めよ。
97:132人目の素数さん
20/03/22 21:37:38.27 liILqu/N.net
>>92
> PCR(54,0)
mean lower upper
0.0282436750 0.0000053869 0.0839653803
98:132人目の素数さん
20/03/23 03:28:01.22 uvHIelYA.net
日本人の平均身長を推測するのにその値は1~2mの間であるという弱情報事前分布は合理的。
現時点での新型コロナの有病率は0.1未満の一様分布という弱情報事前分布として
【富山県最強伝説】新型コロナウイルスPCR検査件数 54人 陽性0人
ある集団から54人を無作為に選んでPCR検査したら陽性0であった。感度0.7 特異度0.9としてこの集団の有病率の期待値と9信頼区間を推測する。
事前分布のパラメータを変えるとstanだとコンパイルが必要になるのでjagsでプログラムを組んでみた。
# 感度SEN, 特異度SPCの検査でN人中X人が陽性であったときの推定有病率prevalence
# 弱情報事前分布はprevalence ~ dunif(0,UL) , UL:上限
library(rjags)
PCRj <- function(N,X,UL=1,SEN=0.7,SPC=0.9,verbose=TRUE){ # UL:upper limit of dunif(0,UL)
modelstring=paste0('
model
{
x ~ dbin(p,n)
p <- prev*sen + (1-prev)*(1-spc)
prev ~ dunif(0,',UL,')
}
')
if(verbose & UL!=1) cat(modelstring)
writeLines(modelstring,'TEMPmodel.txt')
dataList=list(n=N,x=X,sen=SEN,spc=SPC)
jagsModel = jags.model( file="TEMPmodel.txt" ,data=dataList,quiet=TRUE)
update(jagsModel)
codaSamples = coda.samples( jagsModel ,
variable=c("prev","p"), n.iter=1e6, thin=10)
js=as.matrix(codaSamples)
BEST::plotPost(js[,'prev'],xlab='prevalence',showMode = TRUE) ; lines(density(js[,'prev']),col='skyblue')
round(c(mean=mean(js[,'prev']),HDInterval::hdi(js[,'prev'])[1:2]),10)
}
実行結果
> PCRj(54,0,UL=0.1)
model
{
x ~ dbin(p,n)
p <- prev*sen + (1-prev)*(1-spc)
prev ~ dunif(0,0.1)
}
|**************************************************| 100%
mean lower upper
0.0245104429 0.0000003782 0.0703606657
99:132人目の素数さん
20/03/23 09:55:27 LegnQVLy.net
統計のことははぜんぜんわからんが、確率論的には
検査陽性率の期待値=有病率×感度+(1-有病率)×(1ー特異度)
なんだから、
有病率<<1なら、検査陽性率の期待値≒有病率×感度+(1- 特異度)
っちゅうことで、感度70%で特異度が100%なら検査陽性率の3割増し
が有病率だとみなせばよろしいんでしょ。
一方、特異度が90%しかなかったりすると陽性率の期待値が10%も
水増しされちゃうから、有病率の推定が大幅に困難になる。
つまり、特異度がよほど高くなければ(99%とかね)、有病率が数%以下
の状況でランダム検査しても偽陽性が真陽性を上回って混乱をきたす。
(見かけ上致死率は下がる、ってのが南朝鮮の状況か)
100:132人目の素数さん
20/03/23 11:25:41 eyONmTBV.net
>>95
同感。
有病率0.1%で特異度0.9なら偽陽性だらけになるんだよね。
101:132人目の素数さん
20/03/23 14:59:33 9TP9mpqz.net
>>95-96
疾病率=陽性率、すなわち無作為抽出ならその算段も成立するかもね。
現行の制度では推定される市中感染率が1/10000程度で陽性率が5%ほどらしいから現行制度下での検査はうまく行ってますね。
ただ感度+特異度が1.7位あるので検査数増やした方が有益である事は間違いがない。
102:132人目の素数さん
20/03/23 19:12:11 O5lTfF0I.net
>>95
誤解や、誤解を引き起こしかねない内容があったので、勝手に補足させていただきます。
>>っちゅうことで、感度70%で特異度が100%なら検査陽性率の3割増し
>>が有病率だとみなせばよろしいんでしょ。
この場合、検査陽性率の期待値=有病率×感度 なのだから、1/(0.7)=1.4285...
3割増しではなく、4割強増しと言うべき。
>>一方、特異度が90%しかなかったりすると陽性率の期待値が10%も
>>水増しされちゃうから、有病率の推定が大幅に困難になる。
感度70、特異度100で、有病率 1%、0.1%、0.01%、0.001%の時、10万人に検査したときの
検査陽性者数は700,70,7,0.7人です。 有病率に比例して増減し、これらは、はっきり見極めできます。
一方、感度70、特異度90で、同じ事をすると、それぞれ、10600,10060,10006,10000.6人です。
有病率に応じて差はあるのですが、常に偽陽性が約10000人いて、誤差を考えると、見極めは困難です。
「陽性率の期待値が10%も水増しされちゃうから」と書かれていますが、これは、
検査対象者10万人に対する約10%=1万人が常に水増しされているのであって、
陽性と判定される人の数が10%水増しされていると誤解しないよう、補足しておきます。
103:132人目の素数さん
20/03/24 01:26:59.63 EUfp1x4d.net
>>98
心配性だねw
3割が4割でもたいして変わらんがな。
陽性率の期待値が10%水増しってのも、期待値の10%じゃなくて、
期待値そのものが10%上昇するんだってことくらい、元の式見りゃ
自明だしね。
そもそも感度や特異度の具体的な数値はよくわかんないんだから、
具体的な数字にこだわってもしょうがない。
有病率が低いときの、有病率と陽性率と特異度の関係がわかれば
よろし。
104:132人目の素数さん
20/03/24 01:40:13 TnHQvRcs.net
こういう計算が必要になる。
事前分布を選択する(例. 有病率は高々10%として(0.0.1]の一様分布とする)、
陽性確率は真陽性確率と偽陽性確率の和、
陽性数はこの確率で二項分布、
以上を実際に得られた検査数と陽性数から最尤値となるパラメータとして有病率の分布を出して期待値や信頼区間を出す。
手計算では無理。
105:132人目の素数さん
20/03/24 08:34:40.14 TnHQvRcs.net
感度0.7 特異度0.9でコンピュータに計算させると
陽性率が30%なら有病率の推測値は
> PCRj(100,30)
mean lower upper
0.3396123 0.1994400 0.4964517
> PCRj(1000,300)
mean lower upper
0.3343576 0.2876152 0.3832246
> PCRj(10000,3000)
mean lower upper
0.3333387 0.3186388 0.3481056
> PCRj(100000,30000)
mean lower upper
0.3332425 0.3286774 0.3380952
と 期待値も信頼区間もそれらしい値になるが、
陽性率が1%なら有病率の推測値は偽陽性が増えまくるので陽性率と有病率が乖離する、
> PCRs(100,1)
mean lower upper
0.01497496587 0.00004299248 0.04480331292
> PCRs(1000,10)
mean lower upper
0.001642898836 0.000001299175 0.004982795701
> PCRs(10000,100)
mean lower upper
0.000174692810 0.000000052897 0.000560041026
> PCRs(100000,1000)
mean lower upper
0.000016780530193 0.000000002738145 0.000051489256688
陽性率が60%なら、今度は偽陰性が増えまくるので偽陰性が増えまくるので陽性率と有病率が乖離する。
> PCRs(100,1)
> PCRs(100,60)
mean lower upper
0.8280581 0.6766480 0.9764282
> PCRs(1000,600)
mean lower upper
0.8335023 0.7826355 0.8825766
> PCRs(10000,6000)
mean lower upper
0.8334634 0.8187334 0.8492064
> PCRs(100000,60000)
mean lower upper
0.8332873 0.8289242 0.8379535
106:132人目の素数さん
20/03/24 08:42:02.80 TnHQvRcs.net
Rとstanでベイズ統計ができるなら、以下のコードで実行可能。
パッケージ rstan と BEST(信頼区間グラフ描出用)が必要
library("BEST")
library("rstan")
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
options(scipen = 5)
model.string='
data{
int n; // sample size
int x; // number of positive test
real<lower=0,upper=1> sen; // sensitivity 0.7
real<lower=0,upper=1> spc; // specificity 0.9
real<lower=0,upper=1> ul; // uniform(0,ul)
}
parameters{
real<lower=0,upper=1> prev; // prevalence
}
transformed parameters{
real<lower=0,upper=1> p;
p = prev*sen + (1-prev)*(1-spc) ; // probability of positive test result
}
model{
x ~ binomial(n,p);
prev ~ uniform(0,ul);
}
'
writeLines( model.string , con='model.stan' )
corona1.model=stan_model('model.stan')
# saveRDS(corona1.model,file='corona1.rds')
# corona1.model=readRDS('corona1.rds')
PCRs <- function(N=1000,X=10,UL=1,SEN=0.7,SPC=0.9,verbose=FALSE,...){
data = list(n=N,x=X,sen=SEN,spc=SPC,ul=UL)
fit.corona = sampling(corona1.model, data=data,
seed=1234,control=list(adapt_delta=0.99),...)
if(verbose) print(fit.corona, prob=c(0.025,0.5,0.975),pars=c('prev'),digits=8)
ms=rstan::extract(fit.corona)
BEST::plotPost(ms$prev,showMode = T,xlab='prevalence') ; lines(density(ms$prev),col='skyblue')
c(mean=mean(ms$prev),HDInterval::hdi(ms$prev)[1:2])
}
107:132人目の素数さん
20/03/24 08:46:15.46 TnHQvRcs.net
>>102(補足)
非対称分布の信頼区間計算にパッケージHDIntervalも必要。
108:132人目の素数さん
20/03/24 08:51:43.98 TnHQvRcs.net
>>76
偽陽性率が現時点での有病率を大きく上回るから東大生の言い分が正しい。
有病率が3割程度になれば上先生の言い分が正しい。
109:132人目の素数さん
20/03/24 10:27:48 TnHQvRcs.net
>>104
1000人調べたときの検査陽性率と推定陽性率をグラフにしてみた。
灰色直線は検査陽性率=推定陽性率の直線
検査陽性率が低いときは過小評価、高いときは過大評価する。
URLリンク(i.imgur.com)
110:132人目の素数さん
20/03/24 10:29:17 mBslr8ul.net
何言ってんの?
市中感染率をいくらでも正しく推定できるかなんでしよ?
結論のために問題変えるなよ。
政治の話をここに持ち込むなよ。
111:132人目の素数さん
20/03/24 10:52:19 EUfp1x4d.net
>>104
上の言い分も一理あるし、東大生の言い分も一理あるw
特異度が90%を越える高い値だという前提があればこうなる。
1)検査陽性率が低く(数%以下)、特異度がほぼ100%でないの
なら、市中感染率を推定するのは難しい。とはいえ、上限は
(検査陽性率/感度)程度で抑えられる。つまり10%以下くらいの
ことは言えるが、0%かもしれない。
2)一方、検査陽性率が高ければ(数十%以上)、下限は検査
陽性率程度と見込めるが、市中感染率は感度に依存して大きく
変化する。つまり、数十%以上とは言えるが100%近いかどうか
までは不明。
ってことで、検査陽性率からある程度市中感染率の目安は立つが、
それがどこまで意味があるとみなせるかは疑問。TPO次第か。
112:132人目の素数さん
20/03/24 10:58:55 EUfp1x4d.net
あと、補足すると、感度と特異度が正確にわかっているのなら、
統計学的に市中感染率を推定することはある程度可能だけど、
実際はそうではないから、上様の言い分には意味がない。
113:132人目の素数さん
20/03/24 11:03:07 TnHQvRcs.net
>>105
感度と特異度を変化させて、検査陽性率と推定有病率の関係をグラフにしてみた。
URLリンク(i.imgur.com)
114:132人目の素数さん
20/03/24 11:06:16 mBslr8ul.net
だから自分の言いたい結論が先に決まっててそれに合わせて好き勝手に問題読み替えてるだけじゃん?
そんな考え方しかできないなら理系板に来んなよ。
115:132人目の素数さん
20/03/24 11:09:05 TnHQvRcs.net
>>108
感度も特異度も定数でなく何らかの分布に従うパラメータとしてモデルを組めばいいだけの話だろ。
116:132人目の素数さん
20/03/24 11:48:05 TnHQvRcs.net
感度が最頻値0.7 標準偏差0.05のベータ分布β(58.229, 25.527 )
特異度が最頻値0.9 標準偏差0.05のベータ分布β(36.172, 4.908)
に従うと過程して
model
{
x ~ dbin(p,n)
p <- prev*sen + (1-prev)*(1-spc)
sen ~ dbeta(sn[1],sn[2])
spc ~ dbeta(sp[1],sp[2])
prev ~ dunif(0,ul)
}
こういうモデルでMCMCすれば可能。
実行結果
> PCRj2(1000,10) # 陽性率1%で有病率を推定
mean lower upper
0.001667604624 0.000000053909 0.004956969423
> PCRj2(1000,300) # 陽性率30%で有病率を推定
0.33414 0.28797 0.38253
> PCRj2(1000,600) # 陽性率60%で有病率を推定
mean lower upper
0.8
117:3296 0.78428 0.88496
118:132人目の素数さん
20/03/24 11:56:18 TnHQvRcs.net
β分布のパラメータを出すRスクリプト
# a,b to Mode,mean,variance
ab2Mmv<-function(a,b){
M<-(a-1)/(a+b-2)
m<-a/(a+b)
v<-a*b/((a+b)^2*(a+b+1))
cat('Mode =',M,'mean =',m,'variance =',v,'\n')
invisible(c(Mode=M,mean=m,variance=v))
}
# Mode,kappa to mean,variance
Mk2mvab= function( mode , kappa ) {
# if ( mode < 0 | mode > 1) stop("must have 0 <= mode <= 1")
# if ( kappa <2 ) stop("kappa must be >= 2 for mode parameterization")
a = mode * ( kappa - 2 ) + 1
b = ( 1.0 - mode ) * ( kappa - 2 ) + 1
m=a/(a+b)
v=m*(1-m)/(a+b+1)
return( c( mean=m , variance=v,a=a,b=b ) )
}
# Mode,variance to a,b
Mv2ab = function(mode,vari){
f=function(kappa) Mk2mvab(mode,kappa)[2] - vari
kappa=uniroot(f,c(2,10000))$root
ab=Mk2mvab(mode,kappa)[c('a','b')]
ab2Mmv(ab[1],ab[2])
return(ab)
}
(sn=Mv2ab(0.7,0.05^2))
curve(dbeta(x,sn[1],sn[2]),bty='l')
(sp=Mv2ab(0.9,0.05^2))
curve(dbeta(x,sp[1],sp[2]),bty='l')
119:132人目の素数さん
20/03/24 11:57:10 TnHQvRcs.net
上記の準備をして以下で実行
PCRj2 <- function(
N,X,
UL=1,
SEN=0.7,
SPC=0.9,
SD=0.05,
print=TRUE){
# UL:upper limit of dunif(0,UL)
library(rjags)
library(BEST)
sn=Mv2ab(SEN,SD^2)
sp=Mv2ab(SPC,SD^2)
modelstring=paste0('
model
{
x ~ dbin(p,n)
p <- prev*sen + (1-prev)*(1-spc)
sen ~ dbeta(sn[1],sn[2])
spc ~ dbeta(sp[1],sp[2])
prev ~ dunif(0,ul)
}
')
writeLines(modelstring,'TEMPmodelj.txt')
dataList=list(n=N,x=X,ul=UL,sen=SEN,spc=SPC,sn=sn,sp=sp)
jagsModel = jags.model( file="TEMPmodelj.txt" ,data=dataList, quiet=TRUE)
update(jagsModel)
codaSamples = coda.samples( jagsModel ,
variable=c("prev","p","sen","spc"), n.iter=1e5, thin=5)
js=as.matrix(codaSamples)
if(print){
BEST::plotPost(js[,'prev'],xlab='prevalence',showMode = TRUE)
lines(density(js[,'prev']),col='skyblue')}
re=c(mean=mean(js[,'prev']),HDInterval::hdi(js[,'prev'])[1:2])
return(re)
}
options(digits = 5)
options(scipen = 5)
PCRj2(1000,10) # 陽性率1%で有病率を推定
PCRj2(1000,300) # 陽性率30%で有病率を推定
PCRj2(1000,600) # 陽性率60%で有病率を推定
120:132人目の素数さん
20/03/24 13:15:20 /QqkwKRd.net
>>99
期待値というのは、無次元量ではない。観測値とか物理量と同じように単位をつけて議論できる量。
従って「期待値が10%増える」等という言葉があれば、期待値が1.1倍になるのだろうと感じるのが普通。
そのような性質を持つ期待値に対し、「10%増える」と表現し、
「期待値の値そのものが、0.1増えることを意味している」
と説明しなければならないならば、やはり誤解を招きやすい表現だと思う。
今回の期待値は比率であり、無次元量であったから、「10%」と言うのが、
どちらの意味としても、通用したため発生したとは言えるが、読み手の立場に立った表現を望む。
似た議論に、選挙時の投票率がある。前回の投票率が40%。今回の投票率が50%だとする。
「前回に比べ、今回は10%増えました」
「前回に比べ、今回は25%増えました」
どちらも、言い得る表現。聞き手の混乱を避けるため、前者の意味で使う場合、
「10%ポイント増えました」とコメントするのを最近聞くようになった。
私にはよい傾向と感じるが、中には、違いは何かとか、混乱の源の存在さえ意識していない人もいるようだ。
「3割増も4割強増も大した差ではない」には、「式が違っても結果が誤差範囲なら問題ない」
という考えが背景に見える。そのような方が、混乱を引き起こしかねない表現を用いた。
だから、補足した。果たして本当に杞憂だったのだろうか?
121:132人目の素数さん
20/03/24 14:42:23 TnHQvRcs.net
富山では62人PCR検査して陽性0人(3月22日までの集計)有病率を推定とその信頼区間を推定したい。
URLリンク(www.pref.toyama.jp)
PCR検査の感度は最頻値0.6標準偏差0.1、特異度は最頻値0.9標準偏差0.05のベータ分布(正規分布は負になったり1を超えるので不適)、
有病率は一様分布として、推定される有病率の期待値と95%を計算せよ。
図示するとこんな感じ。
URLリンク(i.imgur.com)
stanのモデルのスクリプトはこれ
sn,spはβ分布のパラメータ、その計算法は既述
data{
int n; // sample size
int x; // positive test result
real<lower=0,upper=1> ul; // uniform(0,ul)
real<lower=0> sn[2]; // sen ~ beta(sn[1],sn[2])
real<lower=0> sp[2]; // spc ~ beta(sp[1],sp[2])
}
parameters{
real<lower=0,upper=1> prev; // prevalence
real<lower=0,upper=1> sen; // sensitivity
real<lower=0,upper=1> spc; // specificity
}
transformed parameters{
real<lower=0,upper=1> p;
p = prev*sen + (1-prev)*(1-spc) ; // probability of positive test result
}
model{
x ~ binomial(n,p);
prev ~ uniform(0,ul);
sen ~ beta(sn[1],sn[2]);
spc ~ beta(sp[1],sp[2]);
}
122:132人目の素数さん
20/03/24 14:54:23 TnHQvRcs.net
>>116
ここで問題
感度特異度の分布はそのベータ分布として
何人陰性が続けば95%信頼区間の上限が0.05を下回るか?
123:132人目の素数さん
20/03/24 15:14:26 mBslr8ul.net
>>117
感度、特異度の分布???
124:132人目の素数さん
20/03/24 15:56:48 TnHQvRcs.net
>>118
何でも確率変数にするのがベイズ推計。
p値の分布すら考えるぞ。
125:132人目の素数さん
20/03/24 16:10:56 TnHQvRcs.net
PCR検査の感度は最頻値0.6標準偏差0.1、特異度は最頻値0.9標準偏差0.05のベータ分布を事前分布にしたけど、
事後分布はstanによるMCMCで
感度は期待値0.57 95%信頼区間は[0.37,0.77]
特異度は期待値0.96 95%信頼区間は[0.91,0.99]
とコンピュータが計算してくれる。
126:132人目の素数さん
20/03/24 17:38:05 TnHQvRcs.net
>116のように弱情報事前分布を設定することで事後分布は次のように描ける。
URLリンク(i.imgur.com)
127:132人目の素数さん
20/03/24 17:43:56 TnHQvRcs.net
>>54
いや、特異度の事前分布を設定することで事後分布をMCMCで求めることができる。
>116の設定での結果が>121
128:132人目の素数さん
20/03/24 18:10:37.68 EUfp1x4d.net
>>122
結局事前分布の設定次第ってことはないの?
129:132人目の素数さん
20/03/24 18:45:58.22 TnHQvRcs.net
>>123
日本人の平均身長を推測するのにその値は1~2mの間であるという弱情報事前分布は合理的。
感度特異度の分布に正規分布を使うのはアホ。
負になったり、1を超えたりするから。
130:132人目の素数さん
20/03/24 19:08:27.01 TnHQvRcs.net
>>123
感度を0.4-0.8の一様分布、特異度を0.8-1.0の一様分布にしても有病率の推定値は
> round(re$mci,5)
mean lower upper
0.02827 0.00000 0.08592
であまり変わらないね。
131:132人目の素数さん
20/03/24 19:20:40 TnHQvRcs.net
sensitivity ~ N(m=0.6,sd=0.1) specificity ~ N(m=0.9, sd=0.05)
にしても推測有病率は平均3%弱で 95%CIは0-8%とあまり分布の形にはよらないね。
mean lower upper
0.026841384 0.000000153 0.081071379
確率だと定義域が0-1で計算しやすいのでβ分布を使うことが多い。
132:132人目の素数さん
20/03/24 21:43:16 mBslr8ul.net
>>119
そう?
統計の推定の理論で推計する母数は確率変数ではないと習ったけど?
133:132人目の素数さん
20/03/25 05:45:40 jmNOx22O.net
>>127
時代は頻度主義統計からベイズ統計だよ。
134:132人目の素数さん
20/03/25 06:04:40 jmNOx22O.net
頻度主義統計でも最尤推定では
データを固定してパラメータを動かすだろ。
135:132人目の素数さん
20/03/25 06:30:36 jmNOx22O.net
>>127
階層ベイズモデルを扱ったことないの?
>112は簡単な実例。
136:132人目の素数さん
20/03/25 07:18:32 yWXBkNWD.net
>>128 >>130
何を持ってベイズ統計っていってんのか知らん。
pcr検査の感度とは被験者が感染者である場合の検査結果が陽性となる条件付き確率でしょ?
条件付き確率の分布ってどういうことよ?
確率がまた確率になるってなんの話してんの?
変数Xの平均とか分散とかは統計学においては推定すべき定数であって確定値。
それの分散なんて数学的に意味不明。
一体どこの統計学の教科書にそんなデタラメ書いてあんの?
137:132人目の素数さん
20/03/25 07:23:50 r1V62jxn.net
まちがえた。
確率の平均がまた確率変数になるってどういうことよ、ね。
式でかけば確率変数Xの平均E(X)の分散ってなんの話ってことになる。
確率変数Xはある標本空間上の関数だけどE(X)は実数だよ?
138:132人目の素数さん
20/03/25 09:57:10 2o2
139:7M3ww.net
140:132人目の素数さん
20/03/25 10:06:13 2o27M3ww.net
>>131
ベータ分布は定義域が[0,1]で二項分布の確率の確率密度関数としてベイズ階層モデルでは頻用されるよ。
ベイズ階層モデルを使わずにこの計算できるならやってみてくれ。
020/3/24 11:00時点で検査人数での陽性率は171/2013であるという。
新型コロナ肺炎のPCR検査の感度は5~7割、特異度は9割前後らしい。幅をもたせた値を使って検査をうけたグループの有病率を計算せよ。
141:132人目の素数さん
20/03/25 11:27:50 82yASlvk.net
>>133
まぁ言わんとする事はもちろんわかるし伝わるけど、疫学だから数学やってる人間がなんとなく伝わるではダメだろ?
数学だけの話ではなく、疫学は実社会とキチンと繋がってるんだから?
統計学ではあくまで検定する母数は定数。
それは確率モデルでは実数値であり、定数。
そして統計データを確率変数に割り当てる。
当然それらの確率変数は一つの測度空間の一つしかない確率変数であり、平均も分散もひとつしかない定数値。
それらをいっぱい考えてどうこう言ってるんだろうとは思うけどそんなの統計学や疫学の一般的な考えにはない。
何故なら現実世界はひとつしかなく、確率変数に対応している統計量も一個しかない。
もちろん母数がめちゃめちゃ大きい統計量で例えば10000個のデータを100こずつ切って100個の統計量を100の世界からとってきたなんて考えが無理クリできなくはないが、そんな考え方は普通しない。
それはあくまで100個ずつに区切られた10000個の一つの世界の確率変数としか扱わない。
そういうオリジナルな考えで捉えたいならそれは勝手だけど、それならそれで話の中で明示しないとダメ。
数学の世界なら言わずもがなの話は言わなくてもエスパーしてもらえても、疫学、統計学の世界では実社会とつながる話だからダメ。
142:132人目の素数さん
20/03/25 12:30:50.29 jmNOx22O.net
>>135
能書きいいから、
ベイズ階層モデルを使わずにこの計算できるならやってみてくれ。
020/3/24 11:00時点で検査人数での陽性率は171/2013であるという。
新型コロナ肺炎のPCR検査の感度は5~7割、特異度は9割前後らしい。幅をもたせた値を使って検査をうけたグループの有病率を計算せよ。
143:132人目の素数さん
20/03/25 12:39:07.64 jmNOx22O.net
>>136
こういう判断が現実には必要。
検査特性を無視して単純な割り算だと検査を受けた人の有病率は8.5%弱になるけどこれは過大評価か過小評価か?
144:132人目の素数さん
20/03/25 14:59:13.04 jmNOx22O.net
検査感度が5-7割、特異度が9割前後なら
検査陽性率=有病率とすると常に過大評価かどうか気になったので陽性数を変化させて計算してみた。
検査感度はmode=0.6,sd=0.1 特異度はmode=0.9,sd=0.05のベータ分布に設定してJAGSでベイズ階層モデルをたてて計算。
URLリンク(i.imgur.com)
陽性率が20%未満のときは過大評価、それ以上のときは過小評価である、という結論になった。
ベイズ統計を理解できている人の検証希望。
145:132人目の素数さん
20/03/25 17:30:41 jmNOx22O.net
>>138
プログラムの練習がてらに、
MCMCのアルゴリズムの異なるstanでベイズ階層モデルを組んで検証。
当然ながら、同様の結果。 検査陽性率が20%を境に過大評価と過小評価が入れ替わる。
URLリンク(i.imgur.com)
146:132人目の素数さん
20/03/25 21:21:08.00 jmNOx22O.net
>>136(自己レス)
今日の都の発表で(171+41)/(2013+89) に検査陽性率が増えたので再計算。
URLリンク(i.imgur.com)
147:132人目の素数さん
20/03/25 21:33:22.10 jmNOx22O.net
>>138 サンプリング回数を増やしてグラフを完成。 https://i.imgur.com/kLjCD2y.png
149:132人目の素数さん
20/03/26 16:25:58 +rQz06p5.net
>>140
89は検査数で検査人数は74という。
計算し直すと
> PCRj2(N,r,SEN=0.6,SD1=0.1,SPC=0.9,SD2=0.05,N.ITER=5e5)
|**************************************************| 100%
mean lower upper
0.05720165 0.00000015 0.1332385
150:132人目の素数さん
20/03/26 16:34:28 +rQz06p5.net
41/74の推測有病率は
mean lower upper
0.8121975 0.5957315 0.9999992
151:132人目の素数さん
20/03/27 11:07:27 sdGiAEI7.net
オリンピック延期発表後の検査陽性率は88/169で52%だが、
PCR検査の感度と特異度がはっきりしないので、検査陽性率をこの集団の有病率とするのは正しくない。
88/169のときの感度・特異度と推定有病率の関係をグラフにしてみた。
URLリンク(i.imgur.com)
感度0.6、特異度0.9のときの推定有病率は85%で陽性率からの憶測は過小評価といえる。
152:132人目の素数さん
20/03/27 18:36:04.97 8rq7DP6B.net
検査陽性率が小さいときには、実際の有病率より過大評価してるし、
検査陽性率が高いときは、過小評価してるだろうってことでしょ。
そのくらいは定性的に理解できる。
153:132人目の素数さん
20/03/27 20:59:08.67 sdGiAEI7.net
>>145
どこが境目かは直感じゃわからんね。
154:132人目の素数さん
20/03/27 22:15:19.84 8rq7DP6B.net
そりゃ感度や特異度次第だからな。
まあ、数%と数十%では違うんだということがわかればいいんじゃね?
境目なんかどうでもいいでしょ。
155:132人目の素数さん
20/03/27 22:22:39.69 sdGiAEI7.net
陽性率が15%でこれを有病率の推測値に使うのは過大評価なのか過小評価がわからんのはまずいね。
156:132人目の素数さん
20/03/27 23:24:27 sdGiAEI7.net
オリンピック延期決定以後の検査数と陽性数
subjects=c(74,95,87)
positives=c(41,47,40)
PCRs3(subjects,positives,iter=10000,warmup=1000)
として、
感度・特異度を考慮した推定有病率は
mean lower upper
0.77417 0.56756 0.99944
>
日々の陽性数が二項分布に従うとして計算。
157:132人目の素数さん
20/03/28 03:24:30.89 NK6wIjWT.net
志村けんみたいな有名人がコロナに感染してることから日本全体のコロナ感染者数を推定してみる。
まず、日本の有名人が1000人いるとしよう。
つぎに、日本でコロナに感染していない確率をxとしよう。
すると、有名人1000人が一人も感染していない確率は、xの1000乗となる。これをyとおこう。すると、有名人が一人でも感染している確率は(1-y)となる、これをzとおこう。
まとめると以下の関係がなりたつ。
・コロナに感染しない確率:x
・有名人が一人もコロナに感染しない確率:y=x^1000
・有名人が一人でもコロナに感染している確率:z = 1-y
158:132人目の素数さん
20/03/28 03:25:09.78 NK6wIjWT.net
志村は感染したわけなので、以下、2つのケースにわける
ケース1: zが10%のとき
z=0.1, 故にy = 1-0.1=0.9
故にx = y^0.001よりx=0.9^0.001=0.999894
これがコロナに感染していない確率なので、
コロナ感染確率は、1-0.999894=0.000106
よって日本のコロナ感染者数は推定
120,000,000*0.000106=12,720人
ケース2:zが50%のとき
ケース1と同様の計算で、
日本のコロナ感染者数は推定
120,000,000*0.000693=83,160人
159:132人目の素数さん
20/03/28 09:37:19.12 uwBdnirU.net
検査が少ないから感染者増が緩やか?数学的に検証してみた
URLリンク(agora-web.jp)
主な関係国について、新型コロナ感染者数の片対数グラフがある。
URLリンク(agora-web.jp)
FT.comより
感染者数の伸びが日本は緩やかと解釈するのが普通だが、検査が少ないからとする解釈もある。本当はどうなのか計算してみる。
結論を先に書くと、検査が多いか少ないかは関係ない。
160:132人目の素数さん
20/03/28 09:40:23.67 QZo3p56d.net
対数をとると係数(感染者の発見率)は定数項になり、今回の片対数グラフ�
161:フ整理法の前提としてキャンセルされる 日本が展開しているのは患者認定の精度上昇であり、医療リソースの効果を最大化して死者数を低く抑えている要因の一つといえる
162:132人目の素数さん
20/03/28 10:39:35.99 BJlezchp.net
キャバクラ客100人から無作為に5人から検体を採取してこの検体を混合攪拌してコロナ検査したところ陽性であった。
(1)100人のキャバクラ客の陽性数の期待値と95%信頼区間を求めよ。
(2)PCR検査の感度0.6、特異度0.9として100人のキャバクラ客の感染数の期待値と95%信頼区間を求めよ。
163:132人目の素数さん
20/03/28 11:58:45.95 BJlezchp.net
>>151
> m=1000 # 有名人の人数
> n=1.268e5 # 日本の人口
> x=0:n # 感染者数:x, 非感染数:n-x
> pmf=1- chooseZ(n-x,m)/chooseZ(n,m) # 1 - (m人全員非感染の確率)
> pdf=pmf/sum(pmf) # 確率密度関数化して
> (E=sum(x*pdf)) # 期待値を計算
Big Rational ('bigq') :
[1] 63590201/1002
> as.numeric(E)
[1] 63463.27
6万3000人と計算された。
164:132人目の素数さん
20/03/28 12:46:07.99 NK6wIjWT.net
>>155
良く分からんが、ありがとう。
こちとら高校レベルの確率の知識しかないもんで。
165:132人目の素数さん
20/03/28 15:42:06.57 BJlezchp.net
>>156
n(=10)人の中にi人の感染者がいるとき無作為にm(=2)人を選ぶ。
選ばれた2人の中に少なくとも一人の感染者がいる確率をP[x]として、
n個からr個選ぶ組み合わせの数をChoose(n,r)で表すと
P[xi]=1- choose(10-x,2)/choose(10,2)
xを0から10まで変化させて、
Σx*P[x]/(ΣP[x])で
期待値が求まる。
166:132人目の素数さん
20/03/28 15:42:43.27 BJlezchp.net
タイプミス修正
P[x]=1- choose(10-x,2)/choose(10,2)
167:132人目の素数さん
20/03/28 16:07:51.35 qsSYTF8t.net
何このアホスレ?
168:132人目の素数さん
20/03/28 16:53:08.32 BJlezchp.net
有名人の数を増やしてみても同様の結果になった。
> # 有名人が感染
> library(gmp)
> m=18200 # 有名人の数(桜を見る会参加人数)
> n=1.268e5 # 日本の人口
> x=0:n # 感染者数:x, 非感染数:n-x
> pmf=1- chooseZ(n-x,m)/chooseZ(n,m) # 1 - (m人全員非感染の確率)
> pdf=pmf/sum(pmf) # 確率密度関数化して
> (E=sum(x*pdf)) # 期待値を計算
Big Rational ('bigq') :
[1] 1154070201/18202
> as.numeric(E) # E=63463.27 (m=1000) , E=1154070201/18202=63403.48(m=1.268e5)
[1] 63403.48
169:132人目の素数さん
20/03/28 19:42:56.69 NK6wIjWT.net
>>160
なんだってー。直感に反するな
170:132人目の素数さん
20/03/29 09:23:06.20 WogCQeQk.net
>>161
総人口100人として有名人の数を1~100人まで変化させて、有名人に感染者がいたときの100人中の感染者の数をグラフにすると
URLリンク(i.imgur.com)
有名人の数を変化さえても期待値にさほどの変化はない。
171:132人目の素数さん
20/03/29 10:18:20.74 2PsxdXJm.net
>>162
感染者が1名以上という条件だと、
有名人の割合が一定以上になると飽和するんだな。
172:132人目の素数さん
20/03/29 10:39:47.55 WogCQeQk.net
Ax: x人の感染者がいる(x=0~n)という事象
B:最低一人の感染陽性判定という事象
Pr[Ax|B]=Pr[B|Ax]Pr[Ax]/Pr[B]
Pr[Ax]:事前確率
Pr[B|Ax]:尤度
Pr[B]:周辺尤度(規格化定数)
求めたい期待値Eは
Σ(x*Pr[Ax|B])/ΣPr[Ax|B] = Σ(x*Pr[B|Ax]Pr[Ax])/Σ(Pr[B|Ax]Pr[Ax])
Pr[Ax]がxにかかわらず定数であれば
E=Σ(x*Pr[B|Ax])/Σ(Pr[B|Ax])
事前確率分布を一様分布と仮定しての計算
つまり、感染者が1人の確率も50人の確率も100人の確率,....も一定という前提での計算。
173:132人目の素数さん
20/03/29 10:47:28.57 WogCQeQk.net
>>163
そうみたいですね。
> data.frame(有名人=1:10,期待値=sapply(1:10,function(x) fn(100,x)$mean))
有名人 期待値
1 1 67.00000
2 2 62.75000
3 3 60.20000
4 4 58.50000
5 5 57.28571
6 6 56.37500
7 7 55.66667
8 8 55.10000
9 9 54.63636
10 10 54.25000
> data.frame(有名人=1:10*10,期待値=sapply(1:10*10,function(x) fn(100,x)$mean))
有名人 期待値
1 10 54.25000
2 20 52.31818
3 30 51.59375
4 40 51.21429
5 50 50.98077
6 60 50.82258
7 70 50.70833
8 80 50.62195
9 90 50.55435
10 100 50.50000
174:132人目の素数さん
20/03/29 10:58:30.70 1Oo79tY3.net
「有名人」を「wikに載ってる人」と定義し
その数を10000人としてそのうち4人(志村、藤浪、長坂、伊藤隼人)
感染したとしても結果は変わらない
175:132人目の素数さん
20/03/29 10:58:36.48 WogCQeQk.net
昨日の東京のコロナ陽性者は87人検査して63人陽性であったという。
検査の感度0.6 特異度0.9と仮定して、87人中に感染者は何人と推定されるか?
真陽性率=感度=0.6
偽陽性率=1-特異度=0.1
87人中の感染者数をxとすると
陽性者数= 感染者数*真陽性率 + 非感染者数*偽陽性率
63=x*0.6+(87-x)*0.1
これを解くとあり得ない答になる。
176:132人目の素数さん
20/03/29 11:48:31.42 WogCQeQk.net
>>166
総人口n人、有名人m人、そのうち感染者k人とすると
n人中の感染者の期待値は
x = 0 ~ nとして 、xCkはx人からk人選ぶ組み合わせの数を表す
Σ(x*(xCk/nCm))/Σ(xCk/nCm) = =Σ(x*(xCk))/Σ(xCk)
となるのでmの値には依存しない。
n
177:132人目の素数さん
20/03/29 14:27:34.63 2PsxdXJm.net
>>168
するとこの計算で出てくる推定感染者数6万人って値は意味ない感じですか?
178:132人目の素数さん
20/03/29 14:33:09.96 WogCQeQk.net
>>167
陽性者数が87人中63人になるような感度と特異度を最小二乗法で求めると。
> (opt=optim(c(0.6,0.9,63),nazo,method='CG'))
$par
[1] 0.916014625 0.779617519 63.002729987
179:132人目の素数さん
20/03/29 14:48:03.55 0jXKnAa1.net
学術の巨大掲示板群 - アルファ・ラボ
URLリンク(x0000.net)
数学 物理学 化学 生物学 天文学 地理地学
IT 電子 工学 言語学 国語 方言 など
180:132人目の素数さん
20/03/29 15:31:10 WogCQeQk.net
>>170
初期値に依存するから意味のないスクリプトであると判明したので撤回します。
181:132人目の素数さん
20/03/29 15:31:33 WogCQeQk.net
>>169
単なる数字の遊びだろうね。
182:132人目の素数さん
20/03/29 15:37:58 WogCQeQk.net
>>169
前提となっているのが、
日本人1億2680万人いるとして
日本人の感染者数が1人である確率も1億人である確率も同じと、一様分布を仮定しているのが現実離れしている。
よって現実的には意味がない。
183:132人目の素数さん
20/03/31 03:21:38.60 5/cy/U/F.net
URLリンク(youtu.be)
専門家会議がモデルを出したから議論してくれ
184:132人目の素数さん
20/03/31 06:08:43.61 2llZ2I8j.net
>>175
Reed Frost モデルかな?
何を使ったかには言及がなかった。
185:132人目の素数さん
20/03/31 06:12:02.74 2llZ2I8j.net
Reed -Frostはパラメータが1個ですむから推定しやすいんだろう。
186:132人目の素数さん
20/03/31 08:54:47.69 2llZ2I8j.net
>>76
54119人という値になった。
計算プログラムは以下の通り。
# width of 99% confidence interval when 1000 subjects are examined
p2w <- function(
prevalence,
subjects=1000,
sensitivity=0.6,
specificity=0.9,
conf.level=0.99){
# prevalence -> width of 99% confidence interval
n=subjects
p=prevalence*sensitivity+(1-prevalence)*(1-specificity) # positive rate=prev*TP+(1-prev)*FP
q=1-p
2*qnorm(1-(1-conf.level))*sqrt(p*q/n) # width of 99%CI
}
p2w=Vectorize(p2w)
prevalence=seq
187:(0,1,by=0.01) plot(prevalence,p2w(prevalence),bty='l',type='l',lwd=2,ylab='99%CI width', main='subjects:1000\nsensitivity:0.6\nspecificity:0.9') optimize(p2w,c(0,1),maximum=TRUE) # sj2w <- function(subjects){ # subjects -> maximum 99%CI width & its prevalence optimize(function(prev) p2w(prev,subjects),c(0,1),maximum = TRUE) } # at how many subjects 99%ci width equals 0.01 uniroot(function(x,u0=0.01) sj2w(x)$objective-u0,c(1000,100000))
188:132人目の素数さん
20/03/31 09:55:37.96 cpD4Fk2x.net
上って、灘校東大理IIIの超秀才のはずなのに、なんで
あんなに頭の悪い発言ばかりしてんの?
変な宗教にでも取り憑かれて理性が狂わされてるのかな?
189:132人目の素数さん
20/03/31 10:07:35.24 2llZ2I8j.net
日本人1億2680万人からX人を無作為に抽出してPCR検査して、感染者数(≠検査陽性者数)を信頼区間99%誤差±1%で検定したい。
PCR検査は感度0.6,特異度0.9とする。
何人を抽出すれば十分といえるか?
54000人程度になったけど、あってる?
190:132人目の素数さん
20/03/31 14:43:06 2llZ2I8j.net
>>179
超秀才は理Iに行くんじゃないの?
191:132人目の素数さん
20/03/31 14:50:29 ncBHjUEo.net
>>180
感染率の程度、感度・特異度の値の精度の言及無しに出された結論に、ほとんど説得力は無い。
192:132人目の素数さん
20/03/31 15:19:09 2llZ2I8j.net
>>182
感度 beta(13.6991,9.4661)でmode 0.6 sd=0.1
特異 beta(36.172,4.908) でmode 0.9 sd=0.05
でベイズの階層モデルを組んでみるかな。
193:132人目の素数さん
20/03/31 15:45:31.45 2llZ2I8j.net
>>183
そのβ分布を弱情報事前分布に設定して、乱数発生させて計算すると
54000人で99%信頼区間の幅の分布は
> summary(s2w(54000))
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.008144 0.009912 0.009981 0.009927 0.010005 0.010011
となるから、まあ、概ねあっていると思うな。
194:132人目の素数さん
20/03/31 17:50:11.35 ncBHjUEo.net
最も重要なファクターは事前感染率。
事前感染率はいくらに設定したの?
偽陽性が調査対象の10パーセント程含まれる。
医者が怪しいと判断した場合にのみ検査をする場合は、真陽性が調査対象の数十パーセントが期待できる。
このような場合は、真陽性は偽陽性より多数であることが期待でき、検査対象の正確な感染率は把握できるが、
「日本人1億2680万人からX人を無作為に抽出」のような方法だと、感染率0.01%(←現在確認できている感染者の
7倍程度が実際の感染者数に相当)辺りが妥当だと思われるが、この場合、五万人調査して、真陽性5人、偽陽性5000人
のような数字が出てくる。感染率0.02%だったとすると、真陽性10人、偽陽性5000人だ。
中央値のみで判断すると、例えば、5005人の陽性が出ると、0.01%で、5010人の陽性者が出ると0.02% のような
データが出てくる。誤差との見極めは困難。
このような数字から、信頼できる感染率が出せるのか?
195:132人目の素数さん
20/04/01 07:44:43.76 xwYPMdxl.net
>>185
一様分布
196:132人目の素数さん
20/04/01 07:48:29.51 xwYPMdxl.net
確率の分布を考えずにスポットで考える思考のやつとは議論にならんな。
ベイズ階層モデルやったことないの?
197:132人目の素数さん
20/04/01 09:12:32 bZbNlxPT.net
0%~100% までの一様分布のようだな。
つまり、事前確率全く不明だから、1/2教の経典に従い、0.5=50%でやったということ。
医者が検査を行った方がよいと判断した集団でも、なかなか有病率50%はいかない。
そのような結果は、無作為抽出で必要なの調査人数はどれくらいか等という議論では使えない。
全住民を対象にした無作為抽出なら、十
198:万人に一人 以上いる(いた)のは確実だった一方、 百人に一人 という程たくさんはいないだろう と見積もれる。0.001%~1% 辺りで行うべき。 ちょっと考えれば判ることを指摘しているに過ぎない。 調査対象の有病率0.01以下の集団に対し、特異度90%の性能の機器で調査しても、ほとんどがエラー。 せめて 有病率 は、 1-特異度 と同じオーダーか、1-特異度 より大きくないと、扱えない。 特異度99.99%の機器を用意するか、でなければ、有病率を10パーセント程度以上に煮詰めてからやれというお話。
199:132人目の素数さん
20/04/01 09:19:12 deMoC1lt.net
>>188
東京都の行政検査では陽性率が50%を越える日があるぞ。
200:132人目の素数さん
20/04/01 09:26:31 deMoC1lt.net
有病率の事前分布を一様分布として
日々の陽性数は二項分布に従うとして
オリンピック延期決定後の検査を受けた集団での有病率をMCMC出だすと
(感度特異度は既述のβ分布を仮定)
> subjects=c(74,95,87,143,244,330)
> positives=c(17,41,47,40,63,68)
> PCRs3(subjects,positives,iter=10000,warmup=1000)
mean lower upper
0.37288732 0.09822213 0.63719043
201:132人目の素数さん
20/04/01 09:31:13 deMoC1lt.net
>>188
別に有病率を(0,0.1)の一様分布にしても計算できるけど
都の行政検査も陽性率が50%を越える日もあったから一様分布でいいと思うね。行政検査に回った集団の話だけど。
感度・特異度も弱情報事前分布が設定できる。
202:132人目の素数さん
20/04/01 09:33:18 deMoC1lt.net
一変数のポイント確率しか計算できない奴との議論は不毛だね。
203:132人目の素数さん
20/04/01 09:43:53 HHJL1yTu.net
結局なんの疫学データにも基づかない、疫学データで追試することもできない、なんの理論的根拠もない統計仮説下のお話なんて統計学、疫学できないな意味なんかないんだよな。
計算機で遊んでる以上の意味なんかない。
204:132人目の素数さん
20/04/01 09:46:36 bZbNlxPT.net
>>189
だからきちんと「なかなかいかない」と書きました。
>>191
目的が「日本人1億2680万人からX人を無作為に抽出してPCR検査して、感染者数(≠検査陽性者数)
を信頼区間99%誤差±1%で検定したい。 」なのだから、あなたの主張は前提を無視ししている。
205:132人目の素数さん
20/04/01 12:55:21.15 xwYPMdxl.net
>>193
計算機で遊ぶこともできずに電卓で計算して必死で書いていて虚しくない?
CTの診断能を検討した論文。
URLリンク(doi.org)
誰でも鑑別できるのか疑問に思った
このペーパのTable 3に3人の読影医の結果が載っている。
TP FP TN FN sen spc PPV NPV accuracy
1 158 13 192 61 0.72 0.94 0.92 0.76 0.83
2 157 24 181 62 0.72 0.88 0.87 0.74 0.80
3 206 156 49 13 0.94 0.24 0.57 0.79 0.60
陽性尤度比、陰性尤度比、Diagnostic Odd Ratio(陽性尤度比/陰性尤度比)を計算して加えると
TP FP TN FN sen spc PPV NPV acc PLR NLR DOR
1 158 13 192 61 0.72 0.94 0.92 0.76 0.83 11.4 0.30 38
2 157 24 181 62 0.72 0.88 0.87 0.74 0.80 6.1 0.32 19
3 206 156 49 13 0.94 0.24 0.57 0.79 0.60 1.2 0.25 5
PPV,accuracy,DORから読影医3が劣っているようにみえる。
PPVで三者を検定してみる。多重比較になるので一番厳しいBonferri法で補正
Pairwise comparisons using Pairwise comparison of proportions
data: TP out of TP + FP
1 2
2 0.4 -
3 1e-15 2e-11
明らかに3が劣っている。
206:132人目の素数さん
20/04/01 12:55:59.40 xwYPMdxl.net
読影医1,2を加算して計算すると
感度72% [67-76]
特異度91% [88-94]
という結果が得られた。
しかし、現実には何でもコロナと診断する傾向のある読影医3も紛れこむからこういう読影医も加算して計算しないと現実的でないね。
問題
3人を統合したときの感度・特異度とその95%信頼区間を述べよ。
207:132人目の素数さん
20/04/01 12:59:59.63 YULTPcko.net
昔パソコンは習うより慣れろ、理屈なんかわからなくても使ってたらわかるってのがあったけど、まさに正反対の方向にダメダメだな。
学問に対するなんの畏敬の念もない。
208:132人目の素数さん
20/04/01 13:08:31.65 xwYPMdxl.net
>>197
>学問に対するなんの畏敬の念
ひょっとしてアホなの?
209:132人目の素数さん
20/04/01 13:14:52.37 xwYPMdxl.net
Housefield数の計算原理がわからなくても
この画像が新型コロナ肺炎かどうか、診断できる方が有用なんだよな。
URLリンク(pubs.rsna.org)
中心極限定理の証明できなくても、学問への畏敬とかなくても、二項分布を正規分布で近似して計算できる。
210:132人目の素数さん
20/04/01 18:34:14.58 zMY/D89k.net
>>168 他皆様
有名人の感染者が増えてきましたが
市中感染率に影響はないという県警でよろしいのでしょうか?
211:132人目の素数さん
20/04/02 06:13:09 +vJJzaTC.net
>>200
サンプルサイズは期待値の信頼区間幅に影響するけど期待値そのものに影響しないってことでは?
212:132人目の素数さん
20/04/02 09:26:20.10 mzm7EAoV.net
市中感染率が増加の時はもちろんそうだが、一定、あるいは、減少傾向であっても、
経過日数が多くなれば、感染者数は多くなる。
例えば、十日に一人有名人の感染が報告されるというのが継続されていたなら、感染率は一定と
考えられるが、それが、一週間に一人 → 五日に一人 → 三日に一人 → ほぼ毎日 →...
のように、報告されるペースに変化があると、感染率も変化していると考えられる。
213:132人目の素数さん
20/04/02 09:30:08.05 mzm7EAoV.net
補足だが、あまりにも、有名人感染の報告頻度が多くなると、ニュースとしての価値が低くなり、
以前であったら報告されていたであろうケースが報告されなくなるということもあるので、
その辺も考慮して考える必要はある。
214:132人目の素数さん
20/04/03 11:52:46 cch/ocoF.net
横浜市立大学データサイエンス学部佐藤彰洋教授のCOVID-19(新型肺炎)の感染拡大抑止に関する研究・検討資料内容を共有するページ
URLリンク(www.fttsus.jp)
矢原 徹一:九州大学理学研究院教授の試算
URLリンク(jbpress.ismedia.jp)
215:132人目の素数さん
20/04/04 11:37:34 ZFu90Xbq.net
SEIR MODEL
dS(t)/dt = mu*(N-S) - b*S(t)*I(t)/N - nu*S(t)
dE(t)/dt = b*S(t)I(t)/N - (mu+sig)*E(t)
dI(t)/dt = sig*E(t) - (mu+g)*I(t)
dR(t)/dt = g*I(t) - mu*R + nu*S(t)
mu:自然死亡率 b:感染率(S->I)
nu:ワクチン有効率(S->R) sig:発症率(E->I),g:回復率(I->R)
の微分方程式の数値解を使ってシミュレーション
対策しない(外出を控えず、マスクもしない)方が患者や死者は増えるけど早く収束するな。
contact_rate と trannsmission_probabilityを変化させてグラフにしてみた。
URLリンク(i.imgur.com)
216:132人目の素数さん
20/04/04 15:28:24 zerwqPau.net
一次�
217:Y業ごと消滅していいならそうかもな
218:132人目の素数さん
20/04/05 09:54:53.42 fV/kgtmE.net
オリンピック延期決定以後の東京都の行政PCR検査での陽性率をグラフにすると
URLリンク(i.imgur.com)
(陽性数より検査件数の公表は2~3日遅れる)
PCR検査は感度60%、特異度90%くらいなので検査を受けた集団の有病率はもっと多いはず。
感度(最頻値0.6 標準偏差0.1)、特異度(最頻値0.9 標準偏差0.05)のベータ分布に設定、有病率は(0,1)の一様分布でMCMCしたみた。
URLリンク(i.imgur.com)
有病率40%くらいありそうだな。
219:132人目の素数さん
20/04/05 23:57:53 fV/kgtmE.net
新型コロナ肺炎に再感染があるとして流行具合をシミュレーションしてみた。
赤が感染者
上:再感染率0%
中:再感染率1%
下:再感染1%に治癒確率を5倍にする治療薬がある場合
URLリンク(i.imgur.com)
220:132人目の素数さん
20/04/06 00:03:24 xOX4/rO7.net
>>208
準拠したモデルはこれ
SEIRS MODEL
dS(t)/dt = mu*(N-S) - b*S(t)*I(t)/N - nu*S(t) + rho*R(t)
dE(t)/dt = b*S(t)I(t)/N - (mu+sig)*E(t)
dI(t)/dt = sig*E(t) - (mu+g)*I(t)
dR(t)/dt = g*I(t) - mu*R(t) + nu*S(t) - rho*R(t)
mu:自然死亡率 b:感染率(S->I)
nu:ワクチン有効率(S->R) sig:発症率(E->I),g:回復率(I->R)
rho:再感染率(R->S)
Rのスクリプトはここに置いた
スレリンク(hosp板:417番)-420
221:132人目の素数さん
20/04/06 01:00:30.67 taqqH9Ce.net
>>207
疑い例だけに絞って検査してるんだから、有病率が高いのは当たり前。
検査を拡大すれば感染者数は増えるが、有病率はぐっと下がるだろう。
222:132人目の素数さん
20/04/06 05:54:28.61 xOX4/rO7.net
>>210
検査を受けた集団の有病率の事前分布を(0,0.2)に設定して実行すると
URLリンク(i.imgur.com)
MCMCで感度や特異度の事後分布が出せるのが面白い。
223:132人目の素数さん
20/04/06 08:06:29.84 xOX4/rO7.net
>>209
基本的に鎖国しているモデルだから、
入国者や出国者による人口増減とその保菌率・感染率を組み込んだモデルにしないと実態にそぐわないな。
昨今、輸入されているコロナの発症率・回復率などが異なるとすると益々、複雑なモデルになっちゃう。
224:132人目の素数さん
20/04/06 16:39:50 taqqH9Ce.net
死者数が少ないのは注目すべきことだけど、それもいつまで続くか、、、
発症から集中治療室に入るまで2週間弱。死ぬまでがさらに一週間。
ってことは、この2週間での感染者激増が死者数に反映するまで、
あと一週間かかるってこと。
重症者の激増もそろそろ始まることだと思われる。
検査数とか感染者数、重症者数etc.を可視化したデータってここくらいか?
URLリンク(toyokeizai.net)
225:132人目の素数さん
20/04/06 17:50:09 taqqH9Ce.net
韓国は46万件の検査で感染者が1万人。
感染かどうかをどう確定したのか知らないけど、有病率が2%だと
すれば、特異度がよほど高くないと、陽性者のかなりの割合が
偽陽性ってことになりそう。
226:132人目の素数さん
20/04/06 22:07:12.60 Sgr4CLg/.net
>>207
東京都は陽性者数は公表しても検査人数を迅速に公表しないのでRのパッケージAmeliaを使って多重代入法による欠測データ処理してみる。
#
subjects=c(74,95,87,143,244,330,41,145,164,469,NA,NA,NA)
positives=c(17,41,47,40,63,68,13,78,66,97,89,117,143)
(dataset=data.frame(subjects,positives))
# 検査実施人数を欠測データとしてAmeliaで推定
library(Amelia)
na.idx=which(is.na(subjects))
M=1000
set.seed(1234) ; a.out=amelia(dataset,m=M)
imp=a.out$imputations
NAsubject=NULL
for(i in 1:M){
NAsubject=rbind(NAsubject,imp[[i]][na.idx,1])
}
missing.data=round(apply(NAsubject,2,mean))
subjects[na.idx]=missing.data
227: PCRs3(subjects,positives,iter=10000,warmup=1000,verbose=TRUE) 行政検査例の有病率と95%信頼区間 mean lower upper 0.3745734 0.1048757 0.6575508 https://i.imgur.com/DxVNXad.png
228:132人目の素数さん
20/04/07 07:24:45.67 TwUCHcsI.net
4月5日の陽性者数が検査人数より多いのは単純ミス?
URLリンク(stopcovid19.metro.tokyo.lg.jp)
Ameliaでの欠測推測を確認しようと思ったのにもとのデータがあてにならない。
229:132人目の素数さん
20/04/07 07:48:06 TwUCHcsI.net
晋型コロナ肺炎に感度0.9,特異度0.9の迅速検査が開発されたと仮定する。
日本人1億2595万人からX人を無作為抽出して有病率を推定したい。
有病率の99%信頼区間幅を1%以内で検定したい。
何人を抽出すれば十分といえるか?
230:132人目の素数さん
20/04/07 08:24:28 TwUCHcsI.net
>>76
感度0.6 特異度0.9として1000人検査したときはの有病率の99%信頼区間幅は何%以内におさまるか?
231:132人目の素数さん
20/04/07 10:09:26.39 /UNISUK0.net
>>216
検査結果の集計日と、対応する検査数の集計日が揃ってないからかもね。
データがあてにならないのはその通りかも。
232:132人目の素数さん
20/04/07 20:47:10.10 /iaYsYaQ.net
URLリンク(imgur.com)
スプレッドシートのgrowthで予測してみると、9月には人類滅亡しそうです。
元データは
URLリンク(www.worldometers.info)
233:132人目の素数さん
20/04/07 21:10:10.27 TwUCHcsI.net
ソース不明な数字だけど
"名前:名無しさん@1周年[sage] 投稿日:2020/04/07(火) 20:29:49.57 ID:GcsEpAfR0 (PC)
他国の総感染者数が今の日本と同程度の時の死者数
イタリア, 感染3089人 死者107人
スペイン, 感染2965人 死者 84人
アメリカ, 感染3499人 死者 64人
イラン , 感染3513人 死者107人
フランス, 感染3661人 死者 79人
日本 , 感染3123人 死者 77人
"
country=c('イタリア','スペイン','アメリカ','イラン','フランス','日本')
infected=c(3089,2965,3499,3513,3661,3123)
dead=c(107,84,64,107,79,77)
df=data.frame(country,dead,infected)
df
re.lm=lm(dead~infected,data=df)
summary(re.lm)
chisq.test(dead,infected)
names(infected)=country
names(dead)=country
pairwise.prop.test(dead,infected)
> chisq.test(dead,infected)
Pearson's Chi-squared test
data: dead and infected
X-squared = 24, df = 20, p-value = 0.2424
Warning message:
In chisq.test(dead, infected) : Chi-squared approximation may be incorrect
> pairwise.prop.test(dead,infected)
Pairwise comparisons using Pairwise comparison of proportions
data: dead out of infected
イタリア スペイン アメリカ イラン フランス
スペイン 1.00000 - - - -
アメリカ 0.00066 0.11005 - - -
イラン 1.00000 1.00000 0.01762 - -
フランス 0.01843 0.79242 1.00000 0.24130 -
日本 0.24724 1.00000 0.79242 1.00000 1.00000
P value adjustment method: holm
234:132人目の素数さん
20/04/08 00:23:04 KlmY0AUa.net
アメリカの致死率高くない?
現時点での死者数と感染者数からだと3%くらいなんだが、これは
見かけ上の数値で、実際には発症から死亡まで2週間以上かかかる
から、10日くらい前の感染者数を分母にしないと見誤る。
で、このタイムラグを考慮にいれると少なくとも10%くらいの致死率
になりそ。日本の2,3倍になる計算だ。
235:132人目の素数さん
20/04/08 00:27:12 KlmY0AUa.net
>>221
その数字もタイムラグを考慮に入れないと致死率を見誤るよ。
感染者の増大率が大きいところでは、死者数を過小評価する。
たとえば2週間後の死者数で比較すべき。
236:132人目の素数さん
20/04/08 19:19:57.90 gX6rmgSD.net
>>223
感染から死亡までの期間は一定でないからどうだろね。
237:132人目の素数さん
20/04/08 19:24:23.22 gX6rmgSD.net
SEIRモデルに再感染(免疫を失ってR からSへの変遷)があるとしてグラフ化してみた。
1万人に1人の感染者と9人の保菌者がいるとして、
再感染率0%のとき
URLリンク(i.imgur.com)
再感染率0.1%(1000人に1人の割合で免疫を失う)場合
URLリンク(i.imgur.com)
再感染率0.1%で外出自粛等で接触が8割減になった場合
URLリンク(i.imgur.com)
238:132人目の素数さん
20/04/09 01:14:30 m5ECJkjO.net
>>224
まあ、発症してすぐ感染が発覚するわけでもないから、
症状がある程度続いてから陽性判定されるとすれば、
10日前後のディレイかな。しかも、どのくらいの分散
かわかんないけど、正規分布で畳み込んだ感じで。
239:132人目の素数さん
20/04/09 02:02:32 m5ECJkjO.net
英国もフランスも致死率高いな。
のきなみ10%越え。
240:132人目の素数さん
20/04/09 02:09:10 NyLklR2N.net
英国とフランスは医療崩壊してるの?
241:132人目の素数さん
20/04/09 09:37:43.05 m5ECJkjO.net
しらんがな
242:132人目の素数さん
20/04/09 11:57:57 cpUiyQfM.net
これこそ、国を守るってことだな。
キューバ政府は2020年3月24日(火)以降、外国人観光客の入国を認めないと発表しました。
キューバに入国できるのは居住者のみで、これには出張で短期または長期滞在する外国人も含まれます。
加えて、キューバに入国する人は全員14日間隔離されます。
商業活動は維持されますが、商船や航空機の乗務員の移動は制限されます。キューバの領空は封鎖されません。
今回の制限措置は30日間有効で、状況に応じて延長される可能性もあります。
243:132人目の素数さん
20/04/09 12:26:32 cpUiyQfM.net
>>216
オリンピック延期決定以降のデータで検討
陽性人数が検査人数を上回るデータは検査人数を欠測データ(NA)として計算する。
dates subjects positives
1 2020-03-24 74 17
2 2020-03-25 95 41
3 2020-03-26 87 47
4 2020-03-27 143 40
5 2020-03-28 244 63
6 2020-03-29 330 68
7 2020-03-30 41 13
8 2020-03-31 145 78
9 2020-04-01 164 66
10 2020-04-02 469 97
11 2020-04-03 551 89
12 2020-04-04 NA 117
13 2020-04-05 NA 143
14 2020-04-06 356 83
15 2020-04-07 271 79
16 2020-04-08 NA 144
その設定で
# LittleのMCAR検定 データが MCAR または MAR であるという帰無仮説のもとで検定
library(BaylorEdPsych)
LittleMCAR(dataset)$p.value
LittleMCAR(dataset)$p.value
this could take a while[1] 0.004907289
有意差がでて出鼻をくじかれた
244:132人目の素数さん
20/04/09 14:54:07 pa+n/S2z.net
URLリンク(twitter.com)
(deleted an unsolicited ad)
245:132人目の素数さん
20/04/10 00:25:10 EjAtBG4Z.net
日本もすでに欧州同様指数関数ラインにのってる模様、
URLリンク(i.imgur.com)
欧州並になるのは文字通り時間の問題かと
日本のコロナ感染発覚者数
3月24日 1193
3月25日 1307
3月26日 1387
3月27日 1499
3月28日 1693
3月29日 1866
3月30日 1953
3月31日 2178
4月 1日 2384
4月 2日 2617
4月 3日 2935
4月 4日 3139
4月 5日 3654
4月 6日 3906
4月 7日 4257
4月 8日 4667
オリンピック延期決定
246:した3月24日をDay1として線形回帰すると https://i.imgur.com/Q45Xo8V.png y=c(1193,1307,1387,1499,1693,1866,1953,2178,2384,2617,2935,3139,3654,3906,4257,4667) z=as.Date("2020-03-24")+ 0:(length(y)-1) plot(z,y,bty='l',pch=19) x=1:length(y) plot(x,y,bty='l') (re.lm=lm(log(y)~x)) b0=re.lm$coef[1] b1=re.lm$coef[2] names(b0)='' f <- function(x) exp(b0)*exp(b1*x) # 1059.063*exp(0.09231543*x) curve(f(x),add=T) n=1:50 plot(n,f(n),bty='l') f1 <- function(y) print(log(y*exp(-b0))/b1 + as.Date("2020-03-24"),quote=F) f1(c(5000,1e4,5e4,1e5,5e5,1e6,1e7,1e8,1.2595e8)) 感染者が5千、1万、5万、10万、50万、百万、1千万、1億、1億2595万人に達する日は > f1(c(5000,1e4,5e4,1e5,5e5,1e6,1e7,1e8,1.2595e8)) [1] 2020-04-09 2020-04-17 2020-05-04 2020-05-12 2020-05-29 [6] 2020-06-06 2020-07-01 2020-07-26 2020-07-28
247:132人目の素数さん
20/04/10 06:57:42.69 EjAtBG4Z.net
>>233
1日ずれていた。
> f1(c(5000,1e4,5e4,1e5,5e5,1e6,1e7,1e8,1.2595e8,77e8))
[1] 2020-04-10 2020-04-17 2020-05-04 2020-05-12 2020-05-29
[6] 2020-06-06 2020-07-01 2020-07-26 2020-07-28 2020-09-10
世界の人口77億人にたっするのは9月10日となった。
248:132人目の素数さん
20/04/10 07:03:11.67 EjAtBG4Z.net
>>216
行政検査人数と医療機関の行った検査も含めて、陽性者数を公表しているだけみたいだな。
これじゃぁ、陽性率も有病率も何の解析もできない。
249:132人目の素数さん
20/04/10 07:18:35.21 EjAtBG4Z.net
URLリンク(toyokeizai.net)
のデータ使って
全国と東京の線形回帰の係数を比べると
全国
Call:
lm(formula = log(y) ~ x)
Coefficients:
(Intercept) x
6.84962 0.09253
東京
Coefficients:
(Intercept) x
5.1707 0.1332
東京の方が感染者増加速度が速いね。
東京都の人口13,951,636(令和2年1月1日現在) に達する日は
> f1(13951636)
2020-06-15
となった。
250:132人目の素数さん
20/04/10 07:23:05.91 EjAtBG4Z.net
今日の東京の感染者予想人数は
> d2i("2020-4-10")
daily total
241.1515 1934.2759
251:132人目の素数さん
20/04/10 11:19:37 oS02oOtl.net
URLリンク(pbs.twimg.com)
このグラフを見ろ
252:132人目の素数さん
20/04/10 11:22:37 wWOTJ0FB.net
>>238
だから何なの?
253:132人目の素数さん
20/04/10 11:47:04.58 EjAtBG4Z.net
感染者数の増加は直線より指数関数の方がフィットするな。
URLリンク(i.imgur.com)
254:132人目の素数さん
20/04/10 12:02:00 K2O0fmDQ.net
>>238
エクセルのグラフって、マークと曲線がちょっとずれるよね。
ってか、凡例までそうなってて、すっごく気持ち悪い。
なんでだろ?
255:132人目の素数さん
20/04/10 12:04:52 K2O0fmDQ.net
有名人の感染がボチボチ出てくるってことは、市中感染率も
それなりに高いんだろうな。0.1%は超えてそう。
256:132人目の素数さん
20/04/10 14:27:39 XjWf1V4z.net
>>242
120万人罹患してるってこと?
発症者が5000人いってないのに?
発症率0.4%しかない雑魚ウィルスならサイトカインストームなんて起こさないでしょ?
これウィルスが強すぎて免疫系が異常に発動して起こるもんらしいから。
257:132人目の素数さん
20/04/10 14:29:22 XjWf1V4z.net
あ、0.1%か。
12万人で発症率4%。
どのみちないだろな。
258:132人目の素数さん
20/04/10 14:58:22 K2O0fmDQ.net
>>244
0.1%ってのはあくまで都内や大阪を想定した推測ね。
259:132人目の素数さん
20/04/10 15:52:37.05 Qxo2/Tvi.net
>>216
これひどいよなぁ
陽性率隠しの一環にしか見えない
そもそも陽性者数と検査数って普通に検査してりゃ一対一で対応するだろ
各検査機関がそれらを普通に出せば足し合わせるだけでいい
260:132人目の素数さん
20/04/11 01:01:52 g/e7gi/j.net
そういうミスは統計資料にはつきもの。それでもデータさえ
多ければ大勢に影響なし。
261:132人目の素数さん
20/04/11 06:07:34 Dikmq831.net
>>8
優秀ですね
262:132人目の素数さん
20/04/11 13:47:45 mKhDXvLi.net
東京都のこのデータを使って
URLリンク(stopcovid19.metro.tokyo.lg.jp)
オリンピック延期決定前後で陽性者の増え方に違いがあるかをグラフにしてみた。
URLリンク(i.imgur.com)
黒塗り丸が延期決定前で黒線がその回帰曲線1日毎に1.078594 倍に増加
赤丸が延期決定後で赤線がその回帰曲線 1日毎に1.136801 倍に増加
青は全体での回帰曲線 1日毎に1.088687 倍に増加
263:132人目の素数さん
20/04/11 15:33:13 mKhDXvLi.net
>>247
集計ミスじゃなくて、集計デザインそのものに問題があるよ。
陽性数は東京都全部の集計なのに、検査数は行政検査だけで指定医療機関での検査を除外している集計している。
これじゃ陽性率から有病率を推測するのも不可能。
264:132人目の素数さん
20/04/11 17:18:33.74 HlGO5aHy.net
風邪ひいたみたい。感染予防には抜かりはないのに、
265:132人目の素数さん
20/04/11 21:40:53.17 vjB1aI8B.net
ランダムに検査して、感染者数を推定してほしいんだけな
266:132人目の素数さん
20/04/12 09:01:25.70 uhjyERyF.net
オーストリアがやったらしいよ。
それによりますと、全体のおよそ0.3%が陽性と判定され、分析の結果、この時期に少なくとも1万200人、最大で6万7400人が感染していた可能性があることがわかったということです。
オーストリアの人口はおよそ890万人で、政府は最も可能性が高い数字として、この時期におよそ2万8500人が感染していたと推定しています。
267:132人目の素数さん
20/04/12 09:56:34.10 uhjyERyF.net
>>253
コピペしそこなったけど、検査は1533人に施したそうだ。
感度100%、特異度100%の完璧な検査だったと仮定して、サンプルサイズ1500と
サンプル比率0.3%で、信頼度95%で推定すると母集団での比率は0.12%~0.7%
くらいだと推定できるから、オーストリア政府発表と一致するね。そういう
単純な計算なのかねぇ?
268:132人目の素数さん
20/04/12 10:02:10.44 uhjyERyF.net
もひとつコピペし忘れた、判明している感染者数は12,200人ね。
特異度が100%って考えにくいから、まあ、結構な捕捉率じゃねーの?
>当時医療機関などで確認されていた感染者はおよそ1万2200人だったことから、
>政府は、実際に感染していた人は公表の2倍以上に上るとみられるとしています。
269:132人目の素数さん
20/04/12 12:34:54.53 QLYmTD4y.net
>>254
1533*0.3/100=4.6だけど陽性は5人?
270:132人目の素数さん
20/04/12 12:40:40.53 +MLUgV1E.net
>>254
95%信頼区間も計算方法がいくつかあるよね。
> binom::binom.confint(5,1533)
method x n mean lower upper
1 agresti-coull 5 1533 0.003261579 0.0011557601 0.007850672
2 asymptotic 5 1533 0.003261579 0.0004073974 0.006115760
3 bayes 5 1533 0.003585398 0.0009655853 0.006611832
4 cloglog 5 1533 0.003261579 0.0012670654 0.007342721
5 exact 5 1533 0.003261579 0.0010598477 0.007594894
6 logit 5 1533 0.003261579 0.0013582045 0.007811457
7 probit 5 1533 0.003261579 0.0013079622 0.007530462
8 profile 5 1533 0.003261579 0.0011715789 0.006997012
9 lrt 5 1533 0.003261579 0.0011411357 0.006992346
10 prop.test
271: 5 1533 0.003261579 0.0012013872 0.008065116 11 wilson 5 1533 0.003261579 0.0013939315 0.007612501
272:132人目の素数さん
20/04/12 13:03:50.97 +MLUgV1E.net
有病率:一様分布
感度:最頻値0.6標準偏差0.1のベータ分布
特異度:最頻値0.9標準偏差0.05のベータ分布
として検査陽性数は有病率*感度+(1-有病率)*(1-特異度)の確率に従う二項分布
というモデルでMCMCすると
URLリンク(i.imgur.com)
有病率
mean lower upper
0.002678327 0.000029215 0.005582318
感染者数
> c(0.002678327, 0.000029215, 0.005582318)*8900000
[1] 23837.11 260.01 49682.63
という結果になった。
検査陽性数が少ないから信頼区間幅が広くなるなぁ。
エントリーに5以下があると正規分布での近似が悪くなると習ったような記憶がある。
273:132人目の素数さん
20/04/12 17:10:06 uhjyERyF.net
>>256
それは書いてないけど、そうなるよね。陽性は4人か5人だろうね。
274:132人目の素数さん
20/04/12 17:23:55 uhjyERyF.net
>>258
意外と特異度が高いんだよねぇ。99%をゆうにこえるとは...。
ってことは、検査数を増やしても偽陽性に惑わされる心配はあまりないってことか。
この調査から言えることは、検査の特異度がかなり高いってことと、有病率の上限
が1%には届かないと見積もれることかな?
275:132人目の素数さん
20/04/12 20:15:06 NHP4603i.net
最終的に陽性と判断されたのは、1500件中、4,5件だと思われるが、
これは、PCR検査機器の性能チェックじゃ無く、有病率の調査なんだから、
陽性とでた場合には、PCR機器をもう一度回したかもしれないし、
その人を呼んで、再度検体を採取してチェックしたのではないか?
具体的には、一次検査では、1500件中150件ほど陽性が出る。
この150人に対し、再調査を指示して、20件ほど陽性が出る。
この20人に対し、三次検査、四次検査を行って、最終的に4件程が陽性だと判断したと。
検査機器自体の特異度が高いのではなく、運用による効果だと考える方が無難では?
276:132人目の素数さん
20/04/12 20:37:37 t0enBOyQ.net
ドメイン知識ゼロ感
277:132人目の素数さん
20/04/12 21:48:49.76 +MLUgV1E.net
>>260
特異度の事前分布を最頻値0.9標準偏差0.05にした影響かと思って
特異度の事前分布を最頻値0.5 標準偏差0.2にしてMCMCしても同様の結果だった。(感度は最頻値0.6標準偏差0.1で同じ設定)
PCRs2(1533,5,SEN=0.6,SPC=0.5,SD1=0.1,SD2=0.2)$fit
Inference for Stan model: model3.
4 chains, each with iter=2000; warmup=1000; thin=1;
post-warmup draws per chain=1000, total post-warmup draws=4000.
mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
prev 0.00273 0.00005 0.00249 0.00010 0.00200 0.00912 2406 1.0011
sen 0.57489 0.00205 0.10047 0.37477 0.57377 0.76629 2403 1.0009
spc 0.99592 0.00004 0.00182 0.99181 0.99613 0.99882 2203 1.0007
p 0.00559 0.00003 0.00185 0.00258 0.00539 0.00975 3191 1.0000
lp__ -72.11408 0.03392 1.22437 -75.33316 -71.80596 -70.68418 1303 1.0048
prev:prevalence 有病率
sen : sensitivity 感度
spc : specificity 特異度
278:132人目の素数さん
20/04/12 23:37:41 NHP4603i.net
以前、>>188 で
>>せめて 有病率 は、 1-特異度 と同じオーダーか、1-特異度 より大きくないと、扱えない。
と指摘したことがあります。>>260や、>>263の結果は、将にこれを実証したといえそうですね。
主語を、特異度に変えると、
「特異度は 1-有病率 と同じオーダーか、1-有病率 より大きくないと、扱えない。」となります。
279:132人目の素数さん
20/04/12 23:39:41 NHP4603i.net
上の 260へのリンクは、>>258の間違いです。
280:132人目の素数さん
20/04/12 23:59:37 NHP4603i.net
>>264 さらに訂正
誤:「特異度は 1-有病率 と同じオーダーか、1-有病率 より大きくないと、扱えない。」
正:「1-特異度は 有病率 と同じオーダーか、有病率 より小さくないと、扱えない。」
281:132人目の素数さん
20/04/13 05:55:48.12 qYuv8ZQx.net
1-特異度=偽陽性率
282:132人目の素数さん
20/04/13 05:58:17.03 qYuv8ZQx.net
HIV検査などは有病率<偽陽性率だが、だからといって扱えないわけではない。
283:132人目の素数さん
20/04/13 06:14:55.51 qYuv8ZQx.net
陽性(陰性)的中率や精度は有病率に左右されるけど、陽性(陰性)尤度比や診断的オッズ比は有病率には左右されない。
ただ、これだけの話。
284:132人目の素数さん
20/04/13 06:56:13.98 qYuv8ZQx.net
臨床の世界では、陽性尤度比が10以上、陰性尤度比が0.1以下なら有用な検査とされる。
5以上、0.5以下なら、参考にはなるという程度だな。
Evidence-Based Physical Diagnosis, 4e だと3以上、0.3以下としている。
285:132人目の素数さん
20/04/13 09:47:52 5wetDrgQ.net
>>268
>>188などの書き込みは、>>180の問いに対する流れで書かれたものであって、
そこで書かれた「扱えない」は、「有病率の調査として誤差が大きすぎて使えない」という意味。
「検査自体が無意味」等というような意味で使っていたかのような誤解を与えかねない268の書き込みは、
話のすり替えであり、不適当であることを指摘しておく。
286:132人目の素数さん
20/04/13 10:35:51.38 LCm4R7zM.net
>>269
思うに、臨床の場では陽性が疑わしい患者をみるから、
「市中での」有病率はあまり関係ないってことなんじゃ
ないかな?
陽性かどうか五分五分という状況で検査すると思えば、
事前確率が50%という仮定(有病率50%に対応)になって、
陽性判定が真の陽性である確率と偽の陽性である確率との
比は感度/(1-特異度)になる。これが陽性尤度比だよね?
287:132人目の素数さん
20/04/13 15:30:20 SYLNYopk.net
感度0.6 特異度0.9の検査で1533人検査して3例陽性であった、有病率を計算せよ
という問題をやってみるといいね。
解答不能になる。
288:132人目の素数さん
20/04/13 15:39:19 SYLNYopk.net
1533例中5例陽性のとき、弱情報事前分布を以下のように
有病率:一様分布
感度:最頻値0.6標準偏差0.1のベータ分布
特異度:最頻値0.5標準偏差0.2のベータ分布
として検査陽性数は有病率*感度+(1-有病率)*(1-特異度)の確率に従う二項分布
としてJAGSでMCMCしてみた。的中率・精度・尤度比・診断的オッズ比の分布も出してみた。
URLリンク(i.imgur.com)
URLリンク(i.imgur.com)
289:132人目の素数さん
20/04/13 17:12:33.46 SYLNYopk.net
>>273
>253の数字だと陽性数は5人なので1533人検査して5人陽性とする。
有病者の人数をxとすると 0.6x+(1-0.9)(1533-x)=5を解くとありえない数値になるね。
感度0.6特異度0.9という前提が間違っていると考えると、どんな値が尤もらしい値なのかを探求したくなる。
MCMCはその手段を与えてくれる。