2011-02-07 12 views
32

que tiene una matriz de esta manera:eficiente gama de construcción Numpy 2D de la matriz 1D

A = array([1,2,3,4,5,6,7,8,9,10]) 

y estoy tratando de conseguir un arreglo de esta manera:

B = array([[1,2,3], 
      [2,3,4], 
      [3,4,5], 
      [4,5,6]]) 

donde cada fila (de un fijo ancho arbitrario) se desplaza en uno. La matriz de A tiene 10k registros y estoy tratando de encontrar una manera eficiente de hacerlo en Numpy. Actualmente estoy usando vstack y un bucle for que es lento. ¿Hay una manera mas rápida?

Editar:

width = 3 # fixed arbitrary width 
length = 10000 # length of A which I wish to use 
B = A[0:length + 1] 
for i in range (1, length): 
    B = np.vstack((B, A[i, i + width + 1])) 
+0

¿Puedes publicar tu solución vstack/loop? – eumiro

+0

@wxbx: ¿Por favor, elabore más lo que pretende hacer? Tenga en cuenta que 'B = array ([1,2,3], [2,3,4], [3,4,5], [4,5,6])' no es válido 'numpy'! – eat

+0

@wxbx: su solución es realmente desafortunada. Usted 'vstack' la matriz 10000 veces! Ver mi respuesta, lo "apilaré" una sola vez. – eumiro

Respuesta

47

En realidad, hay una forma aún más eficiente de hacerlo ... La desventaja de usar vstack, etc., es que está haciendo una copia de la matriz.

Por cierto, este es efectivamente idéntico al @ respuesta de Pablo, pero Quiero poner esto sólo para explicar las cosas en un poco más de detalle ...

Hay una manera de hacer esto con sólo puntos de vista de manera que no memoria está duplicada.

Estoy tomando prestado directamente de Erik Rigtorp's post to numpy-discussion, que a su vez, lo prestó de Bottleneck de Keith Goodman (¡lo cual es bastante útil!).

El truco básico es que manipular directamente (arrays para unidimensionales) strides of the array:

import numpy as np 

def rolling(a, window): 
    shape = (a.size - window + 1, window) 
    strides = (a.itemsize, a.itemsize) 
    return np.lib.stride_tricks.as_strided(a, shape=shape, strides=strides) 

a = np.arange(10) 
print rolling(a, 3) 

Dónde a es su matriz de entrada y window es la longitud de la ventana que desea (3, en su caso).

Este rendimientos:

[[0 1 2] 
[1 2 3] 
[2 3 4] 
[3 4 5] 
[4 5 6] 
[5 6 7] 
[6 7 8] 
[7 8 9]] 

Sin embargo, no hay absolutamente ninguna duplicación de memoria entre el original a y la matriz devuelta. Esto significa que es rápido y escala mucho mejor que otras opciones.

Por ejemplo (usando a = np.arange(100000) y window=3):

%timeit np.vstack([a[i:i-window] for i in xrange(window)]).T 
1000 loops, best of 3: 256 us per loop 

%timeit rolling(a, window) 
100000 loops, best of 3: 12 us per loop 

Si generalizamos esta a una "ventana móvil" a lo largo del último eje para una matriz de N dimensiones, obtenemos la función "ventana móvil" de Erik Rigtorp :

import numpy as np 

def rolling_window(a, window): 
    """ 
    Make an ndarray with a rolling window of the last dimension 

    Parameters 
    ---------- 
    a : array_like 
     Array to add rolling window to 
    window : int 
     Size of rolling window 

    Returns 
    ------- 
    Array that is a view of the original array with a added dimension 
    of size w. 

    Examples 
    -------- 
    >>> x=np.arange(10).reshape((2,5)) 
    >>> rolling_window(x, 3) 
    array([[[0, 1, 2], [1, 2, 3], [2, 3, 4]], 
      [[5, 6, 7], [6, 7, 8], [7, 8, 9]]]) 

    Calculate rolling mean of last dimension: 
    >>> np.mean(rolling_window(x, 3), -1) 
    array([[ 1., 2., 3.], 
      [ 6., 7., 8.]]) 

    """ 
    if window < 1: 
     raise ValueError, "`window` must be at least 1." 
    if window > a.shape[-1]: 
     raise ValueError, "`window` is too long." 
    shape = a.shape[:-1] + (a.shape[-1] - window + 1, window) 
    strides = a.strides + (a.strides[-1],) 
    return np.lib.stride_tricks.as_strided(a, shape=shape, strides=strides) 

lo tanto, vamos a ver en lo que está pasando aquí ... la manipulación puede parecer un poco mágica de strides una matriz, pero una vez que entienda lo que está pasando, no es en absoluto. Los pasos de una matriz numpy describen el tamaño en bytes de los pasos que se deben seguir para incrementar un valor a lo largo de un eje dado. Por lo tanto, en el caso de una matriz unidimensional de flotantes de 64 bits, la longitud de cada elemento es de 8 bytes, y x.strides es (8,).

