ARTICLE DETAIL

建站实战干货

来自一线的建站与推广经验沉淀,每一条都经过真实交付验证。

EWMA指数加权移动平均的MATLAB实现与工程应用

2026/9/13 16:27:57 拓冰建站 浏览量
EWMA指数加权移动平均的MATLAB实现与工程应用 简介一套基于MATLAB的指数加权移动平均EWMA模型估计资料主要面向金融风险管理、波动率建模与时间序列分析领域的初学者和开发者。与简单移动平均相比指数加权移动平均对近期数据赋予更高权重能够更快捕捉市场短期波动处理非平稳序列时更具灵活性该资源正是围绕这一核心思路展开。包内共6个文件包含3个M脚本、2个MAT数据文件和1个说明文档压缩包仅12KB便于下载与快速使用。M脚本覆盖了EWMA方差、波动率、协方差等核心估计函数数据文件提供演示用样本说明文档则帮助用户理解平滑参数λ的作用、起始值设定以及参数配置流程。目前已有1290人学习下载结合这些代码读者既能验证EWMA在波动率预测中的实际效果也能进一步结合正态分布等假设计算风险价值VaR为资产组合的风险量化提供实用手段。资源结构紧凑、示例清晰代码便于二次修改适合作为金融量化课程设计、项目实践或入门自学的参考资料。1. 指数加权移动平均的估计值实时系统里的递推均值指数加权移动平均EWMA的估计值解决的是实时系统里一个非常具体的问题当前均值到底是多少。设备监控、策略回测、产线质检都会遇到这个需求但全样本均值对近期变化反应太慢滑动窗口又要反复解释“窗口为什么取 60 而不是 80”。EWMA 把这个问题压缩成一条递推式s_t α·x_t (1-α)·s_{t-1}新样本权重是 α历史信息按(1-α)^k指数衰减因此每个时刻都有一个综合全部历史加权信息的估计值还不用回头重算整段序列。相比全样本均值它够快相比滑动窗口它少一个边界参数相比卡尔曼滤波它不需要维护协方差矩阵在 MATLAB 里用一条 filter 就能落地。下面从“为什么它算估计值而不是平滑值”讲起把 MATLAB 实现、α 怎么定、控制图和波动率怎么用一次说透。2. EWMA递推结构为什么它是“估计值”而不是历史平均2.1 递推展开与半衰期EWMA 的递推定义是s_t α·x_t (1-α)·s_{t-1}其中 α 必须在 0 到 1 之间。把递推式反复展开能得到输入对输出的显式权重s_t α·x_t α(1-α)·x_{t-1} α(1-α)^2·x_{t-2} ... (1-α)^t·s_0系数α(1-α)^i加起来是1-(1-α)^t再加上初始项(1-α)^t·s_0总权重正好为 1。这说明递推输出是一个加权平均而不是某种带偏的累积量。因为权重按几何级数衰减工程上通常直接算半衰期也就是某个历史样本的权重衰减到一半需要多少步alpha 0.1; kHalf -log(2) / log(1 - alpha); % 权重衰减一半所需步数 nEff (2 - alpha) / alpha; % 等效样本量 fprintf(alpha%.2f: 半衰期%.1f步, 等效样本量%.1f\n, ... alpha, kHalf, nEff);α 取 0.1 时半衰期约 6.6 步等效样本量约 19α 取 0.01 时半衰期约 69 步等效样本量约 199。这里的等效样本量N_eff (2-α)/α是由权重平方和推导出来的直观含义是“这个指数加权平均大约相当于多少个等权样本”。后文控制图的标准误计算、波动率估计的自由度修正都会用到它建议当成基础参数存起来。2.2 为什么递推输出可以被当作状态估计值如果目标只是把曲线画得好看任何窗口平均都能做到。EWMA 真正被当成估计值来用背景通常是一个简单观测模型x_t μ_t ε_t其中 μ_t 是随时间缓慢移动的真实状态ε_t 是零均值噪声。把s_{t-1}当作对 μ_t 的预测再定义一步预测误差e_t x_t - s_{t-1}EWMA 的更新就可以改写成新息修正形式s_t s_{t-1} α·e_t也就是说先按上一拍估计值预测再用预测误差的 α 比例修正。这个形式和稳态卡尔曼滤波完全同构当 μ_t 服从随机游走、观测噪声方差与状态噪声方差比值固定时线性最小均方误差滤波的稳态增益就是一个常数这个常数恰对应某个 α。所以s_t在设计目标上就是对 μ_t 的估计值而不是为了让曲线更光滑而做的后处理。代价也很直接。估计值的方差随 α 增大而增大但如果 μ_t 存在趋势估计会出现滞后滞后量与(1-α)/α成正比。α 大跟随快但方差大α 小平滑好但对真实漂移反应慢。实际项目里比较可解释的做法是先按业务能容忍的滞后步数定半衰期再反推 α而不是凭手感试数。2.3 alpha 与遗忘因子的命名与取值规则网上搜 EWMA 资料时很容易看到另一套写法v_t λ·v_{t-1} (1-λ)·x_t^2这里的 λ 是遗忘因子对应本文的1-α。金融里 RiskMetrics 常用 λ0.94意思就是 α0.06。两套命名混在文档里最容易把参数设反。下表按“α 越大响应越快”对齐业务场景平滑系数 α遗忘因子 λ等效样本量 N_eff工业过程突发漂移0.2 ~ 0.40.6 ~ 0.84 ~ 9经典 EWMA 控制图0.1 ~ 0.20.8 ~ 0.99 ~ 19金融日频波动率0.03 ~ 0.060.94 ~ 0.9732 ~ 65长记忆信号去噪0.02 ~ 0.050.95 ~ 0.9839 ~ 99经验做法是先把业务可容忍的滞后步数换算成半衰期再用半衰期公式反解 α。比如业务说“3 分钟内必须跟上阶跃变化”采样间隔 5 秒那就是 36 步内衰减一半解出来的 α 大约 0.02。这个反解代码建议单独留成函数后面做自动调参时要拿它定初始搜索区间。3. 在MATLAB里实现指数加权移动平均模型循环、filter与函数封装3.1 用循环先把递推关系写清楚最小可运行的 MATLAB 实现是循环逐行对应递推公式alpha 0.1; % 平滑系数0alpha1 x randn(1000, 1) 0.01*(1:1000); % 带趋势的测试数据列向量 s zeros(size(x)); s(1) x(1); % 初始估计值取第一个观测 for t 2:numel(x) s(t) alpha * x(t) (1 - alpha) * s(t-1); end初始值取第一个观测是最常用做法也可以取前若干点的均值这样起始段收敛更快。循环版本适合教学、算法核对和后续移植到 C 或 Simulink 代码生成在 10^6 量级数据上也不算慢但长序列批量处理时应该用 filter 重写。3.2 用 filter() 把估计值计算向量化MATLAB 的 filter 函数直接实现差分方程对应 EWMA 递推可以写成b alpha; % 分子的前向系数对应 x_t a [1, -(1-alpha)]; % 分母的自回归系数对应 s_{t-1} zi x(1) * (1 - alpha); % 初始状态让第一拍输出等于 x(1) s filter(b, a, x, zi);filter(b, a, x, zi)里b决定当前输入怎么进来a决定上一拍输出怎么反馈。关键是zi如果不设置第一拍输出会是α·x_1整条曲线起始段会明显偏低。要让s_1 x_1需要把初始状态设成x_1(1-α)这一项等价于“初始值 s_0 x_1 对首拍输出的贡献”。filter 版本在长序列上比循环快一个数量级以上也方便用 parfor 对多列信号并行处理代价是zi不够直观所以函数里最好保留循环版本作为注释参照。3.3 封装成可复用的 ewmaEst 并处理缺失值实际采集数据几乎都带 NaN。标准 EWMA 遇到 NaN 时不能更新但也不应该把整段序列截断常见做法是“缺失时保持上一拍估计值”function [s, nEff, kHalf] ewmaEst(x, alpha, initMethod) arguments x (:,1) double alpha (1,1) double {mustBeInRange(alpha, 0, 1)} initMethod string first end x x(:); N numel(x); s zeros(N, 1); if initMethod first s(1) x(1); elseif initMethod mean s(1) mean(x(~isnan(x(1:min(10,N)))), omitnan); end for t 2:N if isnan(x(t)) s(t) s(t-1); % 缺失值不更新 else s(t) alpha * x(t) (1 - alpha) * s(t-1); end end nEff (2 - alpha) / alpha; kHalf -log(2) / log(1 - alpha); end提示arguments 块从 MATLAB R2019b 开始支持。老版本直接改成function [s, nEff, kHalf] ewmaEst(x, alpha, initMethod)再把参数校验写在函数体开头即可。函数返回值顺带给出等效样本量和半衰期这样调用端写报表时不需要重复计算。缺失值策略这里用的是“保持”适合缺失段较短的情况如果缺失段很长估计值会一直停在缺失前的水平此时应该改成前向后向分段重置或者在缺失段结束后重新初始化。三种写法的取舍可以按表来写法适用场景主要注意点循环教学、算法核对、代码生成最直观起始段初始值显式可见filter长序列、批量并行初始状态 zi 必须按 3.2 节设置封装函数多脚本复用、报表输出参数校验与缺失值策略要提前约定从这一步开始后面的调参、控制图和波动率代码都只调用ewmaEst不再重复写递推逻辑。4. 指数加权移动平均的参数选择alpha经验表与MATLAB自动搜索4.1 平滑系数 alpha 的经验表先给可用的初始值前面那张场景表可以作为初始值但有一个换算逻辑更实用从旧系统的滚动窗口长度迁移到 EWMA。滚动窗口 N 点平均等效样本量解出来大致是α ≈ 2/(N1)。旧系统用 20 点窗口对应的 α 约 0.095旧系统用 60 点窗口α 约 0.033。迁移老代码时这个换算能让控制限和报警率保持大体一致而不是换了算法之后整个监控行为都变了。另一个判断方法是按业务滞后容忍度反推。控制图场景下如果要求 20 个采样周期内识别出阶跃半衰期取 10 步左右解出来的 α 在 0.07 到 0.1 之间。这样定出来的参数至少方向是对的之后再用数据做精细优化。4.2 用一步预测误差搜索最优 alpha 的MATLAB代码拟合 α 时最常见的错误是用整条序列的拟合误差。s_t本身包含x_t拿x_t - s_t当误差必然偏小而且会把参数推向更大的 α。正确的做法是用s_{t-1}作为x_t的一步预测function rmse rmseOfEWMA(x, alpha) s ewmaEst(x, alpha); % 复用 3.3 节封装 pred [NaN; s(1:end-1)]; % 上拍估计值作本拍预测 e x(2:end) - pred(2:end); % 一步预测误差 rmse sqrt(mean(e.^2, omitnan)); end x load(signal.mat).x; % 替换成自己的数据列向量 alphas linspace(1e-3, 0.5, 40); % 粗网格搜索 err arrayfun((a) rmseOfEWMA(x, a), alphas); [~, idx] min(err); a_lo alphas(max(1, idx-1)); a_hi alphas(min(numel(alphas), idx1)); aBest fminbnd((a) rmseOfEWMA(x, a), a_lo, a_hi);pred整体错开一拍所以x(2:end) - pred(2:end)比较的正是x_t与s_{t-1}。arrayfun把 40 个候选 α 各跑一遍ewmaEst在 10^5 量级数据上开销很小粗网格选邻域再交给fminbnd精修是为了避免优化目标在小范围内出现多峰。跑完后用[rmseOfEWMA(x,aBest), min(err)]对一下如果两者差距明显说明网格太粗或数据里存在长趋势需要把 α 搜索上限放宽到 0.8 再做一轮。4.3 波动率估计值场景下的负对数似然调参如果估计的是方差而不是均值比如对收益率残差r_t做波动率估计RMSE 目标就不合适因为方差是二阶量必须用似然。常见做法是假设残差零均值正态分布然后最小化负对数似然function v ewmaVar(r, alpha) v zeros(size(r)); v(1) var(r, omitnan); % 初始方差用全样本 for t 2:numel(r) v(t) alpha * r(t)^2 (1 - alpha) * v(t-1); end end function nll ewmaVarNLL(r, alpha) v ewmaVar(r, alpha); nll 0.5 * mean(r(2:end).^2 ./ v(2:end) log(v(2:end)), omitnan); end zBest fminbnd((z) ewmaVarNLL(r, 1/(1exp(-z))), -4, 4); aBest 1 / (1 exp(-zBest));似然里每个时刻的贡献是r_t^2 / v_t ln(v_t)这一项同时惩罚低估和高估。这里用 logit 变换α 1/(1exp(-z))而不是直接搜 α方差目标函数在 α 接近 0.95 的长记忆区域非常平直接搜索容易贴到边界上换成 z 之后搜索空间变成整条实数轴收敛更顺滑。ewmaVar的初始方差用的整体var(r)样本量小时建议改成前 20 点方差否则首段估计值会被全局方差带偏。5. EWMA模型的工程落地控制图、波动率估计与残差自相关检查5.1 过程监控里的EWMA控制图怎么写工业质量监控里EWMA 控制图对小幅漂移比 Shewhart 图更敏感这是它最常见的工程落点。中心线取均值控制限用稳态标准误alpha 0.2; L 2.7; % 控制限倍数常用 2.7 或 3 s ewmaEst(x, alpha); cl mean(x, omitnan); sigma 1.4826 * median(abs(x - cl), omitnan); % MAD 估计尺度 se sigma * sqrt(alpha / (2 - alpha)); % EWMA 稳态标准误 ucl cl L * se; lcl cl - L * se; figure; plot(t, s, b-); hold on; yline(cl, k-); yline(ucl, r--); yline(lcl, r--);1.4826 * MAD是对正态数据的稳健尺度估计抗异常点能力比直接 std 好控制图场景建议优先用。sqrt(alpha/(2-alpha))来自 EWMA 稳态方差公式和滑动窗口的sigma/sqrt(N)不是一回事两者只在特定 α 下数值巧合相等迁移代码时不要混用。参数对应关系如下控制图参数常用值作用α0.05 ~ 0.2越小对漂移越迟钝越大噪声越多L2.7 ~ 3控制误报警率2.7 对应约 370 步平均运行长度尺度估计1.4826·MAD比 std 更抗尖峰稳态标准误σ·sqrt(α/(2-α))与样本量 N 无关5.2 金融波动率估计先估计方差再谈均值价格序列通常先转成对数收益率r_t ln(P_t/P_{t-1})再对r_t^2做 EWMA。因为收益率均值接近 0一般不再对r_t做均值估计直接估计方差。长记忆场景的标准设置是遗忘因子 λ0.94也就是 α0.06alpha 0.06; % 对应 RiskMetrics 0.94 遗忘因子 v ewmaVar(r, alpha); vol sqrt(v); % 波动率估计值序列ewmaVar返回的是方差估计值开方后才是波动率。预测下一期波动率直接用vol(end)即可。这个模型最大的坑是 α 要按数据频率分档日频用 0.06月频用 0.03 左右周频介于两者之间。把日频参数直接搬去小时数据估计值会抖得非常厉害需要先用retime聚合再重新定参而不是换数据不换 α。5.3 残差自相关检验与模型升级参数确定后还要检查模型有没有吃干净可预报结构。残差就是一步预测误差resid x(2:end) - s(1:end-1); % 一步预测残差 [acf, lags, bounds] autocorr(resid, NumLags, 20);bounds是白噪声假设下的 95% 置信带。前几个滞后若超出边界说明 EWMA 只覆盖了均值漂移没有覆盖趋势或周期成分。此时应该升级到 Holt 双指数平滑在 EWMA 基础上加一个趋势项s_t α·x_t (1-α)·(s_{t-1} b_{t-1}) b_t β·(s_t - s_{t-1}) (1-β)·b_{t-1}MATLAB 的tsmovavg只做单重平滑趋势场景建议用fit工具箱里的指数平滑模型或直接手写上面两行更新。调参目标也要改成两步预测误差否则趋势项会被一步预测目标带偏。这一步是实际项目里最容易被跳过的环节但残差检验往往比调参更快暴露模型缺陷。6. EWMA模型上值得保留的三个MATLAB技巧6.1 用仿真数据先验一遍估计值调参结果很难直接评估好坏因为真实状态 μ_t 未知。仿真时可以自己造真值再对不同 α 做对比rng(7); N 2000; mu cumsum(randn(N, 1) * 0.02); % 设定的状态随机游走 x mu randn(N, 1) * 0.5; % 观测值叠加噪声 aList [0.02, 0.05, 0.1, 0.2]; for k 1:numel(aList) sk ewmaEst(x, aList(k)); mse(k) mean((mu - sk).^2); % 真值已知才能算 MSE end [~, bestIdx] min(mse);最优 α 对应的曲线和mu叠在一起画能直观看到偏差与方差的取舍。这个验证方式也适合回答“为什么 α 不能取 0.5”仿真里 α0.5 的 MSE 一般明显偏高而 α0.05 到 0.1 之间最接近最优区间比空口解释更有说服力。6.2 注意滤波器方向在线估计不允许未来信息filter天然是因果的按时间正序处理就没有问题。离线分析时有人会顺手用filtfilt做零相位滤波这对平滑合适但用在估计值语义上是有问题的filtfilt双向处理输出里混杂了未来样本信息故障报警会“提前”出现这在回测里尤其危险。要保持 EWMA 的估计语义离线数据也应该按正序filter一次最多把序列分段分别初始化不能用filtfilt替换。如果离线分析确实需要双向平滑就明确把它叫平滑器不要和实时估计值放在同一张图里对比。6.3 异常值污染时改判不更新监控数据常有跳变普通 EWMA 会把跳变当成观测慢慢拉走估计值。工程上常用的做法是拿前一段残差的 MAD 做尺度门限超出门限只标记不更新K max(5, round(20 / alpha)); for t 2:N e x(t) - s(t-1); mad_t 1.4826 * median(abs(x(max(1,t-K):t-1) - s(max(1,t-K):t-2))); if abs(e) 3 * mad_t s(t) alpha * x(t) (1 - alpha) * s(t-1); else s(t) s(t-1); % 异常点不参与更新 end end门限用 3 倍 MAD 而不是 3 倍 std因为 MAD 对尖峰不敏感跳变不会把门限自身拉大。加上这条之后EWMA 控制图的假报警率通常会明显下降而真正的小幅漂移仍然能在一两个半衰期之内追上去。本文还有配套的精品资源点击获取