offline ( Gauss Newtone approach ) Prediction error method for ARMAX Model, Algorithm Implemented but not good results

2 次查看（过去 30 天）
Reda Chtiba
Reda Chtiba 2022-5-12
回答： Aman 2024-2-7
% Simulation of an ARMAX model y(t)-0.7y(t-1)=0.5u(t-1)+e(t)+0.6e(t)
np=251;
a0=[1 -0.7 ]; b0=[0 0.5 ]; c0=[1 0.6]; % Set parameter values
u=randn(np,1)*sqrt(1); % Generate white noise with variance =1
e=randn(np,1)*sqrt(1); % Generate white noise with variance = 1
y=filter(b0,a0,u)+filter(c0,a0,e); % generate output data : y=(B/A)*u+(C/A)*e
z=[y u]; % store our inpute output data
x=rpem(z,[1 1 1 0 0 1],'ff',1); % build-in predction error method
x(end,:);
%% %%
theta=[-0.6 ;0.4 ; 0.6];% initiliazing
grad=[]; % dummy
iter_max=10;
for iter=1:iter_max
theta1=theta(1);
theta2=theta(2);
theta3=theta(3);
pe=filter([1 theta1],[1 theta3],y)-filter([0 theta2],[1 theta3],u);% Predictive errors, for the armax is : pe=(A/C)y(t)-(B/C)u(t)
for i=2:251 % Gradient for armax is given by grad(t)=(1/C)(-y(t-1),u(t-1),pe(t-1)) .
grad(1,:)=[0 0 0];% dummy value at time 1, to be deleted later
grad(i,:)=filter(1,[1 theta3],[-y(i-1) u(i-1) pe(i-1)]);% we started at i=2 so we can calculate using past value for other data
grad(1,:)=[];% delete first dummy row , we get 250 data
end
pe(1,:)=[];% delete first row , to get 250 data so we can compare it with grad
P=(grad'*grad)\(grad'*pe); % Gauss-Newtone Formula
theta=theta+0.8*P;%GaussNewtone Iterations
end

请先登录，再回答此问题。

回答（1 个）

Aman
Aman 2024-2-7
Hi Reda,
As per my understanding, you are using the Gauss-Newton method of optimizing the prediction error, but the result is not converging.
The matrix that has been created is near singular, so in that case the delta is not getting computed properly. In this case, you need to add a damping factor so that the proper delta is calculated. Refer to the below code, where I have added a damping factor of "1e-3".
% Simulation of an ARMAX model y(t)-0.7y(t-1)=0.5u(t-1)+e(t)+0.6e(t)
np=251;
a0=[1 -0.7 ];
b0=[0 0.5 ];
c0=[1 0.6]; % Set parameter values
u=randn(np,1)*sqrt(1); % Generate white noise with variance = 1
e=randn(np,1)*sqrt(1); % Generate white noise with variance = 1
y=filter(b0,a0,u)+filter(c0,a0,e); % generate output data : y=(B/A)*u+(C/A)*e
z=[y u]; % store our input output data
%%
y_prev = [0; y(1:end-1)]; % y(t-1)
u_prev = [0; u(1:end-1)]; % u(t-1)
e_prev = zeros(np, 1);
theta=[-0.6 ;0.4 ; 0.6];% initiliazing
iter_max=10;
error = [];
for iter=1:iter_max
residual = y - (theta(1)*y_prev + theta(2)*u_prev + theta(3)*e_prev);
J = [-y_prev, u_prev, e_prev];
lambda = 1e-3; % Damping factor
delta = ((J'*J + lambda*eye(size(J,2))) \ J') * residual;
theta = theta + delta;
e_prev = [0; residual(1:end-1)];
end
I hope it helps!

类别

在 Help Center 和 File Exchange 中查找有关 Time Series Analysis 的更多信息

标签

产品


版本

R2022a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by