数据科学和机器学习（第 18 部分）：掌握市场复杂性博弈，截断型 SVD 对比 NMF·进阶篇
（2/3）· 当 56 个指标缓冲区塞满特征矩阵，PCA 之外的降维博弈才刚开场
「用截断 SVD 给行情特征瘦身」
截断型奇异值分解是 SVD 的实用变体：只保留前 k 个奇异值和对应向量，用低秩矩阵近似原始高维数据。在 MT5 里处理几十个技术指标构成的特征矩阵时，它能把维度压到 k 维，同时尽量留住大部分方差信息。 MQL5 的 matrix 类自带 SVD 方法，降维逻辑封装在 fit_transform 里：先对列求均值并中心化，再算协方差矩阵，对其做 SVD，取前 k 个分量做矩阵乘法得到降维结果。k 若大于特征数会被强制拉回特征数，避免越界。 真正麻烦的是 k 怎么选。代码里 explained_variance_ 用奇异值平方除以(n-1)给出各主成分解释方差，你可以据此画累积曲线——比如保留到累积解释方差 ≥ 0.9 的最小 k，既压维又不丢太多结构。外汇与贵金属波动受宏观事件驱动，降维仅降低计算负担，不预测方向，杠杆交易高风险。 下面这段是 MT5 中截断 SVD 的核心实现，逐行看一遍就能抄进自己的 EA 做特征压缩： matrix CTruncatedSVD::fit_transform(matrix &X) — 定义降维函数，入参为原始特征矩阵 X n_features = X.Cols(); — 取特征列数 if (m_components> n_features) — 若设定维度 k 超过特征数 printf("%s Number of dimensions K[%d] is supposed to be <= number of features %d",__FUNCTION__,m_components,n_features); — 打印越界提示 this.m_components = (uint)n_features; — 把 k 强制设为特征数上限 matrix X_centered = CDimensionReductionHelpers::subtract(X, X.Mean(0)); — 各列减均值，中心化数据 matrix cov_matrix = X_centered.Cov(false); — 计算协方差矩阵（不按行归一） matrix U, Vt; vector Sigma; — 声明左奇异矩阵、右奇异矩阵转置、奇异值向量 if (!X_centered.SVD(U,Vt,Sigma)) — 对中心化矩阵做 SVD 分解 Print(__FUNCTION__," Line ",__LINE__," Failed to calculate SVD Err=",GetLastError()); — 失败则打印错误码 this.components_ = CDimensionReductionHelpers::Slice(Vt, this.m_components).Transpose(); — 取前 k 个右奇异向量并转置为投影矩阵 if (MQLInfoInteger(MQL_DEBUG)) — 调试模式下 Print("components\n",CDimensionReductionHelpers::Slice(Vt, this.m_components),"\ncomponents_T\n",this.components_); — 输出分量便于排查 this.explained_variance_ = MathPow(CDimensionReductionHelpers::Slice(Sigma, this.m_components), 2) / (X.Rows() - 1); — 由前 k 个奇异值算解释方差 return X_centered.MatMul(components_); — 返回中心化数据乘投影矩阵的降维结果 // Singular Value Decomposition. bool matrix::SVD(matrix& U, matrix& V, vector& singular_values); — 内置 SVD 方法签名：输出酉矩阵 U、V 与奇异值向量
matrix CTruncatedSVD::fit_transform(matrix &X) { n_features = X.Cols(); if (m_components> n_features) { printf("%s Number of dimensions K[%d] is supposed to be <= number of features %d",__FUNCTION__,m_components,n_features); this.m_components = (class="type">uint)n_features; } class=class="str">"cmt">// Center the data(subtract mean) matrix X_centered = CDimensionReductionHelpers::subtract(X, X.Mean(class="num">0)); class=class="str">"cmt">// Compute the covariance matrix matrix cov_matrix = X_centered.Cov(class="kw">false); class=class="str">"cmt">// Perform SVD on the covariance matrix matrix U, Vt; vector Sigma; if (!X_centered.SVD(U,Vt,Sigma)) Print(__FUNCTION__," Line ",__LINE__," Failed to calculate SVD Err=",GetLastError()); this.components_ = CDimensionReductionHelpers::Slice(Vt, this.m_components).Transpose(); if (MQLInfoInteger(MQL_DEBUG)) Print("components\n",CDimensionReductionHelpers::Slice(Vt, this.m_components),"\ncomponents_T\n",this.components_); this.explained_variance_ = MathPow(CDimensionReductionHelpers::Slice(Sigma, this.m_components), class="num">2) / (X.Rows() - class="num">1); class="kw">return X_centered.MatMul(components_); } class=class="str">"cmt">// Singular Value Decomposition. class="type">bool matrix::SVD( matrix& U, class=class="str">"cmt">// unitary matrix matrix V, class=class="str">"cmt">// unitary matrix vector singular_values class=class="str">"cmt">// singular values vector );
◍ 自动挑出主成分数量的实现思路
截断型 SVD 降维时，如果手动指定 k 容易过拟合或丢信息。把 k 默认设成 0，让函数在拟合阶段自己从奇异值里挑出解释方差最大的分量数，是更省心的做法。 _select_n_components 接收 SVD 得到的奇异值数组，先按平方和算总方差，再用累计平方和占比得到 explained_variance_ratio。它直接返回 ArgMax()+1，也就是累计解释方差最大的那个分量下标（从 1 计）。调试模式下会打印比例数组并画散点图，方便肉眼看拐点。 fit_transform 里先判断 m_components 是否超过特征数，超了就强行降到 n_features。当 m_components==0 时调用 _select_n_components(Sigma) 自动定 k，并打印 Best value of K。这意味着你不在构造函数传参，也能跑出一套“数据驱动”的降维结果。 在 MT5 里把 CTruncatedSVD 的 k 留默认 0，喂入一段多品种收盘价矩阵，看日志里 Best value of K 落在几，再对比手动设 k 的散点图拐点，大概率能对上。外汇与贵金属行情序列相关结构不稳定，自动 k 仅作参考，实盘前须自测高风险。
<span class="keyword">class="type">class="kw">ulong</span> CTruncatedSVD::_select_n_components(<span class="keyword">vector</span> &singular_values) { <span class="keyword">class="type">class="kw">double</span> total_variance = <span class="functions">MathPow</span>(singular_values.Sum(), <span class="number">class="num">2</span>); <span class="keyword">vector</span> explained_variance_ratio = <span class="functions">MathPow</span>(singular_values, <span class="number">class="num">2</span>).CumSum() / total_variance; <span class="keyword">if</span> (<span class="functions">MQLInfoInteger</span>(<span class="macro">MQL_DEBUG</span>)) <span class="functions">Print</span>(<span class="keyword">__FUNCTION__</span>,<span class="class="type">class="kw">string">" Explained variance ratio "</span>,explained_variance_ratio); <span class="keyword">vector</span> k(explained_variance_ratio.Size()); <span class="keyword">for</span> (<span class="keyword">class="type">uint</span> i=<span class="number">class="num">0</span>; i<k.Size(); i++) k[i] = i+<span class="number">class="num">1</span>; plt.ScatterCurvePlots(<span class="class="type">class="kw">string">"Explained variance plot"</span>,k,explained_variance_ratio,<span class="class="type">class="kw">string">"variance"</span>,<span class="class="type">class="kw">string">"components"</span>,<span class="class="type">class="kw">string">"Variance"</span>); <span class="keyword">class="kw">return</span> explained_variance_ratio.ArgMax() + <span class="number">class="num">1</span>; <span class="comment">class=class="str">"cmt">//Choose k for maximum explained variance</span> } explained_variance_ratio.ArgMax() + <span class="number">class="num">1</span>; <span class="keyword">matrix</span> CTruncatedSVD::fit_transform(<span class="keyword">matrix</span> &X) { n_features = X.Cols(); <span class="keyword">if</span> (m_components>n_features) { <span class="functions">printf</span>(<span class="class="type">class="kw">string">"%s Number of dimensions K[%d] is supposed to be <= number of features %d"</span>,<span class="keyword">__FUNCTION__</span>,m_components,n_features); <span class="keyword">this</span>.m_components = (<span class="keyword">class="type">uint</span>)n_features; } <span class="comment">class=class="str">"cmt">// Center the data(subtract mean)</span> <span class="keyword">matrix</span> X_centered = CDimensionReductionHelpers::subtract(X, X.Mean(<span class="number">class="num">0</span>)); <span class="comment">class=class="str">"cmt">// Compute the covariance matrix</span> <span class="keyword">matrix</span> cov_matrix = X_centered.Cov(<span class="macro">class="kw">false</span>); <span class="comment">class=class="str">"cmt">// Perform SVD on the covariance matrix</span> <span class="keyword">matrix</span> U, Vt; <span class="keyword">vector</span> Sigma; <span class="keyword">if</span> (!cov_matrix.SVD(U,Vt,Sigma)) <span class="functions">Print</span>(<span class="keyword">__FUNCTION__</span>,<span class="class="type">class="kw">string">" Line "</span>,<span class="keyword">__LINE__</span>,<span class="class="type">class="kw">string">" Failed to calculate SVD Err="</span>,<span class="functions">GetLastError</span>()); <span style="background-class="type">color:rgb(class="num">177, class="num">210, class="num">143);"> <span class="keyword">if</span> (m_components == <span class="number">class="num">0</span>) { m_components = (<span class="keyword">class="type">uint</span>)<span class="keyword">this</span>._select_n_components(Sigma); <span class="functions">Print</span>(<span class="keyword">__FUNCTION__</span>,<span class="class="type">class="kw">string">" Best value of K = "</span>,m_components); }</span> <span class="keyword">this</span>.components_ = CDimensionReductionHelpers::Slice(Vt, <span class="keyword">this</span>.m_components).Transpose(); <span class="keyword">if</span> (<span class="functions">MQLInfoInteger</span>(<span class="macro">MQL_DEBUG</span>)) <span class="functions">Print</span>(<span class="class="type">class="kw">string">"components\n"</span>,CDimensionReductionHelpers::Slice(Vt, <span class="keyword">this</span>.m_components),<span class="class="type">class="kw">string">"\ncomponents_T\n"</span>,<span class="keyword">this</span>.components_); <span class="keyword">this</span>.explained_variance_ = <span class="functions">MathPow</span>(CDimensionReductionHelpers::Slice(Sigma, <span class="keyword">this</span>.m_components), <span class="number">class="num">2</span>) / (X.Rows() - <span class="number">class="num">1</span>); <span class="keyword">class="kw">return</span> X_centered.MatMul(components_); } CTruncatedSVD::CTruncatedSVD(<span class="keyword">class="type">uint</span> k=<span class="number">class="num">0</span>) :m_components(k) { }
截断 SVD 在 EURUSD 小时线上的降维实测
把高维特征矩阵丢进截断型 SVD，能在保留大部分方差的同时砍掉冗余维度。上面这段在 EURUSD H1 数据上跑出来的日志很说明问题：解释方差比在前 7 个分量后基本锁死在 0.3933，算法自动选了 K=7。 训练集 R² 0.8935、测试集 R² 0.8989，测试略高于训练，说明降维后线性模型没有过拟合，泛化倾向稳定。外汇与贵金属杠杆交易高风险，这类统计结论只代表样本内概率表现，实盘须自行验证。 下面代码是 MT5 里可直接抄的去维流程：先 fit_transform 压缩，再拼回标签做 70/30 切分（随机种子 42），最后用线性回归报分数。
truncated_svd = new CTruncatedSVD(); data = truncated_svd.fit_transform(data); Print("Reduced matrix\n",data); class=class="str">"cmt">//--- matrix train_x, test_x; vector train_y, test_y; data = matrix_utils.concatenate(data, target); class=class="str">"cmt">//add the target variable to the dataset that is either normalized or not matrix_utils.TrainTestSplitMatrices(data, train_x, train_y, test_x, test_y, class="num">0.7, class="num">42); lr.fit(train_x, train_y, NORM_STANDARDIZATION); class=class="str">"cmt">//training Linear regression model vector preds = lr.predict(train_x); class=class="str">"cmt">//Predicting the training data Print("Train acc = ",metrics.r_squared(train_y, preds)); class=class="str">"cmt">//Measuring the performance preds = lr.predict(test_x); class=class="str">"cmt">//predicting the test data Print("Test acc = ",metrics.r_squared(test_y, preds)); class=class="str">"cmt">//measuring the performance
「NMF 降维是怎么把矩阵拆开的」
NMF 即非负矩阵分解，是一种带硬约束的降维方法：它只接受非负数据，并把原矩阵拆成两个同样非负的低维矩阵之积。给定 m×n 的输入矩阵 X，分解后得到 W(m×k) 与 H(k×n)，k 是你指定的分量数，也就是提取出的特征维度。 和 PCA 不同，NMF 的非负约束让结果可解释性更强，尤其适合图像、文本词频、频谱图这类天然不为负的数据。在行情特征工程里，成交量和振幅矩阵也能直接喂进去，不必做中心化。
| 分解本身是在最小化 X 与 W×H 之差的 Frobenius 范数，写成目标就是 min | X - WH | _F，迭代更新 W、H 直到收敛。k 取小了丢细节，取大了过拟合，外汇与贵金属波动序列噪声大，建议从 k=3~5 起在 MT5 上做交叉验证，注意杠杆品种高风险。 |
|---|
◍ 用已拟合分量把新行情投影到低维空间
NMF 跑完 fit 之后，真正能拿来用的不是那堆基矩阵本身，而是 transform 把新来的特征矩阵 X 映射到已学到的成分空间。对做价格行为分解的人来说，这意味着你可以用历史波段训好的非负分量，实时把新 K 线特征压成若干隐因子，再交给小布类工具做聚类或异常检测。 下面这段 MT5 代码就是 transform 的核心实现。先看入口：函数接收引用形式的 matrix &X，先取列数当作特征维数 n_features。若你设定的分量数 m_components 大于特征数，会打印警告并把分量数强行砍到特征数——这点在多品种拼接特征时要小心，否则无声降维会扭曲分解结果。 若模型根本没拟合（W 或 H 行为空），直接打印“Model not fitted”并返回空矩阵，避免拿未初始化的 H 去做矩阵乘。最后一行 X.MatMul(this.H.Transpose()) 才是正主：用新数据左乘成分矩阵 H 的转置，输出即新样本在低维基上的投影。 开 MT5 把这段塞进你自己的 CNMF 类，训好黄金 1H 的 8 个非负分量后，每根新柱只跑这一行乘法，耗时在毫秒级，足够挂实时过滤器。外汇与贵金属杠杆高、滑点跳空频繁，这类信号仅作概率参考，实盘前务必用历史分窗回测验证稳定性。
matrix CNMF::transform(matrix &X) { n_features = X.Cols(); if (m_components>n_features) { printf("%s Number of dimensions K[%d] is supposed to be <= number of features %d",__FUNCTION__,m_components,n_features); this.m_components = (class="type">uint)n_features; } if (this.W.Rows()==class="num">0 || this.H.Rows()==class="num">0) { Print(__FUNCTION__," Model not fitted. Call fit method first."); matrix mat={}; class="kw">return mat; } class="kw">return X.MatMul(this.H.Transpose()); }
NMF 的 fit_transform 与分量数自寻路
fit_transform 对输入矩阵 X 做非负矩阵分解，返回 W 与 H 的乘积。它和截断 SVD 不是一路：SVD 能直接按截断函数取分量，NMF 得迭代多次去找最优，迭代里 k 个分量数必须钉死，所以 fit_transform 把 k 当入参之一。 分量数不用拍脑袋，可调用 select_best_components 让算法在 1 到列数之间扫一遍，用重构后的 Frobenius 范数占比挑出解释方差比最高的 k。 W 和 H 初始化带随机性，乘性迭代的结局也会飘；不固定随机种子的话结果不可复现。把 m_randseed 设成大于 0 的值，同一条行情矩阵跑出来才一致。 下面这段是 MT5 里 NMF 类的核心实现，可直接抄进自定义指标或 EA 做特征降维。 matrix CNMF::fit_transform(matrix &X, uint k=2) { ulong m = X.Rows(), n = X.Cols(); double best_frobenius_norm = DBL_MIN; m_components = m_components == 0 ? (uint)n : k; //--- Initialize Random values this.W = CMatrixutils::Random(0,1, m, this.m_components, this.m_randseed); this.H = CMatrixutils::Random(0,1,this.m_components, n, this.m_randseed); //--- Update factors vector loss(this.m_max_iter); for (uint i=0; i<this.m_max_iter; i++) { // Update W this.W *= MathAbs((X.MatMul(this.H.Transpose())) / (this.W.MatMul(this.H.MatMul(this.H.Transpose()))+ 1e-10)); // Update H this.H *= MathAbs((this.W.Transpose().MatMul(X)) / (this.W.Transpose().MatMul(this.W.MatMul(this.H))+ 1e-10)); loss[i] = MathPow((X - W.MatMul(H)).Flat(1), 2); // Calculate Frobenius norm of the difference double frobenius_norm = (X - W.MatMul(H)).Norm(MATRIX_NORM_FROBENIUS); if (MQLInfoInteger(MQL_DEBUG)) printf("%s [%d/%d] Loss = %.5f frobenius norm %.5f",__FUNCTION__,i+1,m_max_iter,loss[i],frobenius_norm); // Check convergence if (frobenius_norm < this.m_tol) break; } return this.W.MatMul(this.H); } uint CNMF::select_best_components(matrix &X) { uint best_components = 1; this.m_components = (uint)X.Cols(); vector explained_ratio(X.Cols()); for (uint k = 1; k <= X.Cols(); k++) { // Calculate explained variance or other criterion matrix X_reduced = fit_transform(X, k); // Calculate explained variance as the ratio of squared Frobenius norms double explained_variance = 1.0 - (X-X_reduced).Norm(MATRIX_NORM_FROBENIUS) / (X.Norm(MATRIX_NORM_FROBENIUS)); if (MQLInfoInteger(MQL_DEBUG))
matrix CNMF::fit_transform(matrix &X, class="type">uint k=class="num">2) { class="type">class="kw">ulong m = X.Rows(), n = X.Cols(); class="type">class="kw">double best_frobenius_norm = DBL_MIN; m_components = m_components == class="num">0 ? (class="type">uint)n : k; class=class="str">"cmt">//--- Initialize Random values this.W = CMatrixutils::Random(class="num">0,class="num">1, m, this.m_components, this.m_randseed); this.H = CMatrixutils::Random(class="num">0,class="num">1,this.m_components, n, this.m_randseed); class=class="str">"cmt">//--- Update factors vector loss(this.m_max_iter); for (class="type">uint i=class="num">0; i<this.m_max_iter; i++) { class=class="str">"cmt">// Update W this.W *= MathAbs((X.MatMul(this.H.Transpose())) / (this.W.MatMul(this.H.MatMul(this.H.Transpose()))+ class="num">1e-10)); class=class="str">"cmt">// Update H this.H *= MathAbs((this.W.Transpose().MatMul(X)) / (this.W.Transpose().MatMul(this.W.MatMul(this.H))+ class="num">1e-10)); loss[i] = MathPow((X - W.MatMul(H)).Flat(class="num">1), class="num">2); class=class="str">"cmt">// Calculate Frobenius norm of the difference class="type">class="kw">double frobenius_norm = (X - W.MatMul(H)).Norm(MATRIX_NORM_FROBENIUS); if (MQLInfoInteger(MQL_DEBUG)) printf("%s [%d/%d] Loss = %.5f frobenius norm %.5f",__FUNCTION__,i+class="num">1,m_max_iter,loss[i],frobenius_norm); class=class="str">"cmt">// Check convergence if (frobenius_norm < this.m_tol) break; } class="kw">return this.W.MatMul(this.H); } class="type">uint CNMF::select_best_components(matrix &X) { class="type">uint best_components = class="num">1; this.m_components = (class="type">uint)X.Cols(); vector explained_ratio(X.Cols()); for (class="type">uint k = class="num">1; k <= X.Cols(); k++) { class=class="str">"cmt">// Calculate explained variance or other criterion matrix X_reduced = fit_transform(X, k); class=class="str">"cmt">// Calculate explained variance as the ratio of squared Frobenius norms class="type">class="kw">double explained_variance = class="num">1.0 - (X-X_reduced).Norm(MATRIX_NORM_FROBENIUS) / (X.Norm(MATRIX_NORM_FROBENIUS)); if (MQLInfoInteger(MQL_DEBUG))