Как рассчитать скользящее среднее, используя NumPy?

Кажется, не существует функции, которая просто вычисляет скользящее среднее по numpy / scipy, что приводит к запутанным решениям.

У меня вопрос двоякий:

  • Какой самый простой способ (правильно) реализовать скользящую среднюю с помощью numpy?
  • Поскольку это кажется нетривиальным и подверженным ошибкам, есть ли веская причина не включать батареи в этом случае?

22 ответа

Если вам просто нужна простая невзвешенная скользящая средняя, ​​вы можете легко реализовать ее с помощью np.cumsum, что может быть быстрее, чем методы на основе FFT:

РЕДАКТИРОВАТЬ Исправлено ошибочное индексирование, обнаруженное Бином в коде. РЕДАКТИРОВАТЬ

def moving_average(a, n=3) :
    ret = np.cumsum(a, dtype=float)
    ret[n:] = ret[n:] - ret[:-n]
    return ret[n - 1:] / n

>>> a = np.arange(20)
>>> moving_average(a)
array([  1.,   2.,   3.,   4.,   5.,   6.,   7.,   8.,   9.,  10.,  11.,
        12.,  13.,  14.,  15.,  16.,  17.,  18.])
>>> moving_average(a, n=4)
array([  1.5,   2.5,   3.5,   4.5,   5.5,   6.5,   7.5,   8.5,   9.5,
        10.5,  11.5,  12.5,  13.5,  14.5,  15.5,  16.5,  17.5])

Поэтому я думаю, что ответ таков: это действительно легко реализовать, и, возможно, numpy уже немного раздут со специализированным функционалом.

Простой способ достичь этого с помощью np.convolve, Идея, стоящая за этим, заключается в том, чтобы воспользоваться преимуществами способа вычисления дискретной свертки и использовать его для возврата скользящего среднего. Это может быть сделано путем свертки с последовательностью np.ones длины, равной sliding window Длина, которую мы хотим.

Для этого мы можем определить следующую функцию:

def moving_average(x, w):
    return np.convolve(x, np.ones(w), 'valid') / w

Эта функция будет принимать свертку сигнала x и последовательность единиц длины w, Обратите внимание, что выбран mode является valid так что произведение свертки дается только для точек, где сигналы полностью перекрываются.


Случай использования

Давайте посмотрим на некоторые примеры:

x = np.array([5,3,8,10,2,1,5,1,0,2])

Для скользящей средней с окном длины 2 мы бы хотели иметь:

moving_average(x, 2)
# array([4. , 5.5, 9. , 6. , 1.5, 3. , 3. , 0.5, 1. ])

И для окна длины 4:

moving_average(x, 4)
# array([6.5 , 5.75, 5.25, 4.5 , 2.25, 1.75, 2.  ])

подробности

Давайте более подробно рассмотрим способ вычисления дискретной свертки. Следующая функция призвана повторить способ np.convolve вычисляет выходные значения:

def mov_avg(x, w):
    for m in range(len(x)-(w-1)):
        yield sum(np.ones(w) * x[m:m+w]) / w 

Который, для того же примера выше, также даст:

list(mov_avg(x, 2))
# [4.0, 5.5, 9.0, 6.0, 1.5, 3.0, 3.0, 0.5, 1.0]

Итак, что делается на каждом шаге - это поместить внутренний продукт между массивом единиц и текущим окном. В этом случае умножение на np.ones(w) Излишне, учитывая, что мы принимаем sum последовательности.

Ниже приведен пример того, как первые результаты вычисляются так, чтобы он был немного понятнее. Предположим, мы хотим окно w=4:

[1,1,1,1]
[1,2,3,4,5,6,7,8...]
= (1*1 + 1*2 + 1*3 + 1*4) / w = 2.5

И следующий результат будет вычислен как:

  [1,1,1,1]
[1,2,3,4,5,6,7,8...]
= (1*2 + 1*3 + 1*4 + 1*5) / w = 3.5

И так далее, уступая, как упоминалось ранее, скользящее среднее x,

