254

是否有用于 Python 的 SciPy 函数或 NumPy 函数或模块来计算给定特定窗口的一维数组的运行平均值?

4

30 回答 30

342

更新:已经提出了更有效的解决方案,uniform_filter1d可能scipy是“标准”第 3 方库中最好的,并且还提供了一些更新或专门的库。


您可以np.convolve为此使用:

np.convolve(x, np.ones(N)/N, mode='valid')

解释

运行均值是卷积数学运算的一种情况。对于运行平均值,您沿着输入滑动一个窗口并计算窗口内容的平均值。对于离散的一维信号,卷积是一样的,除了计算任意线性组合而不是平均值,即,将每个元素乘以相应的系数并将结果相加。这些系数,一个对应于窗口中每个位置的系数,有时称为卷积。N 个值的算术平均值是(x_1 + x_2 + ... + x_N) / N,所以对应的核是(1/N, 1/N, ..., 1/N),这正是我们使用 得到的np.ones(N)/N

边缘

mode参数np.convolve指定如何处理边缘。我在这里选择了这种valid模式,因为我认为这是大多数人期望运行均值工作的方式,但您可能有其他优先事项。这是一个说明模式之间差异的图:

import numpy as np
import matplotlib.pyplot as plt
modes = ['full', 'same', 'valid']
for m in modes:
    plt.plot(np.convolve(np.ones(200), np.ones(50)/50, mode=m));
plt.axis([-10, 251, -.1, 1.1]);
plt.legend(modes, loc='lower center');
plt.show()

运行均值卷积模式

于 2014-03-24T22:01:33.470 回答
180

高效的解决方案

卷积比简单的方法要好得多,但是(我猜)它使用 FFT,因此速度很慢。然而,专门用于计算运行意味着以下方法可以正常工作

def running_mean(x, N):
    cumsum = numpy.cumsum(numpy.insert(x, 0, 0)) 
    return (cumsum[N:] - cumsum[:-N]) / float(N)

要检查的代码

In[3]: x = numpy.random.random(100000)
In[4]: N = 1000
In[5]: %timeit result1 = numpy.convolve(x, numpy.ones((N,))/N, mode='valid')
10 loops, best of 3: 41.4 ms per loop
In[6]: %timeit result2 = running_mean(x, N)
1000 loops, best of 3: 1.04 ms per loop

请注意,两种方法是等效的numpy.allclose(result1, result2)TrueN越大,时间差越大。

警告:虽然 cumsum 更快,但浮点错误会增加,这可能导致您的结果无效/不正确/不可接受

评论在这里指出了这个浮点错误问题,但我在答案中让它更加明显。.

# demonstrate loss of precision with only 100,000 points
np.random.seed(42)
x = np.random.randn(100000)+1e6
y1 = running_mean_convolve(x, 10)
y2 = running_mean_cumsum(x, 10)
assert np.allclose(y1, y2, rtol=1e-12, atol=0)
  • 您累积的点越多,浮点误差就越大(因此 1e5 点很明显,1e6 点更重要,超过 1e6 并且您可能想要重置累加器)
  • 您可以通过使用作弊,np.longdouble但是对于相对大量的点,您的浮点错误仍然会变得很重要(大约> 1e5,但取决于您的数据)
  • 您可以绘制错误并看到它相对快速地增加
  • 卷积解决方案较慢,但没有这种浮点精度损失
  • uniform_filter1d 解决方案比这个 cumsum 解决方案更快并且没有这种浮点精度损失
于 2014-12-28T22:50:57.413 回答
92

更新:下面的示例显示了pandas.rolling_mean已在最新版本的 pandas 中删除的旧功能。该函数调用的现代等价物将使用 pandas.Series.rolling

In [8]: pd.Series(x).rolling(window=N).mean().iloc[N-1:].values
Out[8]: 
array([ 0.49815397,  0.49844183,  0.49840518, ...,  0.49488191,
        0.49456679,  0.49427121])

pandas比 NumPy 或 SciPy 更适合这个。它的功能rolling_mean可以方便地完成这项工作。当输入是一个数组时,它还会返回一个 NumPy 数组。

任何定制的纯 Python 实现都很难rolling_mean在性能上被击败。以下是针对两个建议解决方案的示例性能:

In [1]: import numpy as np

In [2]: import pandas as pd

In [3]: def running_mean(x, N):
   ...:     cumsum = np.cumsum(np.insert(x, 0, 0)) 
   ...:     return (cumsum[N:] - cumsum[:-N]) / N
   ...:

In [4]: x = np.random.random(100000)

In [5]: N = 1000

In [6]: %timeit np.convolve(x, np.ones((N,))/N, mode='valid')
10 loops, best of 3: 172 ms per loop

