我在尋找優化 Python 中的三重回圈的方法時遇到了麻煩。我將直接給出代碼,以便更好、更簡單地表示我必須計算的內容:
給定兩個名為樣本(M x N) 和D (N x N)的二維陣列以及輸出結果 (NxN):
for sigma in range(M):
for i in range(N):
for j in range(N):
results[i, j] = (1/N) * (samples[sigma, i]*samples[sigma, j]
- samples[sigma, i]*D[j, i]
- samples[sigma, j]*D[i, j])
return results
它可以完成這項作業,但在 python 中根本無效。我試圖解開 for i.. for j.. 回圈,但我無法用 sigma 正確計算它。
有人知道如何優化這幾行嗎?歡迎提出任何建議,例如 numpy、numexpr 等...
uj5u.com熱心網友回復:
我發現改進代碼(即減少回圈次數)的一種方法是使用np.meshgrid.
這是我發現的改進。這需要一些擺弄,但它提供與您的三重回圈代碼相同的輸出。我保持了相同的代碼結構,所以你可以看到哪些部分對應于哪些部分。我希望這對你有用!
for sigma in range(M):
xx, yy = np.meshgrid(samples[sigma], samples[sigma])
results = (1/N) * (xx * yy
- yy * D.T
- xx * D)
print(results) # or return results
.
編輯:這是一個驗證結果是否符合預期的小腳本:
import numpy as np
M, N = 3, 4
rng = np.random.default_rng(seed=42)
samples = rng.random((M, N))
D = rng.random((N, N))
results = rng.random((N, N))
results_old = results.copy()
results_new = results.copy()
for sigma in range(M):
for i in range(N):
for j in range(N):
results_old[i, j] = (1/N) * (samples[sigma, i]*samples[sigma, j]
- samples[sigma, i]*D[j, i]
- samples[sigma, j]*D[i, j])
print('\n\nresults_old', results_old, sep='\n')
for sigma in range(M):
xx, yy = np.meshgrid(samples[sigma], samples[sigma])
results_new = (1/N) * (xx * yy
- yy * D.T
- xx * D)
print('\n\nresults_new', results_new, sep='\n')
編輯 2:完全擺脫回圈:它有點令人費解,但它基本上做同樣的事情。
M, N = samples.shape
xxx, yyy = np.meshgrid(samples, samples)
split_x = np.array(np.hsplit(np.vsplit(xxx, M)[0], M))
split_y = np.array(np.vsplit(np.hsplit(yyy, M)[0], M))
results = np.sum(
(1/N) * (split_x*split_y
- split_y*D.T
- split_x*D), axis=0)
print(results) # or return results
uj5u.com熱心網友回復:
我發現將問題分解為更小的步驟并對其進行處理更容易,直到我們有一個單一的等式。
從您的原始配方開始:
for sigma in range(M):
for i in range(N):
for j in range(N):
results[i, j] = (1/N) * (samples[sigma, i]*samples[sigma, j]
- samples[sigma, i]*D[j, i]
- samples[sigma, j]*D[i, j])
首先是消除j最內層回圈中的索引。為此,我們開始使用向量而不是單個元素:
for sigma in range(M):
for i in range(N):
results[i, :] = (1/N) * (samples[sigma, i]*samples[sigma, :] - samples[sigma, i]*D[:, i] - samples[sigma, :]*D[i, :])
然后,我們消除第二個回圈,即帶有i索引的回圈。在這一步中,我們開始在矩陣中思考。因此,每個回圈都是“sigma 矩陣”的直接求和。
for sigma in range(M):
results = (1/N) * (samples[sigma, :, np.newaxis] * samples[sigma] - samples[sigma, :, np.newaxis] * D.T - samples[sigma, :] * D)
我強烈建議使用此步驟作為解決方案,因為對于較大的M. 但是,只是為了知識...
將矩陣視為 3 維物件。我們在索引零的末尾進行計算和求和,如下所示:
results = (1/N) * (samples[:, :, np.newaxis] * samples[:,np.newaxis] - samples[:, :, np.newaxis] * D.T - samples[:, np.newaxis, :] * D).sum(axis=0)
uj5u.com熱心網友回復:
為了向量化for回圈,我們可以利用廣播,然后沿輸出陣列未反映的任何軸減少。為此,我們可以為每個for回圈索引“分配”一個軸(作為慣例)。對于您的示例,這意味著所有輸入陣列都可以重新整形為維度 3(即len(a.shape) == 3);軸sigma, i, j分別對應。然后我們可以對廣播的陣列執行所有操作,最后沿sigma軸對結果進行歸約(求和)(因為僅i, j反映在結果中):
# Ordering of axes: (sigma, i, j)
samples_i = samples[:, :, np.newaxis]
samples_j = samples[:, np.newaxis, :]
D_ij = D[np.newaxis, :, :]
D_ji = D.T[np.newaxis, :, :]
return (samples_i*samples_j - samples_i*D_ji - samples_j*D_ij).sum(axis=0) / N
下面是一個完整的例子,將參考代碼(使用for回圈)與上面的版本進行比較;請注意,我洗掉了該1/N部分以保持整數域中的計算,從而使陣列相等性測驗準確。
import time
import numpy as np
def timeit(func):
def wrapper(*args):
t_start = time.process_time()
res = func(*args)
t_total = time.process_time() - t_start
print(f'{func.__name__}: {t_total:.3f} seconds')
return res
return wrapper
rng = np.random.default_rng()
M, N = 100, 200
samples = rng.integers(0, 100, size=(M, N))
D = rng.integers(0, 100, size=(N, N))
@timeit
def reference(samples, D):
results = np.zeros(shape=(N, N))
for sigma in range(M):
for i in range(N):
for j in range(N):
results[i, j] = (samples[sigma, i]*samples[sigma, j]
- samples[sigma, i]*D[j, i]
- samples[sigma, j]*D[i, j])
return results
@timeit
def new(samples, D):
# Ordering of axes: (sigma, i, j)
samples_i = samples[:, :, np.newaxis]
samples_j = samples[:, np.newaxis, :]
D_ij = D[np.newaxis, :, :]
D_ji = D.T[np.newaxis, :, :]
return (samples_i*samples_j - samples_i*D_ji - samples_j*D_ij).sum(axis=0)
assert np.array_equal(reference(samples, D), new(samples, D))
這給了我以下基準測驗結果:
reference: 6.465 seconds
new: 0.133 seconds
轉載請註明出處,本文鏈接:https://www.uj5u.com/gongcheng/369817.html
