How to operate elementwise on a matrix of type scipy.sparse.csr_matrix?

Question

In numpy if you want to calculate the sinus of each entry of a matrix (elementise) then

a = numpy.arange(0,27,3).reshape(3,3)
numpy.sin(a)

will get the job done! If you want the power let's say to 2 of each entry

a**2

will do it.

But if you have a sparse matrix things seem more difficult. At least I haven't figured a way to do that besides iterating over each entry of a lil_matrix format and operate on it.

I've found this question on SO and tried to adapt this answer but I was not succesful.

The Goal is to calculate elementwise the squareroot (or the power to 1/2) of a scipy.sparse matrix of CSR format.

What would you suggest?

Schoolmistress

possible duplicate of Element-wise power of scipy.sparse matrix – Shrimp

Answer 1

The following trick works for any operation which maps zero to zero, and only for those operations, because it only touches the non-zero elements. I.e., it will work for sin and sqrt but not for cos.

Let X be some CSR matrix...

>>> from scipy.sparse import csr_matrix
>>> X = csr_matrix(np.arange(10).reshape(2, 5), dtype=np.float)
>>> X.A
array([[ 0.,  1.,  2.,  3.,  4.],
       [ 5.,  6.,  7.,  8.,  9.]])

The non-zero elements' values are X.data:

>>> X.data
array([ 1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9.])

which you can update in-place:

>>> X.data[:] = np.sqrt(X.data)
>>> X.A
array([[ 0.        ,  1.        ,  1.41421356,  1.73205081,  2.        ],
       [ 2.23606798,  2.44948974,  2.64575131,  2.82842712,  3.        ]])

Update In recent versions of SciPy, you can do things like X.sqrt() where X is a sparse matrix to get a new copy with the square roots of elements in X.

Shrimp

@Aufwind: I got so sick and tired of operating on .data that I submitted a patch to Scipy implementing sqrt, sin and tan on CSR and CSC matrices. – Shrimp

This is great, I am sure many people will be verythankful! :-) – Schoolmistress

There is a reason not to support these operations on matrices directly; exponentiation of a matrix already has a standard and quite different meaning from the elementwise one proposed here. Personally, I think sparse matrices ought to respect linear algebra semantics by default; if what you really want is a sparse array, then maybe there should be a subclass for that, analogous with their dense counterparts. – Molest

@EelcoHoogendoorn SciPy sparse matrices effectively follow linear algebra conventions, which is a pain because NumPy arrays don't. SciPy should really have sparse arrays, but every attempt to implement them so far has failed. – Shrimp

Apparently not, as per the latest scipy, since sparse.sqrt is an elementwise squareroot. Personally, ive often used sparse matrices, and never had a need for sparse arrays. So I am ok with sparse matrices (mostly) acting as matrices. The base class is called spmatrix; not sparray. That said, an sparray could be neat; but then youd also want an sp-nd-array generalization, and its not directly clear how to go forward with that. But as for low hanging fruit, I think spmatrix should be castable to a 2d sparray, in order for array semantics to be enabled, and spmatrix should have pure LA semantics. – Molest

Is there any way to access the index of elements with this? That is, performing A[i,j] = f(A[i,j], i, j) as opposed to A[i,j] = f(A[i,j])? – Uncoil