Получить левый, правый, вверх, вниз ненулевой сосед из SciPy разреженной матрицы

Question

Скажите, что у меня 2D-матрица разреженности SciPy:

import numpy as np
from scipy.sparse import csc_matrix

arr = np.array([[0, 0, 1, 0, 1],
                [1, 0, 0, 1, 0],
                [0, 1, 1, 0, 0],
                [1, 0, 0, 1, 0],
                [0, 1, 0, 0, 0],
               ])

csc = csc_matrix(arr)

Для каждого ненулевого элемента в матрице я хотел бы создать четыре новых разреженных матрицы, которые содержат индекс, соответствующий следующему ближайшему ненулевому соседу Left, Right, Up и Down. Элементы на концах могут иметь соседей, которые обернуты вокруг (представьте круговой двусвязный список в горизонтальном и вертикальном направлениях или тороидальный). В случае, когда элемент является единственным ненулевым элементом в своей строке / столбце, соответствующий индекс будет указывать на себя. Кроме того, поскольку индексы могут иметь нулевое значение (при обращении к первой строке или столбцу) и быть неотличимыми от элементов с естественным нулем, мы устанавливаем эти нулевые индексы в -1, чтобы отсеять реальный индекс от нулевых элементов.

Для вышеприведенной матрицы плотные матрицы Left и Down будут выглядеть так:

left = np.array([[0, 0, 4,  0, 2],
                 [3, 0, 0, -1, 0],
                 [0, 2, 1,  0, 0],
                 [3, 0, 0, -1, 0],
                 [0, 1, 0,  0, 0],
                ])

down = np.array([[0, 0,  2, 0, -1],
                 [3, 0,  0, 3,  0],
                 [0, 4, -1, 0,  0],
                 [1, 0,  0, 1,  0],
                 [0, 2,  0, 0,  0],
                ])

Помните, что элементы со значением индекса -1 на самом деле являются ссылками на нулевой индекс. Конечно, мне нужно иметь эти матрицы в разреженной матричной форме, поскольку мои реальные матрицы слишком велики и разрежены, чтобы поместиться в память.

Answer 1

Более строгий подход:

csc = csc_matrix(arr)
inds = (csc.indices,csc.indptr)
irows = np.split(*inds)[1:-1]

down = csc_matrix((np.hstack([np.roll(row,-1) for row in irows]),*inds))
up = csc_matrix((np.hstack([np.roll(row,1) for row in irows]),*inds))

Проверка:

>>> down.A 
array([[0, 0, 2, 0, 0],
       [3, 0, 0, 3, 0],
       [0, 4, 0, 0, 0],
       [1, 0, 0, 1, 0],
       [0, 2, 0, 0, 0]], dtype=int32)

Слева и справа можно получить с помощью представления CSR.

Я не думаю, что кодирование 0 на -1 - это хорошая идея, так как if сломает все редкие улучшения вычислений. необходимо посетить только места, спроектированные csc.nonzeros().

Answer 2

Один возможный ответ (плотная форма):

ix, iy = csc.nonzero()
w = np.where(np.insert(np.diff(ix), 0,1) != 0)[0]
iy2 = np.concatenate([np.roll(_, 1) for _ in np.split(iy,w)])
iy2[iy2==0] = -1

left = csc_matrix(arr.shape)
left[ix, iy] = iy2

ix, iy = csc.transpose().nonzero()
w = np.where(np.insert(np.diff(ix), 0,1) != 0)[0]
iy2 = np.concatenate([np.roll(_, 1) for _ in np.split(iy,w)])
iy2[iy2==0] = -1

down = csc_matrix(arr.T.shape)
down[ix, iy] = iy2
down = down.transpose()
print(left.todense(), '\n', down.todense())


 >> [[ 0  0  4  0  2]
 [ 3  0  0 -1  0]
 [ 0  2  1  0  0]
 [ 3  0  0 -1  0]
 [ 0  1  0  0  0]]

[[ 0  0  2  0 -1]
 [ 3  0  0  3  0]
 [ 0  4 -1  0  0]
 [ 1  0  0  1  0]
 [ 0  2  0  0  0]]

