如何在2D图上绘制$\alpha $置信区域？

Question

关于绘制置信区间有很多答案.

我正在阅读Lourme A.等人(2016年)的论文,我想从论文中得出90%的置信区间和10%的特殊点,如图2所示:.

我不能使用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)在风险管理框架中测试Gaussian和Student's t copulas.经济模型.

题.任何人都可以帮助我并解释变量v=(v_1,...,v_d)和G(v_1),..., G(v_d)等式吗？

我认为v是非随机矩阵,维度应该是d=2(维度)$ k ^ 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) 的任何分布中抽取的简单随机向量。（这里我使用均匀分布）。由于您想要 2D 样本，因此nu可以使用单个样本nu = runif(2)。鉴于上述解释，G是一个均值为 0 且协方差矩阵 的高斯 pdf Rg。（Rg 的二维尺寸为 2x2）。

现在该段落说的是：如果您有一个随机样本，并且希望从给定的维数和置信水平nu中抽取它，那么您需要计算以下统计量并检查其是否低于 Chi^2 分布的 pdf和。Gammadalpha(G(nu) %*% Rg^-1) %*% G(nu)dalpha

例如：

# 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- 现在在上面的代码中已更正（并更新了数字）。