Scipy редкие... массивы?

Итак, я делаю классификацию Kmeans, используя массивы numpy, которые довольно разрежены - много и много нулей. Я решил, что я использую scipy "разреженный" пакет, чтобы уменьшить накладные расходы на хранение, но я немного смущен тем, как создавать массивы, а не матрицы.

Я прочитал этот урок о том, как создавать разреженные матрицы: http://www.scipy.org/SciPy_Tutorial#head-c60163f2fd2bab79edd94be43682414f18b90df7

Чтобы подражать массиву, я просто создаю матрицу 1xN, но, как вы можете догадаться, Asp.dot(Bsp) работает не так, потому что вы не можете умножить две матрицы 1xN. Мне пришлось бы перенести каждый массив на Nx1, и это было бы довольно хромым, так как я буду делать это для каждого вычисления точечного продукта.

Затем я попытался создать матрицу NxN, где столбец 1 == строка 1 (чтобы вы могли умножить две матрицы и просто взять верхний левый угол в качестве точечного продукта), но это оказалось действительно неэффективным.

Я бы хотел использовать scipy sparse package как волшебную замену для numpy array(), но пока я не уверен, что делать.

Любые советы?

Ответы

Ответ 1

Используйте формат scipy.sparse, основанный на строках или столбцах: csc_matrix и csr_matrix.

Они используют эффективные реализации C под капотом (включая умножение), а транспонирование - это не-op (esp., если вы вызываете transpose(copy=False)), точно так же, как с массивами numpy.

EDIT: некоторые тайминги через ipython:

import numpy, scipy.sparse
n = 100000
x = (numpy.random.rand(n) * 2).astype(int).astype(float) # 50% sparse vector
x_csr = scipy.sparse.csr_matrix(x)
x_dok = scipy.sparse.dok_matrix(x.reshape(x_csr.shape))

Теперь x_csr и x_dok на 50% разрежены:

print repr(x_csr)
<1x100000 sparse matrix of type '<type 'numpy.float64'>'
        with 49757 stored elements in Compressed Sparse Row format>

И тайминги:

timeit numpy.dot(x, x)
10000 loops, best of 3: 123 us per loop

timeit x_dok * x_dok.T
1 loops, best of 3: 1.73 s per loop

timeit x_csr.multiply(x_csr).sum()
1000 loops, best of 3: 1.64 ms per loop

timeit x_csr * x_csr.T
100 loops, best of 3: 3.62 ms per loop

Итак, похоже, что я солгал. Транспонирование очень дешево, но нет эффективной реализации C csr * csc (в последнем scipy 0.9.0). В каждом вызове создается новый объект csr: - (

В качестве взлома (хотя в наши дни scipy относительно стабилен), вы можете делать точечный продукт непосредственно на разреженных данных:

timeit numpy.dot(x_csr.data, x_csr.data)
10000 loops, best of 3: 62.9 us per loop

Обратите внимание, что этот последний подход повторяет многократное плотное умножение. Разреженность составляет 50%, поэтому она быстрее, чем dot(x, x), в 2 раза.

Ответ 2

Вы можете создать подкласс одного из существующих 2d разреженных массивов

from scipy.sparse import dok_matrix

class sparse1d(dok_matrix):
    def __init__(self, v):
        dok_matrix.__init__(self, (v,))
    def dot(self, other):
        return dok_matrix.dot(self, other.transpose())[0,0]

a=sparse1d((1,2,3))
b=sparse1d((4,5,6))
print a.dot(b)

Ответ 3

Я не уверен, что это действительно намного лучше или быстрее, но вы можете сделать это, чтобы избежать использования транспонирования:

Asp.multiply(Bsp).sum()

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

Однако, возможно, проще их трансформировать:

Asp*Bsp.T

Это не так много, но вы также можете сделать подкласс и изменить метод mul().