In [7]: %timeit running_mean(x, N)
100 loops, best of 3: 6.72 ms per loop

In [8]: %timeit pd.rolling_mean(x, N)[N-1:]
100 loops, best of 3: 4.74 ms per loop

In [9]: np.allclose(pd.rolling_mean(x, N)[N-1:], running_mean(x, N))
Out[9]: True

关于如何处理边缘值也有很好的选择。

于 2015-05-09T14:51:13.493 回答
73

您可以使用scipy.ndimage.filters.uniform_filter1d

import numpy as np
from scipy.ndimage.filters import uniform_filter1d
N = 1000
x = np.random.random(100000)
y = uniform_filter1d(x, size=N)

uniform_filter1d

  • 给出具有相同 numpy 形状的输出(即点数)
  • 允许多种方式来处理'reflect'默认的边界,但在我的情况下,我更想要'nearest'

它也相当快(比上面给出的 cumsum 方法快近 50 倍,快np.convolve2-5 倍):

%timeit y1 = np.convolve(x, np.ones((N,))/N, mode='same')
100 loops, best of 3: 9.28 ms per loop

%timeit y2 = uniform_filter1d(x, size=N)
10000 loops, best of 3: 191 µs per loop

这里有 3 个函数可以让您比较不同实现的错误/速度:

from __future__ import division
import numpy as np
import scipy.ndimage.filters as ndif
def running_mean_convolve(x, N):
    return np.convolve(x, np.ones(N) / float(N), 'valid')
def running_mean_cumsum(x, N):
    cumsum = np.cumsum(np.insert(x, 0, 0))
    return (cumsum[N:] - cumsum[:-N]) / float(N)
