2015-10-21 42 views
1

我正在尋找一種向量化方法來應用將二維數組返回給二維數組的每一行併產生的函數一個三維數組。NumPy:應用將矩陣返回矩陣的每一行的函數的一般向量化方法

更具體地說,我有一個函數,它需要一個長度爲p的向量並返回一個二維數組(m乘n)。以下是我的函數的程式化版本:

import numpy as np 
def test_func(x, m, n): 
    # this function is just an example and does not do anything useful. 
    # but, the dimensions of input and output is what I want to convey. 
    np.random.seed(x.sum()) 
    return np.random.randint(5, size=(m, n)) 

我有一件T由P 2維輸入數據:

t = 5 
p = 6 
input_data = np.arange(t*p).reshape(t, p) 
input_data 
Out[403]: 
array([[ 0, 1, 2, 3, 4, 5], 
     [ 6, 7, 8, 9, 10, 11], 
     [12, 13, 14, 15, 16, 17], 
     [18, 19, 20, 21, 22, 23], 
     [24, 25, 26, 27, 28, 29]]) 

我想申請test_func到input_data的每一行。由於test_func返回一個矩陣,我期望創建一個3維(t乘m乘n)的數組。我可以用下面的代碼產生我想要的結果:

output_data = np.array([test_func(x, m=3, n=2) for x in input_data]) 
output_data 
Out[405]: 
array([[[0, 4], 
     [0, 4], 
     [3, 3], 
     [1, 0]], 

     [[1, 0], 
     [1, 0], 
     [4, 1], 
     [2, 4]], 

     [[3, 3], 
     [3, 0], 
     [1, 4], 
     [0, 2]], 

     [[2, 4], 
     [2, 1], 
     [3, 2], 
     [3, 1]], 

     [[3, 4], 
     [4, 3], 
     [0, 3], 
     [3, 0]]]) 

但是,這段代碼似乎並不是最優的代碼。它有一個明確的降低速度,並使用不必要地分配額外內存的中間列表。所以,我喜歡找到一個矢量化的解決方案。我最好的猜測是以下代碼,但它不起作用。

output = np.apply_along_axis(test_func, m=3, n=2, axis=1, arr=input_data) 
Traceback (most recent call last): 

    File "<ipython-input-406-5bef44da348f>", line 1, in <module> 
    output = np.apply_along_axis(test_func, m=3, n=2, axis=1, arr=input_data) 

    File "C:\Anaconda\lib\site-packages\numpy\lib\shape_base.py", line 117, in apply_along_axis 
    outarr[tuple(i.tolist())] = res 

ValueError: could not broadcast input array from shape (3,2) into shape (3) 

請您建議一種有效的方法來解決這個問題。

UPDATE

下面是我想申請的實際功能。它執行多維古典縮放。問題的目的不是優化函數的內部運作,而是找到一個向量化函數apply的泛化方法。但是,本着全面披露的精神,我將實際功能放在這裏。請注意,此功能僅當p == M *(M-1)/ 2

def mds_classical_scaling(v, m, n):  

    # create a symmetric distance matrix from the elements in vector v 
    D = np.zeros((m, m)) 
    D[np.triu_indices(4, k=1)] = v 
    D = (D + D.T) 

    # Transform the symmetric matrix 
    A = -0.5 * (D**2) 
    # Create centering matrix  
    H = np.eye(m) - np.ones((m, m))/m 
    # Doubly center A and store in B 
    B = H*A*H 

    # B should be positive definite otherwise the function 
    # would not work. 
    mu, V = eig(B) 

    #index of largest eigen values 
    ndx = (-mu).argsort() 

    # calculate the point configuration from largest eigen values 
    # and corresponding eigen vectors 
    Mu1 = diag(mu[ndx][:n]) 
    V1 = V[:, ndx[:n]] 
    X = V1*sqrt(Mu1)  

    return X 

任何性能提升,我從量化得到的是微不足道的比較實際的功能。主要原因是學習:)

+5

你可以使用'np.vectorize'或'np。apply_along_axis'使任意的Python函數以「矢量化」的方式運行,但是這些通用的解決方案與標準的Python for循環相比,性能優勢可以忽略不計。爲了獲得任何有意義的性能改進,您需要具體說明要矢量化的實際功能。 –

+0

謝謝ali_m。你的評論和奧利弗的回答提供了我正在尋找的答案。 – Sina

回答

2

ali_m的評論是專注:對於嚴重的速度增益,您應該更具體地瞭解該功能的功能。

話雖這麼說,如果你仍然想使用np.apply_along_axis得到一個(可能)小速度提升,再考慮(重讀that function's docstring後),您可以輕鬆地

  1. 換你的函數產生一維數組,
  2. 使用np.apply_along_axis與包裝和
  3. 重塑所得陣列:

    def test_func_wrapper(*args, **kwargs): 
        return test_func(*args, **kwargs).ravel() 
    
    output = np.apply_along_axis(test_func_wrapper, m=3, n=2, axis=1, arr=input_data) 
    np.allclose(output.reshape(5,3, -1), output_data) 
    # output: True 
    

請注意,這是一個通用方式來加速這樣的循環。如果您使用更具體到實際問題的功能,您可能會獲得更好的性能。

相關問題