开头先交代一下我为什么会写这篇。前阵子给某个预测类项目做方案，数据是常规的多输入单输出回归，样本量不算大，大概几百条记录，特征七八个。一开始我图省事，直接用Matlab自带的fitrgp跑高斯过程回归，默认参数一把梭。结果测试集上预测曲线和真实值之间那条误差带简直像心电图，完全没法交付。后来认真把GPR的超参数调了一遍，精度立刻上了一个台阶，我才意识到：GPR这个模型对超参数极其敏感，而“怎么找超参数”才是真正值得花力气的地方。那段时间我试了网格搜索、贝叶斯优化、粒子群，最后在一篇新算法论文里看到蜣螂算法（DBO），自己动手把DBO和高斯过程回归拼在一起，做了个Matlab版的DBO-GPR多变量回归预测方案，又叠了一层交叉验证来做可靠评估。这套流程跑通之后，预测精度和稳定性都比原来默认参数强太多了。这篇文章就围绕这套方案，把原理、代码骨架、交叉验证的正确姿势和踩过的坑一次性讲透。
1. 这套方案到底在解决什么问题——以及为什么标准GPR不够用
1.1 GPR不是拿来就能用，超参数才是灵魂
高斯过程回归的核心思想是：不去直接学一个f(x)，而是对函数的分布做推断。也就是说，我们认为所有可能的函数都服从一个高斯过程先验，然后根据观测数据更新出后验，预测的时候给出均值作为预测值，同时给出方差作为不确定性区间。
这个思路在样本量小、噪声水平不明确的时候特别香。因为它天然自带不确定性估计，你不仅能告诉用户“预测值是多少”，还能告诉他“这个预测值有多可信”。而且GPR是非参数方法，表达能力强，不需要像神经网络那样精心设计结构。
但代价也很明显：结果好不好，几乎完全取决于核函数和超参数的设置。用得最多的平方指数核（RBF核）长这样：
k(x_i, x_j) = sigma_f^2 * exp( -0.5 * ||x_i - x_j||^2 / sigma_l^2 )其中sigma_l是长度尺度，控制函数光滑程度；sigma_f是信号标准差，控制函数整体的振幅范围。再加一个观测噪声项sigma_n，就是GPR里最常碰到的三个超参数。
问题就出在这三个参数上。sigma_l给太小，模型会把噪声也当信号学进去，预测曲线剧烈抖动；sigma_l给太大，模型又过度光滑，连真实的趋势变化都被磨平了。sigma_f和sigma_n的比例失调，同样会导致模型要么过度自信要么过度保守。如果是多维输入特征，还可以换用自动相关性决定（ARD）核，每一个特征有一个独立的长度尺度，参数数量直接变成特征数加2，手动调的噩梦级别。
1.2 标准fitrgp直接跑会遇到的具体困境
Matlab里fitrgp函数其实提供了OptimizeHyperparameters选项，可以自动优化超参数。但它默认用的是贝叶斯优化，在小样本、低维参数空间上效果尚可，一旦遇到特征维度稍高、数据分布比较复杂的情况，就容易出现两个问题：一是优化迭代次数受限，搜索不充分，结果随初始点波动很大；二是贝叶斯优化的采集函数在非凸的似然面上容易陷入局部极值。
手动调参更不靠谱。三个参数互相耦合，你调一个固定另外两个，永远只能看到局部切片，看不到全局地形。而且GPR的对数边际似然函数并不是简单的高斯形状，它可能是多峰的，甚至存在若干个数值相近但性质完全不同的参数组合。
我当时拿某组多变量数据做了个对比实验：默认参数的测试集RMSE大约是2.6，手动调参后降到1.9，但再用优化算法搜一圈，能到1.3。这个差距在工程上基本就是“能不能用”和“看起来能用”的区别。
1.3 为什么选蜣螂算法而不是粒子群、遗传算法
既然要优化GPR的超参数，那选择哪个优化算法就很重要。我对比过粒子群（PSO）、遗传算法（GA）、灰狼算法（GWO）、樽海鞘算法（SSA），最后选了DBO，原因有三点：
第一，收敛速度快。DBO的种群更新机制里，大部分个体承担的是局部开发任务，少数个体负责全局探索，分工明确，前期不需要像PSO那样反复调节惯性权重和学习因子。
第二，实现简单且参数少。DBO没有遗传算法里交叉率、变异率那一堆设置，也没有粒子群里的速度更新公式，整个算法的核心就是几种位置更新规则，Matlab实现起来非常干净。
第三，不容易早熟。虽然它的局部开发能力强，但偷窃行为这一支始终保留了对外部区域的搜索通道，相当于一直有一个“侦察兵”在外部转悠，降低了种群整体陷入局部最优的概率。
不过要强调一句：不选PSO不代表PSO不行，只是在这个特定场景下，DBO的性价比更高。如果你手头的GPR超参数搜索范围已经缩得很小，PSO也完全够用。
2. DBO优化的核心：蜣螂的六种行为到底在搜什么
2.1 从滚粪球到参数寻优的映射逻辑
蜣螂算法是2022年提出的一种群体智能优化算法。它模拟的是蜣螂的六种生存行为：滚球、跳舞、繁殖、觅食、偷窃和飞行。我第一次看到这个算法时想着估计也就是蹭个热点，读完原始论文之后发现它的行为建模其实挺有意思。
在自然界的蜣螂种群中，滚粪球的主力个体负责把粪球从一个地方推到另一个地方，它们需要根据太阳、风向来直线推进，遇到障碍物就得停下来调整方向重新规划路线。跳舞行为则发生在滚球个体无法确定前进方向的时候，它们通过旋转身体重新定向。繁殖行为是雌性蜣螂把粪球滚到安全区域附近，在周围产卵孵化。觅食行为是成年蜣螂在局部范围内搜索食物资源。偷窃行为是最有意思的——一部分蜣螂不自己滚球，而是专门去偷别的个体藏好的粪球。
映射到优化问题上，每一种行为都代表一种搜索策略：
- 滚球：全局探索，对应大步长的位置更新
- 跳舞：当探索方向受阻时，重新随机化一个方向再探索
- 繁殖产卵：在局部最优区域周围划定边界进行精细搜索
- 觅食：在已知较好区域内做小范围的局部搜索
- 偷窃：从全局最优位置出发向外探索，防止种群停滞
这几种行为在每次迭代中按比例分配给不同的个体，相当于同时开了几条搜索链路。
2.2 滚球与跳舞的Matlab更新规则
滚球行为的核心位置更新公式是这样的：
X_new = X + alpha * k * X_old + beta * |X - X_best| + lambda * theta * (X - X_best)这里alpha是随机扰动权重，k是偏置系数，beta和lambda控制向全局最优位置靠拢的强度，theta是方向向量。简化实现时，很多代码会把最后一项合并到第二项里，只保留两个主要寻优分量，实际效果也够用。
跳舞行为是在随机条件下触发的，当随机数小于某个阈值时，滚球个体不再沿用原来的方向向量，而是重新生成一个随机方向再执行更新，相当于在全局探索中加入了一次重定向。
繁殖行为的位置更新围绕一个动态变化的“产卵区”进行，这个区域的上下界会随着迭代进程逐渐向当前全局最优位置收缩。落在产卵区内的个体在小范围内随机扰动，实现对局部区域的加密搜索。
觅食行为和繁殖类似，但它是以全局最优位置为圆心，在一定半径内做随机游走，而且这个半径也会随迭代次数变化。偷窃行为则是在最优位置的邻域内进行一次较大的随机跳跃，保证部分个体不会永久停留在某个局部区域。
2.3 GPR的优化变量是什么——从3维到D+2维
讲完算法侧，回到GPR侧。DBO要优化的变量就是GPR的核函数超参数。以Matlab中fitrgp使用的平方指数核为例，优化的变量有三个：
| 参数 | 含义 | 典型范围 |
|---|---|---|
| SigmaL | 长度尺度，控制函数光滑程度 | [0.01, 10] |
| SigmaF | 信号标准差，控制函数幅值变化范围 | [0.01, 10] |
| SigmaN | 观测噪声标准差 | [0.001, 1] |
如果换用ARD平方指数核ardsquaredexponential，那么每个输入维度都有一个独立的SigmaL，优化变量数量变成特征数加2。比如你有7个特征，优化变量就是9维。维度一上去，网格搜索直接失效，手动调参更是无从下手，这也是为什么DBO这类群体优化算法在这里特别有优势。
需要注意的是，SigmaL、SigmaF、SigmaN这三个参数在fitrgp中都被限制为正值，所以DBO种群初始化时就要加上下界约束，在更新位置后还要加一个边界截断处理，否则会出现负值导致预测结果直接报错。
3. 把DBO和GPR焊在一起的关键：fitness函数与超参数编码设计
3.1 fitness函数的两种设计路线，以及我最终选的那种
DBO优化GPR，核心在于怎么定义“适应度”。我在实际实现中试过两条路线。
路线一：交叉验证均方根误差。每次迭代中，每个蜣螂个体携带一组超参数，把这组参数传给fitrgp，在训练集上训练模型，然后用验证集预测，计算均方根误差作为适应度值。交叉验证的评估结果更稳健，不容易被某一次划分影响。代价是训练多个GPR模型，计算量很大。一个个体一次迭代要训练K个GPR模型，整个算法要跑几十次迭代，如果数据量稍大，这个方案会很慢。
路线二：负对数边际似然作为适应度。GPR在训练过程中本身就能计算出对数边际似然值，这个值衡量的是“在当前超参数下，观测数据出现的概率”。对数边际似然越大，说明超参数与数据的契合程度越高。所以适应度可以直接取负对数边际似然（NLMLL）。好处是不需要额外汇总多次模型预测误差，训练一次GPR就能得到一个评估值；坏处是边际似然优化容易过拟合训练数据，对噪声数据的适应性不如交叉验证。
两条路线我最后都实现了，结论是：在数据量不太大（比如一千条以内）时，优先选路线一；数据量大、时间紧张时，选路线二。实际上很多公开发表的研究文章里，写的就是路线一，因为交叉验证结果看起来更“稳妥”，审稿人也不会挑毛病。如果你是自己做工程，我更推荐一种折中方案：用路线二做DBO搜索的大方向，最后用一次K折交叉验证来确认最优超参数的泛化效果。
3.2 超参数编码与种群结构
设计fitness函数的另一个关键点是把蜣螂个体的位置向量和GPR超参数对应起来。以RBF核为例，个体位置是一个3维向量：
个体位置 [SigmaL, SigmaF, SigmaN]种群是一个PopulationSize * 3的矩阵。初始化时：
lb = [0.01, 0.01, 0.001]; % 下界 ub = [10, 10, 1]; % 上界 population = lb + rand(N, 3) .* (ub - lb);每一轮迭代，先对种群中每个个体调用fitness函数，得到适应度数组，然后根据适应度找到当前全局最优位置X_best。之后按照DBO的行为规则依次更新种群，再次计算适应度，直到满足终止条件（达到最大迭代次数或误差低于阈值）。
这里很容易踩一个坑：fitrgp要求传入的超参数是正数且有有限值，但DBO更新过程中，有些个体位置可能被更新成负值或无穷大。我的做法是在每次位置更新之后立刻做一个边界截断：
population(i, :) = max(population(i, :), lb); population(i, :) = min(population(i, :), ub);边界截断之后还要检查数值是否为有限值，把NaN向量替换成随机初始化值，避免fitrgp内部报错。
3.3 用Matlab实现fitness函数的关键写法
fitness函数的输入是一个超参数向量，输出是适应度值。假设使用验证集误差作为适应度，核心代码大致如下：
function rmse = fitnessFunc(theta, Xtrain, Ytrain, Xval, Yval) try mdl = fitrgp(Xtrain, Ytrain, ... 'KernelFunction', 'squaredexponential', ... 'KernelParameters', [theta(1), theta(2)], ... 'Sigma', theta(3), ... 'Standardize', true, ... 'OptimizeHyperparameters', 'none'); Ypred = predict(mdl, Xval); rmse = sqrt(mean((Yval - Ypred).^2)); catch rmse = 1e10; % 防止部分超参数组合导致异常，给予极大的适应度惩罚 end end注意KernelParameters这里传入的是一个1×2向量，包含SigmaL和SigmaF，而Sigma单独传第三个参数，因为fitrgp的Sigma选项既可以传数值也可以传优化选项，传数值时就是固定噪声标准差。
3.4 不确定性也要一起评估：不止看RMSE
GPR的独特优势是能输出预测方差，所以评价这套方案时，我除了看RMSE和R²之外，还额外统计了预测区间的覆盖率。就是把测试集每个样本的真实值落在预测均值 ± 1.96倍标准差这个区间内的比例。覆盖率太差，说明模型的不确定性估计失真，哪怕RMSE好看也不能用。
加这个评价维度还有个实际原因：很多时候我们做回归预测，最终业务上需要的不只是单一数值，而是需要知道一个“波动范围”，比如库存备货量的上下限。GPR若不优化好超参数，预测方差会整体偏大或者偏小，给业务决策带来误导。
4. 交叉验证在DBO-GPR里的正确嵌套方式——这一步最容易翻车
4.1 嵌套到哪一层，结果完全不一样
标题里带“交叉验证”，但交叉验证放在哪里，区别非常大。我见过不少代码把K折交叉验证直接放在DBO每代适应度评估里，也就是每个个体都做一遍5折CV，来算出平均RMSE作为适应度。这种做法确实能最真实地反映超参数的泛化能力，但计算开销很大。
更常见、也更推荐的做法是两层结构：
- 外层：把原始数据划分为训练集和测试集。测试集从头到尾不参与超参数搜索，只用于最终评估。
- 内层：DBO搜索过程中，对训练集内部做K折交叉验证，用K折平均RMSE作为适应度值。
这样做的好处是避免超参数选择直接看到测试集信息，防止“数据泄漏”。数据泄漏这词听起来高大上，其实就是你把答案偷看了一眼。如果你用整个数据集（包括测试集）去搜超参数，那交叉验证评估出来的性能会虚高，等真正上线换一批数据，效果立刻现原形。
4.2 时间序列数据要警惕的切分陷阱
如果你的数据有明确的时间顺序，比如传感器时序、销量序列，那交叉验证的切分方式要格外小心。常规的cvpartition随机打乱再切分，会把时间上相邻的样本搅混在一起，训练集和验证集之间可能存在强烈的时序相关性，导致验证集预测误差被低估。
我在某次处理某波动数据集时就吃过这个亏。随机切分的5折CV，平均RMSE只有0.8，看着很棒；结果换成按时间顺序前80%训练、后20%测试，RMSE直接涨到2.1。原因就是相邻时间点的样本高度相关，随机划分相当于让模型“见过了”验证集附近的数据点。
正确做法是按时间顺序留出连续的验证集，或者在切分时要求验证集样本与训练集之间保留一定的时间间隔。代码上可以用Matlab的cvpartition中的holdout，或者直接手动按索引切分：
idxTrain = 1:floor(0.8 * size(data, 1)); idxTest = floor(0.8 * size(data, 1)) + 1:size(data, 1);如果必须做K折且数据是时间序列，推荐用“前折训练、后折验证”的扩展式切分，而不是随机打乱的K折。
4.3 K值怎么选
交叉验证中K的选择对DBO-GPR的影响，很多时候被低估了。K=3折叠，训练集比例只有67%，GPR在小训练集上拟合出的函数形态可能有较大偏差，导致超参数搜索偏向某个并不真实的区域。K=10折叠，训练集比例有90%，GPR能学到更稳定的函数形态，但每次要训练10个模型，配合DBO的种群迭代，计算时间可能翻好几倍。
我自己的经验是：样本量在300到1000条时，用5折比较平衡；超过1000条时可以用3折，因为训练集比例已经足够高；低于200条时，干脆不要做K折，直接用留一法或者70/30固定划分。留一法在GPR上有天然优势——GPR本身依赖于协方差矩阵求逆，留一法训练K个模型的成本极高，但如果数据量只有几十条，这个成本可以接受。
另外还有一点：K折CV作为适应度函数时，折数越多，适应度曲线的噪声越小，DBO的搜索过程越稳定。但K越大，单次适应度评估耗时越长，整个优化过程的实时性会变差。做工程的话，我一般先用3折快速摸一遍参数空间，找到大致最优区域后，再用5折或10折精搜。效果类似先粗调再微调。
5. Matlab代码骨架：麻雀虽小五脏俱全
5.1 整体流程与文件组织
我的实现拆成了四个文件：主脚本、DBO优化主循环函数、fitness函数、结果可视化脚本。主脚本负责读数据、归一化、划分训练测试集、调用DBO优化、评估最终模型。
%% 主脚本：DBO-GPR多变量回归预测 clear; clc; rng(42); % 固定随机种子，保证结果可复现 % 1. 加载数据，假定data是N行D列的表格，最后一列为输出 % data = readmatrix('yourdata.csv'); X = data(:, 1:end-1); Y = data(:, end); % 2. 归一化 [Xnorm, Xps] = mapminmax(X', 0, 1); Xnorm = Xnorm'; [Ynorm, Yps] = mapminmax(Y', 0, 1); Ynorm = Ynorm'; % 3. 划分训练集与测试集 numSamples = size(X, 1); idxTrain = 1:floor(0.8 * numSamples); idxTest = floor(0.8 * numSamples) + 1:numSamples; XTrain = Xnorm(idxTrain, :); YTrain = Ynorm(idxTrain, :); XTest = Xnorm(idxTest, :); YTest = Ynorm(idxTest, :); % 4. DBO参数设置 N = 30; % 种群规模 MaxIter = 60; % 迭代次数 dim = 3; % 优化变量数量 lb = [0.01, 0.01, 0.001]; ub = [10, 10, 1]; % 5. 调用DBO优化 [bestTheta, bestFitness, convergence] = DBOGPR(...);5.2 DBO主循环代码核心片段
完整的DBO实现代码不算长，核心结构如下：
function [bestPos, bestFit, curve] = DBOGPR(...) population = lb + rand(N, dim) .* (ub - lb); fitness = arrayfun(@(i) fitnessFunc(population(i,:), ...), 1:N); [bestFit, bestIdx] = min(fitness); bestPos = population(bestIdx, :); for t = 1:MaxIter % 滚球更新 for i = 1:N alpha = rand; if rand < 0.3 % 跳舞行为：方向紊乱时重新定向 theta = rand(1, dim) - 0.5; population(i,:) = population(i,:) + alpha * abs(bestPos - population(i,:)) .* theta; else % 正常滚球，向全局最优靠拢 k = 0.1; population(i,:) = population(i,:) + alpha * k * population(i,:) ... + rand * (bestPos - population(i,:)); end end % 繁殖、觅食、偷窃行为更新（此处省略细节实现） % 每次更新后都要做边界截断 population = max(population, lb); population = min(population, ub); % 重新计算适应度并更新全局最优 % 记录收敛曲线 curve(t) = bestFit; end end实际完整版会多几块更新逻辑。繁殖行为要在当前最优附近生成一个小边界的搜索区域，只更新落在该区域的个体；觅食行为是从学习因子TempX出发做局部随机扰动；偷窃行为是在全局最优附近的位置进行随机跳跃。这几种行为组合起来，每轮迭代后的种群多样性比重会明显好过单纯全部个体向最优靠拢。
5.3 fitrgp训练与预测的可视化
得到最优超参数后，用整个训练集重新训练一次最终模型，再对测试集做预测。预测结果的可视化我一般画三张图：第一张是测试集预测值与真实值的拟合对比曲线；第二张是误差带图，用patch函数画95%置信区间；第三张是收敛曲线。
finalMdl = fitrgp(XTrain, YTrain, ... 'KernelFunction', 'squaredexponential', ... 'KernelParameters', [bestTheta(1), bestTheta(2)], ... 'Sigma', bestTheta(3), ... 'Standardize', true); [Ypred, Ysd] = predict(finalMdl, XTest); Ypred = mapminmax('reverse', Ypred', Yps)'; Ysd = Ysd'; % 测试点预测标准差预测方差Ysd这个输出很容易被忽略，但它是GPR区别于其他回归模型的最有价值信息。画图时把Ypred ± 1.96 * Ysd画成阴影带，直接能在图上看到模型在哪段区域自信、在哪段区域犹豫。
6. 实测踩坑复盘：我在这套流程里交过的学费
6.1 坑一：维度不一致导致fitrgp直接报错
第一次跑通DBO-GPR时，我在fitness函数里传参少传了一个维度。种群个体是3维，但KernelParameters里我却填了[theta(1), theta(2), theta(3)]三个值，结果fitrgp的RBF核只接受两个参数，直接报错“KernelParameters must be a vector of length 2”。这类问题粗看很低级，但在代码耦合度高的优化流程里非常隐蔽，因为报错信息不会直接告诉你“是你fitness函数的参数结构错了”，它只是告诉你fitrgp内部参数不对。
排查方法其实很朴素：先固定一组可行的超参数，直接调用fitnessFunc验证能否正常返回一个数值，能跑通再去跑DBO主循环。最好在fitness函数开头加一段断言：
assert(numel(theta) == 3, 'theta dimension must be 3');这种防御式编程在多层嵌套流程里能省好几小时的排查时间。
6.2 坑二：噪声标准差被优化成负值，正定矩阵直接崩
GPR内部要计算协方差矩阵的Cholesky分解，这个分解要求矩阵正定。如果噪声标准差SigmaN被更新成负值，或者SigmaL小到接近0，协方差矩阵可能变得病态甚至不正定，fitrgp就会报错。表面上看是数值问题，深层次原因是DBO搜索空间没有完全约束住。
解决方式有两个层面。第一层是边界截断，我没收住的所有情况下都要执行；第二层是给SigmaN设一个下限，不能只依赖初始下界。我最终把lb的第三个分量从0.001改成了0.01，并且在fitness函数的catch块里对所有异常参数组合都返回一个极大惩罚值。这一步看似简单，对最终收敛稳定性的影响非常大。
6.3 坑三：时序数据上的交叉验证形同虚设
前面提到的时间序列切分问题，我在某个传感器数据集上实际吃过亏。那次是用默认的随机5折CV作为DBO的适应度函数，搜出来的最优参数在测试集上表现尚可但不如CV结果那么亮眼。后来我特意去检查了适应度评估阶段用的CV划分，发现随机划分把相邻时间点的样本揉在同一个折里，等于模型通过相邻样本“间接看到”了验证集的模式。
修正方式是用分组时间序列交叉验证。简单说就是先把样本按时间顺序分成M组，训练时永远只使用靠前的组，验证时使用下一组，这样每个验证集都严格晚于训练集。DBO-GPR里实现这一点，只需要在fitness函数中手动维护一个基于索引的划分逻辑，别用cvpartition默认的随机划分。
6.4 坑四：蜣螂种群早熟，所有个体挤进同一个局部最优
DBO虽然局部开发能力强，但在某些地形复杂的适应度面上，仍然会出现种群多样度下降过快的情况。直观表现是收敛曲线前20代快速下降，后40代几乎不动，测试集上的误差明显高于预期。我曾在一个特征强相关的数据集上遇到这个问题，有效特征只有两个但维度有八个，DBO在迭代后期所有个体的三个参数都在非常窄的范围内振荡。
针对早熟，我加了一个种群多样性检测。当所有个体与全局最优位置的距离都很小时，抽取一部分个体重新随机初始化，相当于给种群注入新的多样性。实际效果是测试集RMSE从2.1降到了1.6，收敛后期曲线又继续下探了一段。这个策略和遗传算法里的变异思想很像，但只针对全局最优个体周围的局部拥挤群体做，计算开销很可控。
7. 收敛性与结果解读：怎么判断优化真的跑好了
7.1 收敛曲线怎么看科学
收敛曲线是最直接的工具，但很多人只盯着“最终值有多小”，不看形态。我总结出三个判断要点：
- 曲线初始阶段是否快速下降。如果前5代几乎没有下降，大概率是种群初始范围设计不合理，导致所有个体离最优区域太远。
- 曲线中后期是否还在缓慢下降。后期即便下降幅度小，只要不是完全水平，说明算法还在做有效探索，没有完全停滞。
- 曲线是否出现突然的“大跳水”。如果某代适应度突然从2.0降到0.8，要警惕是不是种群更新逻辑有bug，比如某个个体的位置越界后恰好撞到了边界上的最优——出现这种情况需要检查更新公式。
7.2 种群规模与迭代次数的经验值
我试过几组参数对比。种群规模20、迭代40次的方案，在某些数据上能收敛，但在复杂数据上不够充分；改成种群30、迭代60次，最终RMSE大约能降低10%到20%；再加到种群50、迭代100次，效果提升就不明显了，但计算时间几乎翻倍。
结论是：先小规模跑通，再根据收敛曲线形态决定要不要加大预算。我一般先用N=20、MaxIter=30做一次快速扫描，观察最优适应度在迭代末段是否还在明显下降。如果还在降，就加大规模继续跑；如果已经平缓，就直接用当前参数。
7.3 结果评估指标别只看R²
测试集上最终评估，我最常用的一组指标包括RMSE、MAE、R²和95%预测区间覆盖率。R²在这个场景下容易“虚高”，因为如果数据本身就存在强趋势，即使预测偏了也会因为方差分解显得很好看。RMSE更直观，MAE可以反映绝对误差的平均水平。
我还习惯把预测误差按测试集各个样本分别画出来，观察是否集中在某些区域。如果误差在数据变化剧烈的区段突然变大，说明模型训练数据在该区域覆盖不足，可能需要补采数据而不是换模型。
7.4 参数与性能的对应关系
从多次实验中我观察到几个规律性的现象：SigmaL偏大时，模型光滑，对趋势的拟合滞后明显；SigmaL偏小在训练集上表现很好但在测试集上变差，过拟合信号；SigmaN如果被优化到接近下界，往往意味着数据噪声被低估，预测区间会偏窄。所以跑完DBO之后，我会把最优超参数和默认超参数各训练一遍，对比两者的测试集误差带宽度。如果优化后区间明显更合理（覆盖率更接近95%），那这套方案才算真正起作用。
一点个人体会
DBO-GPR这套流程，技术上不能算什么新东西——群体优化算法套GPR超参数搜索，方法论很成熟了。但这个组合在Matlab里落地时，真正决定成败的是细节：适应度函数的稳定性、交叉验证的切分方式、边界处理与异常防御、种群多样性的保持。我在实际操作中最大的体会是：优化算法只是搜索工具，GPR本身才是建模核心，所以一定不要只顾着调DBO的参数，而忽视GPR模型自身的合理性。
后续想继续扩展的话，有两个方向比较推荐。一是把核函数从RBF换成ARD平方指数核，让每个输入维度拥有独立长度尺度，在特征相关性复杂的场景下能明显提升拟合精度；二是对比贝叶斯优化与DBO在不同维度超参数空间上的表现，两者的耗时与最终精度的权衡在不同数据规模下差异很大。这套Matlab代码本身也预留了接口，换核函数、换数据集基本只需要改几行参数，有需要的话我可以再写一篇针对ARD核和多元输出的扩展实现。
