基于阈值创建间隔斜坡数组 - Python / NumPy答案

【问题标题】：Creating intervaled ramp array based on a threshold - Python / NumPy基于阈值创建间隔斜坡数组 - Python / NumPy
【发布时间】：2017-09-13 03:18:05
【问题描述】：

我想测量满足某些条件（如停止时钟）的子数组的长度，但一旦条件不再满足，该值应重置为零。所以，结果数组应该告诉我，有多少值满足某些条件（例如值 > 1）：

[0, 0, 2, 2, 2, 2, 0, 3, 3, 0]

结果应该是以下数组：

[0, 0, 1, 2, 3, 4, 0, 1, 2, 0]

在python中可以很方便的定义一个函数，返回对应的numy数组：

def StopClock(signal, threshold=1):

    clock = []
    current_time = 0
    for item in signal:
        if item > threshold:
            current_time += 1
        else:
            current_time = 0
        clock.append(current_time)
    return np.array(clock)

StopClock([0, 0, 2, 2, 2, 2, 0, 3, 3, 0])

但是，我真的不喜欢这个 for 循环，特别是因为这个计数器应该在更长的数据集上运行。我想到了一些np.cumsum 结合np.diff 的解决方案，但是我没有通过重置部分。有人知道上述问题的更优雅的 numpy 风格的解决方案吗？

【问题讨论】：

标签： python arrays numpy

【解决方案1】：

此解决方案使用 pandas 执行groupby：

s = pd.Series([0, 0, 2, 2, 2, 2, 0, 3, 3, 0])
threshold = 0
>>> np.where(
        s > threshold, 
        s
        .to_frame()  # Convert series to dataframe.
        .assign(_dummy_=1)  # Add column of ones.
        .groupby((s.gt(threshold) != s.gt(threshold).shift()).cumsum())['_dummy_']  # shift-cumsum pattern
        .transform(lambda x: x.cumsum()), # Cumsum the ones per group.
        0)  # Fill value with zero where threshold not exceeded.
array([0, 0, 1, 2, 3, 4, 0, 1, 2, 0])

【讨论】：

【解决方案2】：

是的，我们可以使用diff-styled differentiation 和cumsum 以矢量化方式创建这样的间隔斜坡，这对于大型输入数组来说应该非常有效。重置部分通过在每个间隔结束时分配适当的值来处理，使用 cum-summing 的想法在每个间隔结束时重置数字。

这是完成所有这些的一种实现 -

def intervaled_ramp(a, thresh=1):
    mask = a>thresh

    # Get start, stop indices
    mask_ext = np.concatenate(([False], mask, [False] ))
    idx = np.flatnonzero(mask_ext[1:] != mask_ext[:-1])
    s0,s1 = idx[::2], idx[1::2]

    out = mask.astype(int)
    valid_stop = s1[s1<len(a)]
    out[valid_stop] = s0[:len(valid_stop)] - valid_stop
    return out.cumsum()

示例运行 -

Input (a) : 
[5 3 1 4 5 0 0 2 2 2 2 0 3 3 0 1 1 2 0 3 5 4 3 0 1]
Output (intervaled_ramp(a, thresh=1)) : 
[1 2 0 1 2 0 0 1 2 3 4 0 1 2 0 0 0 1 0 1 2 3 4 0 0]

Input (a) : 
[1 1 1 4 5 0 0 2 2 2 2 0 3 3 0 1 1 2 0 3 5 4 3 0 1]
Output (intervaled_ramp(a, thresh=1)) : 
[0 0 0 1 2 0 0 1 2 3 4 0 1 2 0 0 0 1 0 1 2 3 4 0 0]

Input (a) : 
[1 1 1 4 5 0 0 2 2 2 2 0 3 3 0 1 1 2 0 3 5 4 3 0 5]
Output (intervaled_ramp(a, thresh=1)) : 
[0 0 0 1 2 0 0 1 2 3 4 0 1 2 0 0 0 1 0 1 2 3 4 0 1]

Input (a) : 
[1 1 1 4 5 0 0 2 2 2 2 0 3 3 0 1 1 2 0 3 5 4 3 0 5]
Output (intervaled_ramp(a, thresh=0)) : 
[1 2 3 4 5 0 0 1 2 3 4 0 1 2 0 1 2 3 0 1 2 3 4 0 1]

运行时测试

进行公平基准测试的一种方法是在问题中使用发布的样本，并多次平铺并将其用作输入数组。有了这个设置，这是时间安排 -

In [841]: a = np.array([0, 0, 2, 2, 2, 2, 0, 3, 3, 0])

In [842]: a = np.tile(a,10000)

# @Alexander's soln
In [843]: %timeit pandas_app(a, threshold=1)
1 loop, best of 3: 3.93 s per loop

# @Psidom 's soln
In [844]: %timeit stop_clock(a, threshold=1)
10 loops, best of 3: 119 ms per loop

# Proposed in this post
In [845]: %timeit intervaled_ramp(a, thresh=1)
1000 loops, best of 3: 527 µs per loop

【讨论】：

虽然 Alexander 的解决方案很优雅，而 Psidom 的解决方案是最易读的，但由于速度，这个解决方案是完美的。谢谢大家！

【解决方案3】：

另一个 numpy 解决方案：

import numpy as np
a = np.array([0, 0, 2, 2, 2, 2, 0, 3, 3, 0])

def stop_clock(signal, threshold=1):
    mask = signal > threshold
    indices = np.flatnonzero(np.diff(mask)) + 1
    return np.concatenate(list(map(np.cumsum, np.array_split(mask, indices))))

stop_clock(a)
# array([0, 0, 1, 2, 3, 4, 0, 1, 2, 0])

【讨论】：