AR自回归模型的完整构建流程涉及多个关键环节，从数据预处理到最终预测输出，每一步都对模型效果产生重要影响。本文将基于MATLAB平台，系统性地演示AR模型的建模全链路，并提供可直接运行的代码与详细注释，方便读者快速上手实践。

MATLAB AR模型建模与预测完整工作流程

%% 1. 数据准备：生成/加载并可视化数据
% 假设您有自己的时间序列数据 `myData`
% 这里以模拟的平稳时间序列为例
rng(0); % 设置随机种子确保可重复性
T = 200; % 数据长度
y = 0.7*randn(T,1); % 生成白噪声
for t = 3:T
    y(t) = 0.5*y(t-1) - 0.3*y(t-2) + y(t); % 生成AR(2)过程
end
figure;
plot(y, 'LineWidth', 1.5);
title('原始时间序列数据');
xlabel('时间点'); ylabel('数值');
grid on;

%% 2. 平稳性检验与处理
% 进行ADF单位根检验（需要Econometric Toolbox的adftest函数）
% 如果检验结果为0（非平稳），则需要进行差分
try
    [h, pValue] = adftest(y, 'Model', 'ARD');
    if ~h
        disp('数据非平稳，进行一阶差分...');
        y_diff = diff(y);
        y = y_diff; % 用差分后的数据继续分析
    else
        disp('数据通过平稳性检验。');
    end
catch
    disp('未找到adftest函数，请确保已安装Econometric Toolbox。');
    % 作为备选，可以可视化观察或手动差分
    % y = diff(y); 
end
% 去中心化：减去均值（对平稳序列很重要）
y = y - mean(y);

%% 3. 模型识别：确定AR阶数(p)
% 方法一：绘制自相关(ACF)和偏自相关(PACF)图
figure;
subplot(2,1,1);
autocorr(y); title('样本自相关函数 (ACF)');
subplot(2,1,2);
parcorr(y); title('样本偏自相关函数 (PACF)');
% PACF在滞后p阶后“截尾”，即为建议的阶数
% 方法二：使用AIC/BIC信息准则自动定阶（更客观）
maxP = 10; % 设定最大搜索阶数
[aic, bic] = aicbic_custom(y, maxP); % 需自定义aicbic_custom函数，下方提供
[~, bestP_aic] = min(aic);
[~, bestP_bic] = min(bic);
fprintf('AIC建议最佳阶数 p = %d', bestP_aic-1); % 阶数从0开始
fprintf('BIC建议最佳阶数 p = %d', bestP_bic-1);

%% 4. 模型估计与拟合
p = bestP_bic - 1; % 选择BIC建议的阶数（通常更简洁）
% 使用aryule函数（Yule-Walker方法）估计AR参数
[ar_coeff, noise_var] = aryule(y, p);
ar_coeff = -ar_coeff(2:end); % 提取真正的自回归系数（去掉首项1）
% 构建AR模型对象（如果使用System Identification Toolbox）
if exist('ar', 'file') == 2
    sys = ar(y, p, 'yw');
    compare(sys, y, 10); % 对比模型拟合效果
end

%% 5. 模型诊断：检验残差是否为白噪声
% 计算残差
y_pred = filter([0; ar_coeff], 1, y); % 单步预测
residuals = y(p+1:end) - y_pred(p+1:end); % 预测残差
% 残差的自相关检验
figure;
subplot(2,1,1);
autocorr(residuals);
title('残差序列的自相关图'); % 如果大部分自相关系数在置信区间内，则通过
subplot(2,1,2);
histogram(residuals, 20);
title('残差直方图'); xlabel('残差'); ylabel('频数');
% Ljung-Box Q检验（更严格的统计检验）
[h_lbq, p_lbq] = lbqtest(residuals, 'Lags', [10, 15]);
if h_lbq
    fprintf('警告：残差可能不是白噪声 (p-value = %.4f)。考虑增加模型阶数或检查模型类型。', p_lbq);
else
    fprintf('残差通过白噪声检验 (p-value = %.4f)。', p_lbq);
end