x = np.arange(9) 
print x.strides 

Ahora bien, si formar de nuevo este en una 2D, matriz de 3x3, los pasos serán (3 * 8, 8), ya que habría que saltar 24 bytes para incrementar un paso a lo largo del primer eje, y 8 bytes para incrementar un paso a lo largo del segundo eje.

y = x.reshape(3,3) 
print y.strides 

mismo modo una transposición es lo mismo que simplemente invirtiendo los pasos de una matriz:

print y 
y.strides = y.strides[::-1] 
print y 

Claramente, los pasos de una matriz y la forma de una matriz están íntimamente ligados. Si cambiamos uno, tenemos que cambiar el otro en consecuencia, de lo contrario no tendremos una descripción válida del búfer de memoria que realmente contiene los valores de la matriz.

Por lo tanto, si desea cambiar tanto la forma y el tamaño de una matriz al mismo tiempo, no puede hacerlo sólo mediante el establecimiento de x.strides y x.shape, incluso si los nuevos avances y la forma son compatibles.

Aquí es donde entra numpy.lib.as_strided. En realidad, es una función muy simple que simplemente establece los pasos y la forma de una matriz al mismo tiempo.

Comprueba que los dos son compatibles, pero no que los pasos anteriores y la nueva forma sean compatibles, como sucedería si establece los dos de forma independiente. (De hecho, lo hace a través de numpy's __array_interface__, que permite que clases arbitrarias describan un búfer de memoria como una matriz numpy.)

Así que, todo lo que hemos hecho se hace de modo que pasos un artículo adelante (8 bytes en el caso de una matriz de 64 bits) a lo largo de un eje, pero también solo pasos 8 bytes hacia adelante a lo largo del otro eje.

En otras palabras, en caso de un tamaño de "ventana" de 3, la matriz tiene una forma de (whatever, 3), pero en lugar de pisar un completo 3 * x.itemsize para la segunda dimensión, se sólo unos pasos un elemento hacia adelante, haciendo efectivamente las filas de la nueva matriz una vista de "ventana móvil" en la matriz original.

(Esto también significa que x.shape[0] * x.shape[1] no será el mismo que x.size para su nueva matriz.)

En cualquier caso, es de esperar que hace las cosas un poco más claro ..

+0

Kinggton: Realmente admiro tu respuesta, pero ¿no crees que es demasiado exagerado para la pregunta de OP? ;-). Gracias – eat

+4

@eat - ¡Lo es! :) Definitivamente es excesivo para una matriz corta (y la matriz de 10K elementos de OP es bastante corta), pero aún es útil saberlo. Honestamente, creo que a veces me gusta escribir respuestas demasiado largas ... –

+1

Kingston: Gracias por una respuesta realmente detallada, aprendí mucho allí. ¡También guardé tu código contra la respuesta de @ eumiro y tu respuesta continua me dio 60 veces más de velocidad! Teniendo en cuenta que planeo usar esto en una matriz mucho más grande, la aceleración es increíblemente útil. :) – wxbx

2

¿Qué método se utiliza?

import numpy as np 
A = np.array([1,2,3,4,5,6,7,8,9,10]) 
width = 3 

np.vstack([A[i:i-len(A)+width] for i in xrange(len(A)-width)]) 
# needs 26.3µs 

np.vstack([A[i:i-width] for i in xrange(width)]).T 
# needs 13.2µs 

Si su anchura es relativamente baja (3) y tiene un gran A (10.000 elementos), entonces la diferencia es aún más importante: 32.4ms para la primera y 44μs para la segunda.

+0

gracias! esto es justo lo que necesitaba! y sí acabo de romper numpy hoy tan lentamente aprendiendo. – wxbx

1

Creo que esto podría ser más rápido que un bucle, cuando la anchura se fija en un número bajo ...

import numpy 
a = numpy.array([1,2,3,4,5,6]) 
b = numpy.reshape(a, (numpy.shape(a)[0],1)) 
b = numpy.concatenate((b, numpy.roll(b,-1,0), numpy.roll(b,-2,0)), 1) 
b = b[0:(numpy.shape(a)[0]/2) + 1,:] 

EDITAR Claramente, las soluciones utilizando pasos son superiores a este, con la única desventaja importante que aún no están bien documentados ...

9

Esta solución no se implementa de manera eficiente mediante un bucle python, ya que viene con todo tipo de comprobación de tipos que se evita mejor cuando se trabaja con matrices numpy. Si la matriz es excepcionalmente alto, usted notará una velocidad grande con esto:

newshape = (4,3) 
newstrides = (A.itemsize, A.itemsize) 
B = numpy.lib.stride_tricks.as_strided(A, shape=newshape, strides=newstrides) 

