<p>要将函数应用于数组的每一行,可以使用:</p>
<pre><code>np.apply_along_axis(mahalanobis_sqdist, 1, d1, mean1, Sig1)
</code></pre>
<p>不过,在这种情况下,还有更好的办法。不必对每一行应用函数。相反,您可以对整个<code>d1</code>数组应用NumPy操作来计算相同的结果。<a href="http://docs.scipy.org/doc/numpy/reference/generated/numpy.einsum.html" rel="noreferrer">np.einsum</a>可以替换<code>for-loop</code>和对<code>np.dot</code>的两个调用:</p>
<hr/>
<pre><code>def mahalanobis_sqdist2(d, mean, Sigma):
Sigma_inv = np.linalg.inv(Sigma)
xdiff = d - mean
return np.einsum('ij,im,mj->i', xdiff, xdiff, Sigma_inv)
</code></pre>
<hr/>
<p>以下是一些基准:</p>
<pre><code>import numpy as np
np.random.seed(1)
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
def mahalanobis_sqdist2(d, mean, Sigma):
Sigma_inv = np.linalg.inv(Sigma)
xdiff = d - mean
return np.einsum('ij,im,mj->i', xdiff, xdiff, Sigma_inv)
def using_loop(d1, mean, Sigma):
expected = []
for r in d1:
expected.append(mahalanobis_sqdist(r[0:4], mean1, Sig1))
return np.array(expected)
d1 = np.random.random((25,4))
mean1 = np.array([ 5.028, 3.48 , 1.46 , 0.248])
Sig1 = np.cov(d1[0:25, 0:4].T)
expected = using_loop(d1, mean1, Sig1)
result = np.apply_along_axis(mahalanobis_sqdist, 1, d1, mean1, Sig1)
result2 = mahalanobis_sqdist2(d1, mean1, Sig1)
assert np.allclose(expected, result)
assert np.allclose(expected, result2)
</code></pre>
<hr/>
<pre><code>In [92]: %timeit mahalanobis_sqdist2(d1, mean1, Sig1)
10000 loops, best of 3: 31.1 µs per loop
In [94]: %timeit using_loop(d1, mean1, Sig1)
1000 loops, best of 3: 569 µs per loop
In [91]: %timeit np.apply_along_axis(mahalanobis_sqdist, 1, d1, mean1, Sig1)
1000 loops, best of 3: 806 µs per loop
</code></pre>
<p>因此<code>mahalanobis_sqdist2</code>比<code>for-loop</code>快18倍,比使用<code>np.apply_along_axis</code>快26倍。</p>
<hr/>
<p>注意,<code>np.apply_along_axis</code>,<code>np.vectorize</code>,<code>np.frompyfunc</code>是Python实用程序函数。在引擎盖下他们使用<code>for-</code>或<code>while-loop</code>s。这里没有真正的“矢量化”。它们可以提供语法帮助,但不要期望它们能使代码的性能比您自己编写的代码更好。</p>