快速累积和和幂算子答案

【问题标题】：Fast cumulative sum and power operator快速累积和和幂算子
【发布时间】：2019-04-18 11:14:23
【问题描述】：

我有一个预测算法，它使用以下代码处理时间序列的趋势直到给定范围：

import numpy as np
horizon = 91
phi = 0.2
trend = -0.004
trend_up_to_horizon = np.cumsum(phi ** np.arange(horizon) + 1) * self.trend

在此示例中，前两个 trend_up_horizon 值是：

array([-0.008 , -0.0128])

有没有一种计算速度更快的方法来实现这一点？目前这需要很长时间，因为我猜使用np.cumsum 方法和** 运算符很昂贵。

感谢您的帮助

【问题讨论】：

您的实际数据集的horizon 值是多少？
@Divakar 地平线值为horizon = 91。
什么是“长时间”？我得到：%timeit np.cumsum(0.2 ** np.arange(91) + 1) 是 7.2 微秒，没有 cumsum 需要 3.9 微秒。也许迭代器可以用更少的乘法更快地生成数组？

标签： python numpy cumulative-sum

【解决方案1】：

你可以使用 Cython 让它稍微快一点，但这并不多

在基本的 np.cumsum(phi ** np.arange(horizon) + 1) * trend 上运行 %timeit 表示在我的笔记本电脑上需要 17.5µs，这并不多

具有相同功能的 Cython 版本是：

import numpy as np
cimport numpy as np
cimport cython

@cython.boundscheck(False)
def do_cumsum(size_t horizon, double phi, double trend):
    cdef np.ndarray[double, ndim=1] out = np.empty(horizon, dtype=np.float)
    cdef double csum = 0
    cdef int i

    for i in range(horizon):
        csum += phi ** i + 1
        out[i] = csum * trend

    return out

这将do_cumsum(horizon, phi, trend) 的时间减少到 6.9µs，而如果我切换到单精度/32 位浮点数，这将下降到 4.5µs

也就是说，微秒并不多，您最好将精力集中在其他地方

【讨论】：

【解决方案2】：

您可以更快地执行此操作。正如您已经假设（不必要的）电源运算符是这里的主要问题。

除此之外，Numpy 没有 power(float64,int64) 的特殊实现，其中指数是一个小的正整数。相反，Numpy 总是计算 power(float64,float64)，这是一项复杂得多的任务。

Numba 和 Numexpr 对简单案例 power(float64,int64) 有一个特殊的实现，所以让我们首先尝试一下。

第一种方法

import numpy as np
import numba as nb

horizon = 91
phi = 0.2
trend = -0.004

@nb.njit()
def pow_cumsum(horizon,phi,trend):
    out=np.empty(horizon)
    csum=0.
    for i in range(horizon):
        csum+=phi**i+1
        out[i]=csum*trend
    return out

如前所述，不需要直接计算幂，可以重写算法来完全避免这种情况。

第二种方法

@nb.njit()
def pow_cumsum_2(horizon,phi,trend):
    out=np.empty(horizon)

    out[0]=2.*trend
    TMP=2.
    val=phi
    for i in range(horizon-1):
        TMP=(val+1+TMP)
        out[i+1]=TMP*trend
        val*=phi
    return out

时间

%timeit np.cumsum(phi ** np.arange(horizon) + 1) * trend
7.44 µs ± 89.1 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
%timeit pow_cumsum(horizon,phi,trend)
882 ns ± 4.91 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)
%timeit pow_cumsum_2(horizon,phi,trend)
559 ns ± 3.8 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)

【讨论】：