0

假设一个无噪声的 AR(1) 过程y(t)= a*y(t-1)。我有以下概念性问题，很高兴得到澄清。

Q1 - 数学公式和实现之间的差异 - AR 模型的数学公式是y(t) = - summmation over i=1 to p[a*y(t-p)] + eta(t)其中 p=模型阶数的形式，并且eta(t)是高斯白噪声。但是，当使用任何方法（如arburg()最小二乘法）估计系数时，我们只需调用该函数。我不知道是否隐式添加了白高斯噪声。然后，当我们用估计的系数求解 AR 方程时，我看到没有考虑负号，也没有添加噪声项。

AR 模型的正确表示是什么？当我只有一个包含 1000 个数据点的样本时，如何找到 k 次试验的平均系数？

Q2 - 如何为 k 次试验模拟fitted_data 中的编码问题，然后找到残差 - 我拟合了从未知系统生成的数据“数据”并通过以下方式获得系数

load('data.txt');

for trials = 1:10

    model = ar(data,1,'ls');
    original_data=data;

    fitted_data(i)=coeff1*data(i-1); %  **OR**
    data(i)=coeff1*data(i-1); 

    fitted_data=data;

    residual= original_data - fitted_data;
    plot(original_data,'r'); hold on; plot(fitted_data);

end

在计算残差时，是否通过用获得的系数解析 AR 方程获得了上述 fit_data？Matlab 有这样做的功能，但我想自己做。那么，在从原始数据中找到系数后，我该如何解决？上面的编码不正确。附上原始数据和fitted_data的图。原始数据与拟合数据的图

4

2 回答 2

2

如果您的模型只是y(n)= a*y(n-1)使用 scalar a，那么这就是解决方案。

y = randn(10, 1);
a = y(1 : end - 1) \ y(2 : end);

y_estim = y * a;
residual = y - y_estim;

当然，您应该将数据分成训练测试，并应用于a测试数据。您可以将此方法推广到y(n)= a*y(n-1) + b*y(n-2)等。

注意\表示mldivide()函数：mldivide

编辑： 

% model: y[n] = c + a*y(n-1) + b*y(n-2) +...+z*y(n-n_order)
n_order = 3;
allow_offset = true; % alows c in the model

% train
y_train = randn(20,1); % from your data
[y_in, y_out] = shifted_input(y_train, n_order, allow_offset);
a = y_in \ y_out;

% now test
y_test = randn(20,1); % from your data
[y_in, y_out] = shifted_input(y_test, n_order, allow_offset);
y_estim = y_in * a; % same a
residual = y_out - y_estim;

这里是shifted_input()：

function [y_in, y_out] = shifted_input(y, n_order, allow_offset)
y_out = y(n_order + 1 : end);
n_rows = size(y, 1) - n_order;
y_in = nan(n_rows, n_order);
for k = 1 : n_order    
    y_in(:, k) = y(1 : n_rows);
    y = circshift(y, -1);
end
if allow_offset
    y_in = [y_in, ones(n_rows, 1)];
end
return
于 2013-07-30T16:26:14.170 回答
1

AR 类型的模型可以用于多种用途，包括线性预测、线性预测编码、过滤噪声。eta(t) 不是我们有兴趣保留的东西，而是算法的一部分是通过在数据中寻找持久模式来尽可能地消除它们的影响。

我有教科书，在线性预测的背景下，不包括在总和之前的表达式中包含的负号。另一方面，Matlab 的功能lpc是：

Xp(n) = -A(2)*X(n-1) - A(3)*X(n-2) - ... - A(N+1)*X(n-N)

lpc如果您还没有，我建议您查看函数，并查看文档中的示例，例如：

randn('state',0);
noise = randn(50000,1);  % Normalized white Gaussian noise
x = filter(1,[1 1/2 1/3 1/4],noise);
x = x(45904:50000);
% Compute the predictor coefficients, estimated signal, prediction error, and autocorrelation sequence of the prediction error: 
p = lpc(x,3);
est_x = filter([0 -p(2:end)],1,x);    % Estimated signal
e = x - est_x;                        % Prediction error
[acs,lags] = xcorr(e,'coeff');        % ACS of prediction error

估计x值计算为est_x。请注意示例如何使用filter. 再次引用 matlab 文档，filter(b,a,x)“是标准差分方程的“直接形式 II 转置”实现：

a(1)*y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb)
                      - a(2)*y(n-1) - ... - a(na+1)*y(n-na)

这意味着在前面的示例est_x(n)中计算为

  est_x(n) = -p(2)*x(n-1) -p(3)*x(n-2) -p(4)*x(n-3)

这就是你所期望的！

编辑：

关于函数ar，matlab文档解释了输出系数与上面讨论的 lp 场景中的含义相同。

评估 AR 模型输出的正确方法是计算

data_armod(i)= -coeff(2)*data(i-1) -coeff(3)*data(i-2) -coeff(4)*data(i-3)

其中 coeff 是返回的系数矩阵

 model = ar(data,3,'ls');
 coeff = model.a;
于 2013-07-30T15:10:20.193 回答