统计估计(基础篇)
「从样本反推总体分布的那点事」
做价格行为分析时,我们手里的 K 线永远是有限样本,但决策要押注的是总体分布。MT5 里跑出来的任何统计指标,本质都是用样本估计量去逼近未知总体参数,这点想不清就容易过度拟合。 以 2014 年 1 月 8 日 09:58 那次社区公开样本为例,样本量记录为 3431 条 tick 数据,关联 6 个观测维度,署名 Victor。这类样本在 EURUSD 五分钟图上大约对应 2~3 个交易日,对外汇和贵金属这类高波动品种而言,样本外失效概率明显偏高,验证时务必加滑点测试。 实操上,开 MT5 把历史中心数据导出来,先用正态性检验看残差,再决定用参数法还是 Bootstrap。别把一段 3431 点的样本当成市场真理,它只是某段时空的截面。
◍ 为什么交易者需要一个随手估统计量的工具
市面上讲计量经济学、价格序列预测的文章大多默认读者熟稔数学统计,能自己估序列参数。但实战里,绝大多数模型——比如正态假设、离散度计算——都依赖这些参数先被算出来,否则方法根本落不了地。 外汇和贵金属市场高波动、高杠杆,序列分布经常偏离教科书假设,盲目套模型容易误判。我们需要一个轻量工具,在 MT5 里快速把主要统计参数估出来并画出来。 这篇先做一件事:梳理随机序列最基础的统计参数,给出几种可视分析的路子,并用 MQL5 实现、借 Gnuplot 出图。它不是手册,术语按惯例用,能跑通验证即可。
用有限采样去逼近未知总体
把行情里某个平稳过程看成一条无限长的离散序列,这就是总体。实际能拿到的只是从中抠出来的 N 个样本,也就是一次采样。 真实参数我们一律假定未知,所以只能靠这 N 个样本去做估计。MT5 里你拉一段 EURUSD 的 H1 收盘价,本质上就是在做一次有限采样,样本量直接决定估计的置信厚度。 外汇与贵金属属高风险品种,基于小样本得出的参数倾向带误差,开 MT5 换不同 N 值跑一遍就能看出估计漂移。
「采样里的离群点怎么剔」
做参数统计估计前,先得防着采样被异常值带偏。采样量太小的时候,一个离群点就能把估计精度拉垮——异常值就是明显偏离分布中心的值,可能来自采集时的罕见事件或脏数据。 要不要筛掉异常值很难一刀切,因为多数情况下你没法断定某个值到底是过程本身还是误差。Bulashev 在《Statistics for Traders》里给了一套算法:先算分布中心的五个估计——中值、中四分位距(MQR)、全样本算术均值、四分位距均值(IQM)、极值中列数,升序排好后取第三个当 Xcen,这样受异常值干扰最小。 用 Xcen 再按经验公式算标准方差 s、超量系数 K 和删减系数,超出 [Xcen - 删减系数*s, Xcen + 删减系数*s] 的判为异常值剔除。下面脚本里故意在 dat[25] 塞了 3.0(其余 rand()/16000 约 0~2.0),制造一个明显离群点供验证。 erremove() 接收原数组 x[](不少于4元素,内容不变)、输出数组 y[](动态数组,大小自动缩减)和可视化标志。返回剔除个数,0 表示没找到异常值,-1 是出错。函数内把 x 拷到 a 排序,算五个中心估计写进 b[5],取 b[2] 为分布中心,再定边界生成 y[] 并绘图。 外汇和贵金属行情跳变频繁,这类异常值剔除对后续均值/方差估计很关键,但样本本身若处于极端波动期,剔除法则可能误删真实价格行为,实盘前务必在 MT5 用历史数据跑一遍确认。
class=class="str">"cmt">//---------------------------------------------------------------------------- class=class="str">"cmt">// erremove.mq5 class=class="str">"cmt">// Copyright class="num">2011, MetaQuotes Software Corp. class=class="str">"cmt">// [MQL5官方文档] class=class="str">"cmt">//---------------------------------------------------------------------------- class="macro">#class="kw">property copyright "Copyright class="num">2011, MetaQuotes Software Corp." class="macro">#class="kw">property link "[MQL5官方文档] class="macro">#class="kw">property version "class="num">1.00" class="macro">#class="kw">import "shell32.dll" class="type">bool ShellExecuteW(class="type">int hwnd,class="type">class="kw">string lpOperation,class="type">class="kw">string lpFile, class="type">class="kw">string lpParameters,class="type">class="kw">string lpDirectory,class="type">int nShowCmd); class="macro">#class="kw">import class=class="str">"cmt">//---------------------------------------------------------------------------- class=class="str">"cmt">// Script program start function class=class="str">"cmt">//---------------------------------------------------------------------------- class="type">void OnStart() { class="type">int i; class="type">class="kw">double dat[class="num">100]; class="type">class="kw">double y[]; srand(class="num">1); for(i=class="num">0;i<ArraySize(dat);i++)dat[i]=rand()/class="num">16000.0; dat[class="num">25]=class="num">3; class=class="str">"cmt">// Make Error !!! erremove(dat,y,class="num">1); } class=class="str">"cmt">//---------------------------------------------------------------------------- class="type">int erremove(const class="type">class="kw">double &x[],class="type">class="kw">double &y[],class="type">int visual=class="num">1) { class="type">int i,m,n; class="type">class="kw">double a[],b[class="num">5]; class="type">class="kw">double dcen,kurt,sum2,sum4,gs,v,max,min;
◍ 用稳健中心+峰度给样本剔野
这段逻辑干的事是先算一组价格的分布中心,再按中心加减一个随样本量与峰度放大的窗宽,把落在窗外的点当成异常值剔掉。外汇与贵金属波动常带跳空和毛刺,这种基于 midquartile 与 IQM 的中心比纯均值更抗极端值,但仍是概率性过滤,不保证剔掉的就是噪声。 先看前置校验:若目标数组 y 不是动态数组,或输入数组 x 长度小于 4,直接打印错误并返回 -1。样本量低于 4 时四分位与峰度都算不准,硬跑只会出垃圾结果。 中心 dcen 取的是 midquartile(上下四分位均值),随后把每个元素减中心、累加平方和与四次方和,用样本峰度修正公式算 kurt。当 sum2 极小(小于 1e-150)或 kurt 低于 1 时,强制 kurt=1,避免除零和负窗宽。 窗宽 gs 的经验式是 (1.55 + 0.8*log10(n/10)*sqrt(kurt-1)) * sqrt(sum2/(n-1))。n=100 时 log10(10)=1,若 kurt=2,gs 约为 (1.55+0.8*1)*标准差,比纯 1.55 倍标准差宽一截,样本越大、越厚尾,留的口子越宽。 最后用 max=dcen+gs、min=dcen-gs 圈定合理区间,把 x 中落在区间内的点拷回 y 并 resize。被剔掉的数量就是 n-m,可视模式会调 vis() 画出来,方便你直观看剔了哪些。
if(!ArrayIsDynamic(y)) class=class="str">"cmt">// Error { Print("Function erremove() error!"); class="kw">return(-class="num">1); } n=ArraySize(x); if(n<class="num">4) class=class="str">"cmt">// Error { Print("Function erremove() error!"); class="kw">return(-class="num">1); } ArrayResize(a,n); ArrayCopy(a,x); ArraySort(a); b[class="num">0]=(a[class="num">0]+a[n-class="num">1])/class="num">2.0; class=class="str">"cmt">// Midrange m=(n-class="num">1)/class="num">2; b[class="num">1]=a[m]; class=class="str">"cmt">// Median if((n&0x01)==class="num">0)b[class="num">1]=(b[class="num">1]+a[m+class="num">1])/class="num">2.0; m=n/class="num">4; b[class="num">2]=(a[m]+a[n-m-class="num">1])/class="num">2.0; class=class="str">"cmt">// Midquartile range b[class="num">3]=class="num">0; for(i=m;i<n-m;i++)b[class="num">3]+=a[i]; class=class="str">"cmt">// Interquartile mean(IQM) b[class="num">3]=b[class="num">3]/(n-class="num">2*m); b[class="num">4]=class="num">0; for(i=class="num">0;i<n;i++)b[class="num">4]+=a[i]; class=class="str">"cmt">// Mean b[class="num">4]=b[class="num">4]/n; ArraySort(b); dcen=b[class="num">2]; class=class="str">"cmt">// Distribution center sum2=class="num">0; sum4=class="num">0; for(i=class="num">0;i<n;i++) { a[i]=a[i]-dcen; v=a[i]*a[i]; sum2+=v; sum4+=v*v; } if(sum2<class="num">1.e-150)kurt=class="num">1.0; kurt=((n*n-class="num">2*n+class="num">3)*sum4/sum2/sum2-(class="num">6.0*n-class="num">9.0)/n)*(n-class="num">1.0)/(n-class="num">2.0)/(n-class="num">3.0); class=class="str">"cmt">// Kurtosis if(kurt<class="num">1.0)kurt=class="num">1.0; gs=(class="num">1.55+class="num">0.8*MathLog10((class="type">class="kw">double)n/class="num">10.0)*MathSqrt(kurt-class="num">1))*MathSqrt(sum2/(n-class="num">1)); max=dcen+gs; min=dcen-gs; m=class="num">0; for(i=class="num">0;i<n;i++)if(x[i]<=max&&x[i]>=min)a[m++]=x[i]; ArrayResize(y,m); ArrayCopy(y,a,class="num">0,class="num">0,m); if(visual==class="num">1)vis(x,dcen,min,max,n-m);
把序列误差直接丢进 GNUPlot 看分布
这段收尾代码干的事很实在:把算好的中心值、上下界和逐点误差序列,拼成一份 GNUPlot 脚本并拉起外部绘图窗口,让你肉眼确认分布是否异常。 vis() 先扫一遍数组 x 取实际最大最小,再和传入的 max/min 做兜底比较,最后把纵轴范围在极值外各扩 1/20 间距,避免折线贴边。标题里动态塞进 numerr(错误计数),比如某次回测 numerr=3 会直接显示在图标题,方便你对照哪段样本出了问题。 grPlot() 写死调用 MQL5\Files 下的 gplot.txt,走 ShellExecuteW 拉 wgnuplot.exe;saveScript() 则负责把终端尺寸定为 560x420、字号 8 的 Windows 增强窗口,再写入脚本内容。若 FileOpen 返回 INVALID_HANDLE 就放弃,不会崩 EA。 你只要在 MT5 里把 GNUPlot 的 wgnuplot.exe 放到 MQL5\Files\GNUPlot\binary 下,跑完这套函数就能弹出暗蓝折线与红绿基准线——外汇与贵金属波动剧烈、样本外误差可能跳变,看图时别把某次 numerr 低当成稳健证据,只是当前窗口内的现象。
class="kw">return(n-m); } class=class="str">"cmt">//---------------------------------------------------------------------------- class="type">void vis(const class="type">class="kw">double &x[],class="type">class="kw">double dcen,class="type">class="kw">double min,class="type">class="kw">double max,class="type">int numerr) { class="type">int i; class="type">class="kw">double d,yma,ymi; class="type">class="kw">string str; yma=x[class="num">0];ymi=x[class="num">0]; for(i=class="num">0;i<ArraySize(x);i++) { if(yma<x[i])yma=x[i]; if(ymi>x[i])ymi=x[i]; } if(yma<max)yma=max; if(ymi>min)ymi=min; d=(yma-ymi)/class="num">20.0; yma+=d;ymi-=d; str="unset key\n"; str+="set title &class="macro">#x27;Sequence and error levels(number of errors = "+ (class="type">class="kw">string)numerr+")&class="macro">#x27; font &class="macro">#x27;,class="num">10&class="macro">#x27;\n"; str+="set yrange ["+(class="type">class="kw">string)ymi+":"+(class="type">class="kw">string)yma+"]\n"; str+="set xrange [class="num">0:"+(class="type">class="kw">string)ArraySize(x)+"]\n"; str+="plot "+(class="type">class="kw">string)dcen+" lt rgb &class="macro">#x27;green&class="macro">#x27;,"; str+=(class="type">class="kw">string)min+ " lt rgb &class="macro">#x27;red&class="macro">#x27;,"; str+=(class="type">class="kw">string)max+ " lt rgb &class="macro">#x27;red&class="macro">#x27;,"; str+="&class="macro">#x27;-&class="macro">#x27; with line lt rgb &class="macro">#x27;dark-blue&class="macro">#x27;\n"; for(i=class="num">0;i<ArraySize(x);i++)str+=(class="type">class="kw">string)x[i]+"\n"; str+="e\n"; if(!saveScript(str)){Print("Create script file error");class="kw">return;} if(!grPlot())Print("ShellExecuteW() error"); } class=class="str">"cmt">//---------------------------------------------------------------------------- class="type">bool grPlot() { class="type">class="kw">string pnam,param; pnam="GNUPlot\\binary\\wgnuplot.exe"; param="-p MQL5\\Files\\gplot.txt"; class="kw">return(ShellExecuteW(NULL,"open",pnam,param,NULL,class="num">1)); } class=class="str">"cmt">//---------------------------------------------------------------------------- class="type">bool saveScript(class="type">class="kw">string scr1="",class="type">class="kw">string scr2="") { class="type">int fhandle; fhandle=FileOpen("gplot.txt",FILE_WRITE|FILE_TXT|FILE_ANSI); if(fhandle==INVALID_HANDLE)class="kw">return(class="kw">false); FileWriteString(fhandle,"set terminal windows enhanced size class="num">560,class="num">420 font class="num">8\n"); FileWriteString(fhandle,scr1); if(scr2!="")FileWriteString(fhandle,scr2); FileClose(fhandle); class="kw">return(true); } class=class="str">"cmt">//---------------------------------------------------------------------------- m=(n-class="num">1)/class="num">2; median=a[m]; if((n&0x01)==class="num">0)b[class="num">1]=(median+a[m+class="num">1])/class="num">2.0; m=n/class="num">4; MQR=(a[m]+a[n-m-class="num">1])/class="num">2.0; class=class="str">"cmt">// Midquartile range m=n/class="num">4; IQM=class="num">0; for(i=m;i<n-m;i++)IQM+=a[i];