pdist para tensor theano

Tengo una matriz simbólica theana.

x = T.fmatrix('input') 

x adelante, x estará poblado por n vectores de dim d (a la hora del tren).

Me gustaría tener el equivalente pdist de pdist ( scipy.spatial.distance.pdist de pdist ), algo así como

 D = theano.pdist( x ) 

¿Cómo puedo conseguir esto?

Llamar a scipy.spatial.distance.pdist en x directamente no funciona, ya que x en esta etapa solo es simbólico …

Actualización: Me gustaría mucho poder imitar el pdist “compacto” de pdist : es decir, calcular solo ~ pdist de las entradas n x n de la matriz de distancia.

pdist de scipy es una colección de diferentes funciones, no existe un equivalente de Theano para todas ellas a la vez. Sin embargo, cada distancia específica, al ser una expresión matemática de forma cerrada, se puede escribir en Theano como tal y luego comstackr.

Tomemos como ejemplo la distancia de la norma minkowski p (copia + pasteable):

 import theano import theano.tensor as T X = T.fmatrix('X') Y = T.fmatrix('Y') P = T.scalar('P') translation_vectors = X.reshape((X.shape[0], 1, -1)) - Y.reshape((1, Y.shape[0], -1)) minkowski_distances = (abs(translation_vectors) ** P).sum(2) ** (1. / P) f_minkowski = theano.function([X, Y, P], minkowski_distances) 

Tenga en cuenta que abs llama al __abs__ , por lo que abs también es una función theano. Ahora podemos comparar esto con pdist :

 import numpy as np from scipy.spatial.distance import pdist rng = np.random.RandomState(42) d = 20 # dimension nX = 10 nY = 30 x = rng.randn(nX, d).astype(np.float32) y = rng.randn(nY, d).astype(np.float32) ps = [1., 3., 2.] for p in ps: d_theano = f_minkowski(x, x, p)[np.triu_indices(nX, 1)] d_scipy = pdist(x, p=p, metric='minkowski') print "Testing p=%1.2f, discrepancy %1.3e" % (p, np.sqrt(((d_theano - d_scipy) ** 2).sum())) 

Esto produce

 Testing p=1.00, discrepancy 1.322e-06 Testing p=3.00, discrepancy 4.277e-07 Testing p=2.00, discrepancy 4.789e-07 

Como puede ver, la correspondencia está ahí, pero la función f_minkowski es ligeramente más general, ya que compara las líneas de dos arreglos posiblemente diferentes. Si se pasa dos veces la misma matriz como entrada, f_minkowski devuelve una matriz, mientras que pdist devuelve una lista sin redundancia. Si se desea este comportamiento, también se puede implementar de forma totalmente dinámica, pero me limitaré al caso general aquí.

Sin embargo, debe observarse una posibilidad de especialización: en el caso de p=2 , los cálculos se vuelven más simples a través de la fórmula binomial, y esto se puede usar para ahorrar un valioso espacio en la memoria: mientras que la distancia general de Minkowski, como se implementó anteriormente, crea una Matriz 3D (debido a la evitación de bucles for y la sum acumulativa), que es prohibitiva, dependiendo de la dimensión d (y nX, nY ), para p=2 podemos escribir

 squared_euclidean_distances = (X ** 2).sum(1).reshape((X.shape[0], 1)) + (Y ** 2).sum(1).reshape((1, Y.shape[0])) - 2 * X.dot(YT) f_euclidean = theano.function([X, Y], T.sqrt(squared_euclidean_distances)) 

que solo usa el espacio O(nX * nY) lugar de O(nX * nY * d) Verificamos la correspondencia, esta vez en el problema general:

 d_eucl = f_euclidean(x, y) d_minkowski2 = f_minkowski(x, y, 2.) print "Comparing f_minkowski, p=2 and f_euclidean: l2-discrepancy %1.3e" % ((d_eucl - d_minkowski2) ** 2).sum() 

flexible

 Comparing f_minkowski, p=2 and f_euclidean: l2-discrepancy 1.464e-11 

No he trabajado con Theano antes, pero aquí hay una solución basada en funciones puras de Numpy (tal vez la conviertas a las funciones equivalentes de theano. Ten en cuenta que estoy usando la transmisión automática en la expresión de abajo, así que es posible que tengas que volver a escribir eso). explícitamente si Theano no lo soporta):

 # X is an m-by-n matrix (rows are examples, columns are dimensions) # D is an m-by-m symmetric matrix of pairwise Euclidean distances a = np.sum(X**2, axis=1) D = np.sqrt((a + a[np.newaxis].T) - 2*np.dot(X, XT)) 

Se basa en el hecho de que: ||uv||^2 = ||u||^2 + ||v||^2 - 2*uv . (Mostré esto en respuestas anteriores mías usando MATLAB)

Aquí hay una comparación con las funciones existentes de Scipy:

 import numpy as np from scipy.spatial.distance import pdist, squareform def my_pdist(X): a = np.sum(X**2, axis=1) D = np.sqrt((a + a[np.newaxis].T) - 2*np.dot(X, XT)) return D def scipy_pdist(X): D = squareform(pdist(X, metric='euclidean')) return DX = np.random.rand(5, 3) D1 = my_pdist(X) D2 = scipy_pdist(X) 

La diferencia debe ser despreciable, cerca de la máquina epsilon ( np.spacing(1) ):

 >>> np.linalg.norm(D1-D2) 8.5368137554718277e-16 

HTH


EDITAR:

Aquí hay otra implementación con un solo bucle:

 def my_pdist_compact(X): D = np.empty(shape=[0,0], dtype=X.dtype) for i in range(X.shape[0]-1): D = np.append(D, np.sqrt(np.sum((X[i,] - X[i+1:,])**2, axis=1))) return D 

Código MATLAB algo equivalente:

 function D = my_pdist_compact(X) n = size(X,1); D = cell(n-1,1); for i=1:n-1 D{i} = sqrt(sum(bsxfun(@minus, X(i,:), X(i+1:end,:)).^2, 2)); end D = vertcat(D{:}); end 

Esto devuelve las distancias de pares en forma compacta (parte triangular superior de la matriz simétrica). Esta es la misma salida que pdist . Usa squareform para convertirlo a matriz completa.

 >>> d1 = my_pdist_compact(X) >>> d2 = pdist(X) # from scipy.spatial.distance >>> (d1 == d2).all() True 

Te lo dejo a ti para ver si es posible escribir el bucle equivalente usando Theano (ver theano.scan ).