提高重复检查效率 - Python

Question

我是一名研究人员，使用Python来处理气候模型输出，以查找某些类型的风暴。我有8个大的numpy阵列（尺寸是109574 x 52 x 57）。这些数组填充1表示当天有风暴（第一维是时间），0表示没有风暴。另外两个维度是经度和纬度。

我要消除这些阵列背到回天。例如，如果第1天和第2天有暴风雨，我只想计1次风暴。如果第一天，第二天和第三天有暴风雨，我想只计算1和3共计两次风暴，第1-4天会有两次风暴，依此类推。我在最后使用np.sum在时间轴上统计了阵列中的1，结果发现了风暴＃。

我运行下面的代码来实现这一点，但我面对的问题，这是非常缓慢的。因为我将不得不为其他数据集重复此过程，所以我想知道是否有办法加快此过程以提高效率。我有我的代码在下面，我会很乐意澄清任何事情。

# If there is a storm that overlaps two two-day periods, only count it once 
print("Eliminating doubles...") 
for i in range(52): 
    for j in range(57): 
     print(i,j) 
     for k in range(109573): 
      if((storms1[k,i,j]) == 1 and (storms1[k+1,i,j] == 1)): 
       storms1[k,i,j] = 0 
      if((storms2[k,i,j]) == 1 and (storms2[k+1,i,j] == 1)): 
       storms2[k,i,j] = 0 
      if((storms3[k,i,j]) == 1 and (storms3[k+1,i,j] == 1)): 
       storms3[k,i,j] = 0 
      if((storms4[k,i,j]) == 1 and (storms4[k+1,i,j] == 1)): 
       storms4[k,i,j] = 0 
      if((storms5[k,i,j]) == 1 and (storms5[k+1,i,j] == 1)): 
       storms5[k,i,j] = 0 
      if((storms6[k,i,j]) == 1 and (storms6[k+1,i,j] == 1)): 
       storms6[k,i,j] = 0 
      if((storms7[k,i,j]) == 1 and (storms7[k+1,i,j] == 1)): 
       storms7[k,i,j] = 0 
      if((storms8[k,i,j]) == 1 and (storms8[k+1,i,j] == 1)): 
       storms8[k,i,j] = 0 

在有人建议用循环遍历数组之前，为了提出这个问题，我改变了变量名以简化它们。

感谢您的帮助。

Answer 1

这里是一个矢量化功能，可替换的最内层循环：

def do(KK): 
    # find stretches of ones 
    switch_points = np.where(np.diff(np.r_[0, KK, 0]))[0] 
    switch_points.shape = -1, 2 
    # isolate stretches starting on odd days and create mask 
    odd_starters = switch_points[switch_points[:, 0] % 2 == 1, :] 
    odd_mask = np.zeros((KK.shape[0] + 1,), dtype=KK.dtype) 
    odd_mask[odd_starters] = 1, -1 
    odd_mask = np.add.accumulate(odd_mask[:-1]) 
    # apply global 1,0,1,0,1,0,... mask 
    KK[1::2] = 0 
    # invert stretches starting on odd days 
    KK ^= odd_mask 

呼叫它从外部对环​​的（i和j）中：

do(storms1[:, i, j]) 
do(storms2[:, i, j]) 
etc. 

它将改变阵列到位。

这应该是比循环（两个外循环不会有所作为）快得多。

工作原理：

它发现的起点和那些块的端点。我们知道在每个这样的块中，每一个块都必须是零。 使用全局1,0,1,0,1,0，...掩码算法每隔一天就会清零。

产生

• 正确的结果在块即开始甚至几天
• 块
• 与所述正确图案的块补上十多天启动外没有变化

算法的最后一步是反转这些奇怪的起始块。

Answer 2

使用一维数组，模拟的第一轴的一个例子。首先，找到1的组开始。接下来，找到每个组的长度。最后，计算出活动的基础上你的逻辑数量：

import numpy 

a = numpy.random.randint(0,2,20) 

# Add an initial 0 
a1 = numpy.r_[0, a] 

# Mark the start of each group of 1's 
d1 = numpy.diff(a1) > 0 

# Indices of the start of groups of 1's 
w1 = numpy.arange(len(d1))[d1] 

# Length of each group 
cs = numpy.cumsum(a) 
c = numpy.diff(numpy.r_[cs[w1], cs[-1]+1]) 

# Apply the counting logic 
storms = c - c//2 

print(a) 
>>> array([0, 1, 0, 1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1]) 
print(c) 
>>> array([1, 2, 4, 1, 3]) 
print(storms) 
>>> array([1, 1, 2, 1, 2]) 

可以为您节省更多的内存比我在这里展示通过重用变量名不再需要后，他们等。

Answer 3

所以我想你想：

storms_in[:,i,j] = [0,0,1,1,0,1,1,1,0,1,0,1,1,1,1,0] 
storms_out[:,i,j]= [0,0,1,0,0,1,0,1,0,1,0,1,0,0,1,0] 

这是不你的代码示例是干什么的，但你说你想在你的第二段做。

要做到这一点，你需要两个步骤

def storms_disc(storms): # put the whole array here, boolean-safe 
    z = np.zeros((1,) + storms.shape[1:]) # zero-pads for the ends 
    changes = np.r_[storms.astype('int8') ,z] - np.r_[z, storms.astype('int8')] #find where the weather changes 
    changes=((changes[:-1] == 1) | (changes[1:] == -1)).astype('int8') # reduce dimension 
    return ((np.r_[changes, z] - np.r_[z, changes])[:-1] == 1).astype(storms.dtype) #find the first of successive changes 

该向量化的全过程，而且你只需要调用它的8倍。该astype电话是因为减去布尔导致一个错误，尽管它们的价值是1和0

测试：

storms=np.random.randint(0,2,90).reshape(10,3,3) 
storms.T 

array([[[1, 0, 0, 1, 1, 1, 1, 1, 1, 0], 
     [0, 0, 1, 1, 0, 1, 1, 0, 0, 1], 
     [1, 0, 1, 0, 1, 0, 1, 0, 1, 0]], 

     [[0, 0, 0, 1, 0, 1, 0, 0, 0, 0], 
     [0, 1, 0, 0, 1, 1, 1, 0, 0, 0], 
     [0, 1, 0, 0, 1, 0, 1, 0, 1, 1]], 

     [[0, 1, 0, 1, 0, 1, 1, 0, 0, 0], 
     [0, 1, 0, 1, 0, 1, 0, 0, 1, 1], 
     [0, 0, 0, 1, 1, 1, 0, 0, 1, 0]]], dtype=int8) 

storms_disc(storms).T 

array([[[1, 0, 0, 1, 0, 0, 0, 0, 1, 0], 
     [0, 0, 1, 0, 0, 1, 0, 0, 0, 1], 
     [1, 0, 1, 0, 1, 0, 1, 0, 1, 0]], 

     [[0, 0, 0, 1, 0, 1, 0, 0, 0, 0], 
     [0, 1, 0, 0, 1, 0, 1, 0, 0, 0], 
     [0, 1, 0, 0, 1, 0, 1, 0, 1, 0]], 

     [[0, 1, 0, 1, 0, 1, 0, 0, 0, 0], 
     [0, 1, 0, 1, 0, 1, 0, 0, 1, 0], 
     [0, 0, 0, 1, 0, 1, 0, 0, 1, 0]]], dtype=int8)