優化與numpy的sparsey陣列

Question

我有一個緩慢的numpy的操作掙扎，使用Python 3

的操作，我有以下操作：

np.sum(np.log(X.T * b + a).T, 1) 

其中

(30000,1000) = X.shape 
(1000,1) = b.shape 
(1000,1) = a.shape 

我的問題這個操作很慢（大約1.5秒），並且它在一個循環內，所以它重複了大約100次，這使得我的代碼的運行時間非常長。

我想知道是否有更快的實現這個功能。

也許有用的事實：X非常稀疏（只有0.08％的條目非零），但是一個NumPy數組。

Answer 1

我們可以優化這似乎是瓶頸對數運算和被超越的一個函數可以用numexpr module加速，然後sum-reduce用NumPy加速，因爲NumPy做得好得多，因此給了我們一個混合版本，就像這樣 -

import numexpr as ne 

def numexpr_app(X, a, b): 
    XT = X.T 
    return ne.evaluate('log(XT * b + a)').sum(0) 

仔細觀察廣播業務：XT * b + a，我們看到廣播有兩個階段，我們可以進一步優化。目的是要看看是否可以減少到一個階段，這在某些部門似乎是可能的。這給我們一個稍微修改後的版本，如下圖所示 -

def numexpr_app2(X, a, b): 
    ab = (a/b) 
    XT = X.T 
    return np.log(b).sum() + ne.evaluate('log(ab + XT)').sum(0) 

---

運行測試和驗證

原始的方法 -

def numpy_app(X, a, b): 
    return np.sum(np.log(X.T * b + a).T, 1) 

計時 -

In [111]: # Setup inputs 
    ...: density = 0.08/100 # 0.08 % sparse 
    ...: m,n = 30000, 1000 
    ...: X = scipy.sparse.rand(m,n,density=density,format="csr").toarray() 
    ...: a = np.random.rand(n,1) 
    ...: b = np.random.rand(n,1) 
    ...: 

In [112]: out0 = numpy_app(X, a, b) 
    ...: out1 = numexpr_app(X, a, b) 
    ...: out2 = numexpr_app2(X, a, b) 
    ...: print np.allclose(out0, out1) 
    ...: print np.allclose(out0, out2) 
    ...: 
True 
True 

In [114]: %timeit numpy_app(X, a, b) 
1 loop, best of 3: 691 ms per loop 

In [115]: %timeit numexpr_app(X, a, b) 
10 loops, best of 3: 153 ms per loop 

In [116]: %timeit numexpr_app2(X, a, b) 
10 loops, best of 3: 149 ms per loop 

只是爲了證明在啓動log部分與原NumPy的方法的瓶頸規定的觀察，這裏是它的時機 - 在其改善是顯著

In [44]: %timeit np.log(X.T * b + a) 
1 loop, best of 3: 682 ms per loop 

-

In [120]: XT = X.T 

In [121]: %timeit ne.evaluate('log(XT * b + a)') 
10 loops, best of 3: 142 ms per loop 

Answer 2

有點不清楚你爲什麼要做np.sum(your_array.T, axis=1)而不是np.sum(your_array, axis=0)。

您可以使用scipy sparse matrix：（對於X使用列壓縮格式，使XT的行壓縮，因爲你用b具有XT的一個排的形狀乘）

X_sparse = scipy.sparse.csc_matrx(X) 

並更換XT * b通過：

X_sparse.T.multiply(b) 

但是，如果一個不稀疏，它不會幫助你儘可能多。

這些都是速度提升我獲得這個操作：

In [16]: %timeit X_sparse.T.multiply(b) 
The slowest run took 10.80 times longer than the fastest. This could mean that an intermediate result is being cached. 
1000 loops, best of 3: 374 µs per loop 

In [17]: %timeit X.T * b 
10 loops, best of 3: 44.5 ms per loop 

有：

import numpy as np 
from scipy import sparse 

X = np.random.randn(30000, 1000) 
a = np.random.randn(1000, 1) 
b = np.random.randn(1000, 1) 
X[X < 3] = 0 
print(np.sum(X != 0)) 
X_sparse = sparse.csc_matrix(X)