def running_mean_uniform_filter1d(x, N):
    return ndif.uniform_filter1d(x, N, mode='constant', origin=-(N//2))[:-(N-1)]
于 2017-04-04T07:06:08.567 回答
52

您可以使用以下方法计算运行平均值:

import numpy as np

def runningMean(x, N):
    y = np.zeros((len(x),))
    for ctr in range(len(x)):
         y[ctr] = np.sum(x[ctr:(ctr+N)])
    return y/N

但它很慢。

幸运的是,numpy 包含一个卷积函数,我们可以使用它来加快速度。运行均值相当于x与一个很长的向量进行卷积N,所有成员都等于1/N。convolve 的 numpy 实现包括起始瞬态,因此您必须删除前 N-1 个点:

def runningMeanFast(x, N):
    return np.convolve(x, np.ones((N,))/N)[(N-1):]

在我的机器上,快速版本要快 20-30 倍,具体取决于输入向量的长度和平均窗口的大小。

请注意,convolve 确实包含一个'same'似乎应该解决起始瞬态问题的模式,但它在开始和结束之间拆分了它。

于 2012-12-05T21:21:38.550 回答
29

对于一个简短、快速的解决方案,它在一个循环中完成整个事情,没有依赖关系,下面的代码效果很好。

mylist = [1, 2, 3, 4, 5, 6, 7]
N = 3
cumsum, moving_aves = [0], []

for i, x in enumerate(mylist, 1):
    cumsum.append(cumsum[i-1] + x)
    if i>=N:
        moving_ave = (cumsum[i] - cumsum[i-N])/N
        #can do stuff with moving_ave here
        moving_aves.append(moving_ave)
于 2017-06-28T08:29:40.840 回答
27

或用于计算的python模块

在我在 Tradewave.net 的测试中,TA-lib 总是获胜:

import talib as ta
import numpy as np
import pandas as pd
import scipy
from scipy import signal
import time as t

PAIR = info.primary_pair
PERIOD = 30

def initialize():
    storage.reset()
    storage.elapsed = storage.get('elapsed', [0,0,0,0,0,0])

def cumsum_sma(array, period):
    ret = np.cumsum(array, dtype=float)
    ret[period:] = ret[period:] - ret[:-period]
    return ret[period - 1:] / period

def pandas_sma(array, period):
    return pd.rolling_mean(array, period)

def api_sma(array, period):
    # this method is native to Tradewave and does NOT return an array
    return (data[PAIR].ma(PERIOD))

def talib_sma(array, period):
    return ta.MA(array, period)

def convolve_sma(array, period):
    return np.convolve(array, np.ones((period,))/period, mode='valid')

def fftconvolve_sma(array, period):    
    return scipy.signal.fftconvolve(
        array, np.ones((period,))/period, mode='valid')    

def tick():

    close = data[PAIR].warmup_period('close')

    t1 = t.time()
    sma_api = api_sma(close, PERIOD)
    t2 = t.time()
    sma_cumsum = cumsum_sma(close, PERIOD)
    t3 = t.time()
    sma_pandas = pandas_sma(close, PERIOD)
    t4 = t.time()
    sma_talib = talib_sma(close, PERIOD)
    t5 = t.time()
    sma_convolve = convolve_sma(close, PERIOD)
    t6 = t.time()
    sma_fftconvolve = fftconvolve_sma(close, PERIOD)
    t7 = t.time()

    storage.elapsed[-1] = storage.elapsed[-1] + t2-t1
    storage.elapsed[-2] = storage.elapsed[-2] + t3-t2
    storage.elapsed[-3] = storage.elapsed[-3] + t4-t3
    storage.elapsed[-4] = storage.elapsed[-4] + t5-t4
    storage.elapsed[-5] = storage.elapsed[-5] + t6-t5    
    storage.elapsed[-6] = storage.elapsed[-6] + t7-t6        

    plot('sma_api', sma_api)  
    plot('sma_cumsum', sma_cumsum[-5])
    plot('sma_pandas', sma_pandas[-10])
    plot('sma_talib', sma_talib[-15])
    plot('sma_convolve', sma_convolve[-20])    
    plot('sma_fftconvolve', sma_fftconvolve[-25])

def stop():

    log('ticks....: %s' % info.max_ticks)

    log('api......: %.5f' % storage.elapsed[-1])
    log('cumsum...: %.5f' % storage.elapsed[-2])
    log('pandas...: %.5f' % storage.elapsed[-3])
    log('talib....: %.5f' % storage.elapsed[-4])
    log('convolve.: %.5f' % storage.elapsed[-5])    
    log('fft......: %.5f' % storage.elapsed[-6])

结果:

[2015-01-31 23:00:00] ticks....: 744
[2015-01-31 23:00:00] api......: 0.16445
[2015-01-31 23:00:00] cumsum...: 0.03189
[2015-01-31 23:00:00] pandas...: 0.03677
[2015-01-31 23:00:00] talib....: 0.00700  # <<< Winner!
[2015-01-31 23:00:00] convolve.: 0.04871
[2015-01-31 23:00:00] fft......: 0.22306

在此处输入图像描述

于 2015-10-10T15:21:39.173 回答
23

有关即用型解决方案,请参阅https://scipy-cookbook.readthedocs.io/items/SignalSmooth.html。它提供flat窗口类型的运行平均值。请注意,这比简单的自己动手的卷积方法要复杂一些,因为它试图通过反映数据来处理数据开头和结尾的问题(这在您的情况下可能有效,也可能无效。 ..)。

首先,您可以尝试:

a = np.random.random(100)
plt.plot(a)
b = smooth(a, window='flat')
plt.plot(b)
于 2012-12-05T19:23:12.210 回答
17

我知道这是一个老问题,但这是一个不使用任何额外数据结构或库的解决方案。它与输入列表的元素数量呈线性关系,我想不出任何其他方法来提高效率(实际上,如果有人知道分配结果的更好方法,请告诉我)。

注意:使用 numpy 数组而不是列表会快得多,但我想消除所有依赖项。也可以通过多线程执行来提高性能

该函数假设输入列表是一维的,所以要小心。

### Running mean/Moving average
def running_mean(l, N):
    sum = 0
    result = list( 0 for x in l)

    for i in range( 0, N ):
        sum = sum + l[i]
        result[i] = sum / (i+1)

    for i in range( N, len(l) ):
        sum = sum - l[i-N] + l[i]
        result[i] = sum / N

    return result

例子

假设我们有一个列表data = [ 1, 2, 3, 4, 5, 6 ],我们想要计算周期为 3 的滚动平均值,并且您还想要一个与输入列表大小相同的输出列表(这是最常见的情况)。

第一个元素的索引为 0,所以滚动平均值应该在索引 -2、-1 和 0 的元素上计算。显然我们没有 data[-2] 和 data[-1](除非你想使用特殊的边界条件),所以我们假设这些元素是 0。这相当于对列表进行零填充,除了我们实际上不填充它,只跟踪需要填充的索引(从 0 到 N-1)。

所以,对于前 N 个元素,我们只是不断地将元素累加到一个累加器中。

result[0] = (0 + 0 + 1) / 3  = 0.333    ==   (sum + 1) / 3
result[1] = (0 + 1 + 2) / 3  = 1        ==   (sum + 2) / 3
result[2] = (1 + 2 + 3) / 3  = 2        ==   (sum + 3) / 3

从元素 N+1 向前简单的累加是行不通的。我们期望result[3] = (2 + 3 + 4)/3 = 3,但这与(sum + 4)/3 = 3.333.

计算正确值的方法是减去data[0] = 1sum+4从而给出 sum + 4 - 1 = 9

发生这种情况是因为 current sum = data[0] + data[1] + data[2],但对于每个i >= N因为在减法之前sumis也是如此data[i-N] + ... + data[i-2] + data[i-1]

于 2016-11-16T03:57:50.270 回答
15

我觉得这可以使用瓶颈优雅地解决

请参阅下面的基本示例:

import numpy as np
import bottleneck as bn

a = np.random.randint(4, 1000, size=100)
mm = bn.move_mean(a, window=5, min_count=1)
  • “mm”是“a”的移动平均值。

  • “窗口”是要考虑的移动平均值的最大条目数。

  • “min_count”是移动均值考虑的最小条目数(例如,对于前几个元素或数组是否具有 nan 值)。

好的部分是瓶颈有助于处理 nan 值,它也非常有效。

于 2019-01-17T11:47:28.870 回答
9

我还没有检查这有多快,但你可以尝试:

from collections import deque

cache = deque() # keep track of seen values
n = 10          # window size
A = xrange(100) # some dummy iterable
cum_sum = 0     # initialize cumulative sum

for t, val in enumerate(A, 1):
    cache.append(val)
    cum_sum += val
    if t < n:
        avg = cum_sum / float(t)
    else:                           # if window is saturated,
        cum_sum -= cache.popleft()  # subtract oldest value
        avg = cum_sum / float(n)
于 2015-03-09T15:42:17.503 回答
7

聚会有点晚了,但我做了我自己的小函数,它不会用零环绕末端或填充物,然后也可以用来找到平均值。作为进一步的处理,它还在线性间隔的点处重新采样信号。随意自定义代码以获得其他功能。

该方法是具有归一化高斯核的简单矩阵乘法。

def running_mean(y_in, x_in, N_out=101, sigma=1):
    '''
    Returns running mean as a Bell-curve weighted average at evenly spaced
    points. Does NOT wrap signal around, or pad with zeros.
    
    Arguments:
    y_in -- y values, the values to be smoothed and re-sampled
    x_in -- x values for array
    
    Keyword arguments:
    N_out -- NoOf elements in resampled array.
    sigma -- 'Width' of Bell-curve in units of param x .
    '''
    import numpy as np
    N_in = len(y_in)

    # Gaussian kernel
    x_out = np.linspace(np.min(x_in), np.max(x_in), N_out)
    x_in_mesh, x_out_mesh = np.meshgrid(x_in, x_out)
    gauss_kernel = np.exp(-np.square(x_in_mesh - x_out_mesh) / (2 * sigma**2))
    # Normalize kernel, such that the sum is one along axis 1
    normalization = np.tile(np.reshape(np.sum(gauss_kernel, axis=1), (N_out, 1)), (1, N_in))
    gauss_kernel_normalized = gauss_kernel / normalization
    # Perform running average as a linear operation
    y_out = gauss_kernel_normalized @ y_in

    return y_out, x_out

添加正态分布噪声的正弦信号的简单用法: 在此处输入图像描述

于 2016-10-17T12:21:30.433 回答
7

我建议 pandas 更快地执行此操作,而不是 numpy 或 scipy:

df['data'].rolling(3).mean()

这采用“数据”列的 3 个周期的移动平均线 (MA)。您还可以计算移位的版本,例如排除当前单元格的版本(向后移位)可以很容易地计算为:

df['data'].shift(periods=1).rolling(3).mean()
于 2018-06-16T13:29:31.827 回答
7

Python标准库解决方案

这个生成器函数接受一个可迭代的和一个窗口大小N ,并产生窗口内当前值的平均值。它使用 a deque,这是一种类似于列表的数据结构,但针对两个端点pop的快速修改 ( , append)进行了优化。

from collections import deque
from itertools import islice

def sliding_avg(iterable, N):        
    it = iter(iterable)
    window = deque(islice(it, N))        
    num_vals = len(window)

    if num_vals < N:
        msg = 'window size {} exceeds total number of values {}'
        raise ValueError(msg.format(N, num_vals))

    N = float(N) # force floating point division if using Python 2
    s = sum(window)
    
    while True:
        yield s/N
        try:
            nxt = next(it)
        except StopIteration:
            break
        s = s - window.popleft() + nxt
        window.append(nxt)
        

这是实际执行的功能:

>>> values = range(100)
>>> N = 5
>>> window_avg = sliding_avg(values, N)
>>> 
>>> next(window_avg) # (0 + 1 + 2 + 3 + 4)/5
>>> 2.0
>>> next(window_avg) # (1 + 2 + 3 + 4 + 5)/5
>>> 3.0
>>> next(window_avg) # (2 + 3 + 4 + 5 + 6)/5
>>> 4.0
于 2018-11-04T15:05:04.517 回答
5

另一种使用numpy或找到移动平均线的方法pandas

import itertools
sample = [2, 6, 10, 8, 11, 10]
list(itertools.starmap(
    lambda a,b: b/a, 
    enumerate(itertools.accumulate(sample), 1))
)

将打印[2.0, 4.0, 6.0, 6.5, 7.4, 7.833333333333333]

  • 2.0 = (2)/1
  • 4.0 = (2 + 6) / 2
  • 6.0 = (2 + 6 + 10) / 3
  • ...
于 2016-06-03T00:00:18.677 回答
5

上面有很多关于计算运行平均值的答案。我的回答增加了两个额外的功能:

  1. 忽略 nan 值
  2. 计算 N 个相邻值的平均值,不包括感兴趣的值本身

第二个特征对于确定哪些值与总体趋势相差一定量特别有用。

我使用 numpy.cumsum 因为它是最省时的方法(参见上面 Alleo 的回答)。

N=10 # number of points to test on each side of point of interest, best if even
padded_x = np.insert(np.insert( np.insert(x, len(x), np.empty(int(N/2))*np.nan), 0, np.empty(int(N/2))*np.nan ),0,0)
n_nan = np.cumsum(np.isnan(padded_x))
cumsum = np.nancumsum(padded_x) 
window_sum = cumsum[N+1:] - cumsum[:-(N+1)] - x # subtract value of interest from sum of all values within window
window_n_nan = n_nan[N+1:] - n_nan[:-(N+1)] - np.isnan(x)
window_n_values = (N - window_n_nan)
movavg = (window_sum) / (window_n_values)

此代码仅适用于 Ns。可以通过更改 padded_x 和 n_nan 的 np.insert 来调整奇数。

示例输出(黑色为 raw,蓝色为 movavg): 原始数据(黑色)和每个值周围 10 个点的移动平均值(蓝色),不包括该值。 nan 值被忽略。

此代码可以很容易地修改为删除从小于 cutoff = 3 个非 nan 值计算的所有移动平均值。

window_n_values = (N - window_n_nan).astype(float) # dtype must be float to set some values to nan
cutoff = 3
window_n_values[window_n_values<cutoff] = np.nan
movavg = (window_sum) / (window_n_values)

原始数据(黑色)和移动平均值(蓝色),同时忽略任何少于 3 个非 nan 值的窗口

于 2018-07-25T18:52:17.637 回答
5

mab有一条评论埋在上面有这种方法 的答案之一中。bottleneckmove_mean一个简单的移动平均线:

import numpy as np
import bottleneck as bn

a = np.arange(10) + np.random.random(10)

mva = bn.move_mean(a, window=2, min_count=1)

min_count是一个方便的参数,它基本上会将移动平均线带到数组中的那个点。如果您不设置min_count,它将等于window,并且所有window分数都将是nan

于 2018-08-05T21:21:13.447 回答
5

使用@Aikude 的变量,我写了一条线。

import numpy as np

mylist = [1, 2, 3, 4, 5, 6, 7]
N = 3

mean = [np.mean(mylist[x:x+N]) for x in range(len(mylist)-N+1)]
print(mean)

>>> [2.0, 3.0, 4.0, 5.0, 6.0]
于 2019-10-29T08:24:54.543 回答
5

所有上述解决方案都很差,因为它们缺乏

  • 由于原生 python 而不是 numpy 矢量化实现,速度快,
  • 由于使用不当导致的数值稳定性numpy.cumsum,或
  • 由于O(len(x) * w)实现为卷积的速度。

给定

import numpy
m = 10000
x = numpy.random.rand(m)
w = 1000

请注意,x_[:w].sum()等于x[:w-1].sum()。因此,对于第一个平均值,numpy.cumsum(...)添加x[w] / w(通过x_[w+1] / w)和减去0(从x_[0] / w)。这导致x[0:w].mean()

通过 cumsum,您将通过另外x[w+1] / w加减来更新第二个平均值x[0] / w,从而得到x[1:w+1].mean().

这种情况一直持续到x[-w:].mean()达到。

x_ = numpy.insert(x, 0, 0)
sliding_average = x_[:w].sum() / w + numpy.cumsum(x_[w:] - x_[:-w]) / w

该解决方案是矢量化的O(m)、可读的和数值稳定的。

于 2020-03-13T15:30:18.847 回答
4

这个问题现在比 NeXuS 上个月写的时候还要老,但我喜欢他的代码如何处理边缘情况。然而,因为它是一个“简单的移动平均线”,它的结果落后于它们所应用的数据。我认为以比 NumPy 的模式更令人满意的方式处理边缘情况validsame并且full可以通过对convolution()基于方法应用类似的方法来实现。

我的贡献使用了一个中央运行平均值来将其结果与他们的数据保持一致。当可用于全尺寸窗口的点太少时,将根据阵列边缘的连续较小窗口计算运行平均值。[实际上,从越来越大的窗口开始,但这是一个实现细节。]

import numpy as np

def running_mean(l, N):
    # Also works for the(strictly invalid) cases when N is even.
    if (N//2)*2 == N:
        N = N - 1
    front = np.zeros(N//2)
    back = np.zeros(N//2)

    for i in range(1, (N//2)*2, 2):
        front[i//2] = np.convolve(l[:i], np.ones((i,))/i, mode = 'valid')
    for i in range(1, (N//2)*2, 2):
        back[i//2] = np.convolve(l[-i:], np.ones((i,))/i, mode = 'valid')
    return np.concatenate([front, np.convolve(l, np.ones((N,))/N, mode = 'valid'), back[::-1]])

它相对较慢,因为它使用convolve(). 并且可能会被真正的 Pythonista 修饰很多,但是,我相信这个想法是成立的。

于 2017-01-02T00:28:53.450 回答
3

仅使用 Python 标准库(内存高效)

只需给出deque仅使用标准库的另一个版本。令我惊讶的是,大多数答案都使用pandasor numpy

def moving_average(iterable, n=3):
    d = deque(maxlen=n)
    for i in iterable:
        d.append(i)
        if len(d) == n:
            yield sum(d)/n

r = moving_average([40, 30, 50, 46, 39, 44])
assert list(r) == [40.0, 42.0, 45.0, 43.0]

实际上我在 python 文档中找到了另一个实现

def moving_average(iterable, n=3):
    # moving_average([40, 30, 50, 46, 39, 44]) --> 40.0 42.0 45.0 43.0
    # http://en.wikipedia.org/wiki/Moving_average
    it = iter(iterable)
    d = deque(itertools.islice(it, n-1))
    d.appendleft(0)
    s = sum(d)
    for elem in it:
        s += elem - d.popleft()
        d.append(elem)
        yield s / n

然而,在我看来,实现比它应该的要复杂一些。但它必须在标准 python 文档中是有原因的,有人可以评论我的实现和标准文档吗?

于 2018-01-27T02:52:25.213 回答
3

通过阅读其他答案,我认为这不是问题所要求的,但我来到这里需要保持一个不断增长的值列表的运行平均值。

因此,如果您想保留从某处(站点、测量设备等)获取的值的列表并n更新最后值的平均值,您可以使用下面的代码,这样可以最大限度地减少添加新值的工作量要素:

class Running_Average(object):
    def __init__(self, buffer_size=10):
        """
        Create a new Running_Average object.

        This object allows the efficient calculation of the average of the last
        `buffer_size` numbers added to it.

        Examples
        --------
        >>> a = Running_Average(2)
        >>> a.add(1)
        >>> a.get()
        1.0
        >>> a.add(1)  # there are two 1 in buffer
        >>> a.get()
        1.0
        >>> a.add(2)  # there's a 1 and a 2 in the buffer
        >>> a.get()
        1.5
        >>> a.add(2)
        >>> a.get()  # now there's only two 2 in the buffer
        2.0
        """
        self._buffer_size = int(buffer_size)  # make sure it's an int
        self.reset()

    def add(self, new):
        """
        Add a new number to the buffer, or replaces the oldest one there.
        """
        new = float(new)  # make sure it's a float
        n = len(self._buffer)
        if n < self.buffer_size:  # still have to had numbers to the buffer.
            self._buffer.append(new)
            if self._average != self._average:  # ~ if isNaN().
                self._average = new  # no previous numbers, so it's new.
            else:
                self._average *= n  # so it's only the sum of numbers.
                self._average += new  # add new number.
                self._average /= (n+1)  # divide by new number of numbers.
        else:  # buffer full, replace oldest value.
            old = self._buffer[self._index]  # the previous oldest number.
            self._buffer[self._index] = new  # replace with new one.
            self._index += 1  # update the index and make sure it's...
            self._index %= self.buffer_size  # ... smaller than buffer_size.
            self._average -= old/self.buffer_size  # remove old one...
            self._average += new/self.buffer_size  # ...and add new one...
            # ... weighted by the number of elements.

    def __call__(self):
        """
        Return the moving average value, for the lazy ones who don't want
        to write .get .
        """
        return self._average

    def get(self):
        """
        Return the moving average value.
        """
        return self()

    def reset(self):
        """
        Reset the moving average.

        If for some reason you don't want to just create a new one.
        """
        self._buffer = []  # could use np.empty(self.buffer_size)...
        self._index = 0  # and use this to keep track of how many numbers.
        self._average = float('nan')  # could use np.NaN .

    def get_buffer_size(self):
        """
        Return current buffer_size.
        """
        return self._buffer_size

    def set_buffer_size(self, buffer_size):
        """
        >>> a = Running_Average(10)
        >>> for i in range(15):
        ...     a.add(i)
        ...
        >>> a()
        9.5
        >>> a._buffer  # should not access this!!
        [10.0, 11.0, 12.0, 13.0, 14.0, 5.0, 6.0, 7.0, 8.0, 9.0]

        Decreasing buffer size:
        >>> a.buffer_size = 6
        >>> a._buffer  # should not access this!!
        [9.0, 10.0, 11.0, 12.0, 13.0, 14.0]
        >>> a.buffer_size = 2
        >>> a._buffer
        [13.0, 14.0]

        Increasing buffer size:
        >>> a.buffer_size = 5
        Warning: no older data available!
        >>> a._buffer
        [13.0, 14.0]

        Keeping buffer size:
        >>> a = Running_Average(10)
        >>> for i in range(15):
        ...     a.add(i)
        ...
        >>> a()
        9.5
        >>> a._buffer  # should not access this!!
        [10.0, 11.0, 12.0, 13.0, 14.0, 5.0, 6.0, 7.0, 8.0, 9.0]
        >>> a.buffer_size = 10  # reorders buffer!
        >>> a._buffer
        [5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0]
        """
        buffer_size = int(buffer_size)
        # order the buffer so index is zero again:
        new_buffer = self._buffer[self._index:]
        new_buffer.extend(self._buffer[:self._index])
        self._index = 0
        if self._buffer_size < buffer_size:
            print('Warning: no older data available!')  # should use Warnings!
        else:
            diff = self._buffer_size - buffer_size
            print(diff)
            new_buffer = new_buffer[diff:]
        self._buffer_size = buffer_size
        self._buffer = new_buffer

    buffer_size = property(get_buffer_size, set_buffer_size)

您可以使用以下方法对其进行测试:

def graph_test(N=200):
    import matplotlib.pyplot as plt
    values = list(range(N))
    values_average_calculator = Running_Average(N/2)
    values_averages = []
    for value in values:
        values_average_calculator.add(value)
        values_averages.append(values_average_calculator())
    fig, ax = plt.subplots(1, 1)
    ax.plot(values, label='values')
    ax.plot(values_averages, label='averages')
    ax.grid()
    ax.set_xlim(0, N)
    ax.set_ylim(0, N)
    fig.show()

这使:

值及其平均值作为值的函数 #

于 2018-02-07T22:42:45.460 回答
3

出于教育目的,让我再添加两个 Numpy 解决方案(比 cumsum 解决方案慢):

import numpy as np
from numpy.lib.stride_tricks import as_strided

def ra_strides(arr, window):
    ''' Running average using as_strided'''
    n = arr.shape[0] - window + 1
    arr_strided = as_strided(arr, shape=[n, window], strides=2*arr.strides)
    return arr_strided.mean(axis=1)

def ra_add(arr, window):
    ''' Running average using add.reduceat'''
    n = arr.shape[0] - window + 1
    indices = np.array([0, window]*n) + np.repeat(np.arange(n), 2)
    arr = np.append(arr, 0)
    return np.add.reduceat(arr, indices )[::2]/window

使用的函数:as_stridedadd.reduceat

于 2018-11-01T09:36:12.287 回答
2

虽然这里有这个问题的解决方案,但请看看我的解决方案。它非常简单并且运行良好。

import numpy as np
dataset = np.asarray([1, 2, 3, 4, 5, 6, 7])
ma = list()
window = 3
for t in range(0, len(dataset)):
    if t+window <= len(dataset):
        indices = range(t, t+window)
        ma.append(np.average(np.take(dataset, indices)))
else:
    ma = np.asarray(ma)
于 2018-01-12T11:33:24.813 回答
2

移动平均滤波器怎么样?它也是单线的,并且具有优势,如果您需要矩形以外的其他东西,您可以轻松地操作窗口类型,即。数组 a 的 N 长简单移动平均线:

lfilter(np.ones(N)/N, [1], a)[N:]

并应用了三角形窗口:

lfilter(np.ones(N)*scipy.signal.triang(N)/N, [1], a)[N:]

注意:我通常将前 N 个样本作为虚假样本丢弃,因此[N:]最后,但这不是必需的,只是个人选择的问题。

于 2018-07-16T19:48:21.043 回答
2

一个新的convolve配方被合并到 Python 3.10 中。

给定


import collections, operator

from itertools import chain, repeat


size = 3 + 1
kernel = [1/size] * size                                              

代码

def convolve(signal, kernel):
    # See:  https://betterexplained.com/articles/intuitive-convolution/
    # convolve(data, [0.25, 0.25, 0.25, 0.25]) --> Moving average (blur)
    # convolve(data, [1, -1]) --> 1st finite difference (1st derivative)
    # convolve(data, [1, -2, 1]) --> 2nd finite difference (2nd derivative)
    kernel = list(reversed(kernel))
    n = len(kernel)
    window = collections.deque([0] * n, maxlen=n)
    for x in chain(signal, repeat(0, n-1)):
        window.append(x)
        yield sum(map(operator.mul, kernel, window))

演示

list(convolve(range(1, 6), kernel))
# [0.25, 0.75, 1.5, 2.5, 3.5, 3.0, 2.25, 1.25]

细节

卷积是一种可以应用于移动平均的通用数学运算。这个想法是,给定一些数据,您将数据子集(窗口)作为“掩码”或“内核”滑过数据,在每个窗口上执行特定的数学运算。在移动平均线的情况下,内核是平均值:

在此处输入图像描述

您现在可以通过more_itertools.convolve. more_itertools是流行的第三方包;通过安装> pip install more_itertools

于 2021-01-05T14:02:00.150 回答
1

另一个仅使用标准库和双端队列的解决方案:

from collections import deque
import itertools

def moving_average(iterable, n=3):
    # http://en.wikipedia.org/wiki/Moving_average
    it = iter(iterable) 
    # create an iterable object from input argument
    d = deque(itertools.islice(it, n-1))  
    # create deque object by slicing iterable
    d.appendleft(0)
    s = sum(d)
    for elem in it:
        s += elem - d.popleft()
        d.append(elem)
        yield s / n

# example on how to use it
for i in  moving_average([40, 30, 50, 46, 39, 44]):
    print(i)

# 40.0
# 42.0
# 45.0
# 43.0
于 2018-09-11T10:14:27.090 回答
1

我的解决方案基于维基百科的“简单移动平均线”。

from numba import jit
@jit
def sma(x, N):
    s = np.zeros_like(x)
    k = 1 / N
    s[0] = x[0] * k
    for i in range(1, N + 1):
        s[i] = s[i - 1] + x[i] * k
    for i in range(N, x.shape[0]):
        s[i] = s[i - 1] + (x[i] - x[i - N]) * k
    s = s[N - 1:]
    return s

与之前建议的解决方案的比较表明,它比 scipy 最快的解决方案“uniform_filter1d”快两倍,并且具有相同的错误顺序。速度测试:

import numpy as np    
x = np.random.random(10000000)
N = 1000

from scipy.ndimage.filters import uniform_filter1d
%timeit uniform_filter1d(x, size=N)
95.7 ms ± 9.34 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
%timeit sma(x, N)
47.3 ms ± 3.42 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

错误比较:

np.max(np.abs(np.convolve(x, np.ones((N,))/N, mode='valid') - uniform_filter1d(x, size=N, mode='constant', origin=-(N//2))[:-(N-1)]))
8.604228440844963e-14
np.max(np.abs(np.convolve(x, np.ones((N,))/N, mode='valid') - sma(x, N)))
1.41886502547095e-13
于 2021-11-02T10:22:53.923 回答
0

如果您必须对非常小的数组(少于大约 200 个元素)重复执行此操作,我发现仅使用线性代数即可获得最快的结果。最慢的部分是设置乘法矩阵 y,您只需要做一次,但之后可能会更快。

import numpy as np
import random 

N = 100      # window size
size =200     # array length

x = np.random.random(size)
y = np.eye(size, dtype=float)

# prepare matrix
for i in range(size):
  y[i,i:i+N] = 1./N
  
# calculate running mean
z = np.inner(x,y.T)[N-1:]

于 2021-09-23T12:24:16.623 回答
-7

如果您确实选择滚动自己,而不是使用现有库,请注意浮点错误并尽量减少其影响:

class SumAccumulator:
    def __init__(self):
        self.values = [0]
        self.count = 0

    def add( self, val ):
        self.values.append( val )
        self.count = self.count + 1
        i = self.count
        while i & 0x01:
            i = i >> 1
            v0 = self.values.pop()
            v1 = self.values.pop()
            self.values.append( v0 + v1 )

    def get_total(self):
        return sum( reversed(self.values) )

    def get_size( self ):
        return self.count

如果您的所有值的数量级大致相同,那么这将有助于通过始终添加大致相似数量级的值来保持精度。

于 2014-03-24T22:29:47.567 回答