Appliquer une fonction à chaque ligne d'un ndarray
J'ai cette fonction pour calculer le carré de la distance de Mahalanobis vecteur x:
def mahalanobis_sqdist(x, mean, Sigma):
'''
Calculates squared Mahalanobis Distance of vector x
to distibutions' mean
'''
Sigma_inv = np.linalg.inv(Sigma)
xdiff = x - mean
sqmdist = np.dot(np.dot(xdiff, Sigma_inv), xdiff)
return sqmdist
J'ai un tableau numpy qui a une forme de (25, 4)
. Donc, je veux appliquer cette fonction à l'ensemble des 25 lignes de mon tableau sans une boucle for. Donc, en gros, comment puis-je écrire le vectorisé forme de cette boucle:
for r in d1:
mahalanobis_sqdist(r[0:4], mean1, Sig1)
où mean1
et Sig1
sont :
>>> mean1
array([ 5.028, 3.48 , 1.46 , 0.248])
>>> Sig1 = np.cov(d1[0:25, 0:4].T)
>>> Sig1
array([[ 0.16043333, 0.11808333, 0.02408333, 0.01943333],
[ 0.11808333, 0.13583333, 0.00625 , 0.02225 ],
[ 0.02408333, 0.00625 , 0.03916667, 0.00658333],
[ 0.01943333, 0.02225 , 0.00658333, 0.01093333]])
J'ai essayé ce qui suit, mais il ne fonctionne pas:
>>> vecdist = np.vectorize(mahalanobis_sqdist)
>>> vecdist(d1, mean1, Sig1)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/lib/python2.7/dist-packages/numpy/lib/function_base.py", line 1862, in __call__
theout = self.thefunc(*newargs)
File "<stdin>", line 6, in mahalanobis_sqdist
File "/usr/lib/python2.7/dist-packages/numpy/linalg/linalg.py", line 445, in inv
return wrap(solve(a, identity(a.shape[0], dtype=a.dtype)))
IndexError: tuple index out of range
source d'informationauteur Vahid Mir
Vous devez vous connecter pour publier un commentaire.
Appliquer une fonction à chaque ligne d'un tableau, vous pouvez utiliser:
Dans ce cas, cependant, il ya une meilleure façon. Vous n'avez pas à appliquer une fonction à chaque ligne. Au lieu de cela, vous pouvez appliquer NumPy des opérations pour l'ensemble de la
d1
tableau pour calculer le même résultat. np.einsum pouvez remplacer lefor-loop
et les deux appels ànp.dot
:Voici quelques repères:
Ainsi
mahalanobis_sqdist2
est d'environ 18 fois plus rapide qu'unfor-loop
et 26x plus rapide que l'utilisationnp.apply_along_axis
.Noter que
np.apply_along_axis
np.vectorize
np.frompyfunc
Python sont les fonctions de l'utilitaire. Sous le capot, ils utilisentfor-
ouwhile-loop
s. Il n'y a pas de véritable "vectorisation" se passe ici. Ils peuvent fournir syntaxique de l'aide, mais ne vous attendez pas à faire de votre code faire mieux qu'unfor-loop
vous écrivez vous-même.La réponse par @unutbu fonctionne très bien pour l'application de n'importe quelle fonction pour les lignes d'un tableau.
Dans ce cas particulier, il y a quelques mathématique symétries vous pouvez utiliser qui permettra d'accélérer considérablement les choses si vous travaillez avec de grands tableaux.
Ici est une version modifiée de votre fonction:
Si vous vous retrouvez à l'aide de toute sorte de grand
Sigma
je vous recommande de cacheSigma_inv
et passer que comme un argument à votre fonction à la place.Depuis qu'il est en 4x4 dans cet exemple, ce n'est pas grave.
Je vais vous montrer comment faire face à de gros
Sigma
de toute façon pour quelqu'un d'autre qui vient dans ce.Si vous n'allez pas être en utilisant le même
Sigma
à plusieurs reprises, vous ne serez pas en mesure de mettre en cache, de sorte que, au lieu de l'inversion de la matrice, vous pouvez utiliser une autre méthode pour résoudre le système linéaire.Ici, je vais utiliser la décomposition LU intégrée dans SciPy.
Cette seulement améliore le temps si le nombre de colonnes de
x
est grande par rapport à son nombre de lignes.Voici une fonction qui montre que l'approche:
Voici quelques temps.
Je vais inclure la version avec
einsum
comme mentionné dans la réponse à faire.:
Toutefois, la modification de la taille des matrices de changements le moment des résultats.
Par exemple, laisser
x = np.random.rand(2500, 4)
, les horaires sont les suivants:Et laisser
x = np.random.rand(1000, 1000)
Sigma1 = np.random.rand(1000, 1000)
etmean1 = np.random.rand(1000)
les horaires sont les suivants:Modifier: j'ai remarqué que l'une des autres réponses utilisé la décomposition de Cholesky.
Étant donné que
Sigma
est symétrique et définie positive, nous pouvons réellement faire mieux que mes résultats ci-dessus.Il y a quelques bonnes routines de BLAS et LAPACK disponibles par le biais de SciPy qui peuvent travailler de façon symétrique définie positive des matrices.
Voici deux versions plus rapide.
Le premier toujours inverse de Sigma.
Si vous pré-calculer l'inverse et de le réutiliser, il est beaucoup plus rapide (le 1000x1000 cas prend de 35,6 ms sur ma machine avec le pré-calculé l'inverse).
J'ai aussi utilisé einsum à prendre le produit, la somme, le long de la dernière de l'axe.
Cela a fini par être légèrement plus rapide que de faire quelque chose comme
(A * B).sum(axis=-1)
.Ces deux fonctions donner le minutage suivant:
Premier cas de test:
Deuxième cas de test:
Troisième cas de test:
Viens de voir un très beau commentaire sur reddit qui pourrait accélérer les choses encore un peu plus:
Le problème est que
np.vectorize
vectorizes plus de tous les arguments, mais vous avez besoin de vectoriser seulement sur le premier. Vous avez besoin d'utiliserexcluded
argument mot-clé àvectorize
: