J'ai un tableau de coordonnées x, y, z de plusieurs points (~ 10 ^ 10) (seulement 5 montré ici) Je veux faire un nouveau tableau avec seulement ceux points qui sont au moins une distance cela fonctionne, mais il prend vraiment longtemps pour exécuter ce calcul. J'ai lu quelque part que Je me demandais si quelqu'un a des suggestions sur la manière d'accélérer ce calcul. EDIT: Pendant mon objectif de trouver les particules qui sont au moins Une certaine distance de tous les autres reste la même chose, je viens de réaliser qu'il y a une forte défaute dans mon code, disons que j'ai 3 particules, mon code fait ce qui suit pour la première itération de d à partir de tous les autres points de la liste. J'ai écrit un code à l'aide de pendant boucle, pendant que Les boucles sont très lentes. i , il calcule les distances 1-> 2 , 1-> 3 , disons 1-> 2 est inférieur à la distance de seuil D , donc le code lance la particule 1 . Pour la prochaine itération de i , il ne fait que 2-> 3 , et disons qu'il trouve qu'il est supérieur à d , donc il Garde les particules 2 , mais c'est faux! Depuis 2 devrait également être jeté avec des particules 1 . La solution de @svohara est la bonne!
4 Réponses :
déposer l'annexe, il doit être vraiment lent. Vous pouvez avoir un vecteur statique de distances et utiliser [] pour mettre le numéro dans la bonne position. P> li>
Utilisez min au lieu de tous. Vous n'avez besoin que de vérifier si la distance minimale est plus grande que x. p> li>
En fait, vous pouvez casser votre annexe au moment où vous trouverez une distance plus petite que votre limite, puis vous pouvez abandonner les deux points. De cette façon, vous n'avez même pas besoin de sauver une distance (sauf si vous en avez besoin plus tard). P>
de votre commentaire, je crois que cela ferait, si vous n'avez pas de points répétés. P>
selected_points = []
for p1 in a:
save_point = True
for p2 in a:
if p1!=p2 and distance.euclidean(p1,p2)<d:
save_point = False
break
if save_point:
selected_points.append(p1)
return selected_points
Merci! Je n'ai pas besoin des distances, mais je suis confus sur la façon de casser l'itération en J dès que je rencontre que la distance est inférieure à D, puis passez à la prochaine itération de I?
La ligne P1! = P2 donne une erreur "ValueError: la valeur de vérité d'un tableau avec plus d'un élément est ambiguë. Utilisez A.ANY () ou A.Al ()" car p1 et P2 sont des matrices [x y z].
Avez-vous vérifié votre dactylographie? Cela a bien fonctionné ici et nous ne demandons pas la valeur de vérité de la matrice. J'ai oublié qu'il pointit (distance.euclidien) parce que je l'ai simulé.
J'en copie en collé votre code, il est d'abord dit "" Retour "Fonction extérieure", puis lorsque j'ai changé de retour à l'impression, je reçois "ValueError: la valeur de vérité d'un tableau avec plus d'un élément est ambigu. Utilisez A.ANY ( ) ou a.all () ", utilisez-vous Python 2.7 ou une autre version?
Voici une approche vectorisée utilisant pour un énorme jeu de données comme distance.pdist - 10e10 , nous devrions peut-être effectuer les opérations en morceaux Basé sur la mémoire système disponible.
La sortie contient tous les points au lieu de seuls ceux qui sont distants de tous les autres.
@Hushu pas sûr si je vous suis. Il donne le même O / P comme avec votre code dans la question.
Désolé mon mauvais. Je viens de réaliser qu'il y a un bogue dans mon code d'origine. S'il vous plaît vérifier la modification.
Votre algorithme est quadratique (10 ^ 20 opérations), voici une approche linéaire si la distribution est presque aléatoire.
Split votre espace dans des cases de taille Puis pour chaque boîte,
S'il n'y a qu'un point, il suffit de calculer la distance avec des points dans un petit quartier.
aire, il n'y a rien à faire.
d / sqrt (3) ^ 3 . Mettre chaque point dans sa boîte.
Pour les grands ensembles de données et les points à faible dimension (tels que vos données à 3 dimensions), il existe parfois un grand avantage pour utiliser une méthode d'indexation spatiale. L'arborescence K-D est un choix populaire pour les données à faible dimension.
La stratégie consiste à indexer le jeu de données. Demandez ensuite à l'index en utilisant le même ensemble de données, pour renvoyer les voisins 2 les plus proches de chaque point. Le premier voisin le plus proche est toujours le point lui-même (avec dist = 0), nous voulons donc vraiment savoir à quelle distance le prochain point le plus proche est (2e voisin le plus proche). Pour ces points où le seuil de 2 nn est>, vous avez le résultat. P>
from scipy.spatial import cKDTree as KDTree
import numpy as np
#a is the big data as numpy array N rows by 3 cols
a = np.random.randn(10**8, 3).astype('float32')
# This will create the index, prepare to wait...
# NOTE: took 7 minutes on my mac laptop with 10^8 rand 3-d numbers
# there are some parameters that could be tweaked for faster indexing,
# and there are implementations (not in scipy) that can construct
# the kd-tree using parallel computing strategies (GPUs, e.g.)
k = KDTree(a)
#ask for the 2-nearest neighbors by querying the index with the
# same points
(dists, idxs) = k.query(a, 2)
# (dists, idxs) = k.query(a, 2, n_jobs=4) # to use more CPUs on query...
#Note: 9 minutes for query on my laptop, 2 minutes with n_jobs=6
# So less than 10 minutes total for 10^8 points.
# If the second NN is > thresh distance, then there is no other point
# in the data set closer.
thresh_d = 0.1 #some threshold, equiv to 'd' in O.P.'s code
d_slice = dists[:, 1] #distances to second NN for each point
res = np.flatnonzero( d_slice >= thresh_d )
Une note sur la complexité de la requête. Chaque requête est O (journal (n)), avec n échantillons, la complexité totale de temps pour compléter la requête de toutes les points est en moyenne O (n log (n)).
Combien de temps cela prend-il ?
J'ai couru pendant la nuit ~ 7 heures et ça marche toujours.
En tant que suggestion rapide, vous ne pouvez pas continuer à calculer des distances si l'on est plus grand que celui
d . Il réduira une autre exécution à travers la matrice dans la clausetout (DIS> D pour DI dans l'interdiste)Pouvez-vous fournir un script fonctionnel? Il n'y a pas de définition de
halos etsélectionné_halos dans le vôtreDésolé pour une idée incomplète, mais je me souviens de lire un moment où une approche du regroupement K-signifie (qui nécessite également des calculs de distance euclidiens) qui utilise l'inégalité du triangle pour réduire le nombre de calculs requis
@max: merci! Je suppose donc que je peux vérifier si la distance est supérieure à D à chaque itération de J et dès qu'il est faux, je brise l'itération et allez à la prochaine itération. Comment rompez-vous l'itération au milieu?
pause est un mot-clé pour casser la boucle. Mais vous devez être conscient des distances de points avec moins d'index quei + 1 car vous ne les prenez pas en considération qui n'est pas juste, je suppose.