scipy.spatial.distance.pdist
restituisce una matrice a distanza condensata. Dalla documentazione :
Restituisce una matrice della distanza condensata Y. Per ogni e (dove), la dist metrica (u = X [i], v = X [j]) viene calcasting e memorizzata nella voce ij.
- Traccia un piano basato su un vettore normale e un punto in Matlab o matplotlib
- Come installare numpy e scipy per Ironpython27? Il vecchio metodo non funziona
- ImportError durante l'importazione di alcuni moduli da SciPY
- Modo vettorizzato per calcolare le matrici di due punti del prodotto riga-punto con Scipy
- Come posso installare qualcosa su Travis CI senza un timeout?
Pensavo che ij
intendesse i*j
. Ma penso che potrei sbagliarmi. Prendere in considerazione
X = array([[1,2], [1,2], [3,4]]) dist_matrix = pdist(X)
allora la documentazione dice che dist(X[0], X[2])
dovrebbe essere dist_matrix[0*2]
. Tuttavia, dist_matrix[0*2]
è 0 – non 2.8 come dovrebbe essere.
Qual è la formula che dovrei usare per accedere alla somiglianza di due vettori, dati i
e j
?
Puoi guardarlo in questo modo: Supponiamo che x
sia m per n. Le possibili coppie di m
righe, scelte due alla volta, sono itertools.combinations(range(m), 2)
, ad es. Per m=3
:
>>> import itertools >>> list(combinations(range(3),2)) [(0, 1), (0, 2), (1, 2)]
Quindi se d = pdist(x)
, la k
th tuple in combinations(range(m), 2))
fornisce gli indici delle righe di x
associate a d[k]
.
Esempio:
>>> x = array([[0,10],[10,10],[20,20]]) >>> pdist(x) array([ 10. , 22.36067977, 14.14213562])
Il primo elemento è dist(x[0], x[1])
, il secondo è dist(x[0], x[2])
e il terzo è dist(x[1], x[2])
.
Oppure puoi vederlo come gli elementi nella parte triangular superiore della matrice della distanza quadrata, uniti insieme in una matrice 1D.
Per esempio
>>> squareform(pdist(x)) array([[ 0. , 10. , 22.361], [ 10. , 0. , 14.142], [ 22.361, 14.142, 0. ]]) >>> y = array([[0,10],[10,10],[20,20],[10,0]]) >>> squareform(pdist(y)) array([[ 0. , 10. , 22.361, 14.142], [ 10. , 0. , 14.142, 10. ], [ 22.361, 14.142, 0. , 22.361], [ 14.142, 10. , 22.361, 0. ]]) >>> pdist(y) array([ 10. , 22.361, 14.142, 14.142, 10. , 22.361])
Una matrice a distanza condensata restituita da pdist può essere convertita in una matrice a distanza intera utilizzando scipy.spatial.distance.squareform
:
>>> import numpy as np >>> from scipy.spatial.distance import pdist, squareform >>> points = np.array([[0,1],[1,1],[3,5], [15, 5]]) >>> dist_condensed = pdist(points) >>> dist_condensed array([ 1. , 5. , 15.5241747 , 4.47213595, 14.56021978, 12. ])
Usa la forma squareform
per convertire in full matrix:
>>> dist = squareform(dist_condensed) array([[ 0. , 1. , 5. , 15.5241747 ], [ 1. , 0. , 4.47213595, 14.56021978], [ 5. , 4.47213595, 0. , 12. ], [ 15.5241747 , 14.56021978, 12. , 0. ]])
La distanza tra il punto i, j è memorizzata in dist [i, j]:
>>> dist[2, 0] 5.0 >>> np.linalg.norm(points[2] - points[0]) 5.0
Si possono convertire gli indici usati per accedere agli elementi della matrice quadrata all’indice nella matrice condensata:
def square_to_condensed(i, j, n): assert i != j, "no diagonal elements in condensed matrix" if i < j: i, j = j, i return n*j - j*(j+1)/2 + i - 1 - j
Esempio:
>>> square_to_condensed(1, 2, len(points)) 3 >>> dist_condensed[3] 4.4721359549995796 >>> dist[1,2] 4.4721359549995796
Anche l'altra direzione è ansible senza sqaureform, che è migliore in termini di tempo di esecuzione e consumo di memoria:
import math def calc_row_idx(k, n): return int(math.ceil((1/2.) * (- (-8*k + 4 *n**2 -4*n - 7)**0.5 + 2*n -1) - 1)) def elem_in_i_rows(i, n): return i * (n - 1 - i) + (i*(i + 1))/2 def calc_col_idx(k, i, n): return int(n - elem_in_i_rows(i + 1, n) + k) def condensed_to_square(k, n): i = calc_row_idx(k, n) j = calc_col_idx(k, i, n) return i, j
Esempio:
>>> condensed_to_square(3, 4) (1.0, 2.0)
>>> import numpy as np >>> from scipy.spatial.distance import pdist, squareform >>> points = np.random.random((10**4,3)) >>> %timeit dist_condensed = pdist(points) 1 loops, best of 3: 555 ms per loop
La creazione di sqaureform risulta essere molto lenta:
>>> dist_condensed = pdist(points) >>> %timeit dist = squareform(dist_condensed) 1 loops, best of 3: 2.25 s per loop
Se stiamo cercando due punti con la massima distanza non è sorprendente che la ricerca in full matrix sia O (n) mentre in forma condensata solo O (n / 2):
>>> dist = squareform(dist_condensed) >>> %timeit dist_condensed.argmax() 10 loops, best of 3: 45.2 ms per loop >>> %timeit dist.argmax() 10 loops, best of 3: 93.3 ms per loop
Ottenere gli elementi necessari per i due punti non richiede quasi tempo in entrambi i casi, ma ovviamente c'è un sovraccarico per il calcolo dell'indice condensato:
>>> idx_flat = dist.argmax() >>> idx_condensed = dist.argmax() >>> %timeit list(np.unravel_index(idx_flat, dist.shape)) 100000 loops, best of 3: 2.28 µs per loop >>> %timeit condensed_to_square(idx_condensed, len(points)) 100000 loops, best of 3: 14.7 µs per loop
Il vettore della matrice compressa corrisponde alla regione triangular inferiore della matrice quadrata. Per convertire per un punto in quella regione triangular, è necessario calcolare il numero di punti a sinistra nel triangolo e il numero sopra nella colonna.
È ansible utilizzare la seguente funzione per convertire:
q = lambda i,j,n: n*j - j*(j+1)/2 + i - 1 - j
Dai un’occhiata:
import numpy as np from scipy.spatial.distance import pdist, squareform x = np.random.uniform( size = 100 ).reshape( ( 50, 2 ) ) d = pdist( x ) ds = squareform( d ) for i in xrange( 1, 50 ): for j in xrange( i ): assert ds[ i, j ] == d[ q( i, j, 50 ) ]
Se vuoi accedere all’elemento di pdist
corrispondente pdist
(i, j) -th della matrice di distanza quadrata, la matematica è la seguente: Assumi i < j
(altrimenti inverti gli indici) se i == j
, la risposta è 0.
X = random((N,m)) dist_matrix = pdist(X)
Quindi l'elemento (i, j) -th è dist_matrix [ind] dove
ind = (N - array(range(1,i+1))).sum() + (j - 1 - i).
Questa è la versione del triangolo superiore ( i
condensed_idx = lambda i,j,n: i*n + j - i*(i+1)/2 - i - 1
Questo è molto facile da capire:
i*n + j
vai nella posizione nella matrice quadrata; - i*(i+1)/2
rimuovi il triangolo inferiore (inclusa la diagonale) in tutte le linee prima di i; - i
rimuovi le posizioni nella riga i prima della diagonale; - 1
rimuovi le posizioni nella riga i sulla diagonale. Dai un’occhiata:
import scipy from scipy.spatial.distance import pdist, squareform condensed_idx = lambda i,j,n: i*n + j - i*(i+1)/2 - i - 1 n = 50 dim = 2 x = scipy.random.uniform(size = n*dim).reshape((n, dim)) d = pdist(x) ds = squareform(d) for i in xrange(1, n-1): for j in xrange(i+1, n): assert ds[i, j] == d[condensed_idx(i, j, n)]
Se qualcuno sta cercando una trasformazione inversa (cioè dato un indice di elemento idx
, capire a quale elemento (i, j)
corrisponde ad esso), ecco una soluzione vettoriale risonabile:
def actual_indices(idx, n): n_row_elems = np.cumsum(np.arange(1, n)[::-1]) ii = (n_row_elems[:, None] - 1 < idx[None, :]).sum(axis=0) shifts = np.concatenate([[0], n_row_elems]) jj = np.arange(1, n)[ii] + idx - shifts[ii] return ii, jj n = 5 k = 10 idx = np.random.randint(0, n, k) a = pdist(np.random.rand(n, n)) b = squareform(a) ii, jj = actual_indices(idx, n)] assert np.allclose(b[ii, jj, a[idx])
L'ho usato per capire gli indici delle righe più vicine in una matrice.
m = 3 # how many closest lowest_dist_idx = np.argpartition(-a, -m)[-m:] ii, jj = actual_indices(lowest_dist_idx, n) # rows ii[0] and jj[0] are closest