在q中给出它的分位数，求解qbeta（）中的shape1和shape2 [英] Solving for shape1 and shape2 in qbeta() given its quantiles in R

问题描述

背景：

qbeta（p，shape1，shape2）是内置基地R功能。在这个函数中， shape1 和 shape2 总是> 0（我的目的总是> 1）。此外， p 是概率，因此0 <= p <= 1。

问题：<假设我知道 qbeta（p = c（.025，.975），shape1 = x，shape2 = y）= c （.3，.8），那么，是否有办法找出 x 和 y 可以合理地（即，在一些误差）？

注意：我相信，有一些相关的评论和代码在 HERE 即可。因此，如果其中一位同事能够理解或使用这些代码和评论并将其应用于此，我将深表感谢。

最后，我的最终目标是用 x & y 找到过程。

这是一个强力方法。

  d = expand.grid（shape1 = 1）我希望有人能想出更明智的方法来做到这一点。 ：100/10，shape2 = 1：100/10）
 
 for（i in 1：NROW（d））{
x = round（qbeta（p = c（0.025，0.975）如果（x [1]> 0.28& x [1] <0.32& x [2]，则形状1 = d [i，1]，形状2 = d [i，2] > 0.78& x [2]< 0.82）{
 print（d [i，]）
 break 
} 
} 
 
＃ shape1 shape2 
＃5368 6.8 5.4 
 
 qbeta（p = c（.025，.975），shape1 = 6.8，shape2 = 5.4）
＃[1] 0.2860268 0.8109389 
  

我正在提供 whuber 的解决方案来自此处

  ################################ 
 #REQUIRED功能
＃Beta CDF的Logistic转换。 
 f.beta函数（alpha，beta，x，lower = 0，upper = 1）{
p <-pbeta（（x-lower）/（upper-lower），alpha， beta）
 log（p /（1-p））
} 
 
＃平方和。 
 delta<  - 函数（fit，actual）{
 sum（（fit-actual）^ 2）
} 
 
＃目标函数处理转换后的参数`theta`和
＃使用`f.beta`和`delta`来适应这些值并测量它们的差异。 
目标<  - 函数（theta，x，prob，...）{
 ab <-exp（θ）＃参数是alpha和beta 
 fit< ;  -  f.beta（ab [1]，ab [2]，x，...）
 return（delta（fit，prob））
} 
 ###### ########################## 
 
 #INPUT 
x<  -  c（.3，.8 ）
 xp < - （函数（p）log（p /（1-p）））（c（0.025,0.975））
 
 #SOLVE 
 start< （1，1，1，1，2，...，1），其中，x（x，x）= 0，upper = 1，
 typsize = c（1， 1），fscale = 1e-12，gradtol = 1e-12）
 parms <-exp（sol $估计）
 
 #ANSWER 
 parms 
＃[ 1] 7.597429 6.016240 
 
 #TEST 
 qbeta（p = c（.025，.975），parms [1]，parms [2]）
＃[1] 0.3 0.8 
  

Background:

qbeta(p, shape1, shape2) is an built-in base R function. In this function, shape1 and shape2 are always > 0 (for my purposes always > 1). Also, p is probability thus 0 <= p <= 1.

Question:

Suppose I know qbeta(p = c(.025, .975), shape1 = x, shape2 = y ) = c(.3, .8 ), THEN, is there a way I can find out what x and y could reasonably be (i.e., within some margin of error)?

Note: There are, I believe, some relevant comments and codes provided in HERE. So, in case one of the colleagues can makes sense or use of those codes and comments and apply them here, I'll deeply appreciate it.

Finally, my ultimate goal is to wrap a function out of this x & y finding process.

解决方案

Here is a brute force method. I'm hoping someone else can come up with a smarter way to do this.

d = expand.grid(shape1 = 1:100/10, shape2 = 1:100/10)

for (i in 1:NROW(d)){
    x = round(qbeta(p = c(0.025, 0.975), shape1 = d[i,1], shape2 = d[i,2]),2)
    if(x[1] > 0.28 & x[1] < 0.32 & x[2] > 0.78 & x[2] < 0.82){
        print(d[i,])
        break
    }
}

#     shape1 shape2
#5368    6.8    5.4

qbeta(p = c(.025, .975), shape1 = 6.8, shape2 = 5.4)
#[1] 0.2860268 0.8109389

I am providing whuber's solution from here below

################################
#REQUIRED FUNCTIONS
# Logistic transformation of the Beta CDF.
f.beta <- function(alpha, beta, x, lower=0, upper=1) {
    p <- pbeta((x-lower)/(upper-lower), alpha, beta)
    log(p/(1-p))
}

# Sums of squares.
delta <- function(fit, actual){
    sum((fit-actual)^2)
}

# The objective function handles the transformed parameters `theta` and
# uses `f.beta` and `delta` to fit the values and measure their discrepancies.
objective <- function(theta, x, prob, ...) {
    ab <- exp(theta) # Parameters are the *logs* of alpha and beta
    fit <- f.beta(ab[1], ab[2], x, ...)
    return (delta(fit, prob))
}
################################

#INPUT
x <- c(.3, .8)
x.p <- (function(p) log(p/(1-p)))(c(0.025, 0.975))

#SOLVE
start <- log(c(1e1, 1e1))
sol <- nlm(objective, start, x=x, prob=x.p, lower=0, upper=1,
                        typsize=c(1,1), fscale=1e-12, gradtol=1e-12)
parms <- exp(sol$estimate)

#ANSWER
parms
#[1] 7.597429 6.016240

#TEST
qbeta(p = c(.025, .975), parms[1], parms[2])
#[1] 0.3 0.8

这篇关于在q中给出它的分位数，求解qbeta（）中的shape1和shape2的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！