如何在二维图上绘制 $\alpha$ 置信区域？

Question

有很多关于 plotting confidence intervals 的答案。

我正在阅读 Lourme A. et al (2016) and I'd like to draw the 90% confidence boundary and the 10% exceptional points like in the Fig. 2 from the paper: 的论文。

我无法使用 LaTeX 插入带有置信区域定义的图片：

library("MASS")
library(copula)
set.seed(612)

n <- 1000 # length of sample
d <- 2    # dimension

# random vector with uniform margins on (0,1)
u1 <- runif(n, min = 0, max = 1)
u2 <- runif(n, min = 0, max = 1)

u = matrix(c(u1, u2), ncol=d)

Rg  <- cor(u)   # d-by-d correlation matrix
Rg1 <- ginv(Rg) # inv. matrix 

# round(Rg %*% Rg1, 8) # check

# the multivariate c.d.f of u is a Gaussian copula 
# with parameter Rg[1,2]=0.02876654

normal.cop = normalCopula(Rg[1,2], dim=d)
fit.cop    = fitCopula(normal.cop, u, method="itau") #fitting
# Rg.hat     = fit.cop@estimate[1]
# [1] 0.03097071
sim        = rCopula(n, normal.cop) # in (0,1)

# Taking the quantile function of N1(0, 1)

y1 <- qnorm(sim[,1], mean = 0, sd = 1)
y2 <- qnorm(sim[,2], mean = 0, sd = 1)

par(mfrow=c(2,2))

plot(y1, y2, col="red");  abline(v=mean(y1), h=mean(y2))
plot(sim[,1], sim[,2], col="blue")
hist(y1); hist(y2)

参考。 Lourme, A., F. Maurer (2016) 在风险管理框架中测试高斯和学生 t copula。经济建模。

问题。谁能帮我解释一下方程中的变量v=(v_1,...,v_d)和G(v_1),..., G(v_d)？

我认为v是非随机矩阵，维度应该是$k^2$（网格点）乘以d=2（维度）。例如，

axis_x <- seq(0, 1, 0.1) # 11 grid points
axis_y <- seq(0, 1, 0.1) # 11 grid points
v <- expand.grid(axis_x, axis_y)
plot(v,  type = "p")

Answer 1

所以，你的问题是关于向量 nu 和对应的 G(nu)。

nu 是从具有定义域 (0,1) 的 any 分布中提取的简单随机向量。（这里我使用均匀分布）。由于您想要 2D 样本，因此单个 nu 可以是 nu = runif(2)。鉴于上面的解释，G 是均值为 0 和协方差矩阵 Rg 的高斯 pdf。（Rg 在 2D 中的尺寸为 2x2）。

现在这一段是怎么说的：如果你有一个随机样本 nu 并且你希望它是从 Gamma 给定的维数 d 和置信水平 alpha 那么您需要计算以下统计数据 (G(nu) %*% Rg^-1) %*% G(nu) 并检查它是否低于 d 和 alpha.

的 Chi^2 分布的 pdf

例如：

# This is the copula parameter
Rg <- matrix(c(1,runif(2),1), ncol = 2)
# But we need to compute the inverse for sampling
Rginv <- MASS::ginv(Rg)

sampleResult <- replicate(10000, {
  # we draw our nu from uniform, but others that map to (0,1), e.g. beta, are possible, too
  nu <- runif(2)
  # we compute G(nu) which is a gaussian cdf on the sample
  Gnu <- qnorm(nu, mean = 0, sd = 1)
  # for this we compute the statistic as given in formula
  stat <- (Gnu %*% Rginv) %*% Gnu
  # and return the result
  list(nu = nu, Gnu = Gnu, stat = stat)
})

theSamples <- sapply(sampleResult["nu",], identity)

# this is the critical value of the Chi^2 with alpha = 0.95 and df = number of dimensions
# old and buggy threshold <- pchisq(0.95, df = 2)
# new and awesome - we are looking for the statistic at alpha = .95 quantile
threshold <- qchisq(0.95, df = 2)
# we can accept samples given the threshold (like in equation)
inArea <- sapply(sampleResult["stat",], identity) < threshold

plot(t(theSamples), col = as.integer(inArea)+1)

红色的点是你要保留的点（我在这里绘制了所有点）。

至于绘制决策边界，我认为它有点复杂，因为您需要计算 nu 的确切对，以便 (Gnu %*% Rginv) %*% Gnu == pchisq(alpha, df = 2)。这是一个线性系统，您可以求解 Gnu，然后应用逆运算在决策边界处得到 nu。

编辑： 再次阅读该段，我注意到，Gnu 的参数没有改变，只是 Gnu <- qnorm(nu, mean = 0, sd = 1).

编辑： 有一个错误：对于阈值，您需要使用分位数函数 qchisq 而不是分布函数 pchisq - 现在已更正上面的代码（并更新了数字）。

Answer 2

这有两个部分：首先，计算作为 X 和 Y 函数的 copula 值；然后，绘制给出 copula 超过阈值的边界的曲线。

计算值基本上是@drey 回答的线性代数。这是重写的版本，因此 copula 由函数给出。

cop1 <- function(x)
{
    Gnu <- qnorm(x)
    Gnu %*% Rginv %*% Gnu
}

copula <- function(x)
{
    apply(x, 1, cop1)
}

可以使用与相同的方法绘制边界曲线（这又是教科书 Modern Applied Stats with S 和 Elements of Stat Learning 使用的方法）。创建一个值网格，并使用插值法找到给定高度的等高线。

Rg <- matrix(c(1,runif(2),1), ncol = 2)
Rginv <- MASS::ginv(Rg)

# draw the contour line where value == threshold
# define a grid of values first: avoid x and y = 0 and 1, where infinities exist
xlim <- 1e-3
delta <- 1e-3
xseq <- seq(xlim, 1-xlim, by=delta)
grid <- expand.grid(x=xseq, y=xseq)
prob.grid <- copula(grid)
threshold <- qchisq(0.95, df=2)

contour(x=xseq, y=xseq, z=matrix(prob.grid, nrow=length(xseq)), levels=threshold,
        col="grey", drawlabels=FALSE, lwd=2)

# add some points
data <- data.frame(x=runif(1000), y=runif(1000))
points(data, col=ifelse(copula(data) < threshold, "red", "black"))

如何在二维图上绘制 $\alpha$ 置信区域？

How to draw an $\alpha$ confidence areas on a 2D-plot?

plot

r

contour

confidence-interval