2013-05-23 33 views
11

Sto implementando un algoritmo che richiede di guardare sottomodelli consecutivi non sovrapposti all'interno di un array numpy (strettamente bidimensionale). ad esempio, per il 12 da 12Indici di sottomatrici di dimensioni fisse di array numpy

>>> a = np.random.randint(20, size=(12, 12)); a 
array([[ 4, 0, 12, 14, 3, 8, 14, 12, 11, 18, 6, 6], 
     [15, 13, 2, 18, 15, 15, 16, 2, 9, 16, 6, 4], 
     [18, 18, 3, 8, 1, 15, 14, 13, 13, 13, 7, 0], 
     [ 1, 9, 3, 6, 0, 4, 3, 15, 0, 9, 11, 12], 
     [ 5, 15, 5, 6, 4, 4, 18, 13, 10, 17, 11, 8], 
     [13, 17, 8, 15, 17, 12, 7, 1, 13, 15, 0, 18], 
     [ 2, 1, 11, 12, 3, 16, 11, 9, 10, 15, 4, 16], 
     [19, 11, 10, 7, 10, 19, 7, 13, 11, 9, 17, 8], 
     [14, 14, 17, 0, 0, 0, 11, 1, 10, 14, 2, 7], 
     [ 6, 15, 6, 7, 15, 19, 2, 4, 6, 16, 0, 3], 
     [ 5, 10, 7, 5, 0, 8, 5, 8, 9, 14, 4, 3], 
     [17, 2, 0, 3, 15, 10, 14, 1, 0, 7, 16, 2]]) 

e guardando sottomatrici 3x3, vorrei la prima sottomatrice 3x3 ad essere dall'angolo in alto a sinistra:

>>> a[0:3, 0:3] 
array([[ 4, 0, 12], 
     [15, 13, 2], 
     [18, 18, 3]]) 

Il prossimo lungo da dare da a[0:3, 3:6] e così via. Non importa se l'ultima serie di indici di ciascuna riga o colonna viene eseguita alla fine dell'array - il comportamento di numpy di dare semplicemente la porzione all'interno della slice esistente è sufficiente.

Desidero un modo per generare questi indici di sezione in modo programmatico per matrici e sottomatrici di dimensioni arbitrarie. Attualmente ho questo:

size = 3 
x_max = a.shape[0] 
xcoords = range(0, x_max, size) 
xcoords = zip(xcoords, xcoords[1:]) 

e similmente per generare y_coords, in modo che la serie di indici è dato da itertools.product(xcoords, ycoords).

La mia domanda è: c'è un modo più diretto per farlo, forse usando numpy.mgrid o qualche altra tecnica numpy?

risposta

6

Come gli indici

Ecco un modo rapido per ottenere una specifica size x size blocco:

base = np.arange(size) # Just the base set of indexes 
row = 1    # Which block you want 
col = 0     
block = a[base[:, np.newaxis] + row * size, base + col * size] 

Se si voleva si poteva costruire matrici simili al vostro xcoords come:

y, x = np.mgrid[0:a.shape[0]/size, 0:a.shape[1]/size] 
y_coords = y[..., np.newaxis] * size + base 
x_coords = x[..., np.newaxis] * size + base 

Quindi è possibile accedere a un blocco come questo:

block = a[y_coords[row, col][:, np.newaxis], x_coords[row, col]] 

Come i blocchi direttamente

Se si desidera solo per ottenere i blocchi (e non gli indici delle voci blocco), userei np.split (due volte):

blocks = map(lambda x : np.split(x, a.shape[1]/size, 1), # Split the columns 
         np.split(a, a.shape[0]/size, 0)) # Split the rows 

allora avete una lista 2D di size x size blocchi:

>>> blocks[0][0] 
array([[ 4, 0, 12], 
     [15, 13, 2], 
     [18, 18, 3]]) 

>>> blocks[1][0] 
array([[ 1, 9, 3], 
     [ 5, 15, 5], 
     [13, 17, 8]]) 

si potrebbe quindi rendere questo un allineamento NumPy e utilizzare lo stesso stile di indicizzazione come sopra:

>>> blocks = np.array(blocks) 
>>> blocks.shape 
(4, 4, 3, 3) 
+1

FYI, ho appena controllato e [la risposta di Saulio] (http://stackoverflow.com/a/16715845/310165) è circa 3-4 volte più veloce del mio metodo 'map', se tutto ciò di cui hai bisogno sono i blocchi . – Geoff

5

È possibile utilizzare l'one-liner:

r = 3 
c = 3 
lenr = a.shape[0]/r 
lenc = a.shape[1]/c 
np.array([a[i*r:(i+1)*r,j*c:(j+1)*c] for (i,j) in np.ndindex(lenr,lenc)]).reshape(lenr,lenc,r,c) 
+2

+1 - Non sapevo di 'ndindex'. Sembra abbastanza utile. – Geoff

2

Sto aggiungendo questa risposta ad una vecchia questione in quanto una modifica ha ha risposto a questa domanda. Ecco un modo alternativo per calcolare i blocchi:

size = 3 
lenr, lenc = int(a.shape[0]/size), int(a.shape[1]/size) 

t = a.reshape(lenr,size,lenc,size).transpose(0, 2, 1, 3) 

Il profilo mostra che questo è il più veloce. Profilazione fatta con python 3.5, e i risultati della mappa passati a array() per la compatibilità, dal momento che in 3.5 mappa restituisce un iteratore.

reshape/transpose: 643 ns per loop 
reshape/index:  45.8 µs per loop 
Map/split:   10.3 µs per loop 

È interessante notare che la versione iteratore della mappa è più veloce. In ogni caso, l'uso di reshape e transpose è più veloce.

+0

Questo è senza dubbio la strada da percorrere per questo problema; nessuna domanda al riguardo –

Problemi correlati