r - R中随机指数相关矩阵的Cholesky分解

Question

我有一组使用以下代码创建的指数相关矩阵。

x=runif(n)
R=matrix(0,n,n)  
for (j in 1:n)
{
  for(k in 1:n)
  {
    R[j,k]=exp(-(x[j]-x[k])^2);
  }
}

现在我想得到他们的 Cholesky 分解。但其中许多是负确定的。我该如何解决这个问题？

Answer 1

空间或时间建模中使用的指数相关矩阵有一个alpha控制衰减速度的因素：

exp(- alpha * (x[i] - x[j]) ^ 2))

您已将此类因子固定为 1。但实际上，此类因子是根据数据估算的。

请注意，这alpha是确保数值正定性所必需的。该矩阵原则上是正定的，但如果alpha不足以快速衰减，则在数值上不是。

鉴于此， 和 之间的x <- runif(n, 0, 1)距离聚集在一个短距离内。看到相关性衰减的范围并不大，也许您想尝试一下。x[i]x[j][0, 1]alpha = 10000

或者，如果您想留在alpha = 1，则需要使距离更加分散。试试x <- runif(n, 0, 100)。衰减非常快，即使alpha = 1.

所以我们看到了距离和之间的对偶alpha。这也是为什么这种相关矩阵可以稳定地用于统计建模的原因。什么时候alpha估计，可以自适应距离，使得相关矩阵总是正定的。

---

例子：

f <- function (xi, xj, alpha) exp(- alpha * (xi - xj) ^ 2)

n <- 100

# large alpha, small distance
x <- runif(n, 0, 1)
A <- outer(x, x, f, alpha = 10000)
R <- chol(A)

# small alpha, large distance
x <- runif(n, 0, 100)
A <- outer(x, x, f, alpha = 1)
R <- chol(A)

Answer 2

尝试使用它来构造正定矩阵

A<-matrix(runif(n^2),n,n)
dim(A)
A<-A%*%t(A)
chol(A)