cf.m

对平稳时间序列中的参数估计方法（包括矩估计

% ARMA(p,q)模型差分方程预测方法
clear;
p=input('输入自回归阶数：');
q=input('输入滑动平均阶数：');
N=input('输入已知数据的个数：');
l=input('输入做几期预测：');
if p>0
for t=1:p
    a(t)=input('输入自回归参数：');
end
else
    a=0;
end
if q>0
for t=1:q
    b(t)=input('输入滑动平均参数：');
end
else 
    b=0;
end
for t=1:N
    x(t)=input('输入数据：');
end
sum1=0;
sum2=0;
res=zeros(N,1);        % 对残差赋初值为0
for t=p+1:N            % 对残差进行迭代计算
    for k=1:p
        sum1=sum1+a(k)*x(t-k);
    end
    for r=1:q
        sum2=sum2+b(r)*x(t-r);
    end
    res(t)=x(t)-sum1+sum2;
    sum1=0;
    sum2=0;
end
sum3=0;
sum4=0;
%  进行递推预测
for h=1:l
    for j=1:p
        sum3=sum3+a(j)*x(N-j+h);
    end
if h<q
    for j=h:q
        sum4=sum4+b(j)*res(N-j+h);
    end
else
    sum4=0;
end
    x(N+h)=sum3-sum4;               % 这里的x(N+h)表示预测的值
    sum3=0;
    sum4=0;
end
%  显示预测的值
plot(x);
hold on;
for j=1:l
    disp(x(N+j));
    plot(N+j,x(N+j),'or');
end 
grid on;

% 得到新值进行修正预报
new=input('输入新值:');
max=input('输入要计算权的个数:');
G=Green(a,b,p,q,max);   % 计算G权
y=x;
CC=new-x(N+1);
for j=1:l-1
    y(N+1+j)=x(N+1+j)+G(j)*CC          % 得到的y的值就是修正后得到的值
    hold on;
    plot(N+1+j,y(N+1+j),'+g');
end