代码之家  ›  专栏  ›  技术社区  ›  MiniQuark

是否有一个函数可以将NumPy的广播规则应用于形状列表并返回最终形状?

  •  0
  • MiniQuark  · 技术社区  · 7 年前

    这是 关于广播工作原理的问题(即,它不是 these questions ).

    我只想找到一个函数,可以将NumPy的广播规则应用于形状列表,并返回最终形状,例如:

    >>> broadcast_shapes([6], [4, 2, 3, 1], [2, 1, 1])
    [4, 2, 3, 6]
    

    谢谢

    3 回复  |  直到 7 年前
        1
  •  2
  •   Paul Panzer    7 年前

    下面是另一个直接实现,它恰好在示例中击败了其他实现。值得一提的是@hpaulj和@Warren Weckesser的hack,它几乎同样快,而且更简洁:

    def bs_pp(*shapes):
        ml = max(shapes, key=len)
        out = list(ml)
        for l in shapes:
            if l is ml:
                continue
            for i, x in enumerate(l, -len(l)):
                if x != 1 and x != out[i]:
                    if out[i] != 1:
                        raise ValueError
                    out[i] = x
        return (*out,)
    
    def bs_mq1(*shapes):
        max_rank = max([len(shape) for shape in shapes])
        shapes = [[1] * (max_rank - len(shape)) + shape for shape in shapes]
        final_shape = [1] * max_rank
        for shape in shapes:
            for dim, size in enumerate(shape):
                if size != 1:
                    final_size = final_shape[dim]
                    if final_size == 1:
                        final_shape[dim] = size
                    elif final_size != size:
                        raise ValueError("Cannot broadcast these shapes")
        return (*final_shape,)
    
    import numpy as np
    
    def bs_mq2(*shapes):
        max_rank = max([len(shape) for shape in shapes])
        shapes = np.array([[1] * (max_rank - len(shape)) + shape
                          for shape in shapes])
        shapes[shapes==1] = -1
        final_shape = shapes.max(axis=0)
        final_shape[final_shape==-1] = 1
        return (*final_shape,)
    
    def bs_hp_ww(*shapes):
        return np.broadcast(*[np.empty(shape + [0,], int) for shape in shapes]).shape[:-1]
    
    L = [6], [4, 2, 3, 1], [2, 1, 1]
    
    from timeit import timeit
    
    print('pp:       ', timeit(lambda: bs_pp(*L), number=10_000)/10)
    print('mq 1:     ', timeit(lambda: bs_mq1(*L), number=10_000)/10)
    print('mq 2:     ', timeit(lambda: bs_mq2(*L), number=10_000)/10)
    print('hpaulj/ww:', timeit(lambda: bs_hp_ww(*L), number=10_000)/10)
    
    assert bs_pp(*L) == bs_mq1(*L) and bs_pp(*L) == bs_mq2(*L) and bs_pp(*L) == bs_hp_ww(*L)
    

    样本运行:

    pp:        0.0021552839782088993
    mq 1:      0.00398325570859015
    mq 2:      0.01497043427079916
    hpaulj/ww: 0.003267909213900566
    
        2
  •  1
  •   user2357112    7 年前

    我在NumPy文档中没有看到任何与此相关的内容。可以向每个目标形状广播单个0维数组,然后相互广播所有结果:

    def broadcast_shapes(*shapes):
        base = numpy.array(0)
        broadcast1 = [numpy.broadcast_to(base, shape) for shape in shapes]
        return numpy.broadcast(*broadcast1).shape
    

    这样可以避免为大型形状分配大量内存。不过,需要创建阵列感觉有点傻。

        3
  •  1
  •   MiniQuark    7 年前

    下面是一个简单的实现,以防有人需要它(它可能有助于理解广播)。不过我更喜欢使用NumPy函数。

    def broadcast_shapes(*shapes):
        max_rank = max([len(shape) for shape in shapes])
        shapes = [[1] * (max_rank - len(shape)) + shape for shape in shapes]
        final_shape = [1] * max_rank
        for shape in shapes:
            for dim, size in enumerate(shape):
                if size != 1:
                    final_size = final_shape[dim]
                    if final_size == 1:
                        final_shape[dim] = size
                    elif final_size != size:
                        raise ValueError("Cannot broadcast these shapes")
        return final_shape
    

    编辑

    我将这个函数与其他几个答案进行了比较,结果证明它是最快的( 编辑 ,Paul Panzer编写了一个更快的函数,请参见他的答案,我将其添加到下面的列表中):

    %timeit bs_pp(*shapes) # Peter Panzer's answer
    2.33 µs ± 10.2 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
    
    %timeit broadcast_shapes1(*shapes)  # this answer
    4.21 µs ± 11.2 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
    
    %timeit broadcast_shapes2(*shapes) # my other answer with shapes.max(axis=0)
    12.8 µs ± 67.3 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
    
    %timeit broadcast_shapes3(*shapes) # user2357112's answer
    18 µs ± 26.9 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
    
    %timeit broadcast_shapes4(*shapes) # hpaulj's answer
    18.1 µs ± 263 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
    
        4
  •  0
  •   hpaulj    7 年前
    In [120]: shapes = [6], [4, 2, 3, 1], [2, 1, 1]                                 
    In [121]: arrs = np.broadcast_arrays(*[np.empty(shape,int) for shape in shapes])
         ...:                                                                       
    In [122]: [a.shape for a in arrs]                                               
    Out[122]: [(4, 2, 3, 6), (4, 2, 3, 6), (4, 2, 3, 6)]
    
    In [124]: np.lib.stride_tricks._broadcast_shape(*[np.empty(shape,int) for shape 
         ...: in shapes])                                                           
    Out[124]: (4, 2, 3, 6)
    
    In [131]: np.broadcast(*[np.empty(shape,int) for shape in shapes]).shape        
    Out[131]: (4, 2, 3, 6)
    

    第二个速度快了很多,4.79秒对42.4秒。第三个速度快了一点。

    正如我第一次评论的那样,我从 broadcast_arrays ,然后查看代码。那是给我的吗 _broadcast_shape ,然后 np.broadcast .

        5
  •  0
  •   MiniQuark    7 年前

    假设这些形状实际上可以被广播,那么这是可行的:

    def broadcast_shapes(*shapes):
        max_rank = max([len(shape) for shape in shapes])
        shapes = np.array([[1] * (max_rank - len(shape)) + shape
                          for shape in shapes])
        shapes[shapes==1] = -1
        final_shape = shapes.max(axis=0)
        final_shape[final_shape==-1] = 1
        return final_shape
    

    如果假设没有空维度,那么 -1 黑客是没有必要的:

    def broadcast_shapes(*shapes):
        max_rank = max([len(shape) for shape in shapes])
        shapes = np.array([[1] * (max_rank - len(shape)) + shape
                          for shape in shapes])
        return shapes.max(axis=0)
    
    推荐文章