R中多元正态分布的数值积分

Question

我正在尝试找到满足以下条件的 t1 and t2 的解决方案

# x1 and x2 are bivariate normal distributed with 
mu1=mu2=0
cov=matrix(c(1,1/2,1/2,1),2,2)

# trying to find the solution satisfying the equality with  t1=t2
# p(x1>t1) + p(x1<t2,x2>t2) = 0.05

我想知道是否有人可以帮助我编写 R.THX!

中的代码

Answer 1

p(x1>t) + p(x1<t,x2>t) = 0.05 等同于 p(x1<t,x2<t) = 0.95。因此：

library(mvtnorm)
qmvnorm(0.95, mean = rep(0, 2), sigma = cov)
#$quantile
#[1] 1.916399
#
#$f.quantile
#[1] 2.600961e-06
#
#attr(,"message")
#[1] "Normal Completion"

检查：

pmvnorm(c(1.916399, -Inf), c(Inf, Inf), , mean = rep(0, 2), sigma = cov) +
  pmvnorm(c(-Inf, 1.916399), c(1.916399, Inf), , mean = rep(0, 2), sigma = cov)
#[1] 0.04999262
#attr(,"error")
#[1] 2e-16
#attr(,"msg")
#[1] "Normal Completion"

Answer 2

在 R 中，您感兴趣的概率是

  pnorm(t, lower.tail = FALSE) + 
    pmvnorm(c(-Inf,t), c(t,Inf), sigma = Sigma) 



library(mvtnorm)
Sigma <- rbind(c(1,0.5),c(0.5,1))
f <- function(t){
  pnorm(t, lower.tail = FALSE) + 
    pmvnorm(c(-Inf,t), c(t,Inf), sigma = Sigma) - 0.05
}

uniroot(f, c(-6,6))
# $root
# [1] 1.916334

R中多元正态分布的数值积分

numerical integration of multivariate normal distribution in R

r

numerical-integration