3

我对 uniroot 命令有疑问。我无法以合乎逻辑的方式提出我的问题,因为我不知道为什么在下面示例的第二种情况下每次结果都不同。在第一种情况下,我的 f 函数的结果总是相同的:

library(mvtnorm)
f <- function(y1_upr,y2_upr = -0.05453663,target = 25e-4,df=3) {
    pmvt(upper = c(y1_upr,y2_upr),df = df) - target
}
uniroot(f,c(-10000,10000))$root

但是我不知道为什么当我在同一个函数中添加另一个参数时,结果每次都会改变(见下文):

g <- function(y1_upr,
              y2_upr = -0.05453663,
              y3_upr = -0.06236616,
              target = 25e-4,
              df = 3) {
    pmvt(upper = c(y1_upr,y2_upr,y3_upr),df = df) - target
} 
uniroot(g,c(-10000,10000))$root

您会看到,每次使用 g 函数应用 uniroot 命令时,都会得到不同的结果。有人对此有想法吗?我能找到解决结果的方法吗?

4

2 回答 2

3

该算法可能涉及一些随机性度量以选择起点。尝试:

set.seed(1) 
uniroot(g,c(-10000,10000))$root
于 2013-01-15T23:19:27.660 回答
2

@joran 是正确的,因为它的默认GenzBretz算法pmvt涉及随机元素。有一种替代算法,TVPACK它给出了一致的结果,但可能不太普遍适用:

replicate(10,pmvt(upper=c(0,-0.05453663,-0.06236616),df=3))
 [1] 0.1145384 0.1145367 0.1145365 0.1145336 0.1145381 0.1145354 0.1145385
 [8] 0.1145369 0.1145384 0.1145385
replicate(10,pmvt(upper=c(0,-0.05453663,-0.06236616),df=3,algorithm=TVPACK()))
 [1] 0.1145364 0.1145364 0.1145364 0.1145364 0.1145364 0.1145364 0.1145364
 [8] 0.1145364 0.1145364 0.1145364
于 2013-01-16T07:56:56.953 回答