2017-07-05 35 views
0

我寫了一個python程序,我嘗試cythonize。 有什麼建議如何讓for循環更高效,因爲這需要99%的時間?使我的代碼更有效

這是for循環:

for i in range(l): 
     b1[i] = np.nanargmin(locator[i,:]) # Closer point 
     locator[i, b1[i]] = NAN # Do not consider Closer point 
     b2[i] = np.nanargmin(locator[i,:]) # 2nd Closer point 
     Adjacents[i,0] = np.array((Existed_Pips[b1[i]]), dtype=np.double) 
     Adjacents[i,1] = np.array((Existed_Pips[b2[i]]), dtype=np.double) 

這是代碼的其餘部分:

import numpy as np 
cimport numpy as np 
from libc.math cimport NAN #, isnan 

def PIPs(np.ndarray[np.double_t, ndim=1, mode='c'] ys, unsigned int nofPIPs, unsigned int typeofdist): 
    cdef: 
     unsigned int currentstate, j, i 
     np.ndarray[np.double_t, ndim=1, mode="c"] D 
     np.ndarray[np.int64_t, ndim=1, mode="c"] Existed_Pips 
     np.ndarray[np.int_t, ndim=1, mode="c"] xs 
     np.ndarray[np.double_t, ndim=2] Adjacents, locator, Adjy, Adjx, Raw_Fire_PIPs, Raw_Fem_PIPs 
     np.ndarray[np.int_t, ndim=2, mode="c"] PIP_points, b1, b2 

    cdef unsigned int l = len(ys) 
    xs = np.arange(0,l, dtype=np.int) # Column vector with xs 
    PIP_points = np.zeros((l,1), dtype=np.int) # Binary indexation 
    PIP_points[0] = 1 # One indicate the PIP points.The first two PIPs are the first and the last observation. 
    PIP_points[-1] = 1 
    Adjacents = np.zeros((l,2), dtype=np.double) 
    currentstate = 2 # Initial PIPs 

    while currentstate <= nofPIPs: # for eachPIPs in range(nofPIPs) 
     Existed_Pips = np.flatnonzero(PIP_points) 
     currentstate = len(Existed_Pips) 
     locator = np.full((l,currentstate), NAN, dtype=np.double) #np.int* 
     for j in range(currentstate): 
      locator[:,j] = np.absolute(xs-Existed_Pips[j]) 
     b1 = np.zeros((l,1), dtype=np.int) 
     b2 = np.zeros((l,1), dtype=np.int) 
     for i in range(l): 
      b1[i] = np.nanargmin(locator[i,:]) # Closer point 
      locator[i, b1[i]] = NAN # Do not consider Closer point 
      b2[i] = np.nanargmin(locator[i,:]) # 2nd Closer point 
      Adjacents[i,0] = np.array((Existed_Pips[b1[i]]), dtype=np.double) 
      Adjacents[i,1] = np.array((Existed_Pips[b2[i]]), dtype=np.double) 

     ##Calculate Distance 
     Adjx = Adjacents   
     Adjy = np.array([ys[np.array(Adjacents[:,0], dtype=np.int)], ys[np.array(Adjacents[:,1], dtype=np.int)]]).transpose() 
     Adjx[Existed_Pips,:] = NAN # Existed PIPs are not candidates for new PIP. 
     Adjy[Existed_Pips,:] = NAN 

     if typeofdist == 1: #Euclidean Distance 
      ##[D] = EDist(ys,xs,Adjx,Adjy) 
      ED = np.power(np.power((Adjx[:,1]-xs),2) + np.power((Adjy[:,1]-ys),2),(0.5)) + np.power(np.power((Adjx[:,0]-xs),2) + np.power((Adjy[:,0]-ys),2),(0.5)) 

     EDmax = np.nanargmax(ED) 
     PIP_points[EDmax]=1 

     currentstate=currentstate+1 

    return np.array([Existed_Pips, ys[Existed_Pips]]).transpose() 

回答

1

一對夫婦的建議:

  1. 採取調用np.nanargmin跳出循環(使用axis參數讓您在整個陣列上進行操作一次。這減少了Python函數調用的次數,你必須做:

    b1 = np.nanargmin(locator,axis=1) 
    locator[np.arange(locator.shape[0]),b1] = np.nan 
    b2 = np.nanargmin(locator,axis=1) 
    
  2. 你分配Adjacents是奇 - 你似乎是首先創建一個長度爲1的陣列的右側。而不是僅僅做

    Adjacents[i,0] = Existed_Pips[b1[i]] 
    # ... 
    

    然而,在這種情況下,你也可以採取兩條線外循環,消除了整個循環:

    Adjacents = np.vstack((Existing_Pips[b1], Existings_Pips[b2])).T 
    

所有這一切都依賴於numpy的,而比Cython,加速,但它可能擊敗你的版本。

+0

非常感謝!在一些大數據集中運行多次時,運行時間從55秒增加到5.8秒,所以這個工作非常好! – Kristin

+0

這比我預期的還要好。我認爲這一定很好。 – DavidW