0
$\begingroup$

I have a sparse symmetric matrix of dimension 1393x1393 (8308 no zero elements), with bandwidth 1380. By Cuthill–McKee algorithm, I could achieve a new matrix with bandwidth 89. If we call this matrix B, then I'm interested in calculate the inverse I-rho*B, where I is a diagonal matrix and rho is a parameter updated in my MCMC algorithm. I need to calculate matrix inverse many times (during MCMC) and someone told me by reducing the matrix bandwidth, I can obtain faster matrix inverse.

I am using R as my programming language. The solve function is used to calculate the inverse of a matrix. I noticed no difference when calculating the inverse of a sparse matrix with bandwidth 1380 or bandwidth 89.

I guess I need some explicit commands/options in order to take advantage of the bandwidth reduction, but I don't know what keywords should I aim at. Can anyone give me some suggestions?

A former post with more information.

| cite | improve this question | |
$\endgroup$
  • 1
    $\begingroup$ Are you timing the solves or are you eyeballing the time it takes to solve? $\endgroup$ – KyleW Jun 29 '17 at 0:09
  • $\begingroup$ I know very little about R but, from looking at the documentation for solve and its underlying classes, it appears that the factorization methods are probably doing their own reordering of the equations for efficiency; so they are essentially discarding the ordering you are providing. Ideally, you would like to call a lower-level factorization routine once on your matrix and then reuse that factorization repeatedly in your algorithm. $\endgroup$ – Bill Greene Jun 29 '17 at 12:29
1
$\begingroup$

Three things immediately come to mind:

  • R might not take advantage of sparsity when using the solve command to compute the inverse of a matrix. Usually, the inverse of a sparse matrix is dense, so the built-in procedure might have been written to just convert it to a dense matrix, figuring that there's rarely anything to be gained from sparsity. To check this, see if you can print out any type info for the inverse matrix, or failing that check the number of non-zero entries in the inverse.
  • Reducing the bandwidth of a sparse matrix does reduce the fill-in of sparse matrix factorization. However, this is an indirect way to get at what you really want: matrix factors that are as sparse as possible. You might instead try using an approximate minimum degree ordering, which could give you better results if the factorization algorithm in R knows how to take advantage of sparsity.
  • Finally, computing the inverse of a matrix is very poorly conditioned from a numerical perspective in a lot of cases. Do you really need to be computing the inverse matrix, or is it enough to be able to multiply several vectors by the inverse matrix? In that case, a sparse LU factorization is enough for your purposes.
| cite | improve this answer | |
$\endgroup$

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.