Отсутствие в NumPy конкретной функции, специфичной для предметной области, возможно, связано с дисциплиной основной команды и верностью основной директиве NumPy: предоставлять N-мерный тип массива, а также функции для создания и индексации этих массивов. Как и многие основополагающие цели, эта задача не маленькая, и NumPy делает это блестяще.

(Гораздо больше) SciPy содержит гораздо большую коллекцию доменных библиотек (так называемые подпакеты от разработчиков SciPy) - например, численная оптимизация (оптимизация), обработка сигналов (сигнал) и интегральное исчисление (интеграция).

Я предполагаю, что функция, за которой вы работаете, находится по крайней мере в одном из подпакетов SciPy (возможно, scipy.signal); Тем не менее, я бы сначала посмотрел в коллекции SciPy scikits, определить соответствующие scikit (ы) и искать там интересующую функцию.

Scikits - это независимо разработанные пакеты, основанные на NumPy / SciPy и ориентированные на конкретную техническую дисциплину (например, scikits-image, scikits-learn и т. Д.). Некоторые из них (в частности, великолепный OpenOpt для численной оптимизации) получили высокую оценку, зрелые проекты задолго до того, как выбрать проживание в относительно новой рубрике. Домашняя страница Scikits понравилась выше перечисляющим около 30 таких scikits, хотя по крайней мере некоторые из них больше не находятся в активной разработке.

Следуя этому совету, вы попадете в серию scikits-time; однако этот пакет больше не находится в активной разработке; По сути, Pandas стала, AFAIK, де-факто библиотекой временных рядов на основе NumPy.

У Панд есть несколько функций, которые можно использовать для вычисления скользящего среднего; самый простой из них, вероятно, Rolling_mean, который вы используете так:

>>> # the recommended syntax to import pandas
>>> import pandas as PD
>>> import numpy as NP

>>> # prepare some fake data:
>>> # the date-time indices:
>>> t = PD.date_range('1/1/2010', '12/31/2012', freq='D')

>>> # the data:
>>> x = NP.arange(0, t.shape[0])

>>> # combine the data & index into a Pandas 'Series' object
>>> D = PD.Series(x, t)

Теперь просто вызовите функцию roll_mean, передающую объект Series и размер окна, который в моем примере ниже равен 10 дням.

>>> d_mva = PD.rolling_mean(D, 10)

>>> # d_mva is the same size as the original Series
>>> d_mva.shape
    (1096,)

>>> # though obviously the first w values are NaN where w is the window size
>>> d_mva[:3]
    2010-01-01         NaN
    2010-01-02         NaN
    2010-01-03         NaN

убедитесь, что это сработало - например, сравнили значения 10 - 15 в исходной серии с новой серией, сглаженной с помощью скользящего среднего

>>> D[10:15]
     2010-01-11    2.041076
     2010-01-12    2.041076
     2010-01-13    2.720585
     2010-01-14    2.720585
     2010-01-15    3.656987
     Freq: D

>>> d_mva[10:20]
      2010-01-11    3.131125
      2010-01-12    3.035232
      2010-01-13    2.923144
      2010-01-14    2.811055
      2010-01-15    2.785824
      Freq: D

Функция roll_mean, а также около дюжины или около того других функций неофициально сгруппированы в документации Pandas под рубрикой "Функции движущегося окна"; вторая связанная группа функций в Pandas называется экспоненциально взвешенными функциями (например, ewma, которая вычисляет экспоненциально скользящее взвешенное среднее). Тот факт, что эта вторая группа не включена в первую (функции движущегося окна), возможно, объясняется тем, что экспоненциально-взвешенные преобразования не полагаются на окно фиксированной длины.

Вот несколько способов сделать это, а также некоторые тесты. Лучшие методы - это версии, использующие оптимизированный код из других библиотек. Вbottleneck.move_meanметод, вероятно, лучший во всех отношениях. Вscipy.convolveподход также очень быстрый, расширяемый, синтаксически и концептуально простой, но плохо масштабируется для очень больших значений окна. Вnumpy.cumsum метод хорош, если вам нужен чистый numpy подходить.

Примечание. Некоторые из них (например,bottleneck.move_mean) не по центру, и ваши данные сместятся.

import numpy as np
import scipy as sci
import scipy.signal as sig
import pandas as pd
import bottleneck as bn
import time as time