Answer 3

In [183]: arr = np.array([[0, 0, 1, 0, 1],
     ...:                 [1, 0, 0, 1, 0],
     ...:                 [0, 1, 1, 0, 0],
     ...:                 [1, 0, 0, 1, 0],
     ...:                 [0, 1, 0, 0, 0],
     ...:                ])
     ...:                
In [184]: from scipy import sparse
In [185]: M = sparse.lil_matrix(arr)
In [186]: M.rows
Out[186]: 
array([list([2, 4]), list([0, 3]), list([1, 2]), list([0, 3]), list([1])],
      dtype=object)

Это та же информация, что и из плотного массива:

In [187]: [np.where(row)[0] for row in arr]
Out[187]: [array([2, 4]), array([0, 3]), array([1, 2]), array([0, 3]), array([1])]

Полагаю, вы уже выяснили, как сгенерировать желаемый left (или right) из плотного массива, поэтому я не буду вдаваться в эти детали (мне лень разбираться с вашими спецификациями упаковки) ).

Для столбцов:

 In [189]: M.T.rows
 Out[189]: 
 array([list([1, 3]), list([2, 4]), list([0, 2]), list([1, 3]), list([0])],
  dtype=object)

Из формата csc вы можете использовать:

In [190]: Mc = sparse.csc_matrix(arr)
In [191]: Mc.indptr
Out[191]: array([0, 2, 4, 6, 8, 9], dtype=int32)
In [192]: Mc.indices
Out[192]: array([1, 3, 2, 4, 0, 2, 1, 3, 0], dtype=int32)
In [193]: for i in range(5):
     ...:     print(Mc.indices[Mc.indptr[i]:Mc.indptr[i+1]])
     ...:     
[1 3]
[2 4]
[0 2]
[1 3]
[0]

В этом примере все строки или столбцы имеют только 1 или 2 ненулевых значения. Я предполагаю, что в более широком, более общем, будет много ненулевых. Также с csc (и csr есть вероятность того, что индексы для каждой строки не сортируются - есть редкий метод для решения этой проблемы.

Что касается построения разреженной матрицы возврата, вы можете изменить атрибут data копии (он будет иметь такую ​​же разреженность).

In [194]: M.data
Out[194]: 
array([list([1, 1]), list([1, 1]), list([1, 1]), list([1, 1]), list([1])],
      dtype=object)
In [195]: Mc.data
Out[195]: array([1, 1, 1, 1, 1, 1, 1, 1, 1], dtype=int64)

Или создать разреженную матрицу из массивов (как обычно для ввода в стиле coo).

---

С моей lil версией tch's решение немного быстрее:

ind = sparse.lil_matrix(M.shape,dtype='int')
for i,row in enumerate(M.rows):
    k = np.array(row)
    ind[i,k] = np.roll(k+1,1)

Еще лучше с моей идеей заменить data:

ind = M.copy()
for row,dat in zip(ind.rows,ind.data):
    k = np.array(row)
    dat[:] = np.roll(k+1,1).tolist()

или с Mr = Mc.tocsr()

ind = Mr.copy()
for i in range(Mr.shape[0]):
    slc = slice(Mr.indptr[i],Mr.indptr[i+1])
    k = Mr.indices[slc]
    ind.data[slc] = np.roll(k+1,1)

Answer 4

Вот возможный способ сделать левого соседа. Это не особенно эффективно, но, вероятно, работает хорошо, если во всей матрице не много ненулевых записей. Вы можете немного оптимизировать его, получая ненулевые записи каждой строки по ходу и вычисляя j[i==row] один раз.

Обратите внимание, что я просто сдвигаю индексы на один, а не устанавливаю 0 на -1.

i,j = csc.nonzero()
ind = sp.sparse.csc_matrix(csc.shape,dtype='int')
for row in range(csc.shape[0]):
    ind[row,j[i==row]] = np.roll(j[i==row]+1,1)

ind.A = array([[0, 0, 5, 0, 3],
   [4, 0, 0, 1, 0],
   [0, 3, 2, 0, 0],
   [4, 0, 0, 1, 0],
   [0, 2, 0, 0, 0]])