4

我在 Python 中有一个函数。我想让它更快吗？有没有人有任何提示？

def garchModel(e2, omega=0.01, beta=0.1, gamma=0.8 ):

    sigma    = np.empty( len( e2 ) )
    sigma[0] = omega

    for i in np.arange(  1, len(e2) ):

        sigma[i] = omega + beta * sigma[ i-1 ] + gamma * e2[ i-1 ]

    return sigma
4

5 回答 5

5

下面的代码可以工作，但是有太多的技巧，我不确定它是否取决于一些最终可能会崩溃的未记录的实现细节：

from numpy.lib.stride_tricks import as_strided
from numpy.core.umath_tests import inner1d

def garch_model(e2, omega=0.01, beta=0.1, gamma=0.8):
    n = len(e2)
    sigma = np.empty((n,))
    sigma[:] = omega
    sigma[1:] += gamma * e2[:-1]
    sigma_view = as_strided(sigma, shape=(n-1, 2), strides=sigma.strides*2)
    inner1d(sigma_view, [beta, 1], out=sigma[1:])
    return sigma

In [75]: e2 = np.random.rand(1e6)

In [76]: np.allclose(garchModel(e2), garch_model(e2))
Out[76]: True

In [77]: %timeit garchModel(e2)
1 loops, best of 3: 6.93 s per loop

In [78]: %timeit garch_model(e2)
100 loops, best of 3: 17.5 ms per loop
于 2013-10-24T22:28:10.580 回答
3

我尝试使用 Numba，它对我的​​数据集提供了 200 倍的改进！

感谢以上所有建议，但我无法让他们给我正确的答案。我将尝试阅读有关线性滤波器的信息，但现在是星期五晚上，我有点太累了，无法再接收信息了。

from numba import autojit
@autojit
def garchModel2(e2, beta=0.1, gamma=0.8, omega=0.01, ):

    sigma    = np.empty( len( e2 ) )
    sigma[0] = omega

    for i in range(  1, len(e2) ):
        sigma[i] = omega + beta * sigma[ i-1 ] + gamma * e2[ i-1 ]

    return sigma
于 2013-10-25T20:18:00.427 回答
2

这是基于@stx2 想法的解决方案。一个潜在的问题是当变大beta**N时可能会导致浮点溢出（与 相同）。Ncumprod

>>> def garchModel2(e2, omega=0.01, beta=0.1, gamma=0.8):
    wt0=cumprod(array([beta,]*(len(e2)-1)))
    wt1=cumsum(hstack((0.,wt0)))+1
    wt2=hstack((wt0[::-1], 1.))*gamma
    wt3=hstack((1, wt0))[::-1]*beta
    pt1=hstack((0.,(array(e2)*wt2)[:-1]))
    pt2=wt1*omega
    return cumsum(pt1)/wt3+pt2

>>> garchModel([1,2,3,4,5])
array([ 0.01    ,  0.811   ,  1.6911  ,  2.57911 ,  3.467911])  
>>> garchModel2([1,2,3,4,5])
array([ 0.01    ,  0.811   ,  1.6911  ,  2.57911 ,  3.467911])
>>> f1=lambda: garchModel2(range(5)) 
>>> f=lambda: garchModel(range(5))
>>> T=timeit.Timer('f()', 'from __main__ import f')
>>> T1=timeit.Timer('f1()', 'from __main__ import f1')
>>> T.timeit(1000)
0.01588106868331031
>>> T1.timeit(1000) #When e2 dimension is samll, garchModel2 is slower
0.04536693909403766
>>> f1=lambda: garchModel2(range(10000))
>>> f=lambda: garchModel(range(10000))
>>> T.timeit(1000)
35.745981961394534
>>> T1.timeit(1000) #When e2 dimension is large, garchModel2 is faster
1.7330512676890066
>>> f1=lambda: garchModel2(range(1000000))
>>> f=lambda: garchModel(range(1000000))
>>> T.timeit(50)
167.33835501439427
>>> T1.timeit(50) #The difference is even bigger.
8.587259274572716

我没有使用beta**N，而是cumprod。**可能会减慢很多。

于 2013-10-25T00:25:25.267 回答
2

您的计算是序列的线性过滤器omega + gamma*e2，因此您可以使用scipy.signal.lfilter. 这是您的计算的一个版本，对初始条件和过滤器的输入进行了适当的调整，以生成与以下相同的输出garchModel：

import numpy as np
from scipy.signal import lfilter

def garch_lfilter(e2, omega=0.01, beta=0.1, gamma=0.8):
    # Linear filter coefficients:
    b = [1]
    a = [1, -beta]

    # Initial condition for the filter:
    zi = np.array([beta*omega])

    # Preallocate the output array, and set the first value to omega:
    sigma = np.empty(len(e2))
    sigma[0] = omega

    # Apply the filter to omega + gamma*e2[:-1]
    sigma[1:], zo = lfilter(b, a, omega + gamma*e2[:-1], zi=zi)

    return sigma

验证它是否给出与@Jaime 函数相同的结果：

In [6]: e2 = np.random.rand(1e6)

In [7]: np.allclose(garch_model(e2), garch_lfilter(e2))
Out[7]: True

它比 快很多garchModel，但不如@Jaime 的函数快。

@Jaime 的时间安排garch_model：

In [8]: %timeit garch_model(e2)
10 loops, best of 3: 21.6 ms per loop

时间garch_lfilter：

In [9]: %timeit garch_lfilter(e2)
10 loops, best of 3: 26.8 ms per loop
于 2013-10-25T19:50:25.643 回答
1

正如@jaime 所示，有一种方法。但是，我不知道是否有办法重写该函数，使其更快，并保持简单。

然后，另一种方法是使用优化“魔法”，例如cython 或 numba。

于 2013-10-25T15:53:58.450 回答