Global Index (short | long) | Local contents | Local Index (short | long)
TH=ivx(z,nn,x,maxsize,Tsamp,p)
IVX Computes instrumental variable estimates for ARX-models. TH = ivx(Z,NN,X) TH: returned as the IV-estimate of the ARX-model A(q) y(t) = B(q) u(t-nk) + v(t) along with relevant structure information. See HELP THETA for the exact structure of TH. Z : the output-input data Z=[y u], with y and u as column vectors. For multi-input systems u=[u1 u2 ... un] NN: NN=[na nb nk] gives the orders and delays associated with the above model. X : is the vector of instrumental variables. This should be of the same size as the y-vector (i.e. the first column of Z). See IV4 for good, automatic choices of instruments. A multioutput variant ig given by IVXMV TH=ivx(Z,NN,X,maxsize,T) allows access to some parameters associated with the algorithm. See HELP AUXVAR for an explanation of these.
This function calls | This function is called by |
---|---|
function TH=ivx(z,nn,x,maxsize,Tsamp,p) % L. Ljung 10-1-86, 4-12-87 % Copyright (c) 1986-98 by The MathWorks, Inc. % $Revision: 2.3 $ $Date: 1997/12/02 03:42:56 $ [Ncap,nz]=size(z); nu=nz-1;[Nc,nx]=size(x);[nnr,nnc]=size(nn); maxsdef=idmsize(Ncap); % % *** Some initial tests on the input arguments *** % if nargin<6, p=1;end if nargin<5, Tsamp=1;end if nargin<4, maxsize=maxsdef;end if maxsize<0,maxsize=maxsdef;end if Tsamp<0,Tsamp=1;end if p<0,p=1;end if isempty(Tsamp),Tsamp=1;end,if isempty(maxsize),maxsize=maxsdef;end if nnr>1,eval('TH=ivxmv(z,nn,x,maxsize,Tsamp,p);'),return,end if Ncap~=Nc | nx~=1,error('The x-vector should be a column vector with the same number of elements as z!'),return,end if length(nn)~=1+2*nu,disp('Incorrect number of orders specified!'), disp('nn should be nn=[na nb nk]'), disp('where nb and nk are row vectors of length equal to the number of inputs'),error('see above'),return,end na=nn(1);nb=nn(2:1+nu);nk=nn(2+nu:1+2*nu);n=na+sum(nb); % % construct regression matrix % nmax=max([na+1 nb+nk])-1; M=floor(maxsize/n); Rxx=zeros(na);Ruu=zeros(sum(nb));Rxu=zeros(na,sum(nb));Rxy=zeros(na); Ruy=zeros(sum(nb),na); F=zeros(n,1); for k=nmax:M:Ncap-1 jj=(k+1:min(Ncap,k+M)); phix=zeros(length(jj),na); phiy=phix; phiu=zeros(length(jj),sum(nb)); for kl=1:na, phiy(:,kl)=-z(jj-kl,1); phix(:,kl)=-x(jj-kl); end ss=0; for ku=1:nu for kl=1:nb(ku), phiu(:,ss+kl)=z(jj-kl-nk(ku)+1,ku+1);end ss=ss+nb(ku); end Rxy=Rxy+phix'*phiy; if nu>0,Ruy=Ruy+phiu'*phiy; Rxu=Rxu+phix'*phiu; Ruu=Ruu+phiu'*phiu; end Rxx=Rxx+phix'*phix; if na>0, F(1:na)=F(1:na)+phix'*z(jj,1);end F(na+1:n)=F(na+1:n)+phiu'*z(jj,1); end clear phiu, clear phix, clear phiy, % % compute estimate % if nu==0,TH=Rxy\F;end if nu>0,TH=[Rxy Rxu;Ruy Ruu]\F;end if p==0, return,end % % proceed to build up THETA-matrix % t=TH;, clear TH; I=[0 Tsamp nu na nb 0 0 zeros(1,nu) nk]; n=na+sum(nb); TH=zeros(n+3,max(length(I),n)); TH(1,1:length(I))=I; ti=fix(clock); ti(1)=ti(1)/100; TH(2,2:6)=ti(1:5); TH(2,7)=11; TH(3,1:n)=t.'; if p==2 e=pe(z,TH); TH(1,1)=e'*e/(length(e)-max(na,sum(nb))); TH(2,1)=TH(1,1)*(1+n/Ncap)/(1-n/Ncap); TH(4:3+n,1:n)=TH(1,1)*inv([Rxx Rxu;Rxu.' Ruu]);end