Панды: увеличить скорость прокручивания окна (применить пользовательскую функцию) - PullRequest
Новая
       11

Панды: увеличить скорость прокручивания окна (применить пользовательскую функцию)

4 голосов
/ 22 апреля 2019

Я использую этот код для применения функции (funcX) к моему фрейму данных с помощью скользящего окна. Основная проблема заключается в том, что размер этого фрейма данных (data) очень велик, и я ищу более быстрый способ выполнения этой задачи. 

import numpy as np

def funcX(x):
    x = np.sort(x)
    xd = np.delete(x, 25)
    med = np.median(xd)
    return (np.abs(x - med)).mean() + med

med_out = data.var1.rolling(window = 51, center = True).apply(funcX, raw = True)

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

1 Ответ

5 голосов
/ 23 апреля 2019

Чтобы быть эффективным, оконный алгоритм должен связать результаты двух наложенных окон.

Здесь, с: med0 медианой, med медианой x \ med0, xl элементов передmed и xg элементов после med в отсортированных элементах, funcX(x) можно рассматривать как: 

<|x-med|> + med = [sum(xg) - sum(xl) - |med0-med|] / windowsize + med

Так что идея состоит в том, чтобы поддерживать буфер, который представляет отсортированное текущее окно, sum(xg) и sum(xl).Используя Numba как раз во время компиляции, здесь возникает очень хорошая производительность.

Первое управление буфером:

init сортирует первое окно и вычисляет влево (xls) и вправо (* 1022)*) суммы.

import numpy as np
import numba
windowsize = 51 #odd, >1
halfsize = windowsize//2

@numba.njit
def init(firstwindow):
    buffer = np.sort(firstwindow)
    xls = buffer[:halfsize].sum()
    xgs = buffer[-halfsize:].sum()   
    return buffer,xls,xgs

shift - линейная часть.Он обновляет буфер, сохраняя его отсортированным.np.searchsorted вычисляет позиции вставки и удаления в O(log(windowsize)).Это технически, поскольку xin<xout и xout<xin не являются симметричными ситуациями. 

@numba.njit
def shift(buffer,xin,xout):
    i_in = np.searchsorted(buffer,xin) 
    i_out = np.searchsorted(buffer,xout)
    if xin <= xout :
        buffer[i_in+1:i_out+1] = buffer[i_in:i_out] 
        buffer[i_in] = xin                        
    else:
        buffer[i_out:i_in-1] = buffer[i_out+1:i_in]                      
        buffer[i_in-1] = xin
    return i_in, i_out

update обновляет буфер и суммы левой и правой частей.Это технически, поскольку xin<xout и xout<xin не являются симметричными ситуациями. 

@numba.njit
def update(buffer,xls,xgs,xin,xout):
    xl,x0,xg = buffer[halfsize-1:halfsize+2]
    i_in,i_out = shift(buffer,xin,xout)

    if i_out < halfsize:
        xls -= xout
        if i_in <= halfsize:
            xls += xin
        else:    
            xls += x0
    elif i_in < halfsize:
        xls += xin - xl

    if i_out > halfsize:
        xgs -= xout
        if i_in > halfsize:
            xgs += xin
        else:    
            xgs += x0
    elif i_in > halfsize+1:
        xgs += xin - xg

    return buffer, xls, xgs

func эквивалентно исходному funcX в буфере.O(1). 

@numba.njit       
def func(buffer,xls,xgs):
    med0 = buffer[halfsize]
    med  = (buffer[halfsize-1] + buffer[halfsize+1])/2
    if med0 > med:
        return (xgs-xls+med0-med) / windowsize + med
    else:               
        return (xgs-xls+med-med0) / windowsize + med

med - глобальная функция.O(data.size * windowsize). 

@numba.njit
def med(data):
    res = np.full_like(data, np.nan)
    state = init(data[:windowsize])
    res[halfsize] = func(*state)
    for i in range(windowsize, data.size):
        xin,xout = data[i], data[i - windowsize]
        state = update(*state, xin, xout)
        res[i-halfsize] = func(*state)
    return res

Производительность: 

import pandas
data=pandas.DataFrame(np.random.rand(10**5))

%time res1=data[0].rolling(window = windowsize, center = True).apply(funcX, raw = True)
Wall time: 10.8 s

res2=med(data[0].values)

np.allclose((res1-res2)[halfsize:-halfsize],0)
Out[112]: True

%timeit res2=med(data[0].values)
40.4 ms ± 462 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

Это ~ 250X быстрее, с размером окна = 51. Час становится 15 секунд.

...