本文介绍了如何在不采用矩阵逆的情况下有效计算diag(X%*%solve(A)%*%t(X))?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我需要以下对角线:

diag(X %*% solve(A) %*% t(X))

其中,A是满秩方阵,X是矩形矩阵. AX都是稀疏的.

where A is a full rank square matrix, X is a rectangular matrix. Both A and X are sparse.

我知道找到矩阵的逆是不好的,除非您真的需要它.但是,我看不到如何重写公式,以便用两个参数将solve(A)替换为solve,这样就可以在不显式求逆的情况下求解线性系统.有可能吗?

I know finding the inverse of a matrix is bad unless you really need it. However, I can't see how to rewrite the formula such that solve(A) is replaced by solve with two arguments, such that a linear system is solved without explicitly inverting. Is that possible?

推荐答案

也许在我真正开始之前,我要提一下,如果你这样做

Perhaps before I really start, I need to mention that if you do

diag(X %*% solve(A, t(X)))

避免矩阵逆. solve(A, B)执行LU分解,并使用所得的三角矩阵因子求解线性系统A x = B.如果未指定B,则默认为对角矩阵,因此将显式计算A的矩阵逆.

matrix inverse is avoided. solve(A, B) performs LU factorization and uses the resulting triangular matrix factors to solve linear system A x = B. If you leave B unspecified, it is default to a diagonal matrix hence you will be explicitly computing matrix inverse of A.

您应该仔细阅读?solve,并可能需要多次阅读提示.它说它基于LAPACK例程 DGESV ,您可以在其中找到数字线性代数.

You should read ?solve carefully, and possibly many times for hints. It says it is based on LAPACK routine DGESV, where you can find the numerical linear algebra behind the scene.

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本身昂贵得多.

The difference between solve(A, t(X)) and solve(A) %*% t(X) is a matter of efficiency. The general matrix multiplication %*% in the latter is much more expensive than solve itself.

但是,即使您使用solve(A, t(X)),您也不会处于最快的轨道,因为您还有另一个%*%.

Yet, even if you use solve(A, t(X)), you are not on the fastest track, as you have another %*%.

此外,由于只希望使用对角线元素,因此在第一次获取完整矩阵时会浪费大量时间. 我在下面的回答将使您处在最快的轨道上.

Also, as you only want diagonal elements, you waste considerable time in first getting the full matrix. My answer below will get you on the fastest track.

我已经用LaTeX重写了所有内容,并且内容也大大扩展了,包括对R实现的引用.最后,如果发现有用的话,还会提供有关QR分解,奇异值分解和Pivoted Cholesky分解的更多资源.

I have rewritten everything in LaTeX, and the content is greatly expanded, too, including reference to R implementation. Extra resources are given on QR factorization, Singular Value decomposition and Pivoted Cholesky factorization in the end, in case you find them useful.

额外资源

如果您对透视式Cholesky分解感兴趣,可以参考通过QR分解,SVD(和Cholesky)计算投影/帽子矩阵因数分解?).在那里,我还将讨论QR因式分解和奇异值分解.

In case you are interested in pivoted Cholesky factorization, you can refer to Compute projection / hat matrix via QR factorization, SVD (and Cholesky factorization?). There I also discuss QR factorization and singular value decomposition.

以上链接是在普通最小二乘回归上下文中设置的.有关加权最小二乘,您可以参考从QR分解中获取帽子矩阵以进行加权最小二乘回归.

The above link is set in ordinary least square regression context. For weighted least square, you can refer to Get hat matrix from QR decomposition for weighted least square regression.

QR因式分解也采用紧凑形式.如果您想进一步了解QR分解的完成方式和存储方式,请参考 QR分解返回的"qraux"是什么.

QR factorization also takes compact form. Should you want to know more on how QR factorization is done and stored, you may refer to What is "qraux" returned by QR decomposition.

这些问题和答案都集中在数值矩阵计算上.下面给出一些统计应用:

These questions and answers all focus on numerical matrix computations. The following gives some statistical application:

  • linear model with lm(): how to get prediction variance of sum of predicted values
  • Sampling from multivariate normal distribution with rank-deficient covariance matrix via Pivoted Cholesky Factorization
  • Boosting lm() by a factor of 3, using Cholesky method rather than QR method

这篇关于如何在不采用矩阵逆的情况下有效计算diag(X%*%solve(A)%*%t(X))?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!

09-09 18:13