SciPy разреженная матрица (COO, CSR): Очистить строку

Question

Для создания scipy sparse matrix у меня есть индексы массива или строки и столбца I и J вместе с массивом данных V. Я использую те построить матрицу в COO format, а затем преобразовать его в CSR,

matrix = sparse.coo_matrix((V, (I, J)), shape=(n, n)) 
matrix = matrix.tocsr() 

У меня есть набор индексов строк, для которых единственная запись должна быть 1.0 по диагонали. До сих пор, я иду через I, найти все индексы, которые нужно протирая, и делать только что:

def find(lst, a): 
    # From <http://stackoverflow.com/a/16685428/353337> 
    return [i for i, x in enumerate(lst) if x in a] 

# wipe_rows = [1, 55, 32, ...] # something something 

indices = find(I, wipe_rows) # takes too long 
I = numpy.delete(I, indices).tolist() 
J = numpy.delete(J, indices).tolist() 
V = numpy.delete(V, indices).tolist() 

# Add entry 1.0 to the diagonal for each wipe row 
I.extend(wipe_rows) 
J.extend(wipe_rows) 
V.extend(numpy.ones(len(wipe_rows))) 

# construct matrix via coo 

Это работает хорошо, но find стремится занять некоторое время.

Любые подсказки о том, как ускорить это? (Возможно, вытирая строки в СОО или формат CSR является лучшей идеей.)

Answer 1

Если вы намерены очистить несколько строк сразу, это

def _wipe_rows_csr(matrix, rows): 
    assert isinstance(matrix, sparse.csr_matrix) 

    # delete rows 
    for i in rows: 
     matrix.data[matrix.indptr[i]:matrix.indptr[i+1]] = 0.0 

    # Set the diagonal 
    d = matrix.diagonal() 
    d[rows] = 1.0 
    matrix.setdiag(d) 

    return 

на сегодняшний день является самым быстрым способом. Он действительно не удаляет строки, но устанавливает все записи в нули, а затем путает с диагональю.

Если записи действительно должны быть удалены, необходимо выполнить некоторые манипуляции с массивами. Это может быть довольно дорогостоящим, но если скорость не проблема: Это

def _wipe_row_csr(A, i): 
    '''Wipes a row of a matrix in CSR format and puts 1.0 on the diagonal. 
    ''' 
    assert isinstance(A, sparse.csr_matrix) 

    n = A.indptr[i+1] - A.indptr[i] 

    assert n > 0 

    A.data[A.indptr[i]+1:-n+1] = A.data[A.indptr[i+1]:] 
    A.data[A.indptr[i]] = 1.0 
    A.data = A.data[:-n+1] 

    A.indices[A.indptr[i]+1:-n+1] = A.indices[A.indptr[i+1]:] 
    A.indices[A.indptr[i]] = i 
    A.indices = A.indices[:-n+1] 

    A.indptr[i+1:] -= n-1 

    return 

заменяет данную строку i матрицы matrix вхождением 1.0 по диагонали.

Answer 2

np.in1d должен быть более быстрый способ нахождения indices:

In [322]: I # from a np.arange(12).reshape(4,3) matrix 
Out[322]: array([0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3], dtype=int32) 

In [323]: indices=[i for i, x in enumerate(I) if x in [1,2]] 

In [324]: indices 
Out[324]: [2, 3, 4, 5, 6, 7] 

In [325]: ind1=np.in1d(I,[1,2]) 

In [326]: ind1 
Out[326]: 
array([False, False, True, True, True, True, True, True, False, 
     False, False], dtype=bool) 

In [327]: np.where(ind1) # same as indices 
Out[327]: (array([2, 3, 4, 5, 6, 7], dtype=int32),) 

In [328]: I[~ind1] # same as the delete 
Out[328]: array([0, 0, 3, 3, 3], dtype=int32) 

Прямая манипуляция coo входов, как этот часто является хорошим способом. Но другой должен использовать математические способности csr. Вы должны быть в состоянии построить диагональную матрицу, которая обнуляет правильные строки, а затем добавляет те обратно в

Вот что я имею в виду:.

In [357]: A=np.arange(16).reshape(4,4) 
In [358]: M=sparse.coo_matrix(A) 
In [359]: M.A 
Out[359]: 
array([[ 0, 1, 2, 3], 
     [ 4, 5, 6, 7], 
     [ 8, 9, 10, 11], 
     [12, 13, 14, 15]]) 

In [360]: d1=sparse.diags([(1,0,0,1)],[0],(4,4)) 
In [361]: d2=sparse.diags([(0,1,1,0)],[0],(4,4)) 

In [362]: (d1*M+d2).A 
Out[362]: 
array([[ 0., 1., 2., 3.], 
     [ 0., 1., 0., 0.], 
     [ 0., 0., 1., 0.], 
     [ 12., 13., 14., 15.]]) 

In [376]: x=np.ones((4,),bool);x[[1,2]]=False 
In [378]: d1=sparse.diags([x],[0],(4,4),dtype=int) 
In [379]: d2=sparse.diags([~x],[0],(4,4),dtype=int) 

Doing это с форматом lil выглядит просто:

In [593]: Ml=M.tolil() 
In [594]: Ml.data[wipe]=[[1]]*len(wipe) 
In [595]: Ml.rows[wipe]=[[i] for i in wipe] 

In [596]: Ml.A 
Out[596]: 
array([[ 0, 1, 2, 3], 
     [ 0, 1, 0, 0], 
     [ 0, 0, 1, 0], 
     [12, 13, 14, 15]], dtype=int32) 

это вроде того, что вы делаете с форматом csr, но это легко заменить каждый список строк с помощью соответствующей кнопки [1] и [I] списка. Но времена конверсии (tolil и т. Д.) Могут повредить время выполнения.