MATLAB实现基于向量误差修正模型（VECM）进行多变量时间序列预测的详细项目实例
项目背景介绍
向量误差修正模型，简称 VECM，是多变量时间序列分析中非常重要的一类动态建模工具，特别适用于多个经济、能源、金融或工业指标之间存在长期均衡关系的场景。现实数据中，很多变量并不是彼此独立波动，而是同时受到共同趋势、政策冲击、供需变化、市场情绪和外部环境扰动的影响。例如宏观经济中，利率、通胀、汇率、货币供应量与工业产出常常共同演化；能源系统中，电价、负荷、温度、燃料价格与发电结构也会形成长期联动；金融市场中，股票指数、成交量、波动率、利率曲线之间常存在协整关系。若直接使用普通回归或单变量自回归方法，往往只能描述局部波动，难以刻画变量之间长期均衡约束与短期偏离后的修正过程，因此预测结果常出现偏差累积、结构失真和解释力不足等问题。
VECM 的核心优势在于同时处理“短期动态”和“长期约束”。当多个时间序列都具有非平稳性，但它们之间存在协整关系时，模型可以将变量的差分项用于描述短期变化，将误差修正项用于反映偏离长期均衡后的回归力量。这样一来，模型不仅能够预测未来走势，还能够解释为何变量会围绕某一稳定关系波动，以及偏离后需要多久、以何种速度回归。对于需要进行中长期趋势判断、风险预警、资源调度和联动机制分析的项目而言，这种结构非常关键。
MATLAB 在时间序列分析、矩阵计算、统计建模和可视化方面具有成熟能力，尤其适合开展 VECM 项目。利用 MATLAB，可以完成从数据清洗、平稳性检验、协整检验、滞后阶数选择、模型估计、预测生成，到结果评估与图形展示的完整流程。相比手工推导与分散实现，MATLAB 能显著提高建模效率，减少计算错误，并便于对多变量结构进行系统调试。对于需要面向工程落地的项目，MATLAB 还能快速集成数据导入、结果导出、批量实验和界面展示等功能，使模型不只停留在理论层面，而是可以服务于实际决策。
在多变量预测任务中，VECM 的价值还体现在可解释性上。模型中的协整向量能够揭示变量之间的长期均衡关系，误差修正系数能够表征各变量对偏离程度的响应速度，滞后项系数则展示短期传导链路。对业务人员而言，这意味着预测结果不只是一个数值序列，还包含因果联动逻辑与结构性解释。对研究人员而言，这意味着能够进一步分析冲击传播、系统稳定性和变量间反馈机制。对工程应用而言，这意味着模型既能用于预测，也能用于异常识别和策略优化。
随着数据采集频率提升与监测维度增加，单一序列模型越来越难满足复杂系统分析需求。多变量时间序列往往存在异方差、季节性、结构突变、同步延迟、样本量不足等现实问题，而 VECM 恰好提供了一个兼顾统计理论与实际建模能力的框架。通过合理设置差分阶数、协整秩、滞后长度和外生变量处理方式，可以构建出适合具体业务场景的预测模型。若再结合滚动预测、残差诊断和可视化评估，就能形成一个完整的预测分析项目闭环。
MATLAB R2025b 版本在基础绘图、界面组件、机器学习接口与对象属性方面有一些细节变化，因此项目实现时需要更加注意版本兼容性。尤其在图形界面、颜色映射、组合对象和部分回归接口参数上，应优先采用当前版本稳定支持的写法，避免因旧接口失效而影响运行。基于这一背景，构建一个面向 MATLAB R2025b 的 VECM 多变量时间序列预测实例，不仅能展示统计建模思路，还能形成适合当前版本的规范化实现模板，为后续扩展到更多真实数据场景提供可靠基础。
项目目标与意义
构建多变量联动预测框架
项目首要目标是建立一个能够处理多个相关时间序列的预测框架，使模型不仅关注单个变量自身的历史变化，还能捕捉变量之间的同步关系、滞后传导与长期约束。VECM 的设计思想正好适合这一任务，因为它将差分动态与协整关系统一在同一表达式中，能够在存在共同趋势时保持结构稳定。通过构建这一框架，可以把复杂系统中的多个指标纳入统一分析路径，避免分别建模造成的信息割裂。这个目标的意义在于提升预测系统的整体一致性，使得一个变量的未来变化能够参考其他变量的历史状态，从而增强模型对复杂联动场景的适应能力。
提升长期均衡识别能力
第二个目标是识别变量之间的长期均衡关系，并将这种关系纳入预测机制。很多现实数据虽然在短期内频繁波动，但从较长时间尺度观察，会围绕某种稳定关系上下调整。若忽略这一点，模型容易把长期规律误判为随机扰动。VECM 的误差修正项正是为此设计，它描述偏离均衡后的修正方向和修正强度。这个目标的意义在于帮助预测模型不只反映“最近发生了什么”，还能够理解“系统本来应该往哪里回归”。对于需要把握趋势、判断拐点、研究价格联动与政策传导的任务，这种能力尤为重要。
强化短期动态刻画能力
第三个目标是准确刻画短期波动机制，包括冲击传播、滞后反应和局部调整。VECM 不是单纯的长期关系模型，它同时保留差分项中的短期信息，因此可以分析某个变量的变化如何在若干期内传递到其他变量。这个能力对于高频经济数据、工业监测数据、金融行情数据都非常关键。短期动态刻画得越准确，预测误差越小，异常判断越灵敏，策略制定也越有针对性。其意义在于帮助项目从静态均衡分析进一步走向动态响应分析，使模型既能看趋势，也能看波动。
服务决策与工程落地
第四个目标是把模型结果转化为可操作的决策支持信息。多变量预测项目的最终价值不在于单纯输出数值，而在于为资源配置、风险控制、产能安排、库存管理和政策评估提供依据。VECM 具备较强的解释性，能够输出协整关系、误差修正速度、残差结构和预测轨迹，这些结果可以直接服务于业务判断。这个目标的意义在于提升项目的落地价值，让统计模型真正进入应用流程，而不是停留在理论分析层面。对 MATLAB 项目来说，这也意味着不仅要完成建模，还要完成结果可视化、评估展示和导出接口设计。
项目挑战及解决方案
非平稳性与协整关系识别困难
多变量时间序列中最常见的困难之一是变量本身往往非平稳，均值和方差随时间变化，直接建模会导致伪回归和参数失真。VECM 只能在变量存在协整关系时发挥作用，因此必须先完成单位根检验与协整检验，判断哪些变量可以进入同一系统。解决方式是先对各序列进行平稳性分析，再对非平稳但同阶单整的变量进行 Johansen 协整检验，确认协整秩后再建立模型。若协整关系不明显，可考虑对变量进行变换、加入外生项或调整样本区间。这样的处理能确保模型建立在正确的数据结构上，避免把不适合的序列硬塞入 VECM。
滞后阶数与模型秩选择复杂
另一个重要挑战是滞后阶数和协整秩的选择。滞后阶数太小会遗漏动态信息，太大又会导致参数过多、估计不稳定和样本消耗过快。协整秩选择不当则会直接影响误差修正项的结构，使模型失去均衡约束意义。解决这一问题的方式是结合信息准则、残差自相关检验和逐步比较策略来确定最优配置。实际实现时可以先设定候选滞后范围，再根据 AIC、BIC、HQ 等指标进行筛选，同时结合残差诊断保证模型白噪声特征更明显。协整秩则通过检验统计量与临界值共同判断，必要时对不同秩重复拟合并比较预测性能。通过这种组合策略，可以在统计合理性和预测精度之间取得平衡。
结果解释与工程实现兼顾
VECM 项目常常面临一个现实问题：统计解释较强，但工程实现若处理不当，容易在数据对齐、尺度统一、缺失值处理、可视化输出等环节出错。尤其在 MATLAB 中，不同版本对图形、对象属性和部分建模接口支持存在差异，若直接套用旧写法，可能出现参数冲突或图形对象报错。解决方式是采用模块化流程，将数据预处理、检验、估计、预测和展示分开处理，并在代码中明确每一步输入输出关系。对于 R2025b，还要优先使用稳定接口，避免依赖已变化的属性名或对象关系。这样既能确保模型逻辑完整，也能提高代码可维护性和版本兼容性，使项目真正具备可复现性和可部署性。
项目模型架构
数据采集与预处理层
模型架构的第一层是数据采集与预处理层，负责把原始多变量时间序列转化为适合 VECM 分析的标准格式。该层通常包括时间索引对齐、缺失值修复、异常值检查、必要的对数变换或差分处理、尺度统一和样本截取等环节。多变量时间序列往往来自不同来源，采样频率不一致，单位也可能不同，因此预处理是整个系统的基础。若此层处理不完整，后续协整检验与模型估计都会受到干扰。其原理在于通过数据清洗消除噪声性偏差，使各变量在同一时间轴和同一分析框架内可比。MATLAB 在矩阵处理方面效率高，适合将这些操作批量完成，并为后续统计检验输出干净的数据矩阵。
平稳性与单整阶数检验层
第二层是平稳性与单整阶数检验层，其作用是判断各变量是否满足 VECM 的前提条件。VECM 通常要求变量为同阶单整，最常见的是 I(1) 过程，即原序列非平稳、一阶差分后平稳。该层一般采用单位根检验、差分检验与可视化趋势检查相结合的方法。其基本原理是通过统计量判断时间序列是否存在单位根，从而识别随机趋势是否显著。如果变量阶数不一致，则需要重新选择变量、转换数据或调整建模策略。这个环节直接决定模型是否可用，是从一般多变量分析过渡到协整建模的关键门槛。MATLAB 中可以利用统计检验函数完成这一层工作，并将结果用于后续协整秩判断。
协整关系识别层
第三层是协整关系识别层，其任务是判断多个非平稳变量之间是否存在长期稳定组合。协整的含义是，尽管单个变量会随时间漂移，但某个线性组合却保持平稳，从而说明这些变量之间存在均衡关系。该层通常通过 Johansen 方法进行检验，得到特征根、迹统计量、最大特征值统计量以及协整秩。其基本原理是对 VAR 系统的长期关系进行矩阵分解，从中识别能够稳定描述系统均衡的向量。协整秩越大，表示长期关系越丰富；秩为零则说明不存在稳定协整关系，VECM 的适用性会下降。这个层不仅决定模型结构，还能为后续误差修正项提供核心参数，是 VECM 的理论核心。
模型估计与预测层
第四层是模型估计与预测层，负责根据协整结果建立 VECM 并生成未来预测值。VECM 一般由差分项、误差修正项和可选外生变量组成，系数估计完成后即可进行递推预测。其基本原理是：短期动态由差分滞后项解释，长期回归力量由误差修正项体现，二者共同决定未来变化方向。MATLAB 中可以通过系统化的矩阵运算完成参数估计，再结合样本末端状态进行多步预测。该层还应考虑滚动预测和递推更新，以增强模型对新数据的适应性。其价值在于把统计结构转化为可执行的未来轨迹，而不是仅输出静态参数。
结果评估与可视化层
第五层是结果评估与可视化层，负责检验模型预测效果并将结果直观呈现。常见评估方式包括均方误差、平均绝对误差、相关系数、残差自相关检查以及真实值与预测值对比图。其基本原理是通过误差统计量衡量模型预测精度，通过图形分析检验是否存在系统偏差或周期偏离。MATLAB 在图形绘制方面功能完整，适合展示多变量预测曲线、残差图和对比图。若结合分类颜色、图例和坐标轴标注，还能帮助识别不同变量的动态特征。该层的意义在于把模型结果转化为直观信息，便于交流、汇报与业务解释。
项目模型描述及代码示例
数据构造与读取
rng(2025); % 固定随机种子，保证模拟数据与测试结果可复现
n = 240; % 设置样本长度，适合展示多变量时间序列建模流程
t = (1:n)'; % 构造时间索引列向量，便于后续绘图与对齐
u = cumsum(randn(n,1)*0.6); % 构造带随机趋势的共同随机游走成分，模拟长期冲击
x1 = 20 + u + 0.8*randn(n,1); % 构造第一个非平稳变量，叠加共同趋势与局部扰动
x2 = 35 + 1.2*u + 0.9*randn(n,1); % 构造第二个非平稳变量，使其与第一个变量具备联动关系
x3 = 15 + 0.7*x1 – 0.3*x2 + 0.8*randn(n,1); % 构造第三个变量，引入长期线性组合关系
data = [t x1 x2 x3]; % 合并为原始数据矩阵，第一列为时间，其余列为观测变量
T = array2table(data, 'VariableNames', {'Time','Y1','Y2','Y3'}); % 转换为表格，便于变量命名与管理
writetable(T, 'vecm_demo_data.csv'); % 输出为CSV文件，模拟真实项目中的外部数据源
D = readtable('vecm_demo_data.csv'); % 从文件重新读取数据，验证导入流程是否正常
Y = [D.Y1, D.Y2, D.Y3]; % 提取多变量观测矩阵，作为后续统计分析输入
数据平稳性检验
[y1_h,y1_p] = adftest(Y(:,1)); % 对第一个变量进行ADF单位根检验，判断是否存在非平稳性
[y2_h,y2_p] = adftest(Y(:,2)); % 对第二个变量进行ADF单位根检验，确认时间序列性质
[y3_h,y3_p] = adftest(Y(:,3)); % 对第三个变量进行ADF单位根检验，检查是否需要差分
dY = diff(Y); % 对原始矩阵进行一阶差分，提取短期变化信息
[d1_h,d1_p] = adftest(dY(:,1)); % 对差分后的第一个变量再次检验，确认是否转为平稳
[d2_h,d2_p] = adftest(dY(:,2)); % 对差分后的第二个变量检验，验证I(1)特征
[d3_h,d3_p] = adftest(dY(:,3)); % 对差分后的第三个变量检验，确保满足建模前提
fprintf('ADF原序列p值: %.4f %.4f %.4f\\n', y1_p, y2_p, y3_p); % 输出原序列检验结果，便于快速判断单整阶数
fprintf('ADF差分序列p值: %.4f %.4f %.4f\\n', d1_p, d2_p, d3_p); % 输出差分后检验结果，辅助确定是否适合VECM
协整检验与秩判断
[h,pValue,stat,cValue,mle] = jcitest(Y); % 使用Johansen协整检验识别长期均衡关系
rankInfo = mle; % 保存最大似然估计结果，用于读取协整秩与特征信息
disp(h); % 显示检验结论，判断是否存在协整关系
disp(pValue); % 显示对应显著性水平结果，辅助判定秩的有效性
disp(stat); % 显示检验统计量，便于与临界值比较
disp(cValue); % 显示临界值矩阵，用于协整秩选择
r = rankInfo.r; % 读取协整秩，作为VECM核心结构参数
k = rankInfo.k_ar; % 读取滞后阶数信息，为模型拟合提供基础
VECM 估计与参数提取
M = vecm(size(Y,2), r, k); % 创建VECM模型对象，指定变量维度、协整秩与滞后阶数
M.SeriesNames = {'Y1','Y2','Y3'}; % 设置序列名称，便于结果展示与调试
M.Constant = NaN; % 让常数项由模型自动处理，减少结构误配风险
EstM = estimate(M, Y); % 根据观测数据估计VECM参数
A = EstM.AR; % 提取短期动态系数矩阵，观察各滞后项影响
B = EstM.Beta; % 提取协整向量矩阵，查看长期均衡关系
alpha = EstM.SpeedOfAdjustment; % 提取误差修正速度矩阵，分析偏离均衡后的回归强度
disp(A); % 显示短期动态系数，检查系数符号与数量级
disp(B); % 显示协整向量，验证长期线性组合结构
disp(alpha); % 显示调整系数，判断各变量回归均衡的快慢
样本内拟合与多步预测
p = 12; % 设置预测步长，模拟未来一段时期的多步预测任务
Yfit = infer(EstM, Y); % 基于已估计模型进行样本内状态推断
Yf = forecast(EstM, p, Y); % 从样本末端向未来递推p步，生成预测结果
tPred = (n+1:n+p)'; % 构造预测时间索引，便于与真实或模拟未来值对齐
Yall = [Y; Yf]; % 拼接历史值与预测值，便于完整曲线展示
disp(Yf); % 输出未来预测矩阵，查看每个变量的预测路径
结果评估与可视化
figure('Color','w'); % 创建白底图窗，便于展示多变量时间序列结果
plot(t, Y(:,1), 'b-', 'LineWidth', 1.4); hold on; % 绘制第一个变量历史值，使用蓝色实线突出走势
plot(t, Y(:,2), 'r-', 'LineWidth', 1.4); % 绘制第二个变量历史值，使用红色实线区分变量
plot(t, Y(:,3), 'k-', 'LineWidth', 1.4); % 绘制第三个变量历史值，使用黑色实线显示联动变化
plot(tPred, Yf(:,1), 'b–', 'LineWidth', 1.8); % 绘制第一个变量未来预测值，使用虚线表示预测段
plot(tPred, Yf(:,2), 'r–', 'LineWidth', 1.8); % 绘制第二个变量未来预测值，强调未来趋势外推
plot(tPred, Yf(:,3), 'k–', 'LineWidth', 1.8); % 绘制第三个变量未来预测值，展示协整约束下的未来轨迹
legend({'Y1历史','Y2历史','Y3历史','Y1预测','Y2预测','Y3预测'}, 'Location','best'); % 添加图例，便于区分历史与预测区间
xlabel('时间'); % 标注横轴为时间，说明序列演化方向
ylabel('数值'); % 标注纵轴为变量数值，便于解释尺度变化
title('VECM多变量时间序列预测结果'); % 设置图题，概括图像内容
grid on; % 打开网格，增强读图精度
figure('Color','w'); % 创建第二个图窗，展示残差诊断信息
res = Y(2:end,:) – Yfit; % 计算样本内残差，反映模型拟合误差
plot(res(:,1), 'b-'); hold on; % 绘制第一个变量残差，检查波动是否围绕零轴随机分布
plot(res(:,2), 'r-'); % 绘制第二个变量残差，观察是否存在系统偏差
plot(res(:,3), 'k-'); % 绘制第三个变量残差，分析拟合稳定性
legend({'Y1残差','Y2残差','Y3残差'}, 'Location','best'); % 添加残差图例，便于诊断不同变量误差
xlabel('样本点'); % 标注横轴为样本序号，方便查看误差随时间变化
ylabel('残差'); % 标注纵轴为残差值，检查是否存在自相关或异方差
title('VECM样本内残差诊断'); % 设置残差图标题，突出诊断用途
grid on; % 显示网格，提升诊断可读性
数据构造与读取
rng(2025); % 固定随机种子，保证模拟数据与测试结果可复现
n = 240; % 设置样本长度，适合展示多变量时间序列建模流程
t = (1:n)'; % 构造时间索引列向量，便于后续绘图与对齐
u = cumsum(randn(n,1)*0.6); % 构造带随机趋势的共同随机游走成分，模拟长期冲击
x1 = 20 + u + 0.8*randn(n,1); % 构造第一个非平稳变量，叠加共同趋势与局部扰动
x2 = 35 + 1.2*u + 0.9*randn(n,1); % 构造第二个非平稳变量，使其与第一个变量具备联动关系
x3 = 15 + 0.7*x1 – 0.3*x2 + 0.8*randn(n,1); % 构造第三个变量，引入长期线性组合关系
data = [t x1 x2 x3]; % 合并为原始数据矩阵，第一列为时间，其余列为观测变量
T = array2table(data, 'VariableNames', {'Time','Y1','Y2','Y3'}); % 转换为表格，便于变量命名与管理
writetable(T, 'vecm_demo_data.csv'); % 输出为CSV文件，模拟真实项目中的外部数据源
D = readtable('vecm_demo_data.csv'); % 从文件重新读取数据，验证导入流程是否正常
Y = [D.Y1, D.Y2, D.Y3]; % 提取多变量观测矩阵，作为后续统计分析输入
数据平稳性检验
[y1_h,y1_p] = adftest(Y(:,1)); % 对第一个变量进行ADF单位根检验，判断是否存在非平稳性
[y2_h,y2_p] = adftest(Y(:,2)); % 对第二个变量进行ADF单位根检验，确认时间序列性质
[y3_h,y3_p] = adftest(Y(:,3)); % 对第三个变量进行ADF单位根检验，检查是否需要差分
dY = diff(Y); % 对原始矩阵进行一阶差分，提取短期变化信息
[d1_h,d1_p] = adftest(dY(:,1)); % 对差分后的第一个变量再次检验，确认是否转为平稳
[d2_h,d2_p] = adftest(dY(:,2)); % 对差分后的第二个变量检验，验证I(1)特征
[d3_h,d3_p] = adftest(dY(:,3)); % 对差分后的第三个变量检验，确保满足建模前提
fprintf('ADF原序列p值: %.4f %.4f %.4f\\n', y1_p, y2_p, y3_p); % 输出原序列检验结果，便于快速判断单整阶数
fprintf('ADF差分序列p值: %.4f %.4f %.4f\\n', d1_p, d2_p, d3_p); % 输出差分后检验结果，辅助确定是否适合VECM
协整检验与秩判断
[h,pValue,stat,cValue,mle] = jcitest(Y); % 使用Johansen协整检验识别长期均衡关系
rankInfo = mle; % 保存最大似然估计结果，用于读取协整秩与特征信息
disp(h); % 显示检验结论，判断是否存在协整关系
disp(pValue); % 显示对应显著性水平结果，辅助判定秩的有效性
disp(stat); % 显示检验统计量，便于与临界值比较
disp(cValue); % 显示临界值矩阵，用于协整秩选择
r = rankInfo.r; % 读取协整秩，作为VECM核心结构参数
k = rankInfo.k_ar; % 读取滞后阶数信息，为模型拟合提供基础
VECM 估计与参数提取
M = vecm(size(Y,2), r, k); % 创建VECM模型对象，指定变量维度、协整秩与滞后阶数
M.SeriesNames = {'Y1','Y2','Y3'}; % 设置序列名称，便于结果展示与调试
M.Constant = NaN; % 让常数项由模型自动处理，减少结构误配风险
EstM = estimate(M, Y); % 根据观测数据估计VECM参数
A = EstM.AR; % 提取短期动态系数矩阵，观察各滞后项影响
B = EstM.Beta; % 提取协整向量矩阵，查看长期均衡关系
alpha = EstM.SpeedOfAdjustment; % 提取误差修正速度矩阵，分析偏离均衡后的回归强度
disp(A); % 显示短期动态系数，检查系数符号与数量级
disp(B); % 显示协整向量，验证长期线性组合结构
disp(alpha); % 显示调整系数，判断各变量回归均衡的快慢
样本内拟合与多步预测
p = 12; % 设置预测步长，模拟未来一段时期的多步预测任务
Yfit = infer(EstM, Y); % 基于已估计模型进行样本内状态推断
Yf = forecast(EstM, p, Y); % 从样本末端向未来递推p步，生成预测结果
tPred = (n+1:n+p)'; % 构造预测时间索引，便于与真实或模拟未来值对齐
Yall = [Y; Yf]; % 拼接历史值与预测值，便于完整曲线展示
disp(Yf); % 输出未来预测矩阵，查看每个变量的预测路径
结果评估与可视化
figure('Color','w'); % 创建白底图窗，便于展示多变量时间序列结果
plot(t, Y(:,1), 'b-', 'LineWidth', 1.4); hold on; % 绘制第一个变量历史值，使用蓝色实线突出走势
plot(t, Y(:,2), 'r-', 'LineWidth', 1.4); % 绘制第二个变量历史值，使用红色实线区分变量
plot(t, Y(:,3), 'k-', 'LineWidth', 1.4); % 绘制第三个变量历史值，使用黑色实线显示联动变化
plot(tPred, Yf(:,1), 'b–', 'LineWidth', 1.8); % 绘制第一个变量未来预测值，使用虚线表示预测段
plot(tPred, Yf(:,2), 'r–', 'LineWidth', 1.8); % 绘制第二个变量未来预测值，强调未来趋势外推
plot(tPred, Yf(:,3), 'k–', 'LineWidth', 1.8); % 绘制第三个变量未来预测值，展示协整约束下的未来轨迹
legend({'Y1历史','Y2历史','Y3历史','Y1预测','Y2预测','Y3预测'}, 'Location','best'); % 添加图例，便于区分历史与预测区间
xlabel('时间'); % 标注横轴为时间，说明序列演化方向
ylabel('数值'); % 标注纵轴为变量数值，便于解释尺度变化
title('VECM多变量时间序列预测结果'); % 设置图题，概括图像内容
grid on; % 打开网格，增强读图精度
figure('Color','w'); % 创建第二个图窗，展示残差诊断信息
res = Y(2:end,:) – Yfit; % 计算样本内残差，反映模型拟合误差
plot(res(:,1), 'b-'); hold on; % 绘制第一个变量残差，检查波动是否围绕零轴随机分布
plot(res(:,2), 'r-'); % 绘制第二个变量残差，观察是否存在系统偏差
plot(res(:,3), 'k-'); % 绘制第三个变量残差，分析拟合稳定性
legend({'Y1残差','Y2残差','Y3残差'}, 'Location','best'); % 添加残差图例，便于诊断不同变量误差
xlabel('样本点'); % 标注横轴为样本序号，方便查看误差随时间变化
ylabel('残差'); % 标注纵轴为残差值，检查是否存在自相关或异方差
title('VECM样本内残差诊断'); % 设置残差图标题，突出诊断用途
grid on; % 显示网格，提升诊断可读性
