2015-07-21 118 views
0

對於200000浮點數據集,此代碼需要半個多小時。可以優化這個動態編程代碼嗎?

import numpy as np 
try: 
    import progressbar 
    pbar = progressbar.ProgressBar(widgets=[progressbar.Percentage(), 
     progressbar.Counter('%5d'), progressbar.Bar(), progressbar.ETA()]) 
except: 
    pbar = list 

block_length = np.loadtxt('bb.txt.gz') # get data file from http://filebin.ca/29LbYfKnsKqJ/bb.txt.gz (2MB, 200000 float numbers) 
N = len(block_length) - 1 

# arrays to store the best configuration 
best = np.zeros(N, dtype=float) 
last = np.zeros(N, dtype=int) 
log = np.log 

# Start with first data cell; add one cell at each iteration 
for R in pbar(range(N)): 
    # Compute fit_vec : fitness of putative last block (end at R) 
    #fit_vec = fitfunc.fitness(
    T_k = block_length[:R + 1] - block_length[R + 1] 
    #N_k = np.cumsum(x[:R + 1][::-1])[::-1] 
    N_k = np.arange(R + 1, 0, -1) 
    fit_vec = N_k * (log(N_k) - log(T_k)) 

    prior = 4 - log(73.53 * 0.05 * ((R+1) ** -0.478)) 
    A_R = fit_vec - prior #fitfunc.prior(R + 1, N) 

    A_R[1:] += best[:R] 

    i_max = np.argmax(A_R) 
    last[R] = i_max 
    best[R] = A_R[i_max] 

# Now find changepoints by iteratively peeling off the last block 
change_points = np.zeros(N, dtype=int) 
i_cp = N 
ind = N 
while True: 
    i_cp -= 1 
    change_points[i_cp] = ind 
    if ind == 0: 
     break 
    ind = last[ind - 1] 
    change_points = change_points[i_cp:] 

print edges[change_points] # show result 

第一個循環非常慢,因爲在每次迭代中數組長度爲R,即增加,導致N^2複雜度。

有什麼辦法可以進一步優化這個代碼,例如:通過預先計算?我也很滿意使用其他編程語言的解決方案。

+1

http://codereview.stackexchange.com/ – Mihai

+0

[跨代發表於代碼評論](http://codereview.stackexchange.com/q/97539/9357) –

+0

這並沒有得到太多的關注,因爲它缺乏一個關鍵標籤,'numpy'。 CR上沒有那麼多'numpy'知識豐富的海報。加CR對於問題格式更挑剔。但我同意CR的意見,這個問題需要更多的解釋。它也應該有一個小的測試數據集。 'progressbar'是一個不必要的複雜因素。 – hpaulj

回答

0

我可以複製A_R(直至fit-prior步驟),爲上三角NxN矩陣與:

def trilog(n): 
    nn = n[:-1,None]-n[None,1:] 
    nn[np.tril_indices_from(nn,-1)]=1 
    return nn 

T_k = trilog(block_length) 
N_k = trilog(-np.arange(N+1)) 
fit_vec = N_k * (np.log(N_k) - np.log(T_k)) 
R = np.arange(N)+1 
prior = 4 - log(73.53 * 0.05 * (R ** -0.478)) 
A_R = fit_vec - prior 
A_R = np.triu(A_R,0) 
print(A_R) 

我還沒有經計算的邏輯和應用best工作。

我只對小數組做過這件事。對於你的完整問題,相應的矩陣對於我的記憶來說太大了。

B=np.ones((200000,200000),float) 

所以,僅僅從內存的考慮,你可能會堅持迭代for R in range(N)