R中加权最小二乘协方差矩阵的加速逆

Question

我需要加速计算R中WLS协方差矩阵的逆，其中矩阵wls.cov.matrix由(完整示例如下)给出:

n = 10000
X = matrix(c(rnorm(n,1,2), sample(c(1,-1), n, replace = TRUE), rnorm(n,2,0.5)), nrow = 1000, ncol = 3)
Q = diag(rnorm(n, 1.5, 0.3))
wls.cov.matrix = solve(t(X)%*%diag(1/diag(Q))%*%X)

有可能加快计算速度吗?

更多信息非常相关的最终目标:这仍然是很少的信息，让我解释更多我的目标，如果有办法加快我的代码会更清楚。
我运行了成千上万次，所以我需要更快的速度。

然而，我每次运行它都使用相同的X，唯一改变的矩阵是Q，它是一个对角矩阵。

如果X是一个与Q相同的方阵，我可以预先计算X^-1和(X^T)^(-1)，

X.inv = solve(X)
X.inv.trans = solve(t(X))

，然后对于每个迭代运行:

Q.inv = diag(1/diag(Q))
wls.cov.matrix = X.inv%*%Q.inv%*%X.inv.trans

但是我的X不是正方形的，所以还有其他的技巧吗?

Answer 1

这里最耗时的部分是t(X)%*%diag(1/diag(Q))%*%X，而不是计算它的逆。

一个不错的技巧是将其计算为

crossprod(X / sqrt(diag(Q)));

确认:

all.equal( (t(X) %*% diag(1/diag(Q)) %*% X) , crossprod(X / sqrt(diag(Q))) );
[1] TRUE

比较计时运行:

Qdiag = diag(Q);
system.time({(t(X) %*% diag(1/Qdiag) %*% X)})
system.time({crossprod(X / sqrt(Qdiag))})

Answer 2

好吧，Q是一个对角矩阵，所以它的逆矩阵就是对角元素的逆矩阵。因此，您可以执行

X = matrix(c(rnorm(n,1,2), sample(c(1,-1), n, replace = TRUE), rnorm(n,2,0.5)), nrow = 1000, ncol = 3)
Qinv = diag(1/rnorm(n, 1.5, 0.3))
wls.cov.matrix = solve(t(X)%*%Qinv%*%X)

事实上，这使速度提高了大约20倍。