重复,但在 numpy 中以可变大小的块

Repeat but in variable sized chunks in numpy

我有一个数组,它是不同块的串联:

a = np.array([0, 1, 2, 10, 11, 20, 21, 22, 23])
#             >     <  >    <  >            <
chunks = np.array([3, 2, 4])
repeats = np.array([1, 3, 2])

上面示例中以新十年开始的每个片段都是一个单独的“块”,我想重复一下。每个块的大小和重复次数都是已知的。我无法在 kronrepeat 之后进行重塑,因为这些块的大小不同。

我想要的结果是

np.array([0, 1, 2, 10, 11, 10, 11, 10, 11, 20, 21, 22, 23, 20, 21, 22, 23])
# repeats:>  1  <  >         3          <  >              2             <

这很容易在循环中完成:

in_offset = np.r_[0, np.cumsum(chunks[:-1])]
out_offset = np.r_[0, np.cumsum(chunks[:-1] * repeats[:-1])]
output = np.zeros((chunks * repeats).sum(), dtype=a.dtype)
for c in range(len(chunks)):
    for r in range(repeats[c]):
        for i in range(chunks[c]):
            output[out_offset[c] + r * chunks[c] + i] = a[in_offset[c] + i]

这导致以下向量化:

regions = chunks * repeats
index = np.arange(regions.sum())

segments = np.repeat(chunks, repeats)
resets = np.cumsum(segments[:-1])
offsets = np.zeros_like(index)
offsets[resets] = segments[:-1]
offsets[np.cumsum(regions[:-1])] -= chunks[:-1]

index -= np.cumsum(offsets)

output = a[index]

有没有更有效的方法来矢量化这个问题?我们很清楚,我不是在要求代码审查。我对这些函数调用如何协同工作感到满意。我想知道是否可以使用完全不同(更有效)的函数调用组合来实现相同的结果。

这个问题的灵感来自 my answer to this question

    for rep, num in zip(repeats, chunks):
        res.extend(list(range(num))*rep)

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

一种更 numpythonic 的方法来完成你的任务(比其他答案)是:

result = np.concatenate([ np.tile(tbl, rpt) for tbl, rpt in
    zip(np.split(a, np.cumsum(chunks[:-1])), repeats) ])

结果是:

array([ 0,  1,  2, 10, 11, 10, 11, 10, 11, 20, 21, 22, 23, 20, 21, 22, 23])

对于范围数组的那些块,我们可以直接处理输入数组,从而避免最后的索引步骤,这应该会有所改善 -

#  @Divakar
def create_ranges(starts, ends, l):
    clens = l.cumsum()
    ids = np.ones(clens[-1],dtype=int)
    ids[0] = starts[0]
    ids[clens[:-1]] = starts[1:] - ends[:-1]+1
    out = ids.cumsum()
    return out

s = np.r_[0,chunks.cumsum()]
starts = a[np.repeat(s[:-1],repeats)]
l = np.repeat(chunks, repeats)
ends = starts+l
out = create_ranges(starts, ends, l)

比其他答案更“numpythonic”的解决方法是 -

np.concatenate(np.repeat(np.split(a, np.cumsum(chunks))[:-1], repeats))
array([ 0,  1,  2, 10, 11, 10, 11, 10, 11, 20, 21, 22, 23, 20, 21, 22, 23])

注意,没有明确的 for-loops。

np.split 有一个隐式循环,正如@Divakar 所指出的那样)。


编辑:基准测试(MacBook pro 13)-

正如@Mad Physicist 在他的 post.

中指出的那样,Divakar 的解决方案更适用于更大的数组、块和重复

出于参考目的,我在此处对工作解决方案进行了基准测试。:

def MadPhysicist1(a, chunks, repeats):
    in_offset = np.r_[0, np.cumsum(chunks[:-1])]
    out_offset = np.r_[0, np.cumsum(chunks[:-1] * repeats[:-1])]
    output = np.zeros((chunks * repeats).sum(), dtype=a.dtype)
    for c in range(len(chunks)):
        for r in range(repeats[c]):
            for i in range(chunks[c]):
                output[out_offset[c] + r * chunks[c] + i] = a[in_offset[c] + i]
    return output

def MadPhysicist2(a, chunks, repeats):
    regions = chunks * repeats
    index = np.arange(regions.sum())

    segments = np.repeat(chunks, repeats)
    resets = np.cumsum(segments[:-1])
    offsets = np.zeros_like(index)
    offsets[resets] = segments[:-1]
    offsets[np.cumsum(regions[:-1])] -= chunks[:-1]

    index -= np.cumsum(offsets)

    output = a[index]
    return output

def create_ranges(starts, ends, l):
    clens = l.cumsum()
    ids = np.ones(clens[-1],dtype=int)
    ids[0] = starts[0]
    ids[clens[:-1]] = starts[1:] - ends[:-1]+1
    out = ids.cumsum()
    return out

def Divakar(a, chunks, repeats):
    s = np.r_[0, chunks.cumsum()]
    starts = a[np.repeat(s[:-1], repeats)]
    l = np.repeat(chunks, repeats)
    ends = starts+l
    return create_ranges(starts, ends, l)

def Valdi_Bo(a, chunks, repeats):
    return np.concatenate([np.tile(tbl, rpt) for tbl, rpt in
                           zip(np.split(a, np.cumsum(chunks[:-1])), repeats)])

def AkshaySehgal(a, chunks, repeats):
    return np.concatenate(np.repeat(np.split(a, np.cumsum(chunks))[:-1], repeats))

我查看了三种输入大小的计时:~100、~1000 和~10k 元素:

np.random.seed(0xA)
chunksA = np.random.randint(1, 10, size=20)   # ~100 elements
repeatsA = np.random.randint(1, 10, size=20)
arrA = np.random.randint(100, size=chunksA.sum())

np.random.seed(0xB)
chunksB = np.random.randint(1, 100, size=20)  # ~1000 elements
repeatsB = np.random.randint(1, 10, size=20)
arrB = np.random.randint(100, size=chunksB.sum())

np.random.seed(0xC)
chunksC = np.random.randint(1, 100, size=200)  # ~10000 elements
repeatsC = np.random.randint(1, 10, size=200)
arrC = np.random.randint(100, size=chunksC.sum())

以下是一些结果:

|               |    A    |    B    |    C    |
+---------------+---------+---------+---------+
| MadPhysicist1 | 1.92 ms |   16 ms |  159 ms |
| MadPhysicist2 | 85.5 µs |  153 µs |  744 µs |
| Divakar       | 75.9 µs | 95.9 µs |  312 µs |
| Valdi_Bo      |  370 µs |  369 µs |  3.4 ms |
| AkshaySehgal  |  163 µs |  165 µs | 1.24 ms |