def rollavg_direct(a,n): 
    'Direct "for" loop'
    assert n%2==1
    b = a*0.0
    for i in range(len(a)) :
        b[i]=a[max(i-n//2,0):min(i+n//2+1,len(a))].mean()
    return b

def rollavg_comprehension(a,n):
    'List comprehension'
    assert n%2==1
    r,N = int(n/2),len(a)
    return np.array([a[max(i-r,0):min(i+r+1,N)].mean() for i in range(N)]) 

def rollavg_convolve(a,n):
    'scipy.convolve'
    assert n%2==1
    return sci.convolve(a,np.ones(n,dtype='float')/n, 'same')[n//2:-n//2+1]  

def rollavg_convolve_edges(a,n):
    'scipy.convolve, edge handling'
    assert n%2==1
    return sci.convolve(a,np.ones(n,dtype='float'), 'same')/sci.convolve(np.ones(len(a)),np.ones(n), 'same')  

def rollavg_cumsum(a,n):
    'numpy.cumsum'
    assert n%2==1
    cumsum_vec = np.cumsum(np.insert(a, 0, 0)) 
    return (cumsum_vec[n:] - cumsum_vec[:-n]) / n

def rollavg_cumsum_edges(a,n):
    'numpy.cumsum, edge handling'
    assert n%2==1
    N = len(a)
    cumsum_vec = np.cumsum(np.insert(np.pad(a,(n-1,n-1),'constant'), 0, 0)) 
    d = np.hstack((np.arange(n//2+1,n),np.ones(N-n)*n,np.arange(n,n//2,-1)))  
    return (cumsum_vec[n+n//2:-n//2+1] - cumsum_vec[n//2:-n-n//2]) / d

def rollavg_roll(a,n):
    'Numpy array rolling'
    assert n%2==1
    N = len(a)
    rolling_idx = np.mod((N-1)*np.arange(n)[:,None] + np.arange(N), N)
    return a[rolling_idx].mean(axis=0)[n-1:] 

def rollavg_roll_edges(a,n):
    # see https://stackru.com/questions/42101082/fast-numpy-roll
    'Numpy array rolling, edge handling'
    assert n%2==1
    a = np.pad(a,(0,n-1-n//2), 'constant')*np.ones(n)[:,None]
    m = a.shape[1]
    idx = np.mod((m-1)*np.arange(n)[:,None] + np.arange(m), m) # Rolling index
    out = a[np.arange(-n//2,n//2)[:,None], idx]
    d = np.hstack((np.arange(1,n),np.ones(m-2*n+1+n//2)*n,np.arange(n,n//2,-1)))
    return (out.sum(axis=0)/d)[n//2:]

def rollavg_pandas(a,n):
    'Pandas rolling average'
    return pd.DataFrame(a).rolling(n, center=True, min_periods=1).mean().to_numpy()

def rollavg_bottlneck(a,n):
    'bottleneck.move_mean'
    return bn.move_mean(a, window=n, min_count=1)

N = 10**6
a = np.random.rand(N)
functions = [rollavg_direct, rollavg_comprehension, rollavg_convolve, 
        rollavg_convolve_edges, rollavg_cumsum, rollavg_cumsum_edges, 
        rollavg_pandas, rollavg_bottlneck, rollavg_roll, rollavg_roll_edges]

print('Small window (n=3)')
%load_ext memory_profiler
for f in functions : 
    print('\n'+f.__doc__+ ' : ')
    %timeit b=f(a,3)

print('\nLarge window (n=1001)')
for f in functions[0:-2] : 
    print('\n'+f.__doc__+ ' : ')
    %timeit b=f(a,1001)

print('\nMemory\n')
print('Small window (n=3)')
N = 10**7
a = np.random.rand(N)
%load_ext memory_profiler
for f in functions[2:] : 
    print('\n'+f.__doc__+ ' : ')
    %memit b=f(a,3)

print('\nLarge window (n=1001)')
for f in functions[2:-2] : 
    print('\n'+f.__doc__+ ' : ')
    %memit b=f(a,1001)

Время, маленькое окно (n=3)

Direct "for" loop : 

4.14 s ± 23.7 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

List comprehension : 
3.96 s ± 27.9 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

scipy.convolve : 
1.07 ms ± 26.7 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

scipy.convolve, edge handling : 
4.68 ms ± 9.69 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

numpy.cumsum : 
5.31 ms ± 5.11 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

numpy.cumsum, edge handling : 
8.52 ms ± 11.1 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

Pandas rolling average : 
9.85 ms ± 9.63 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

bottleneck.move_mean : 
1.3 ms ± 12.2 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

Numpy array rolling : 
31.3 ms ± 91.9 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

Numpy array rolling, edge handling : 
61.1 ms ± 55.9 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

Время, большое окно (n=1001)

Direct "for" loop : 
4.67 s ± 34 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

List comprehension : 
4.46 s ± 14.6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

scipy.convolve : 
103 ms ± 165 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

scipy.convolve, edge handling : 
272 ms ± 1.23 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

numpy.cumsum : 
5.19 ms ± 12.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

numpy.cumsum, edge handling : 
8.7 ms ± 11.5 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

Pandas rolling average : 
9.67 ms ± 199 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

bottleneck.move_mean : 
1.31 ms ± 15.7 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

Память, маленькое окно (n=3)

The memory_profiler extension is already loaded. To reload it, use:
  %reload_ext memory_profiler

scipy.convolve : 
peak memory: 362.66 MiB, increment: 73.61 MiB

scipy.convolve, edge handling : 
peak memory: 510.24 MiB, increment: 221.19 MiB

numpy.cumsum : 
peak memory: 441.81 MiB, increment: 152.76 MiB

numpy.cumsum, edge handling : 
peak memory: 518.14 MiB, increment: 228.84 MiB

Pandas rolling average : 
peak memory: 449.34 MiB, increment: 160.02 MiB

bottleneck.move_mean : 
peak memory: 374.17 MiB, increment: 75.54 MiB

Numpy array rolling : 
peak memory: 661.29 MiB, increment: 362.65 MiB

Numpy array rolling, edge handling : 
peak memory: 1111.25 MiB, increment: 812.61 MiB

Память, большое окно (n=1001)

scipy.convolve : 
peak memory: 370.62 MiB, increment: 71.83 MiB

scipy.convolve, edge handling : 
peak memory: 521.98 MiB, increment: 223.18 MiB

numpy.cumsum : 
peak memory: 451.32 MiB, increment: 152.52 MiB

numpy.cumsum, edge handling : 
peak memory: 527.51 MiB, increment: 228.71 MiB

Pandas rolling average : 
peak memory: 451.25 MiB, increment: 152.50 MiB

bottleneck.move_mean : 
peak memory: 374.64 MiB, increment: 75.85 MiB

Начиная с Numpy 1.20, то предоставляет способ скользить / катиться по окнам элементов. Окна, которые вы можете затем индивидуально усреднить.

Например, для 4-элементное окно:

      from numpy.lib.stride_tricks import sliding_window_view

# values = np.array([5, 3, 8, 10, 2, 1, 5, 1, 0, 2])
np.average(sliding_window_view(values, window_shape = 4), axis=1)
# array([6.5, 5.75, 5.25, 4.5, 2.25, 1.75, 2])

Обратите внимание на промежуточный результат sliding_window_view:

      # values = np.array([5, 3, 8, 10, 2, 1, 5, 1, 0, 2])
sliding_window_view(values, window_shape = 4)
# array([[ 5,  3,  8, 10],
#        [ 3,  8, 10,  2],
#        [ 8, 10,  2,  1],
#        [10,  2,  1,  5],
#        [ 2,  1,  5,  1],
#        [ 1,  5,  1,  0],
#        [ 5,  1,  0,  2]])

Этот ответ с использованием панд адаптирован сверху, так как rolling_mean больше не является частью панд

# the recommended syntax to import pandas
import pandas as pd
import numpy as np

# prepare some fake data:
# the date-time indices:
t = pd.date_range('1/1/2010', '12/31/2012', freq='D')

# the data:
x = np.arange(0, t.shape[0])

# combine the data & index into a Pandas 'Series' object
D = pd.Series(x, t)

Теперь просто вызовите функцию rolling на фрейме данных с размером окна, который в моем примере ниже составляет 10 дней.

d_mva10 = D.rolling(10).mean()

# d_mva is the same size as the original Series
# though obviously the first w values are NaN where w is the window size
d_mva10[:11]

2010-01-01    NaN
2010-01-02    NaN
2010-01-03    NaN
2010-01-04    NaN
2010-01-05    NaN
2010-01-06    NaN
2010-01-07    NaN
2010-01-08    NaN
2010-01-09    NaN
2010-01-10    4.5
2010-01-11    5.5
Freq: D, dtype: float64

Я чувствую, что это может быть легко решено с помощью узкого места

См основной образец ниже:

import numpy as np
import bottleneck as bn

a = np.random.randint(4, 1000, size=(5, 7))
mm = bn.move_mean(a, window=2, min_count=1)

Это дает среднее значение перемещения по каждой оси.

  • "мм" - это скользящее среднее для "а".

  • "Окно" - это максимальное количество записей, которые нужно учитывать для скользящего среднего.

  • "min_count" - это минимальное количество записей, которые нужно учитывать для скользящего среднего (например, для первого элемента или если массив имеет значения nan).

Хорошая часть заключается в том, что Bottleneck помогает справляться со значениями наночастиц, а также очень эффективен.

если кому-то нужно простое решение, вот одно

      def moving_average(a,n):
    N=len(a)
    return np.array([np.mean(a[i:i+n]) for i in np.arange(0,N-n+1)])

вы можете изменить перекрытие между окнами, добавив аргумент шага в np.arange(0,N-n+1,step)

Если вы хотите тщательно позаботиться о краевых условиях (вычислите среднее значение только из доступных элементов по краям), следующая функция поможет.

import numpy as np

def running_mean(x, N):
    out = np.zeros_like(x, dtype=np.float64)
    dim_len = x.shape[0]
    for i in range(dim_len):
        if N%2 == 0:
            a, b = i - (N-1)//2, i + (N-1)//2 + 2
        else:
            a, b = i - (N-1)//2, i + (N-1)//2 + 1

        #cap indices to min and max indices
        a = max(0, a)
        b = min(dim_len, b)
        out[i] = np.mean(x[a:b])
    return out

>>> running_mean(np.array([1,2,3,4]), 2)
array([1.5, 2.5, 3.5, 4. ])

>>> running_mean(np.array([1,2,3,4]), 3)
array([1.5, 2. , 3. , 3.5])

Все ответы, кажется, сосредоточены на случае предварительно вычисленного списка. Для фактического рабочего варианта использования, когда числа приходят один за другим, вот простой класс, который предоставляет услугу усреднения последних N значений:

      import numpy as np
class RunningAverage():
    def __init__(self, stack_size):
        self.stack = [0 for _ in range(stack_size)]
        self.ptr = 0
        self.full_cycle = False
    def add(self,value):
        self.stack[self.ptr] = value
        self.ptr += 1
        if self.ptr == len(self.stack):
            self.full_cycle = True
            self.ptr = 0
    def get_avg(self):
        if self.full_cycle:
            return np.mean(self.stack)
        else:
            return np.mean(self.stack[:self.ptr])

Применение:

      N = 50  # size of the averaging window
run_avg = RunningAverage(N)
for i in range(1000):
    value = <my computation>
    run_avg.add(value)
    if i % 20 ==0: # print once in 20 iters:
        print(f'the average value is {run_avg.get_avg()}')

Вы также можете написать собственное расширение Python C.

Это, конечно, не самый простой способ, но он позволит вам работать быстрее и эффективнее использовать память, чем использованиеnp.cumsumкак строительный блок.

      // moving_average.c
#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION
#include <Python.h>
#include <numpy/arrayobject.h>

static PyObject *moving_average(PyObject *self, PyObject *args) {
    PyObject *input;
    int64_t window_size;
    PyArg_ParseTuple(args, "Ol", &input, &window_size);
    if (PyErr_Occurred()) return NULL;
    if (!PyArray_Check(input) || !PyArray_ISNUMBER((PyArrayObject *)input)) {
        PyErr_SetString(PyExc_TypeError, "First argument must be a numpy array with numeric dtype");
        return NULL;
    }
    
    int64_t input_size = PyObject_Size(input);
    double *input_data;
    if (PyArray_AsCArray(&input, &input_data, (npy_intp[]){ [0] = input_size }, 1, PyArray_DescrFromType(NPY_DOUBLE)) != 0) {
        PyErr_SetString(PyExc_TypeError, "Failed to simulate C array of type double");
        return NULL;
    }
    
    int64_t output_size = input_size - window_size + 1;
    PyObject *output = PyArray_SimpleNew(1, (npy_intp[]){ [0] = output_size }, NPY_DOUBLE);
    double *output_data = PyArray_DATA((PyArrayObject *)output);
    
    double cumsum_before = 0;
    double cumsum_after = 0;
    for (int i = 0; i < window_size; ++i) {
        cumsum_after += input_data[i];
    }
    for (int i = 0; i < output_size - 1; ++i) {
        output_data[i] = (cumsum_after - cumsum_before) / window_size;
        cumsum_after += input_data[i + window_size];
        cumsum_before += input_data[i];
    }
    output_data[output_size - 1] = (cumsum_after - cumsum_before) / window_size;

    return output;
}

static PyMethodDef methods[] = {
    {
        "moving_average", 
        moving_average, 
        METH_VARARGS, 
        "Rolling mean of numpy array with specified window size"
    },
    {NULL, NULL, 0, NULL}
};

static struct PyModuleDef moduledef = {
    PyModuleDef_HEAD_INIT,
    "moving_average",
    "C extension for finding the rolling mean of a numpy array",
    -1,
    methods
};

PyMODINIT_FUNC PyInit_moving_average(void) {
    PyObject *module = PyModule_Create(&moduledef);
    import_array();
    return module;
}
  • METH_VARARGSуказывает, что метод принимает только позиционные аргументы.PyArg_ParseTupleпозволяет вам анализировать эти позиционные аргументы.

  • ИспользуяPyErr_SetStringи возвращая NULL из метода, вы можете сигнализировать интерпретатору Python о возникновении исключения из расширения C.

  • PyArray_AsCArrayпозволяет вашему методу быть полиморфным, когда речь идет о dtype входного массива, выравнивании , является ли массив C-непрерывным (см. «Может ли массив numpy 1d не быть непрерывным?» ) и т. д. без необходимости создавать копию массива. Если вместо этого вы использовалиPyArray_DATA, вам придется разобраться с этим самостоятельно.

  • PyArray_SimpleNewпозволяет создать новый массив numpy. Это похоже на использованиеnp.empty. Массив не будет инициализирован и может содержать недетерминированный мусор, который может вас удивить, если вы забудете перезаписать его.

Создание расширения C

      # setup.py
from setuptools import setup, Extension
import numpy

setup(
  ext_modules=[
    Extension(
      'moving_average',
      ['moving_average.c'],
      include_dirs=[numpy.get_include()]
    )
  ]
)

# python setup.py build_ext --build-lib=.

Ориентиры

      import numpy as np

# Our compiled C extension:
from moving_average import moving_average as moving_average_c

# Answer by Jaime using npcumsum
def moving_average_cumsum(a, n) :
    ret = np.cumsum(a, dtype=float)
    ret[n:] = ret[n:] - ret[:-n]
    return ret[n - 1:] / n

# Answer by yatu using np.convolve
def moving_average_convolve(a, n):
    return np.convolve(a, np.ones(n), 'valid') / n

a = np.random.rand(1_000_000)
print('window_size = 3')
%timeit moving_average_c(a, 3)
%timeit moving_average_cumsum(a, 3)
%timeit moving_average_convolve(a, 3)

print('\nwindow_size = 100')
%timeit moving_average_c(a, 100)
%timeit moving_average_cumsum(a, 100)
%timeit moving_average_convolve(a, 100)
      window_size = 3
958 µs ± 4.68 µs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)
4.52 ms ± 15.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
809 µs ± 463 ns per loop (mean ± std. dev. of 7 runs, 1,000 loops each)

window_size = 100
977 µs ± 937 ns per loop (mean ± std. dev. of 7 runs, 1,000 loops each)
6.16 ms ± 19.1 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
14.2 ms ± 12.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
for i in range(len(Data)):
    Data[i, 1] = Data[i-lookback:i, 0].sum() / lookback

Попробуйте этот фрагмент кода. Я думаю, что это проще и работает. Lookback - это окно скользящей средней.

в Data[i-lookback:i, 0].sum() я положил 0 для ссылки на первый столбец набора данных, но вы можете поместить любой столбец, который вам нравится, если у вас есть несколько столбцов.

Talib содержит простой инструмент скользящего среднего, а также другие аналогичные инструменты усреднения (например, экспоненциальное скользящее среднее). Ниже сравнивается метод с некоторыми другими решениями.


%timeit pd.Series(np.arange(100000)).rolling(3).mean()
2.53 ms ± 40.5 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

%timeit talib.SMA(real = np.arange(100000.), timeperiod = 3)
348 µs ± 3.5 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%timeit moving_average(np.arange(100000))
638 µs ± 45.1 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

Одно предостережение: на самом деле должны быть элементы dtype = float. В противном случае возникает следующая ошибка

Исключение: настоящий не двойной

Я использую либо решение принятого ответа, слегка измененное, чтобы иметь ту же длину для вывода, что и ввод, либоpandas'версия, как указано в комментарии к другому ответу. Я резюмирую их здесь с помощью воспроизводимого примера для использования в будущем:

import numpy as np
import pandas as pd

def moving_average(a, n):
    ret = np.cumsum(a, dtype=float)
    ret[n:] = ret[n:] - ret[:-n]
    return ret / n

def moving_average_centered(a, n):
    return pd.Series(a).rolling(window=n, center=True).mean().to_numpy()

A = [0, 0, 1, 2, 4, 5, 4]
print(moving_average(A, 3))    
# [0.         0.         0.33333333 1.         2.33333333 3.66666667 4.33333333]
print(moving_average_centered(A, 3))
# [nan        0.33333333 1.         2.33333333 3.66666667 4.33333333 nan       ]

скользящая средняя

  • переверните массив в i и просто возьмите среднее значение от i до n.

  • используйте понимание списков для создания мини-массивов на лету.

x = np.random.randint(10, size=20)

def moving_average(arr, n):
    return [ (arr[:i+1][::-1][:n]).mean() for i, ele in enumerate(arr) ]
n = 5

moving_average(x, n)

Сравнивая решение ниже с тем, которое использует cumsum of numpy, это занимает почти половину времени. Это связано с тем, что не нужно проходить через весь массив, чтобы произвести суммирование, а затем выполнять все вычитание. Более того, cumsum может быть "опасным", если массив огромен, а число огромно (возможно переполнение). Конечно, и здесь опасность существует, но по крайней мере суммируются только существенные числа.

def moving_average(array_numbers, n):
    if n > len(array_numbers):
      return []
    temp_sum = sum(array_numbers[:n])
    averages = [temp_sum / float(n)]
    for first_index, item in enumerate(array_numbers[n:]):
        temp_sum += item - array_numbers[first_index]
        averages.append(temp_sum / float(n))
    return averages

Скользящее усреднение просто:

      def avg(x,n=1):
    return avg((x[1:]+x[:-1])/2 ,n-1) if n>1 else (x[1:]+x[:-1])/2 

усреднение numpy:

      def avg(x,n=1):
    return np.convolve(x, np.ones(n):Valid)/n

Дело в том, что до 300 точек данных одноядерный ryzen 5 5900 быстрее, выше этого numpy быстрее ...

До 5000 точек данных с использованием многопроцессорной библиотеки 24 потока (процесса) ryzen 5 5900 быстрее. Выше этого снова numpy быстрее.

CuPy — это другая история с использованием графического процессора:

      import cupy as cp
import numpy as np

def avgnp(x,n=1):
    return np.convolve(x,np.ones(n),mode='valid')/n
def avgcp(x,n=1):
    return cp.convolve(x,cp.ones(n)/n,mode='full')

print('Moving Average for, ')
for i in range(3,8):
   
    s=time.time();
    t_cpu=np.linspace(0,2*np.pi,10**i)
    x_cpu=np.sin(t_cpu)
    avgnp(x_cpu,5)
    print(10**i,'data points via Numpy(CPU):\t',time.time()-s,'*')
    
    time.sleep(1)
    
    s=time.time()
    t=cp.linspace(0,2*np.pi,10**i)
    x=cp.sin(t)
    avgcp(x,5)
    print(10**i,'data points via CuPy(GPU) :\t',time.time()-s)

    time.sleep(1)

Результаты:

      Moving Average for, 
1000 data points via Numpy(CPU):     0.0005373954772949219 *
1000 data points via CuPy(GPU) :     0.00930333137512207
10000 data points via Numpy(CPU):    0.00034356117248535156 *
10000 data points via CuPy(GPU) :    0.0007426738739013672
100000 data points via Numpy(CPU):   0.0014045238494873047 *
100000 data points via CuPy(GPU) :   0.0006258487701416016
1000000 data points via Numpy(CPU):  0.014684438705444336 *
1000000 data points via CuPy(GPU) :  0.0006639957427978516
10000000 data points via Numpy(CPU): 0.1448352336883545 *
10000000 data points via CuPy(GPU) : 0.0006761550903320312

Скользящее среднее для массива с Pandas и Numpy

      import pandas as pd
import numpy as np


def moving_average(arr, n):
    s = pd.Series(a)
    return s.rolling(n).mean().to_numpy()


a = np.array([15,20,34,16,57,89])
print(moving_average(a, 3))

Выход

      array([        nan,         nan, 23.        , 23.33333333, 35.66666667,
       54.        ])

Это самое простое решение, которое я нашел.

Если у вас уже есть массив известного размера

      import numpy as np                                         
M=np.arange(12)
                                                               
avg=[]                                                         
i=0
while i<len(M)-2: #for n point average len(M) - (n-1)
        avg.append((M[i]+M[i+1]+M[i+2])/3) #n is denominator                       
        i+=1     
                                                                                                    
print(avg)

Думаю, самый простой способ - использовать np.mean.
В следующем примере мы вычисляем скользящую среднюю. ma над k последние образцы в массиве a.

      import numpy as np
a = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
k = 3
ma = np.mean(a[-k:])
print(ma)

Вот быстрая реализация с использованием numba (обратите внимание на типы). Обратите внимание, что там со смещением nans есть.

import numpy as np
import numba as nb

@nb.jit(nb.float64[:](nb.float64[:],nb.int64),
        fastmath=True,nopython=True)
def moving_average( array, window ):    
    ret = np.cumsum(array)
    ret[window:] = ret[window:] - ret[:-window]
    ma = ret[window - 1:] / window
    n = np.empty(window-1); n.fill(np.nan)
    return np.concatenate((n.ravel(), ma.ravel())) 

На самом деле мне хотелось немного другого поведения, чем принятый ответ. Я создавал средство извлечения скользящих средних дляsklearnконвейер, поэтому я потребовал, чтобы выходные данные скользящего среднего имели тот же размер, что и входные. Я хочу, чтобы скользящее среднее предполагало, что ряд остается постоянным, то есть скользящее среднее значение[1,2,3,4,5] с окном 2 даст [1.5,2.5,3.5,4.5,5.0].

Для векторов-столбцов (мой вариант использования) мы получаем

def moving_average_col(X, n):
  z2 = np.cumsum(np.pad(X, ((n,0),(0,0)), 'constant', constant_values=0), axis=0)
  z1 = np.cumsum(np.pad(X, ((0,n),(0,0)), 'constant', constant_values=X[-1]), axis=0)
  return (z1-z2)[(n-1):-1]/n

А для массивов

def moving_average_array(X, n):
  z2 = np.cumsum(np.pad(X, (n,0), 'constant', constant_values=0))
  z1 = np.cumsum(np.pad(X, (0,n), 'constant', constant_values=X[-1]))
  return (z1-z2)[(n-1):-1]/n

Конечно, не обязательно принимать постоянные значения для заполнения, но в большинстве случаев этого должно быть достаточно.

Другие вопросы по тегам