Esto da una visión de la matriz A. Si desea una nueva matriz puede editar, hacer lo mismo pero con .copy() al final.

detalles sobre pasos:

La tupla newstrides en este caso será (4,4), ya que la matriz tiene artículos 4 bytes y desea continuar con el paso a través de sus datos en los pasos de un único elemento en el i-dimension. El segundo valor '4' se refiere a las zancadas en la dimensión j (en una matriz 4x4 normal sería 16). Porque en este caso también desea incrementar su lectura desde el búfer en pasos de 4 bytes en la dimensión j.

Joe da una descripción agradable y detallada y deja las cosas claras cuando dice que todo lo que hace este truco es cambiar los pasos y la forma de forma simultánea.

+1

+1 ¡Me ganaste! Estaba en el medio de escribir esto ... Aún publicaré mi respuesta, ya que entra un poco más de detalle. Además, su 'strides = (4,4)' supone que 'A.itemsize' es 4 (es decir, flotantes o enteros de 32 bits).Lo mejor es hacer 'strides = (A.itemsize, A.itemsize)'. –

+0

¿Puede indicarme los documentos para esto? Nunca he visto esta función antes ... – Benjamin

+0

Gracias Joe. ¡Estaba buscando una documentación en línea para vincular pero no mucho por ahí! Esto fue lo mejor que pude encontrar: http://mentat.za.net/numpy/numpy_advanced_slides/ – Paul

2

sólo para ir más allá con la respuesta de @ Joe general

import numpy as np 
def rolling(a, window): 
    step = 2 
    shape = ((a.size-window)/step + 1 , window) 


    strides = (a.itemsize*step, a.itemsize) 

    return np.lib.stride_tricks.as_strided(a, shape=shape, strides=strides) 

a = np.arange(10) 

print rolling(a, 3) 

que emite:

[[0 1 2] 
[2 3 4] 
[4 5 6] 
[6 7 8]] 

Para generalizar más para el caso 2D, es decir, lo utilizan para la extracción de parche desde una imagen

def rolling2d(a,win_h,win_w,step_h,step_w): 

    h,w = a.shape 
    shape = (((h-win_h)/step_h + 1) * ((w-win_w)/step_w + 1) , win_h , win_w) 

    strides = (step_w*a.itemsize, h*a.itemsize,a.itemsize) 


    return np.lib.stride_tricks.as_strided(a, shape=shape, strides=strides) 

a = np.arange(36).reshape(6,6) 
print a 
print rolling2d (a,3,3,2,2) 

que da salida:

[[ 0 1 2 3 4 5] 
[ 6 7 8 9 10 11] 
[12 13 14 15 16 17] 
[18 19 20 21 22 23] 
[24 25 26 27 28 29] 
[30 31 32 33 34 35]] 
[[[ 0 1 2] 
    [ 6 7 8] 
    [12 13 14]] 

[[ 2 3 4] 
    [ 8 9 10] 
    [14 15 16]] 

[[ 4 5 6] 
    [10 11 12] 
    [16 17 18]] 

[[ 6 7 8] 
    [12 13 14] 
    [18 19 20]]] 
+1

Es posible en el ejemplo anterior no obtener resultados que envuelvan el borde derecho de la matriz original. por ejemplo, la tercera salida '[4,5,6; 10,11,12; 16,17,18] ' 'vuelve a enrollarse. Para el procesamiento de imágenes me gustaría evitar esto y simplemente saltar al próximo resultado devuelto. –

0

estoy usando una función más generalizada similar a la de @JustInTime pero aplicable a ndarray

def sliding_window(x, size, overlap=0): 
    step = size - overlap # in npts 
    nwin = (x.shape[-1]-size)//step + 1 
    shape = x.shape[:-1] + (nwin, size) 
    strides = x.strides[:-1] + (step*x.strides[-1], x.strides[-1]) 
    return stride_tricks.as_strided(x, shape=shape, strides=strides) 

Un ejemplo,

x = np.arange(10) 
M.sliding_window(x, 5, 3) 
Out[1]: 
array([[0, 1, 2, 3, 4], 
     [2, 3, 4, 5, 6], 
     [4, 5, 6, 7, 8]]) 


x = np.arange(10).reshape((2,5)) 
M.sliding_window(x, 3, 1) 
Out[2]: 
array([[[0, 1, 2], 
     [2, 3, 4]], 

     [[5, 6, 7], 
     [7, 8, 9]]]) 
1

Eche un vistazo a: view_as_windows.

import numpy as np 
from skimage.util.shape import view_as_windows 
window_shape = (4,) 
aa = np.arange(1000000000) # 1 billion 
bb = view_as_windows(aa, window_shape) 

Alrededor de 1 segundo.

Cuestiones relacionadas