1
$\begingroup$

In MATLAB, i compared elapsed time to invert a Hermitian matrix using inverse(), svd(), and chol(). svd() took the longest. So is there any reason to prefer svd() to the other two methods?

$\endgroup$
1
  • $\begingroup$ in all honesty, this is not a signal processing, but a numerical math question, and might be better off on the math.SE sister site. However, VVT has made the right point: different problems scream for different solutions, and you can't compare things that are applied to different problems. $\endgroup$ – Marcus Müller Dec 28 '20 at 9:18
3
$\begingroup$

The two methods differ, above all, by their applicability to matrix classes.

col (cholesky) decomposes Hermitian, positive-definite rectangular matrices into the product of a lower triangular matrix and its conjugate transpose;

svd (singular value decomposition) factorizes any m×n matrix into the form UΣV*, where U and V are square real or compex unitary matrices, m×m and n×n, respectively, and Σ is an m×n rectangular diagonal matrix with non-negative real numbers on the diagonal.

The method inv internally performs an LU decomposition of the input matrix (or an LDL decomposition if the input matrix is Hermitian), but outputs only the inverse of square matrix only.

Both SVD and Cholesky can be used for computing pseudoinverse of a matrix, provided the matrix satisfies requirement for the method used. The pseudoinverse operation is used to solve linear least squares problems and the other signal processing, image processing, and big data problems.

UPDATE on OP's comment

The matrix can be both Hermitian and not a positive/negative (semi)definite, in which case it is called an Hermitian indefinite matrix. The (generalized) Cholesky decomposes Hermitian positive semidefinite matrices. You can also decompose a negative (semi)definite matrix, say, Anegdef, just call chol()) on -Anegdef, but you cannot compute an inverse of an indefinite matrix with Cholesky, because of unavoidable square roots of negative numbers.

Rank-deficient matrices (having zero singular values) are not invertible, as is the case for most general-case matrices. That is, there exists no matrix A-1 such that A-1A = AA-1 = I, for rank-deficient matrices. Still, many problems which you solve through matrix inversion can be solved for indefinite (and consequently non-invertible) matrices with a generalization of matrix inversion, pseudoinverses of matrices.

I won't bore you with the description of Moore-Penrose pseudoinverse concept, on which subject many articles are readily available online and in textbooks: in short, the matrix A+ is a pseudoinverse of matrix A, if AA+A = A (consequently, A+AA+ = A+). If A has linearly independent columns, the pseudoinverse is a left inverse, because in this case A+A = I, you can check it with a sample low-dimension matrix constructed for this purpose. If A has linearly independent rows, the pseudoinverse is a right inverse, because in this case AA+ = I. In general case, A may have no left/right inverses, but still have pseudoinverse A+ satisfying the equation AA+A = A.

Back to quadratic programming problems: Below, all matrices and vectors are real-valued, so instead of Hermitian conjugation H we have a matrix transpose T here. For an invertible (in a "classic" sense) matrix A, a quadratic form $$ f(x) = {1\over2}x^{\textrm{T}}Ax-x^{\textrm{T}}b $$ has a minimum value at xopt = A-1b: $$ f(x_{opt}) = -{1\over2}b^{\textrm{T}}A^{\textrm{-1}}b $$

To compute an inverse of A, you use function chol() here.

You can solve a minimization problem for a quadratic form with a non-invertible matrix A, provided A is positive semidefinite, even if A has no inverse in this case. For the singular matrix A, the minimum value is $$ f(x_{opt}) = -{1\over2}b^{\textrm{T}}A^{\textrm{+}}b $$ which is reached for all x given by $$ x = A^{\textrm{+}}b + U^{\textrm{T}}{\left(\begin{align}{0\\z}\end{align}\right)}, \hspace{1cm} z \in \mathbf{R}^{\textrm{n-r}} $$ , n is the dimension of x, r is the rank of A.

Matrix U in this solution is the factor of singular value decomposition of A = UΣVT, which you compute with the function svd(). So, maybe not all "the intermediate products, whether they be triangular matrices, eigenvectors, diagonal matrices, etc. are irrelevant when the only things that matter are the computation time(FLOPS)"?

In these examples, all matrices and vectors are real-valued, and I am not sure this example can be readily applied to your beamforming problem. Still, SVD is certainly useful in signal processing, e.g., see Image Reconstruction using SVD.

For the sake of completeness: 1) there exists a generalization of Cholesky decomposition useful for computation of a positive semidefinite matrix, and the pseudoinverse can be computed with the help of Greville's algorithm. But this algorithm is implemented through recursion on matrix columns, and the Cholesky speed can be offset by this recursion; 2) the "stability issues" are addressed by the LDL decomposition.

Any way, these algebraic technicalities are not "too much information" for those who work in signal processing: if, even without a preliminary analysis on positive-definiteness of matrices, you start compute using a chol() function and it turns out that you are getting weird results, you now know what might be the culprit.

$\endgroup$
8
  • $\begingroup$ Thank you for your responses. To be precise, this question pertained to inverting a Hermitian matrix in a signal processing application. I guess I still don't understand why svd() would ever be preferred in such a case, since it takes longer than the other two methods. My understanding is that the intermediate products, whether they be triangular matrices, eigenvectors, diagonal matrices, etc. are irrelevant when the only things that matter are the computation time(FLOPS), memory usage, and accuracy/stability of the resulting inv matrix. I read that the Cholesky method has stability issues... $\endgroup$ – j03y_ Dec 28 '20 at 23:47
  • $\begingroup$ Now that you've made known the context of your study in your question on unconstrained array gain, I can be more specific in my answer on the svd()/chol() inter-comparison, see the answer UPDATE. $\endgroup$ – V.V.T Jan 3 at 18:46
  • $\begingroup$ Thank you - So, are you saying that ... if a matrix is not ... directly invertible, you can still ... get its "pseudo-inverse" by taking its svd()? I may have misunderstood that the svd() enables you to get the literal inverse instead of the "pseudo-inverse." I think i should read this, one of my references, more closely: pillowlab.princeton.edu/teaching/statneuro2018/slides/… $\endgroup$ – j03y_ Jan 3 at 21:43
  • $\begingroup$ Basically, this is what I said, you get it right. Only, before diving deep into the linear algebra, make sure you definitely need pseudoinverses in your project. In your case, examine if your Hermitian matrix in beamforming calculations can be indefinite or rank deficient. If this is not the case, the awareness of the problem is sufficient for success. But, if your matrices can become degenerate, you have no other choice as to master pseudoinverses. Sort of "learning-on-demand". $\endgroup$ – V.V.T Jan 4 at 11:01
  • $\begingroup$ Also, the separation of use cases for application of Cholesky vs SVD was blurred, see the edited text after UPDATE subtitle. $\endgroup$ – V.V.T Jan 4 at 11:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.