Vectorização da operação de índice para uma matriz Scipy.Sparse

Question

O código a seguir funciona muito devagar, mesmo que tudo pareça ser vetorizado.

from numpy import *
from scipy.sparse import *

n = 100000;
i = xrange(n); j = xrange(n);
data = ones(n);

A=csr_matrix((data,(i,j)));

x = A[i,j]

O problema parece ser que a operação de indexação é implementada como uma função python e invocando A[i,j] resulta na seguinte saída de perfil

         500033 function calls in 8.718 CPU seconds

   Ordered by: internal time

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
   100000    7.933    0.000    8.156    0.000 csr.py:265(_get_single_element)
        1    0.271    0.271    8.705    8.705 csr.py:177(__getitem__)
(...)

Ou seja, a função python _get_single_element é chamado 100000 vezes, o que é realmente ineficiente. Por que isso não é implementado no puro C? Alguém conhece uma maneira de contornar essa limitação e acelerar o código acima? Devo estar usando um tipo de matriz esparsa diferente?

Answer 1

Você pode usar A.diagonal() Para recuperar a diagonal muito mais rapidamente (0,0009 segundos vs. 3,8 segundos na minha máquina). No entanto, se você deseja fazer indexação arbitária, essa é uma pergunta mais complicada, porque você não está usando fatias tanto quanto uma lista de índices. A função _Get_Single_Element está sendo chamada 100000 vezes porque apenas iterando os iteradores (i e j) que você passou para ela. Uma fatia seria uma [30: 60,10] ou algo semelhante a isso.

Além disso, eu usaria csr_matrix(eye(n,n)) Para fazer a mesma matriz que você fez com iteradores apenas por simplicidade.

Atualizar:

OK, como sua pergunta é realmente sobre o acesso a muitos elementos aleatórios rapidamente, responderei suas perguntas da melhor maneira possível.

• Por que isso não é implementado no puro C?

A resposta é simples: ninguém chegou a ela. Ainda há muito trabalho a ser feito na área de módulos de matriz esparsa do Scipy pelo que vi. Uma parte que é implementada em C é as conversões entre diferentes formatos de matriz esparsa.

• Alguém conhece uma maneira de contornar essa limitação e acelerar o código acima?

Você pode tentar realmente mergulhar nos módulos de matriz esparsa e tentar acelerar. Fiz isso e consegui reduzir o tempo para menos de um terço do original ao experimentar seu código acima para acessos aleatórios usando matrizes de RSE. Eu tive que acessar diretamente _get_single_element e reduzir significativamente o código para fazer isso, incluindo a retirada de cheques vinculados.

No entanto, ainda era mais rápido usar uma lil_matrix (embora mais lenta para inicializar a matriz), mas eu tive que fazer o acesso com uma compreensão da lista porque as matrizes LIL não estão configuradas para o tipo de indexação que você está fazendo. O uso de uma compreensão da lista para o CSR_Matrix ainda sai do método Lil Matrix muito à frente. Por fim, a matriz LIL é mais rápida para acessar elementos aleatórios porque não é compactada.

Usando o LIL_Matrix em seu formulário original é executado em cerca de um quinto da época do código que você listou acima. Se eu fizer alguns cheques vinculados e ligar diretamente ao método _get1 () da LIL_MATRIX, posso reduzir o tempo mais de 7% da hora original. Para maior clareza, é uma aceleração de 3,4 a 3,8 segundos para cerca de 0,261 segundos.

Por fim, tentei fazer minha própria função que acesse diretamente os dados da matriz LIL e evite as chamadas de função repetidas. O tempo para isso foi de cerca de 0,136 segundos. Isso não aproveitou os dados que estão sendo classificados, que é outra otimização potencial (em particular se você estiver acessando muitos elementos que estão nas mesmas linhas).

Se você quiser mais rápido do que isso, terá que escrever sua própria implementação de matriz escassa de código C, provavelmente.

• Devo estar usando um tipo de matriz esparsa diferente?

Bem, sugiro a matriz Lil se sua intenção é acessar muitos elementos, mas tudo depende do que você precisa fazer. Você também precisa multiplicar matrizes, por exemplo? Lembre -se de que a mudança entre matrizes pode pelo menos às vezes (em determinadas circunstâncias) ser bastante rápida; portanto, não descarte a mudança para um formato de matriz diferente para fazer operações diferentes.

Se você não precisar fazer nenhuma operação de álgebra em sua matriz, talvez você deva apenas usar um padrão ou algo semelhante. O perigo com o DefaultDicts é que, sempre que um elemento é solicitado que não esteja no ditado, ele define esse item para o padrão e o armazena para que possa ser problemático.

Answer 2

Eu acho que _get_single_element é chamado apenas quando o dtype padrão de 'objeto' é usado. Você já tentou fornecer um dtype, como csr_matrix((data, (i,j)), dtype=int32)