我需要以下对角线：
diag(X %*% solve(A) %*% t(X))
其中A 是满秩方阵，X 是矩形矩阵。 A 和 X 都是稀疏的。
我知道除非你真的需要它，否则找到矩阵的逆是不好的。但是，我看不到如何重写公式，以便将 solve(A) 替换为带有两个参数的 solve，从而在不显式反转的情况下求解线性系统。这可能吗？
【问题讨论】：
标签： r algorithm matrix matrix-inverse
也许在我真正开始之前，我需要提一下，如果你这样做了
diag(X %*% solve(A, t(X)))
避免了矩阵求逆。 solve(A, B) 执行 LU 分解并使用生成的三角矩阵因子来求解线性系统 A x = B。如果您未指定B，则默认为对角矩阵，因此您将显式计算A 的逆矩阵。
您应该仔细阅读?solve，并可能多次阅读以获取提示。它说它基于LAPACK 例程DGESV，您可以在其中找到幕后的数值线性代数。
DGESV computes the solution to a real system of linear equations
A * X = B,
where A is an N-by-N matrix and X and B are N-by-N RHS matrices.
The LU decomposition with partial pivoting and row interchanges is
used to factor A as
A = P * L * U,
where P is a permutation matrix, L is unit lower triangular, and U is
upper triangular. The factored form of A is then used to solve the
system of equations A * X = B.
solve(A, t(X)) 和solve(A) %*% t(X) 之间的区别在于效率问题。后者中的一般矩阵乘法%*%比solve本身要贵很多。
然而，即使你使用solve(A, t(X))，你也不是最快的，因为你还有另一个%*%。
此外，由于您只需要对角线元素，因此在第一次获取完整矩阵时会浪费大量时间。 我在下面的回答会让你走上最快的轨道。
我用 LaTeX 重写了所有内容，内容也大大扩展，包括对 R 实现的引用。最后给出了关于 QR 分解、奇异值分解和 Pivoted Cholesky 分解的额外资源，以防您发现它们有用。
额外资源
如果您对旋转 Cholesky 分解感兴趣，可以参考Compute projection / hat matrix via QR factorization, SVD (and Cholesky factorization?)。在那里我还讨论了 QR 分解和奇异值分解。
上面的链接是在普通的最小二乘回归上下文中设置的。加权最小二乘可以参考Get hat matrix from QR decomposition for weighted least square regression。
QR 分解也采用紧凑形式。如果您想了解更多关于 QR 分解是如何完成和存储的，您可以参考What is "qraux" returned by QR decomposition。
这些问题和答案都集中在数值矩阵计算上。下面给出一些统计应用：
【讨论】：