intesa strano comportamento indicizzazione 2d-array booleano in NumPy

Question

Perché questo lavoro:

a=np.random.rand(10,20) 
x_range=np.arange(10) 
y_range=np.arange(20) 

a_tmp=a[x_range<5,:] 
b=a_tmp[:,np.in1d(y_range,[3,4,8])] 

e questo non lo fa:

a=np.random.rand(10,20) 
x_range=np.arange(10) 
y_range=np.arange(20)  

b=a[x_range<5,np.in1d(y_range,[3,4,8])] 

Answer 1

del page on indexing documentazione di riferimento Numpy contiene le risposte, ma richiede un po ' di lettura attenta.

La risposta è che l'indicizzazione con valori booleani è equivalente all'indicizzazione con matrici di interi ottenute trasformando prima gli array booleani con np.nonzero. Pertanto, con array booleani m1, m2

a[m1, m2] == a[m1.nonzero(), m2.nonzero()] 

che (quando riesce, cioè, m1.nonzero().shape == m2.nonzero().shape) è equivalente a:

[a[i, i] for i in range(a.shape[0]) if m1[i] and m2[i]] 

Io non sono sicuro perché è stato progettato per funzionare in questo modo - - Di solito, questo è non quello che vuoi.

per ottenere il risultato più intuitivo, si può invece fare

a[np.ix_(m1, m2)] 

che produce un risultato equivalente a

[[a[i,j] for j in range(a.shape[1]) if m2[j]] for i in range(a.shape[0]) if m1[i]] 

Answer 2

Un'alternativa a np.ix_ è quello di convertire gli array booleani in numero intero array (usando np.nonzero()), quindi utilizzare np.newaxis per creare matrici della forma corretta per sfruttare la trasmissione.

import numpy as np 

a=np.random.rand(10,20) 
x_range=np.arange(10) 
y_range=np.arange(20) 

a_tmp=a[x_range<5,:] 
b_correct=a_tmp[:,np.in1d(y_range,[3,4,8])] 

m1=(x_range<5).nonzero()[0] 
m2=np.in1d(y_range,[3,4,8]).nonzero() 
b=a[m1[:,np.newaxis], m2] 
assert np.allclose(b,b_correct) 

b2=a[np.ix_(x_range<5,np.in1d(y_range,[3,4,8]))] 
assert np.allclose(b2,b_correct) 

np.ix_ tende ad essere più lento del doppio indicizzazione. La soluzione a lungo modulo sembra essere un po 'più veloce:

long-form:

In [83]: %timeit a[(x_range<5).nonzero()[0][:,np.newaxis], (np.in1d(y_range,[3,4,8])).nonzero()[0]] 
10000 loops, best of 3: 131 us per loop 

doppia indicizzazione:

In [85]: %timeit a[x_range<5,:][:,np.in1d(y_range,[3,4,8])] 
10000 loops, best of 3: 144 us per loop 

utilizzando np.ix_:

In [84]: %timeit a[np.ix_(x_range<5,np.in1d(y_range,[3,4,8]))] 
10000 loops, best of 3: 160 us per loop 

Nota: sarebbe una buona idea testare questi tempi sul dispositivo poiché le classifiche potrebbero cambiare in base alla versione di Python, numpy o hardware.