2015-06-27 13 views
8

Per ciascun elemento di una serie casuale di indici 2D (con potenziali duplicati), voglio "+ = 1" sulla griglia corrispondente in un array zero 2D. Tuttavia, non so come ottimizzare il calcolo. Utilizzando lo standard per ciclo, come mostrato qui,Vectorize aggiunta iterativa negli array NumPy

def interadd(): 
    U = 100 
    input = np.random.random(size=(5000,2)) * U 
    idx = np.floor(input).astype(np.int) 

    grids = np.zeros((U,U))  
    for i in range(len(input)): 
     grids[idx[i,0],idx[i,1]] += 1 
    return grids 

Il runtime può essere abbastanza significativo:

>> timeit(interadd, number=5000) 
43.69953393936157 

C'è un modo per vectorize questo processo iterativo?

risposta

8

Si potrebbe accelerarlo un po 'utilizzando np.add.at, che gestisce correttamente il caso degli indici duplicati:

def interadd(U, idx): 
    grids = np.zeros((U,U))  
    for i in range(len(idx)): 
     grids[idx[i,0],idx[i,1]] += 1 
    return grids 

def interadd2(U, idx): 
    grids = np.zeros((U,U)) 
    np.add.at(grids, idx.T.tolist(), 1) 
    return grids 

def interadd3(U, idx): 
    # YXD suggestion 
    grids = np.zeros((U,U)) 
    np.add.at(grids, (idx[:,0], idx[:,1]), 1) 
    return grids 

che dà

>>> U = 100 
>>> idx = np.floor(np.random.random(size=(5000,2))*U).astype(np.int) 
>>> (interadd(U, idx) == interadd2(U, idx)).all() 
True 
>>> %timeit interadd(U, idx) 
100 loops, best of 3: 8.48 ms per loop 
>>> %timeit interadd2(U, idx) 
100 loops, best of 3: 2.62 ms per loop 

E suggerimento di YXD:

>>> (interadd(U, idx) == interadd3(U, idx)).all() 
True 
>>> %timeit interadd3(U, idx) 
1000 loops, best of 3: 1.09 ms per loop 
+4

Stava scrivendo quasi la stessa thing.You può cambiare 'IDX. T.tolist() 'a' (idx [:, 0], idx [:, 1]) ', dovrebbe essere ancora più veloce. – YXD

+1

(errore di battitura appena corretto nel commento sopra) – YXD

5

È possibile convertire gli indici R,C da idx a indici lineari, quindi individuare quelli unici insieme ai loro conteggi e infine memorizzarli nell'output grids come output finale. Ecco l'implementazione per ottenere gli stessi -

# Calculate linear indices corressponding to idx 
lin_idx = idx[:,0]*U + idx[:,1] 

# Get unique linear indices and their counts 
unq_lin_idx,idx_counts = np.unique(lin_idx,return_counts=True) 

# Setup output array and store index counts in raveled/flattened version 
grids = np.zeros((U,U)) 
grids.ravel()[unq_lin_idx] = idx_counts 

test Runtime -

Qui ci sono le prove di runtime che coprono tutti gli approcci (tra cui @DSM's approaches) ed utilizzando le stesse definizioni come elencato in questa soluzione -

In [63]: U = 100 
    ...: idx = np.floor(np.random.random(size=(5000,2))*U).astype(np.int) 
    ...: 

In [64]: %timeit interadd(U, idx) 
100 loops, best of 3: 7.57 ms per loop 

In [65]: %timeit interadd2(U, idx) 
100 loops, best of 3: 2.59 ms per loop 

In [66]: %timeit interadd3(U, idx) 
1000 loops, best of 3: 1.24 ms per loop 

In [67]: def unique_counts(U, idx): 
    ...:  lin_idx = idx[:,0]*U + idx[:,1] 
    ...:  unq_lin_idx,idx_counts = np.unique(lin_idx,return_counts=True) 
    ...:  grids = np.zeros((U,U)) 
    ...:  grids.ravel()[unq_lin_idx] = idx_counts 
    ...:  return grids 
    ...: 

In [68]: %timeit unique_counts(U, idx) 
1000 loops, best of 3: 595 µs per loop 

I runtime suggeriscono che l'approccio basato su np.unique proposto è più del 50% più veloce del secondo approccio più veloce.

+0

'np.unique' usa l'ordinamento sotto il cofano, quindi ha una complessità temporale peggiore di' np.add.at', ma d'altra parte il tuo approccio ha un pattern di accesso alla memoria più veloce per matrice 'grids'. –

+0

@moarningsun Sì, penso che usi 'sorting' e' differentiation' sotto il cofano. Questo ha un senso, credo, su quei runtime più veloci. Sarebbe interessante scoprire cosa sta succedendo con 'add.at'. – Divakar

+0

Questo mi ha fatto pensare ad un approccio interessante: 'grids = np.bincount (lin_idx, minlength = U ** 2) .reshape (U, U)' –

6

risposta di Divakar mi portano a provare la seguente, che sembra essere il metodo più veloce ancora:

lin_idx = idx[:,0]*U + idx[:,1] 
grids = np.bincount(lin_idx, minlength=U**2).reshape(U, U) 

Tempi:

In [184]: U = 100 
    ...: input = np.random.random(size=(5000,2)) * U 
    ...: idx = np.floor(input).astype(np.int) 

In [185]: %timeit interadd3(U, idx) # By DSM/XYD 
1000 loops, best of 3: 1.68 ms per loop 

In [186]: %timeit unique_counts(U, idx) # By Divakar 
1000 loops, best of 3: 676 µs per loop 

In [187]: %%timeit 
    ...: lin_idx = idx[:,0]*U + idx[:,1] 
    ...: grids = np.bincount(lin_idx, minlength=U*U).reshape(U, U) 
    ...: 
10000 loops, best of 3: 97.5 µs per loop 
+0

Un enorme miglioramento sembra! – Divakar