【R言語】統計解析フリーソフトＲ第6章【GNU R】

【R言語】統計解析フリーソフトＲ第6章【GNU R】at MATH

【R言語】統計解析フリーソフトＲ第6章【GNU R】 - 暇つぶし2ch597:信頼区間を求めよ" decken <- function(M, Mmax, N, conf.level=0.95){ if(Mmax < M) return(0) n=Mmax:N pmf=choose(Mmax-1,M-1)/choose(n,M) pdf=pmf/sum (pmf) mean=sum(n*pdf) upr=n[which(cumsum(pdf)>conf.level)[1]] lwr=Mmax c(lower=lwr,mean=mean,upper=upr) } decken(M=5,Mmax=60,N=100) > decken(M=5,Mmax=60,N=100) lower mean upper 60.0000 71.4885 93.0000 これをシミュレーションで確認したい。 # simulation M=5 ; Mmax=60 ; N=100 sub <- function(M,Mmax,N){ n=sample(Mmax:N,1) # n : 参加人数n m.max=max(sample(n,M)) # m.max : n人からM人選んだ最大番号 if(m.max==Mmax) return(n) } sim <- function(){ n=sub(M,Mmax,N) while(is.null(n)){ n=sub(M,Mmax,N) # 最大番号が一致するまで繰り返す } return(n) } runner=replicate(1e4,sim()) summary(runner) ; hist(runner,freq=F,col="lightblue") quantile(runner,prob=c(0,0.95)) cat(names(table(runner)[which.max(table(runner))])) > summary(runner) ; hist(runner,freq=F,col="lightblue") Min. 1st Qu. Median Mean 3rd Qu. Max. 60.00 63.00 68.00 71.43 77.00 100.00 > quantile(runner,prob=c(0,0.95)) 0% 95% 60 93 > cat(names(table(runner)[which.max(table(runner))])) # 最頻値 60 結果は確認できたけど、もっと高速なシミュレーションアルゴリズムはあるだろうか？

次ページ

続きを表示

1を表示