Crea array 2D Numpy dalle coordinate

Question

I punti dati rappresentano le coordinate di un array 2D (matrice). I punti sono regolarmente reticolati, ad eccezione dei punti dati mancanti in alcune posizioni della griglia.

Ad esempio, considerare alcuni dati XYZ che si adattano a una normale griglia di 0.1 con forma (3, 4). Ci sono lacune e punti mancanti, quindi ci sono 5 punti, e non 12:

import numpy as np 
X = np.array([0.4, 0.5, 0.4, 0.4, 0.7]) 
Y = np.array([1.0, 1.0, 1.1, 1.2, 1.2]) 
Z = np.array([3.3, 2.5, 3.6, 3.8, 1.8]) 
# Evaluate the regular grid dimension values 
Xr = np.linspace(X.min(), X.max(), np.round((X.max() - X.min())/np.diff(np.unique(X)).min()) + 1) 
Yr = np.linspace(Y.min(), Y.max(), np.round((Y.max() - Y.min())/np.diff(np.unique(Y)).min()) + 1) 
print('Xr={0}; Yr={1}'.format(Xr, Yr)) 
# Xr=[ 0.4 0.5 0.6 0.7]; Yr=[ 1. 1.1 1.2] 

cosa vorrei vedere è indicato in questa immagine (ambiti: nero index = base-0; grigio = valore della coordinata; colore = valore matrice; bianco = mancante).

Ecco quello che ho, che è intuitiva con un ciclo for:

ar = np.ma.array(np.zeros((len(Yr), len(Xr)), dtype=Z.dtype), mask=True) 
for x, y, z in zip(X, Y, Z): 
    j = (np.abs(Xr - x)).argmin() 
    i = (np.abs(Yr - y)).argmin() 
    ar[i, j] = z 
print(ar) 
# [[3.3 2.5 -- --] 
# [3.6 -- -- --] 
# [3.8 -- -- 1.8]]  

C'è un modo più NumPythonic del vectorising l'approccio per restituire un array 2D ar? O è necessario il ciclo for?

Answer 1

si può fare su una linea con np.histogram2d

data = np.histogram2d(Y, X, bins=[len(Yr),len(Xr)], weights=Z) 
print(data[0]) 
[[ 3.3 2.5 0. 0. ] 
[ 3.6 0. 0. 0. ] 
[ 3.8 0. 0. 1.8]] 

Answer 2

È possibile utilizzare un coo_matrix SciPy. Permette di costruire una matrice sparsa da coordinate e dati. Vedi esempi sul link allegato.

http://docs.scipy.org/doc/scipy-0.15.1/reference/generated/scipy.sparse.coo_matrix.html

Speranza che aiuta.

Answer 3

La matrice sparse è la prima soluzione che è venuto in mente, ma poiché X e Y sono galleggianti, è un po 'disordinato:

In [624]: I=((X-.4)*10).round().astype(int) 
In [625]: J=((Y-1)*10).round().astype(int) 
In [626]: I,J 
Out[626]: (array([0, 1, 0, 0, 3]), array([0, 0, 1, 2, 2])) 

In [627]: sparse.coo_matrix((Z,(J,I))).A 
Out[627]: 
array([[ 3.3, 2.5, 0. , 0. ], 
     [ 3.6, 0. , 0. , 0. ], 
     [ 3.8, 0. , 0. , 1.8]]) 

Occorre ancora, in un modo o nell'altro, per abbinare quelle coordinate con indici [0,1,2 ...]. Il mio trucco veloce era di scalare i valori in modo lineare. Anche così dovevo fare attenzione quando convertivo i float in ints.

sparse.coo_matrix opere perché un modo naturale di definire una matrice sparsa è con (i, j, data) tuple, che ovviamente può essere tradotto I, J, Data elenchi o matrici.

Mi piace la soluzione historgram, anche se non ho avuto occasione di usarlo.

Answer 4

È possibile utilizzare X e Y per creare coordinate X-Y su una griglia 0.1 distanziata estende dalla min to max of X e min to max of Y e quindi inserendo Z's in quelle posizioni specifiche. Questo eviterebbe di utilizzare linspace per ottenere Xr e Yr e come tale deve essere abbastanza efficiente. Ecco l'implementazione -

def indexing_based(X,Y,Z): 
    # Convert X's and Y's to indices on a 0.1 spaced grid 
    X_int = np.round((X*10)).astype(int) 
    Y_int = np.round((Y*10)).astype(int) 
    X_idx = X_int - X_int.min() 
    Y_idx = Y_int - Y_int.min() 

    # Setup output array and index it with X_idx & Y_idx to set those as Z 
    out = np.zeros((Y_idx.max()+1,X_idx.max()+1)) 
    out[Y_idx,X_idx] = Z 

    return out 

test Runtime -

Questa sezione confrontare l'approccio indexing-based contro l'altro np.histogram2d based solution per le prestazioni -

In [132]: # Create unique couples X-Y (as needed to work with histogram2d) 
    ...: data = np.random.randint(0,1000,(5000,2)) 
    ...: data1 = data[np.lexsort(data.T),:] 
    ...: mask = ~np.all(np.diff(data1,axis=0)==0,axis=1) 
    ...: data2 = data1[np.append([True],mask)] 
    ...: 
    ...: X = (data2[:,0]).astype(float)/10 
    ...: Y = (data2[:,1]).astype(float)/10 
    ...: Z = np.random.randint(0,1000,(X.size)) 
    ...: 

In [133]: def histogram_based(X,Y,Z): # From other np.histogram2d based solution 
    ...: Xr = np.linspace(X.min(), X.max(), np.round((X.max() - X.min())/np.diff(np.unique(X)).min()) + 1) 
    ...: Yr = np.linspace(Y.min(), Y.max(), np.round((Y.max() - Y.min())/np.diff(np.unique(Y)).min()) + 1) 
    ...: data = np.histogram2d(Y, X, bins=[len(Yr),len(Xr)], weights=Z) 
    ...: return data[0] 
    ...: 

In [134]: %timeit histogram_based(X,Y,Z) 
10 loops, best of 3: 22.8 ms per loop 

In [135]: %timeit indexing_based(X,Y,Z) 
100 loops, best of 3: 2.11 ms per loop