在命令行中辨识时间序列模型

生成模型并仿真模型输出

通过创建并仿真形式为 $y_k=a_1{y}_{k-1}+a_2{y}_{k-2}+e_k$ 的自回归 (AR) 多项式模型 ts_orig 来生成时间序列数据，其中 $e_k$ 表示高斯随机噪声。该噪声代表了模型中一个未测量的输入。由于该模型是时间序列模型，因此不存在观测输入。

在计算 $e_k$ 之前，请初始化随机数生成器的种子，以确保噪声值具有可重复性。

ts_orig = idpoly([1 -1.75 0.9]);
rng('default')
e = idinput(300,'rgs');

仿真该模型的观测输出 y_obs，并将 y_obs 转换为默认采样时间为一秒的 iddata 对象 y。将模型输出结果与输入噪声一起绘制成图。

y_obs = sim(ts_orig,e); 
y = iddata(y_obs);

plot(e)
hold on
plot(y_obs)
title('Input Noise and Original Model Output')
legend('RGS Noise','Model Output')
hold off

估计模型并比较频谱

函数 etfe 和 spa 提供了两种用于进行谱分析的非参数方法。将 etfe 和 spa 得出的估计功率谱与原始模型进行比较。

ts_etfe = etfe(y);
ts_spa = spa(y);
spectrum(ts_etfe,ts_spa,ts_orig);
legend('ts_{etfe}','ts_{spa}','ts_{orig}')

现在，使用 AR 结构估计一个参数模型。估计一个二阶 AR 模型，并将其频谱与原始模型和 spa 估计值进行比较。

ts_ar = ar(y,2);
spectrum(ts_spa,ts_ar,ts_orig);
legend('ts_{spa}', 'ts_{ar}', 'ts_{orig}')

与非参数模型相比，AR 模型的频谱与原始模型的频谱吻合得更为紧密。

协方差的估计与比较

通过将每个模型的输出与其自身进行卷积，计算原始模型和 AR 模型的协方差函数。

ir_orig = sim(ts_orig,[1;zeros(24,1)]);
Ry_orig = conv(ir_orig,ir_orig(25:-1:1));
ir_ar = sim(ts_ar,[1;zeros(24,1)]);
Ry_ar = conv(ir_ar,ir_ar(25:-1:1));

此外，利用 xcorr 直接从观测到的输出 y 估计协方差 Ry。

Ry = xcorr(y.y,24,'biased');

绘制原始协方差和估计协方差的绘图，并进行比较。

plot((-24:24)',[Ry_orig,Ry_ar,Ry]);
legend('Ry_{orig}', 'Ry_{ar}', 'Ry')

估计的 AR 模型的协方差 Ry_ar 更接近原始协方差 Ry_orig。

预测和比较模型输出结果

使用函数 compare，比较原始模型和 AR 模型的三步预测精度（或拟合百分比）。在此，compare 利用原始模型输出数据 y 计算 ts_orig 和 ts_ar 模型的预测响应，假设未测量的输入 $e_k$ 为零。第四个参量 3 表示预测的步数。

compare(y,ts_orig,ts_ar,3);

图例中的百分比是拟合百分比，代表拟合优度。即使是原始模型，其预测准确率也远未达到 100%，因为预测过程中未将未测量的模型输入 $e_k$ 考虑在内。估计的 AR 模型的拟合值与原始模型接近，这表明该 AR 模型是一个良好的估计。