高校数学の質問スレ（医者・東大卒専用） Part438

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

879: 2025/06/15(日) 08:32:56.90 ID:MIIBNstg(1)調 AA×

[240|320|480|600|100%|JPG|べ|ﾚｽ栞|ﾚｽ消]

879:  [sage] 2025/06/15(日) 08:32:56.90 ID:MIIBNstg

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を求めたりといった解析が可能になります。

http://rio2016.5ch.net/test/read.cgi/math/1723152147/879

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

あと 3 ﾚｽあります
ｽﾚ情報赤ﾚｽ抽出画像ﾚｽ抽出歴の未読ｽﾚ

ぬこの手ぬこTOP 0.044s