[SciPy-dev] Sparse matrix indexing.

[Viral Shah]

Hello,

We are working to convert Circuitscape (a landscape ecology tool) from  
MATLAB to Python.
http://www.nceas.ucsb.edu/~mcrae/circuitscape.html

It exercises a lot of sparse array functionality (indexing,  
assignment, concatenation, and MATLAB sparse()/find(), and sparse  
linear algebra). I am new to Numpy/Scipy, and I'd be delighted to  
figure out if there were a better way to do sparse matrix indexing.  
The general array indexing didn't seem to work for arbitrary index  
sets, at least in the latest release.

I did come across this on scipy-dev, but it seems like it works for  
only CSR ?
http://article.gmane.org/gmane.comp.python.scientific.devel/6702/match=sparse+indexing

I currently use a wrapper that can be used to implement general sparse  
matrix indexing and row/col deletion using sparse matrix  
multiplication. It implements sparse indexing. Deletion, is then  
implemented by simply taking arange(0,len) and removing the indices  
you want to delete, and then calling the indexing code.

This may be a useful way to complete the functionality for data  
structures that may not have full indexing yet. I'm not sure if this  
can be made to handle slices correctly. This is the code snippet I  
use. If its useful, I can help make it robust for inclusion.

     # Implement B = A[I, J]
     def subsref(self,A,I,J):
         nr = A.shape[0]
         nc = A.shape[1]
         nI = len(I)
         nJ = len(J)

         IM = sparse.csc_matrix((ones(nI), c_[arange(0,nI), I].T),  
dims=[nI,nr])
         JM = sparse.csc_matrix((ones(nJ), c_[J, arange(0,nJ)].T),  
dims=[nc,nJ])
         B = IM*A*JM

         return B


-viral