(R)谁能给我一个最高密度区间的算法?

bzzcjhmw  于 2023-09-27  发布在  其他
关注(0)|答案(1)|浏览(101)

在研究贝叶斯统计时,引入了一种称为HPD区间的区间估计,我试图找到beta分布95%的HPD区间。
我尝试了几种方法,但找到一个算法,发现HPD间隔是困难的。
我想用途:

p <- seq(0, 1, length = 100)
density <- dbeta(p, 7, 25)

找到HPD间期
另一个问题是,我知道我们通过计算dbeta(p, 7, 25)得到密度值,但是否存在dbeta的反函数?我试图得到p和相应的密度值。
提前感谢您的帮助!

velaa5lx

velaa5lx1#

看看TeachingDemos::hpd的代码。你需要的是相反的cdf(icdf)。

function(posterior.icdf, conf=0.95, tol=0.00000001,...){
    conf <- min( conf, 1-conf )
    f <- function(x,posterior.icdf,conf,...){
        posterior.icdf(1-conf+x,...) - posterior.icdf(x,...)
    }
    out <- optimize(f, c(0,conf), posterior.icdf = posterior.icdf,
                        conf=conf, tol=tol, ...)
    return( c( posterior.icdf(out$minimum,...), 
               posterior.icdf(1-conf+out$minimum,...) ) )
}

画一张图,你会发现这很容易。

相关问题