Algoritmi veloci per la ricerca della distanza euclidea a coppie

Question

So che MATLAB ha una funzione pdist incorporata che calcolerà le distanze a coppie. Tuttavia, la mia matrice è così grande che i suoi 60000 per 300 e MATLAB esauriscono la memoria.

Questa domanda è un seguito su Matlab euclidean pairwise square distance function.

C'è qualche soluzione per questa inefficienza computazionale. Ho provato a codificare manualmente i calcoli della distanza a coppie e di solito ci vuole un intero giorno per correre (a volte da 6 a 7 ore).

Qualsiasi aiuto è molto apprezzato!

Answer 1

I computer non sono infinitamente grandi o infinitamente veloci. Le persone pensano di avere molta memoria, una CPU veloce, quindi creano solo problemi sempre più grandi e alla fine si chiedono perché il loro problema si verifica lentamente. Il fatto è che questa NON è inefficienza computazionale. È SOLO una CPU sovraccaricata.

Come indicato da Oli in un commento, ci sono qualcosa come 2e9 valori da calcolare, anche assumendo che si calcoli solo la metà superiore o inferiore della matrice della distanza. (6e4^2/2 è approssimativamente 2e9.) Ciò richiederà circa 16 gigabyte di RAM da memorizzare, assumendo che solo una copia dell'array viene creata in memoria. Se il tuo codice è sciatto, potresti facilmente raddoppiarlo o triplicarlo. Non appena entri nella memoria virtuale, le cose si fanno molto più lente.

La ricerca di un grosso problema per correre velocemente non è sufficiente. Per aiutarti davvero, abbiamo bisogno di sapere quanta RAM è disponibile. È un problema di memoria virtuale? Stai usando MATLAB a 64 bit, su una CPU in grado di gestire tutta la RAM necessaria?

Answer 2

Beh, non ho potuto resistere a giocare in giro. Ho creato un Matlab mex C file chiamato pdistc che implementa la distanza euclidea a coppie per precisione singola e doppia. Sulla mia macchina con Matlab R2012b e R2015a è 20 – 25% più veloce di pdist (e la funzione di supporto pdistmex sottostante) per ingressi di grandi dimensioni (ad es. 60.000-per-300).

Come è stato sottolineato, questo problema è fondamentalmente limitato dalla memoria e tu stai chiedendo molto. Il mio codice C mex utilizza una memoria minima oltre a quella necessaria per l'output. Confrontando l'utilizzo della memoria con quello di pdist, sembra che i due siano praticamente gli stessi. In altre parole, pdist non utilizza molta memoria aggiuntiva. Il tuo problema di memoria è probabilmente nella memoria utilizzata prima di chiamare pdist (puoi usare clear per rimuovere qualsiasi array di grandi dimensioni?) O semplicemente perché stai cercando di risolvere un grosso problema computazionale su hardware minuscolo.

Quindi, la mia funzione pdistc probabilmente non sarà in grado di salvarvi la memoria in generale, ma potreste essere in grado di utilizzare un'altra funzione che ho incorporato. Potete calcolare pezzi del vostro vettore di distanza globale a coppie. Qualcosa di simile a questo:

m = 6e3; 
n = 3e2; 
X = rand(m,n); 
sz = m*(m-1)/2; 

for i = 1:m:sz-m 
    D = pdistc(X', i, i+m); % mex C function, X is transposed relative to pdist 
    ...      % Process chunk of pairwise distances 
end 

Questo è notevolmente più lento (10 volte o giù di lì) e questa parte del mio codice C non è ottimizzato bene, ma permetterà molto meno uso di memoria – partendo dal presupposto che non è necessario il intero array alla volta. Nota che potresti fare la stessa cosa in modo molto più efficiente con pdist (o pdistc) creando un loop in cui hai passato direttamente sottoinsiemi di X, piuttosto che tutto.

Se si dispone di un Intel Mac a 64 bit, non è necessario compilare poiché ho incluso il binario .mexmaci64, ma in caso contrario sarà necessario capire come compilare il codice per la macchina. Non posso aiutarti con quello.È possibile che tu non sia in grado di farlo compilare o che ci siano problemi di compatibilità che dovrai risolvere modificando il codice tu stesso. È anche possibile che ci siano errori e il codice si bloccherà Matlab. Inoltre, si noti che è possibile ottenere risultati leggermente diversi rispetto a pdist con differenze tra i due nel campo della macchina epsilon (eps). pdist può o non può fare cose di fantasia per evitare overflow per grandi input e altri problemi numerici, ma sii consapevole che il mio codice non lo fa.

Inoltre, ho creato un semplice pure Matlab implementation. È molto più lento del codice mex, ma ancora più veloce di un'implementazione ingenua o del codice trovato in pdist.

Tutti i file can be found here. L'archivio ZIP include tutti i file. È concesso in licenza BSD. Sentitevi liberi di ottimizzare (ho provato le chiamate BLAS e OpenMP nel codice C senza alcun risultato – forse un po 'di magia del puntatore o GPU/OpenCL potrebbe accelerarlo ulteriormente). Spero che possa essere utile a te o a qualcun altro.

Answer 3

Sul mio sistema il seguente è il più veloce (anche più veloce rispetto al codice C pdistc da @horchler):

function [ mD ] = CalcDistMtx (mX)  
    vSsqX = sum(mX .^ 2); 
    mD = sqrt(bsxfun(@plus, vSsqX.', vSsqX) - (2 * (mX.' * mX)));  
end 

avrete bisogno di un codice C molto ben sintonizzato per battere questo, credo.

P. S.
Utilizzo di MATLAB pdist per il confronto: squareform(pdist(mX.')) è equivalente a CalcDistMtx(mX).
Vale a dire che l'input deve essere trasposto.