在Python中的陣列中的每個元素快速計算特徵向量

Question

我要計算特徵向量數據的陣列（在我的實際情況，我雲多邊形）

爲此我寫了這個功能：

import numpy as np 

def eigen(data): 
    eigenvectors = [] 
    eigenvalues = [] 

    for d in data: 
     # compute covariance for each triangle 
     cov = np.cov(d, ddof=0, rowvar=False) 

     # compute eigen vectors 
     vals, vecs = np.linalg.eig(cov) 
     eigenvalues.append(vals) 
     eigenvectors.append(vecs) 

    return np.array(eigenvalues), np.array(eigenvectors) 

的一些測試數據運行此：

import cProfile 
triangles = np.random.random((10**4,3,3,)) # 10k 3D triangles 
cProfile.run('eigen(triangles)') # 550005 function calls in 0.933 seconds 

工作正常，但它會因爲迭代循環的速度很慢。有沒有更快的方法來計算我所需要的數據而無需迭代數組？如果沒有，任何人都可以建議加快速度的方法嗎？

Answer 1

本事！

嗯，我攻入covariance func definition，並把在規定輸入狀態：ddof=0, rowvar=False和事實證明，一切都降低了，只是三線 -

nC = m.shape[1] # m is the 2D input array 
X = m - m.mean(0) 
out = np.dot(X.T, X)/nC 

把它擴大到我們的3D陣列的情況下，我寫的向下多圈版本與這三條線被重複用於2D陣列的部分從3D輸入陣列，像這樣 -

for i,d in enumerate(m): 

    # Using np.cov : 
    org_cov = np.cov(d, ddof=0, rowvar=False) 

    # Using earlier 2D array hacked version : 
    nC = m[i].shape[0] 
    X = m[i] - m[i].mean(0,keepdims=True) 
    hacked_cov = np.dot(X.T, X)/nC 

升壓它向上

我們需要加快最後三行。在所有迭代X計算可以用做broadcasting -

diffs = data - data.mean(1,keepdims=True) 

接下來，點積運算的各個版本可以與transpose和np.dot來完成，但transpose可能是這樣的多昂貴的事情維數組。存在更好的選擇在np.einsum，像這樣 -

cov3D = np.einsum('ijk,ijl->ikl',diffs,diffs)/data.shape[1] 

使用它！

綜上所述：

for d in data: 
    # compute covariance for each triangle 
    cov = np.cov(d, ddof=0, rowvar=False) 

可以預先計算，如下所示：

diffs = data - data.mean(1,keepdims=True) 
cov3D = np.einsum('ijk,ijl->ikl',diffs,diffs)/data.shape[1] 

這些預計算值可以跨迭代被用來計算像這樣的本徵矢量 -

for i,d in enumerate(data): 
    # Directly use pre-computed covariances for each triangle 
    vals, vecs = np.linalg.eig(cov3D[i]) 

測試它！

這裏有一些運行測試，以評估的前期計算的協方差結果的影響 -

In [148]: def original_app(data): 
    ...:  cov = np.empty(data.shape) 
    ...:  for i,d in enumerate(data):  
    ...:   # compute covariance for each triangle 
    ...:   cov[i] = np.cov(d, ddof=0, rowvar=False) 
    ...:  return cov 
    ...: 
    ...: def vectorized_app(data):    
    ...:  diffs = data - data.mean(1,keepdims=True) 
    ...:  return np.einsum('ijk,ijl->ikl',diffs,diffs)/data.shape[1] 
    ...: 

In [149]: data = np.random.randint(0,10,(1000,3,3)) 

In [150]: np.allclose(original_app(data),vectorized_app(data)) 
Out[150]: True 

In [151]: %timeit original_app(data) 
10 loops, best of 3: 64.4 ms per loop 

In [152]: %timeit vectorized_app(data) 
1000 loops, best of 3: 1.14 ms per loop 

In [153]: data = np.random.randint(0,10,(5000,3,3)) 

In [154]: np.allclose(original_app(data),vectorized_app(data)) 
Out[154]: True 

In [155]: %timeit original_app(data) 
1 loops, best of 3: 324 ms per loop 

In [156]: %timeit vectorized_app(data) 
100 loops, best of 3: 5.67 ms per loop 

Answer 2

我不知道你實際能達到多少加速。

這是一個輕微的修改，可以幫助一點：

%timeit -n 10 values, vectors = \ 
    eigen(triangles) 
10 loops, best of 3: 745 ms per loop 

%timeit values, vectors = \ 
    zip(*(np.linalg.eig(np.cov(d, ddof=0, rowvar=False)) for d in triangles)) 
10 loops, best of 3: 705 ms per loop