Para cada punto de una matriz, busque el punto más cercano en una segunda matriz y genere ese índice

Si tengo dos matrices:

X = np.random.rand(10000,2) Y = np.random.rand(10000,2) 

¿Cómo puedo, para cada punto en X, averiguar qué punto en Y es el más cercano? Así que al final tengo una matriz que muestra:

 x1_index y_index_of_closest 1 7 2 54 3 3 ... ... 

Quiero hacer esto para ambas columnas en X y comparar cada una con cada columna y valor en Y

Esta pregunta es bastante popular. Dado que las preguntas similares se siguen cerrando y enlazando aquí, creo que vale la pena señalar que aunque las respuestas existentes son bastante rápidas para miles de puntos de datos, después de eso comienzan a descomponerse. Mis papas seguras en 10k elementos en cada matriz.

El problema potencial con las otras respuestas es la complejidad algorítmica. Ellos comparan todo en X con todo en Y Para solucionar eso, al menos en promedio, necesitamos una mejor estrategia para descartar algunas de las cosas en Y

En una dimensión, esto es fácil: solo ordena todo y comienza a sacar a los vecinos más cercanos. En dos dimensiones hay una variedad de estrategias, pero los árboles KD son razonablemente populares y ya están implementados en la stack de scipy . En mi máquina, hay un cruce entre los diversos métodos en torno al punto en el que cada uno de X e Y tiene 6k elementos en ellos.

 from scipy.spatial import KDTree tree = KDTree(X) neighbor_dists, neighbor_indices = tree.query(Y) 

El rendimiento extremadamente deficiente de la scipy de scipy de scipy ha sido un punto delicado en mi caso por un tiempo, especialmente con tantas cosas como se han construido sobre él. Probablemente hay conjuntos de datos donde se desempeña bien, pero no he visto ninguno todavía.

Si no le importa una dependencia adicional, puede obtener un aumento de velocidad de 1000x simplemente cambiando su biblioteca KDTree. El paquete pykdtree es pykdtree pip, y casi garantizo que los paquetes conda también funcionan bien. Con este enfoque, mi Chromebook de presupuesto usado puede procesar X e Y con 10 millones de puntos cada uno en apenas 30 segundos. Eso le gana a segfaulting en 10 mil puntos con un tiempo de pared;)

 from pykdtree.kdtree import KDTree tree = KDTree(X) neighbor_dists, neighbor_indices = tree.query(Y) 

Esta tiene que ser la pregunta más frecuente (la he contestado dos veces en la última semana), pero ya que puede expressse de un millón de maneras:

 import numpy as np import scipy.spatial.distance.cdist as cdist def withScipy(X,Y): # faster return np.argmin(cdist(X,Y,'sqeuclidean'),axis=0) def withoutScipy(X,Y): #slower, using broadcasting return np.argmin(np.sum((X[None,:,:]-Y[:,None,:])**2,axis=-1), axis=0) 

También hay un método de solo einsum usando einsum que es más rápido que mi función (pero no cdist ) pero no lo entiendo lo suficiente como para explicarlo.

EDITAR + = 21 meses:

Sin embargo, la mejor manera de hacer esto de forma algorítmica es con KDTree.

 from sklearn.neighbors import KDTree # since the sklearn implementation allows return_distance = False, saving memory y_tree = KDTree(Y) y_index_of_closest = y_tree.query(X, k = 1, return_distance = False) 

@HansMusgrave tiene una aceleración bastante buena para KDTree a continuación.

Y para completar, la respuesta np.einsum , que ahora entiendo:

 np.argmin( # (X - Y) ** 2 np.einsum('ij, ij ->i', X, X)[:, None] + # = X ** 2 \ np.einsum('ij, ij ->i', Y, Y) - # + Y ** 2 \ 2 * X.dot(YT), # - 2 * X * Y axis = 1) 

@Divakar explica bien este método en la página wiki de su paquete eucl_dist