n=input('请输入阶次n=:');
u=input('请输入输入矩阵u=');
y=input('请输入输出矩阵y=');
c=size(u);
%求x(N)
r=[];
for i=1:c(2)-n+1  %行
    for j=1:n     %列
        z(j)=y(n+i-j);
        U(j)=u(n+i-j);
    end
    s=[-z U];
    r(i,:)=s;   
end
X=r';
%求P(N),K(N),cs(N)
cs=zeros(1,2*n)';%初始化参数a,b都为零
I=eye(2*n,2*n);
p=(1.0e+8)*I; %初始化p
for i=1:c(2)-n
    k=p*X(:,i)/(1+(X(:,i))'*p*X(:,i)); %求出k(i)
    p=p-p*X(:,i)/(1+(X(:,i))'*p*X(:,i))*(X(:,i))'*p; %求出p(i)
    cs=cs+k*(y(n+i)-(X(:,i))'*cs); %求出cs(i)
end
a=cs(1:n,:)';
b=cs(n+1:2*n,:)';
a  
b