2011-11-23 66 views
8

我有一個包含標籤的numpy數組。我想根據其大小和邊界框爲每個標籤計算一個數字。我怎樣才能更有效地編寫這個文件,以便在大型數組(〜15000個標籤)上使用它是現實的?如何提高這個numpy循環的效率

A = array([[ 1, 1, 0, 3, 3], 
      [ 1, 1, 0, 0, 0], 
      [ 1, 0, 0, 2, 2], 
      [ 1, 0, 2, 2, 2]]) 

B = zeros(4) 

for label in range(1, 4): 
    # get the bounding box of the label 
    label_points = argwhere(A == label) 
    (y0, x0), (y1, x1) = label_points.min(0), label_points.max(0) + 1 

    # assume I've computed the size of each label in a numpy array size_A 
    B[ label ] = myfunc(y0, x0, y1, x1, size_A[label]) 
+0

真實用例中的A'有多大? –

+1

球場7000x9000 – ajwood

+0

你有沒有做一些分析,看看你的哪些言論是減慢你的言論?也許它是'myfunc'函數,它可能通過將單獨數組中的y0,x0,y1,x1保存在循環中並僅調用一次函數來實現並行化。否則,如果速度真的很重要,你可能需要考慮是否值得做一些C代碼。在使用numpy數組時,我發現cython非常舒適。 –

回答

7

我是不是真的能夠有效地實施這一使用一些NumPy的矢量化功能,所以也許聰明的Python實現會更快。

def first_row(a, labels): 
    d = {} 
    d_setdefault = d.setdefault 
    len_ = len 
    num_labels = len_(labels) 
    for i, row in enumerate(a): 
     for label in row: 
      d_setdefault(label, i) 
     if len_(d) == num_labels: 
      break 
    return d 

這個函數返回一個字典中的每個標籤映射到它出現在第一行的索引。應用的功能AA.TA[::-1]A.T[::-1]也給你的第一列和最後一排,並柱。

如果您希望使用列表而不是字典,可以使用map(d.get, labels)將字典轉換爲列表。或者,您可以從頭開始使用NumPy數組而不是字典,但只要找到所有標籤,您將無法儘早離開循環。

我會感興趣的是(實際上這會加快您的代碼的速度(而且有多少),但我相信它比您的原始解決方案更快。

+0

這不是最好的方法,但它絕對有效。我原來的方式跑了很長時間,我甚至從未讓它完成(20分鐘後放棄)。我剛剛運行了你的方法,並在6分30秒拿到了。 – ajwood

+0

@ajwood:感謝您的反饋。我知道這不太好,但是我能想到的最簡單的解決方案。如果你想更快地做到這一點,我會建議在Cython中實現它。 –

1

性能瓶頸似乎確實是對argmax的呼叫。它可以通過改變如下循環來避免(只計算Y0,Y1,但很容易推廣到X0,X1):

for label in range(1, 4): 
    comp = (A == label) 
    yminind = comp.argmax(0) 
    ymin = comp.max(0) 
    ymaxind = comp.shape[0] - comp[::-1].argmax(0) 
    y0 = yminind[ymin].min() 
    y1 = ymaxind[ymin].max() 

我不知道有關的性能差異的原因,但其中一個原因可能是像==argmaxmax這樣的所有操作都可以直接從輸入數組的形狀中預分配它們的輸出數組,這對argwhere是不可能的。

5

算法:

  1. 改變陣列爲一名維
  2. 獲得通過argsort()
  3. 獲得的維數組排序的版本sorted_A
  4. 使用其中()和排序索引diff()在sorted_A中查找標籤更改位置
  5. 使用更改位置和排序索引來獲取標籤在一個維度中的原始位置。
  6. 從一維位置計算二維位置。

對於像(7000,9000)這樣的大陣列,可以在30s內完成計算。

這裏是代碼:

import numpy as np 

A = np.array([[ 1, 1, 0, 3, 3], 
      [ 1, 1, 0, 0, 0], 
      [ 1, 0, 0, 2, 2], 
      [ 1, 0, 2, 2, 2]]) 

def label_range(A): 
    from itertools import izip_longest 
    h, w = A.shape 
    tmp = A.reshape(-1) 

    index = np.argsort(tmp) 
    sorted_A = tmp[index] 
    pos = np.where(np.diff(sorted_A))[0]+1 
    for p1,p2 in izip_longest(pos,pos[1:]): 
     label_index = index[p1:p2] 
     y = label_index // w 
     x = label_index % w 

     x0 = np.min(x) 
     x1 = np.max(x)+1 
     y0 = np.min(y) 
     y1 = np.max(y)+1 
     label = tmp[label_index[0]] 

     yield label,x0,y0,x1,y1 

for label,x0,y0,x1,y1 in label_range(A): 
    print "%d:(%d,%d)-(%d,%d)" % (label, x0,y0,x1,y1) 

#B = np.random.randint(0, 100, (7000, 9000)) 
#list(label_range(B)) 
+0

我意外地低估了你的帖子,因爲我認爲算法錯了。我不得不做一個虛擬編輯來解鎖投票 - 取而代之的是upvoted。 :) –

5

另一種方法:

使用bincount()來獲取標籤中每行每列數,並保存在行列數陣列的信息。

對於每個標籤,您只需要在行和列中搜索範圍。它比排序更快,在我的電腦上,它可以在幾秒鐘內完成計算。

def label_range2(A): 
    maxlabel = np.max(A)+1 
    h, w = A.shape 
    rows = np.zeros((h, maxlabel), np.bool) 
    for row in xrange(h): 
     rows[row,:] = np.bincount(A[row,:], minlength=maxlabel) > 0 

    cols = np.zeros((w, maxlabel), np.bool) 
    for col in xrange(w): 
     cols[col,:] =np.bincount(A[:,col], minlength=maxlabel) > 0 

    for label in xrange(1, maxlabel): 
     row = rows[:, label] 
     col = cols[:, label] 
     y = np.where(row)[0] 
     x = np.where(col)[0] 
     x0 = np.min(x) 
     x1 = np.max(x)+1 
     y0 = np.min(y) 
     y1 = np.max(y)+1   
     yield label, x0,y0,x1,y1 
+0

這看起來很有希望,我會盡快嘗試。 – ajwood

1

使用PyPy,你可以運行循環,而不用擔心矢量化。它應該很快。