用于日期时间的最快滚动窗口方法
我现在在循环中对时间数据的一个30分钟滚动窗口应用一个有点复杂的函数,结果发现花费的时间太长,已经不太合理。
下面给出一个最小可复现的示例,目前耗时约94 ms,结合数据集的实际规模和循环情况,这个时间显然不可行。
import pandas as pd
import numpy as np
from datetime import datetime,timedelta
import statsmodels.tsa.api as smt
def instrument_variance(ts): #Lernschow et al 2000
ts_acov = smt.acovf(ts)[:3] #Find first few lags of autocovariance
M_0 = ts_acov[0] #Get autocovariance at lag 0
M_fit = np.polyfit(range(3),ts_acov,1) #Get linear fit of the first few lags
M_lim0 = M_fit[0] #Get limit as fit goes to 0
var_inst = M_0-M_lim0 #Get error variance from M11(0)-lim>0(M11)
return var_inst
start = datetime.now()
base = datetime.today()
date_list = [base + timedelta(minutes=x) for x in range(100)]
df = pd.DataFrame(np.random.randint(0,100,size=(100, 3)), index=date_list,columns=list('ABC'))
df_cop = df.copy()
var_inst = df_cop.sort_index().rolling('30min',center=True,min_periods=3).apply(instrument_variance)
df = df.sort_index().rolling('30min',center=True,min_periods=3).var()
df = df.sub(var_inst) #Subtract uncorrelated var from total var to get correlated var
df = df.resample('5min').first()
timed = datetime.now()-start
print(timed)
我的Python知识有限,操作也比较复杂,但有没有什么方法可以让它更快一些?
我也看过类似的问题:
在Python(Numpy或 Pandas)中计算滚动(线性)回归的最有效/最快方法
但我很难找到一个针对基于日期时间的窗口的解决方案。我不能假设采样频率是恒定的,因此我更愿意使用datetime而不是整数索引的方式来实现。
解决方案
not sure if it will be lean enough for your dataset... run_fast_prefix gets about 45x speed-up:
tier 1 optimizations of per-window overhead (37ms → 2.9ms, ~13x). in run_fast:
- 将inline
acovf+polyfit结合使用。每次调用在一个6 元素窗口上都有校验/奇异值分解的开销。用少量的numpy运算和闭式斜率 (M2 - M0)/2取代。 - 将两个滚动遍历合并成一个。原始实现对滚动
.apply和滚动.var的处理都是对窗口的完整遍历。 raw=True跳过为每个窗口构建一个pandas Series。
tier 2优化,灵感来自你分享的另一个SO问题的答案。完全消除每个窗口的调用,在 run_fast_prefix:
- 前缀和:对整列进行一次
cumsum,然后每个窗口的sum只需两个数组查找。所有 ~100个窗口的M0、M1、M2变成了四个向量表达式。逐窗口的Python解释器循环消失,所有工作都在已编译的numpy中完成。
总速度提升:约45x。tier 1是在现有结构上的纯常数因子优化;tier 2是链接的SO问题中所指向的结构性变动。
import numpy as np
import pandas as pd
import statsmodels.tsa.api as smt
# original for reference
def instrument_variance(ts): # Lenschow et al 2000
ts_acov = smt.acovf(ts)[:3]
M_0 = ts_acov[0]
M_fit = np.polyfit(range(3), ts_acov, 1)
M_lim0 = M_fit[0]
var_inst = M_0 - M_lim0
return var_inst
def run_original(df):
var_inst = (df.sort_index()
.rolling('30min', center=True, min_periods=3)
.apply(instrument_variance))
var_tot = (df.sort_index()
.rolling('30min', center=True, min_periods=3)
.var())
return var_tot.sub(var_inst).resample('5min').first()
# combined `apply` path
def correlated_var(a):
"""Rolling correlated variance in one pass.
Biased autocovariance (÷ n) for the instrument term, ddof=1 total var.
Slope of the line through (0,M0),(1,M1),(2,M2) is (M2 - M0)/2 in closed
form, reproducing `np.polyfit(...)[0]` semantics from the original. If
the intent is the lag-0 intercept, use: inst = M0 - ((M0+M1+M2)/3 - slope).
"""
n = a.size
if n < 3:
return np.nan
c = a - a.mean()
ssq = (c * c).sum()
M0 = ssq / n
M1 = (c[:-1] * c[1:]).sum() / n
M2 = (c[:-2] * c[2:]).sum() / n
slope = 0.5 * (M2 - M0)
inst = M0 - slope
var_u = ssq / (n - 1)
return var_u - inst
def run_fast(df):
return (df.sort_index()
.rolling('30min', center=True, min_periods=3)
.apply(correlated_var, raw=True)
.resample('5min').first())
# prefix-sum vectorized path
def _rolling_bounds(df, offset, center=True, closed=None):
"""Exact [start, end) bounds that df.rolling(offset, center=..., closed=...) uses.
Extracted by running pandas' own rolling on a row-index series with a
scalar-returning apply. This is the robust way — pandas' offset-based
centering is implemented as an output-alignment shift, not a bounds shift,
so VariableOffsetWindowIndexer.get_window_bounds(center=True) and any hand-
derived (t - w/2, t + w/2] formula both drift by a row in places.
The bootstrap apply itself is O(N) but cheap per window (scalar return),
so downstream prefix-sum math still gets most of the speedup.
"""
N = len(df)
if N == 0:
return np.zeros(0, dtype=np.int64), np.zeros(0, dtype=np.int64)
idx_ser = pd.Series(np.arange(N, dtype=np.float64), index=df.index)
r = idx_ser.rolling(offset, center=center, closed=closed, min_periods=1)
first = r.apply(lambda w: w[0], raw=True).to_numpy()
last = r.apply(lambda w: w[-1], raw=True).to_numpy()
empty = np.isnan(first)
start = np.where(empty, 0, first).astype(np.int64)
end = np.where(empty, 0, last + 1).astype(np.int64)
return start, end
def _correlated_var_prefix(x, start, end, min_periods=3):
"""Vectorized correlated-variance per window via prefix sums.
Window for output i is x[start[i]:end[i]]. Requires min_periods >= 3
(lag-2 autocov needs at least 3 samples).
"""
assert min_periods >= 3
x = np.asarray(x, dtype=float)
N = x.size
out = np.full(start.shape, np.nan, dtype=float)
if N == 0:
return out
# Padded prefix sums: S[0]=0, S[k]=sum(x[:k]) so sum(x[a:b]) = S[b]-S[a].
S = np.empty(N + 1, dtype=float); S[0] = 0.0
np.cumsum(x, out=S[1:])
Sq = np.empty(N + 1, dtype=float); Sq[0] = 0.0
np.cumsum(x * x, out=Sq[1:])
SP1 = np.zeros(max(N, 1), dtype=float) # lag-1 pair-product cumsum
if N >= 2:
np.cumsum(x[:-1] * x[1:], out=SP1[1:])
SP2 = np.zeros(max(N - 1, 1), dtype=float) # lag-2 pair-product cumsum
if N >= 3:
np.cumsum(x[:-2] * x[2:], out=SP2[1:])
n = (end - start).astype(np.int64)
valid = n >= min_periods
if not valid.any():
return out
s = start[valid]
e = end[valid]
nn = n[valid].astype(float)
Sx = S[e] - S[s]
Sxx = Sq[e] - Sq[s]
mean = Sx / nn
ssq = Sxx - nn * mean * mean
M0 = ssq / nn
var_u = ssq / (nn - 1.0)
# Lag 1: sum over pairs (x[t], x[t+1]) for t in [s, e-1).
sum_P1 = SP1[e - 1] - SP1[s]
sumL1 = Sx - x[e - 1] # sum of left members
sumR1 = Sx - x[s] # sum of right members
n1 = nn - 1.0
M1 = (sum_P1 - mean * (sumL1 + sumR1) + n1 * mean * mean) / nn
# Lag 2: sum over pairs (x[t], x[t+2]) for t in [s, e-2).
sum_P2 = SP2[e - 2] - SP2[s]
sumL2 = Sx - x[e - 2] - x[e - 1]
sumR2 = Sx - x[s] - x[s + 1]
n2 = nn - 2.0
M2 = (sum_P2 - mean * (sumL2 + sumR2) + n2 * mean * mean) / nn
slope = 0.5 * (M2 - M0)
inst = M0 - slope
out[valid] = var_u - inst
return out
def run_fast_prefix(df, offset='30min', min_periods=3, resample_rule='5min'):
df = df.sort_index()
start, end = _rolling_bounds(df, offset, center=True)
arr = df.to_numpy(dtype=float, copy=False)
out = np.empty_like(arr, dtype=float)
for c in range(arr.shape[1]):
out[:, c] = _correlated_var_prefix(arr[:, c], start, end, min_periods)
result = pd.DataFrame(out, index=df.index, columns=df.columns)
return result.resample(resample_rule).first() if resample_rule else result
你可能已经知道:numpy与 pandas在处理NaN时有差异。如果数据中有NaN,请在运行前进行插值或删除。
站内所有文章版权归属LeftHeroAI导航站,无授权禁止任何主体转载、抄袭、复制内容,亦不得私自架设镜像站点。一经侵权,本站将通过法律途径追责。