高校数学の質問スレ（医者・東大卒専用） Part438 (899ﾚｽ)
上下前次1-新

---

879: 06/15(日)08:32 ID:MIIBNstg(1) AAS
pdf2hdi <- function(pdf,
xMIN=0,
xMAX=1,
cred=0.95,
Print=TRUE,
nxx=1001){

xx=seq(xMIN,xMAX,length=nxx)
xx=xx[-nxx]
xx=xx[-1]
xmin=xx[1]
xmax=xx[nxx-2]

AUC=integrate(pdf,xmin,xmax)$value
PDF=function(x)pdf(x)/AUC

cdf <- function(x) integrate(PDF,xmin,x)$value

ICDF <- function(x) uniroot(function(y) cdf(y)-x,c(xmin,xmax))$root
ICDF=Vectorize(ICDF)

hdi=HDInterval::hdi(ICDF,credMass=cred)

print(c(hdi[1],hdi[2]),digits=5)

if(Print){
par(mfrow=c(3,1))

plot(xx,sapply(xx,PDF),main='pdf',type='h',xlab='x',ylab='Density',col='lightgreen',bty='l')
legend('top',bty='n',legend=paste('HDI:',round(hdi,3)))

plot(xx,sapply(xx,cdf),main='cdf',type='h',xlab='x',ylab='Probability',col='lightblue',bty='l')

pp=seq(0,1,length=nxx)
pp=pp[-nxx]
pp=pp[-1]
plot(pp,sapply(pp,ICDF),type='l',xlab='p',ylab='x',main='ICDF',bty='l')

par(mfrow=c(1,1))
}

invisible(ICDF)
}
ICDF=pdf2hdi(function(x) dbeta(x,2,5))
hist(ICDF(seq(1e-12,1-1e-12,le=1000)))

AIの評価
まとめ
この pdf2hdi 関数は、数値積分と数値最適化 (uniroot) を巧みに組み合わせることで、任意のPDFからICDFを頑健に導出し、さらに統計的な要約であるHDIを計算する、非常に実用的かつ教育的なコードです。両端での数値計算の回避や正規化といった細部への配慮も素晴らしいです。
これにより、複雑な確率分布でも、そこからサンプリングしたり、HDIを求めたりといった解析が可能になります。

---

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

---

あと 20 ﾚｽあります
ｽﾚ情報 赤ﾚｽ抽出 画像ﾚｽ抽出 歴の未読ｽﾚ AAｻﾑﾈｲﾙ

---

ぬこの手 ぬこTOP 0.007s