分からない問題はここに書いてね458

[過去ﾛｸﾞ] 分からない問題はここに書いてね458 (1002ﾚｽ)
上下前次1-新
抽出解除必死ﾁｪｯｶｰ(本家) (べ) 自ID ﾚｽ栞あぼーん

このｽﾚｯﾄﾞは過去ﾛｸﾞ倉庫に格納されています｡
次ｽﾚ検索歴削→次ｽﾚ栞削→次ｽﾚ過去ﾛｸﾞﾒﾆｭｰ

ﾘﾛｰﾄﾞ規制です｡10分ほどで解除するので､他のﾌﾞﾗｳｻﾞへ避難してください｡

844(1): １３２人目の素数さん [sage] 2020/03/22(日) 07:50:24.15 ID:liILqu/N(1/6) AAS
コロナの件で疑問に思ったのだけど

無作為に１０００人を抽出してPCR検査を行ったら１０人が陽性であった。
PCR検査の感度０．７、特異度０．９として有病率の期待値と９５％信頼区間は？

有病率の事前確率分布を一様分布と仮定して検査陽性後の有病率をベータ分布で出す。
その有病率と感度・特異度をから陽性的中率を計算。
その陽性的中率を使って有病率を逆算というアルゴリズムで計算したら

> sim()
mean mode lower upper
0.010969194 0.009870096 0.005034295 0.017532847

このアルゴリズムってあってる？

861: １３２人目の素数さん [sage] 2020/03/22(日) 14:02:38.08 ID:liILqu/N(2/6) AAS
>>844(自己レス)
誤謬に気がついたので撤回しますm(_ _)m

872: １３２人目の素数さん [sage] 2020/03/22(日) 23:00:24.46 ID:liILqu/N(3/6) AAS
>>871
n=10のときで総当たりで計算してみた

> rm(list=ls())
> n=10
> dec2nw <- function(num, N, digit = n){
+ r=num%%N
+ q=num%/%N
+ while(q > 0 | digit > 1){
+ r=append(q%%N,r)
+ q=q%/%N
+ digit=digit-1
+ }
+ return(r)
+ }
> d=t(sapply(0:(2^n-1),function(num) dec2nw(num,2,n))) # 10枚のコインの裏表の順列（０を取り除く表とする）
> head(d) ; tail(d)
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,] 0 0 0 0 0 0 0 0 0 0
[2,] 0 0 0 0 0 0 0 0 0 1
[3,] 0 0 0 0 0 0 0 0 1 0
[4,] 0 0 0 0 0 0 0 0 1 1
[5,] 0 0 0 0 0 0 0 1 0 0
[6,] 0 0 0 0 0 0 0 1 0 1
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1019,] 1 1 1 1 1 1 1 0 1 0
[1020,] 1 1 1 1 1 1 1 0 1 1
[1021,] 1 1 1 1 1 1 1 1 0 0
[1022,] 1 1 1 1 1 1 1 1 0 1
[1023,] 1 1 1 1 1 1 1 1 1 0
[1024,] 1 1 1 1 1 1 1 1 1 1
> d1=apply(d,1,sum) #　各行ごとの合計＝裏がでた枚数
> sum((1/2)^d1)/2^n　 # 各行ごとに　(1/2) ^ 裏の枚数を計算して合算する　
[1] 0.056314
> (1-1/4)^n
[1] 0.056314

合致した。

873: １３２人目の素数さん [sage] 2020/03/22(日) 23:20:34.47 ID:liILqu/N(4/6) AAS
n=20での100万回のシミュレーション

> # simulation
> sim <- function(n){
+ (flip1=sample(0:1,n,rep=T)) # 0:1から重複を許してn個選んだ配列をflip1とする
+ (flip2=sample(0:1,sum(flip1),rep=T))　# 0:1からflip1の総和個選んだ配列をflip2とする
+ sum(flip2)==0 # flip2の総和が０か否かを返す
+ }
> mean(replicate(1e7,sim(20))) # 100万回試行して全部取り去った割合を計算
[1] 0.0031899
> (1-1/4)^20
[1] 0.0031712

まあ、近似した。

874: １３２人目の素数さん [sage] 2020/03/22(日) 23:24:01.87 ID:liILqu/N(5/6) AAS
>>867
1枚の硬貨が取り除かれない確率は、裏裏と続くときだから1/4
その余事象（1枚の硬貨が取り除かれる事象）の確率は3/4
どの硬貨の裏表がでるかは独立事象だから、(3/4)^nでいいんじゃないの？

876(1): １３２人目の素数さん [sage] 2020/03/22(日) 23:54:02.20 ID:liILqu/N(6/6) AAS
>>828
シミュレーションプログラムの練習　数理解は賢者にお任せ

[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
En 1.0000 1.1991 1.4378 1.6744 1.9163 2.1544 2.3992 2.6329 2.9352 3.1542 3.4106
Pn 0.3405 0.1997 0.1483 0.1148 0.0917 0.0775 0.0700 0.0608 0.0516 0.0458 0.0416
[,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22]
En 3.6105 3.8742 4.1257 4.3780 4.6665 4.7706 5.0866 5.4512 5.6246 5.8261 6.1474
Pn 0.0375 0.0381 0.0383 0.0308 0.0295 0.0296 0.0269 0.0223 0.0227 0.0221 0.0239
[,23] [,24] [,25] [,26] [,27] [,28] [,29] [,30]
En 6.5768 6.5863 6.8209 7.0221 7.3743 7.652 7.8525 8.1056
Pn 0.0217 0.0204 0.0189 0.0202 0.0183 0.017 0.0166 0.0161

fn <- function(n){
B=c(rep(1,n),0,0) # 1:青玉 0:白玉
flg=3 # drawを初期値
i=0 # 試行の回数カウンター
while(flg==3){ # drawなら繰り返す　1:win 2:lose 3:draw
i=i+1
flg=(1:3)[sum(sample(B,2))+1] # (1:3)[sum(c(0,0))+1] : win
}
c(i=i,flg=flg)
}
sim <- function(n){
k=1e4
re=t(replicate(k,fn(n)))
c(mean(re[,'i']),mean(re[,'flg']==1)) # 回数と勝率を返す
}
n=1:30
re=sapply(n,sim)
En=re[1,]
plot(n,En,bty='l',pch=19)
Pn=re[2,]
plot(n,Pn,bty='l',pch=19)
rownames(re)=c('En','Pn')
re

上下前次1-新書関写板覧索設栞歴

ｽﾚ情報赤ﾚｽ抽出画像ﾚｽ抽出歴の未読ｽﾚ AAｻﾑﾈｲﾙ

ぬこの手ぬこTOP 0.037s