Scipy sparse ... matrices?

Question

Por lo tanto, estoy haciendo una clasificación de Kmeans usando matrices numpy que son bastante dispersas-- montones y montones de ceros. Pensé que usaría el paquete 'disperso' de scipy para reducir la sobrecarga de almacenamiento, pero estoy un poco confundido acerca de cómo crear matrices, no matrices.

He pasado por este tutorial sobre cómo crear matrices dispersas: http://www.scipy.org/SciPy_Tutorial#head-c60163f2fd2bab79edd94be43682414f18b90df7

para imitar una matriz, que simplemente creo una matriz 1xN, pero como pueden suponer, Asp.dot (BSP) doesn' Trabaja bastante porque no puedes multiplicar dos matrices 1xN. Tendría que transponer cada arreglo a Nx1, y eso es bastante cojo, ya que lo haría por cada cálculo de producto de punto.

A continuación, traté de crear una matriz NxN donde la columna 1 == fila 1 (de forma que se pueden multiplicar dos matrices y simplemente tomar la esquina superior izquierda como el producto escalar), pero eso resultó ser realmente ineficiente.

Me encantaría usar el paquete escaso de scipy como un reemplazo mágico para la matriz de numpy(), pero aún no estoy seguro de qué hacer.

¿Algún consejo?

Answer 1

Utilice un formato scipy.sparse basado en filas o columnas: csc_matrix y csr_matrix.

Utilizan implementaciones C eficientes (incluida la multiplicación), y la transposición no es operativa (especialmente si llama al transpose(copy=False)), al igual que con las matrices numpy.

EDIT: algunos tiempos a través de 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)) 

Ahora x_csr y x_dok son 50% escasa:

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

y los tiempos:

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 

lo que parece que me dijo una mentira. La transposición es muy barata, pero no hay una implementación C eficiente de csr * csc (en la última versión scipy 0.9.0).Un nuevo objeto RSE se construye en cada llamada :-(

Como un hack (aunque scipy es relativamente estable en estos días), se puede hacer el producto escalar directamente sobre la escasez de datos:

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

Nota: este El último enfoque hace una multiplicación nude denso de nuevo. La dispersión es 50%, por lo que es más rápido que dot(x, x) por un factor de 2.

Answer 2

Se puede crear una subclase de una de las matrices dispersas existentes en 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) 

Answer 3

no estoy seguro de que es realmente mucho mejor o más rápido, pero se puede hacer esto para evitar el uso de la transpuesta:

Asp.multiply(Bsp).sum() 

Esto solo toma el elemento por elemento producto de las dos matrices y suma los productos. Podría hacer una subclase de cualquier formato de matriz que esté utilizando que tenga la declaración anterior como producto de puntos.

Sin embargo, es probable que sea más fácil de transposición de ellas:

Asp*Bsp.T 

que no parece como mucho que hacer, pero también se puede hacer una subclase y modificar el método mul().