本文介绍一种比 NumPy 原生 sliding_window_view + nanmedian 快 13 倍、比 SciPy 滑动中位数快 10 倍的高性能实现——通过手写 Cython 扩展，利用 C++ std::multiset 维护有序窗口，在 O(n log k) 时间内完成单次滑动更新，显著加速 10,000 条长度为 750 的信号处理任务。
本文介绍一种比 numpy 原生 `sliding_window_view` + `nanmedian` 快 13 倍、比 scipy 滑动中位数快 10 倍的高性能实现——通过手写 cython 扩展，利用 c++ `std::multiset` 维护有序窗口，在 o(n log k) 时间内完成单次滑动更新，显著加速 10,000 条长度为 750 的信号处理任务。
在时间序列分析、生物信号去噪或工业传感器数据预处理中，对数千条中等长度（如 750 点）信号批量计算滑动中位数是常见需求。但标准方法（如 np.lib.stride_tricks.sliding_window_view 配合 np.nanmedian）会产生大量中间数组，内存占用高、缓存不友好，导致性能瓶颈。以原文场景为例：处理 2000 条信号，每条 750 点、窗口大小 150，原方法耗时达 6.11 秒。
✅ 更优解：Cython + std::multiset 实现增量式滑动中位数
核心思想是避免重复排序——用 C++ 标准库的 multiset（底层为红黑树）动态维护当前窗口内所有元素，并始终持有一个指向中位数位置的迭代器。每次滑窗仅需：
- 删除最左侧旧值（O(log k)），
- 插入最右侧新值（O(log k)），
- 调整中位数迭代器（均摊 O(1)），
- 直接读取中位数值（O(1)）。
该策略将单次窗口计算从 O(k log k) 降至 O(log k)，整体复杂度由 O(n·k log k) 优化至 O(n log k)，实测 median.sliding_median 函数仅需 0.20 秒（纯计算部分），较 SciPy 加速 10×，较原始 NumPy 方法加速 30×。
以下是关键 Cython 实现（median.pyx）精简说明：
# cython: language_level=3, boundscheck=False, wraparound=False
import numpy as np
from numpy cimport float64_t
from libcpp.set cimport multiset
from cython.operator cimport dereference as deref, preincrement as inc, predecrement as dec
def sliding_median(float64_t[:] arr, int64_t n):
assert n >= 3 and arr.size >= n
cdef int64_t i, arr_size = arr.size
cdef float64_t[:] out = np.empty(arr_size - n + 1, dtype=np.float64)
# 初始化 multiset 并填充首窗口
cdef multiset[float64_t] sorted_items
for i in range(n):
sorted_items.insert(arr[i])
# 定位中位数迭代器（支持偶数长度双中位平均）
cdef int64_t mid1 = (n - 1) // 2
cdef multiset[float64_t].iterator it = sorted_items.begin()
for i in range(mid1): inc(it)
out[0] = _median_value(it, n)
# 滑动更新：删除左端、插入右端、智能调整迭代器
for i in range(n, arr_size):
cdef float64_t old_val = arr[i - n], new_val = arr[i]
cdef multiset[float64_t].iterator old_it = sorted_items.find(old_val)
# 关键：根据新旧值与当前中位数关系，预判迭代器偏移
cdef bint delayed_inc = False, delayed_dec = False
if old_it == it:
inc(it) # 避免擦除后迭代器失效
sorted_items.erase(old_it)
sorted_items.insert(new_val)
if delayed_inc: inc(it)
if delayed_dec: dec(it)
out[i - n + 1] = _median_value(it, n)
return np.asarray(out)
? 注意：
_median_value()是一个辅助函数，处理奇偶长度逻辑（偶数时取中间两数平均）。完整代码需配合setup.py编译（见原文），并确保系统已安装 Cython 和 C++ 编译器（如 GCC/Clang）。
? 集成到 DataFrame 流水线
为无缝替代原 moving_median，封装加速版函数：
import median # 编译后的 Cython 模块
def moving_median_fast(signal, n=150):
arr = np.asarray(signal, dtype=np.float64)
b = median.sliding_median(arr, n) # 0.2s 核心计算
# 高效补齐前 n-1 个值（可进一步 Cython 化以消除 concat 开销）
b_full = np.concatenate([np.full(n-1, np.nanmedian(arr)), b])
return signal - b_full
# 应用到 DataFrame（推荐先转为 2D 数组以规避 list 转换开销）
df['detrended'] = df['values'].apply(moving_median_fast)
⚠️ 关键优化建议与注意事项
-
数据结构先行：避免
DataFrame存储 Pythonlist。若所有信号等长（如本例均为 750），请直接使用 2D NumPy 数组（shape(2000, 750)）：signals_2d = np.array(column_data) # shape (2000, 750) # 向量化调用（需修改 Cython 支持 axis=1）
-
消除 Python 层开销：
np.concatenate([...])和np.nanmedian()在循环中调用会成为瓶颈。生产环境应：- 预分配输出数组；
- 在 Cython 中直接计算首段均值并写入；
- 使用
prange+nogil实现多线程并行（需显式释放 GIL）。
浮点精度警告：当前实现假设窗口内无重复浮点值精确等于中位数（否则
find()可能定位错误节点）。对高精度场景，建议改用lower_bound+ 计数器增强鲁棒性。编译与部署：首次编译后生成
.so（Linux/macOS）或.pyd（Windows）文件，可像普通模块一样导入。CI/CD 中建议缓存编译产物。
? 性能对比（i5-9600KF，Numpy 2.1.3）
| 方法 | 耗时 | 加速比（vs 原始） |
|---|---|---|
sliding_window_view + nanmedian
|
6.11 s | 1× |
SciPy ndimage.median_filter（串行） |
2.02 s | 3× |
ndfilters（并行，冷启动 8.6 s） |
0.98 s | 6× |
median.sliding_median（纯 Cython） |
0.20 s | 30× |
moving_median_fast（含 Python 封装） |
0.45 s | 13.6× |
✅ 结论：当处理万级中等长度信号时，手写 Cython 是目前最实用的性能解。它平衡了开发成本与执行效率，且易于集成到现有 NumPy/Pandas 生态中。下一步可结合
joblib.Parallel或numba.prange进一步榨干多核性能。