%% 6. 预测与可视化
nForecast = 20; % 预测未来20个点
[y_forecast, y_forecast_mse] = forecast_ar(y, ar_coeff, nForecast); % 需自定义forecast_ar函数
% 计算95%置信区间
ci_bound = 1.96 * sqrt(y_forecast_mse);
lower_bound = y_forecast - ci_bound;
upper_bound = y_forecast + ci_bound;
% 绘制结果
figure;
hold on;
plot(1:length(y), y, 'b-', 'LineWidth', 1.2, 'DisplayName', '历史数据');
forecast_start = length(y) + 1;
forecast_x = forecast_start:forecast_start+nForecast-1;
plot(forecast_x, y_forecast, 'r--', 'LineWidth', 1.5, 'DisplayName', '预测值');
fill([forecast_x, fliplr(forecast_x)], ...
    [lower_bound', fliplr(upper_bound')], ...
    'r', 'FaceAlpha', 0.2, 'EdgeColor', 'none', 'DisplayName', '95%置信区间');
xline(length(y), 'k--', 'HandleVisibility', 'off');
legend('show', 'Location', 'best');
title('AR模型样本外预测');
xlabel('时间点'); ylabel('数值');
grid on;
hold off;

%% 自定义函数1：计算AIC和BIC值
function [aic, bic] = aicbic_custom(y, maxP)
T = length(y);
aic = zeros(maxP+1, 1);
bic = zeros(maxP+1, 1);
for p = 0:maxP
    [coeff, noise_var] = aryule(y, p);
    % 计算对数似然值（假设高斯分布）
    logL = -T/2 * log(2*pi*noise_var) - (T-p)/(2);
    % 计算AIC和BIC
    aic(p+1) = -2*logL + 2*(p+1); % 参数个数为 p+1 (包括方差)
    bic(p+1) = -2*logL + (p+1)*log(T-p);
end
end

%% 自定义函数2：AR模型预测
function [y_forecast, forecast_mse] = forecast_ar(y, ar_coeff, nForecast)
p = length(ar_coeff);
y_extended = [y; zeros(nForecast, 1)]; % 扩展向量以存放预测值
forecast_mse = zeros(nForecast, 1);
noise_var = var(y - filter([0; ar_coeff], 1, y)); % 估计噪声方差
for t = length(y)+1 : length(y)+nForecast
    % 使用过去p个值进行预测
    past_values = y_extended(t-p:t-1);
    y_forecast_local = sum(ar_coeff .* past_values(end:-1:1)); % 注意系数顺序
    y_extended(t) = y_forecast_local;
    % 计算预测误差的均方误差（逐步增大）
    if t == length(y)+1
        forecast_mse(t-length(y)) = noise_var;
    else
        % 对于多步预测，误差会累积
        forecast_mse(t-length(y)) = forecast_mse(t-length(y)-1) + noise_var * sum(ar_coeff(1:t-length(y)-1).^2);
    end
end
y_forecast = y_extended(end-nForecast+1:end);
forecast_mse = forecast_mse(end-nForecast+1:end);
end

AR模型定阶选择与优化关键要点

平稳性检验是AR模型建模的前提条件，也是确保模型有效性的基础约束。若原始数据存在趋势项，可使用detrend函数去除趋势或直接进行差分处理；对于含有季节性成分的时间序列，季节差分是行业标准的处理方式。在模型定阶方面，偏自相关函数（PACF）图的"截尾"特征可作为直观判断依据，但实际工程中更推荐采用AIC（赤池信息准则）或BIC（贝叶斯信息准则）进行客观筛选——其中BIC准则对模型复杂度惩罚更为严格，倾向于选择更加简洁的模型结构。

模型诊断的核心在于残差分析。若残差序列仍存在显著自相关，表明模型尚未充分提取数据中的动态规律，此时可考虑增加模型阶数，或改用ARMA（自回归移动平均模型）进行建模。对于更复杂的应用场景，可以扩展为多变量VAR（向量自回归）模型、引入外生变量的ARX模型，甚至尝试非线性NAR（非线性自回归）神经网络——这些进阶方法适用于不同类型的复杂时间序列预测任务。

AR模型常见问题与解决方案