Documentation of ivstruc

Global Index (short | long) | Local contents | Local Index (short | long)

Function Synopsis

V=ivstruc(z,zv,nn,p,maxsize)

Help text

IVSTRUC Computes the output error fit for families of single-output ARX-models.
   V = IVSTRUC(ZE,ZV,NN)

   V: The first row of V is returned as the output error fit for
   models of the structures, defined by the rows of NN.
   These models are computed by the IV-method. The next rows
   of V are returned as NN'. The last row of V contains the logarithms
   of the conditioning numbers of the IV-matrix for the structures in
   question. The last column of V contains the number of data points in Z

   ZE : the output-input data Z=[y u], with y and u as column vectors.This
   data is used for estimation. For multi-input systems u=[u1 u2 ... un]
   ZV : output-input data used to evaluate the models (could equal ZE).
   NN: Each row of NN is of the format [na nb nk], the orders and
   delays associated with the corresponding model. (See also IV4 and
   STRUC. See HELP SELSTRUC for analysis of V.

   With V = IVSTRUC(Z,NN,p,maxsize) the calculation of conditioning numbers
   can be suppressed (p=0), and the speed/memory trade-off affected.
   See also AUXVAR and ARXSTRUC.

Cross-Reference Information

Listing of function ivstruc

function V=ivstruc(z,zv,nn,p,maxsize)

%   L. Ljung 4-12-87,12-27-91
%   Copyright (c) 1986-98 by The MathWorks, Inc.
%   $Revision: 2.3 $  $Date: 1997/12/02 03:43:01 $

if nargin < 3
   disp('Usage: V = IVSTRUC(EST_DATA,VAL_DATA,ORDERS)')
   return
end

[Ncap,nz]=size(z); nu=nz-1;[nm,nl]=size(nn);
na=nn(:,1);nb=nn(:,2:1+nu);nk=nn(:,2+nu:1+2*nu);
[snbr,snbc]=size(nb);if any(any(nb==zeros(snbr,snbc))'),
error('nb=0 is not feasible for ivstruc!'),end % Mod 911227
nma=max(na);nbkm=max(nb+nk)-ones(1,nu);nkm=min(nk);
n=nma+sum((nbkm-nkm))+nu;nbm=max(nb);
%
% *** Some initial tests on the input arguments ***
%
maxsdef=idmsize(Ncap,n);

if nargin<5,maxsize=maxsdef;,end
if nargin<4, p=1;,end

% *** construct instruments (see (7.111)-(7.112)) ***
tha=arx(z,[nma nbm nkm],maxsize,1,0);
NF=fstab([1 tha(1:nma)']);

     x=zeros(Ncap,1);rs=nma;
     for k=1:nu
         MF=[zeros(1,nkm(k)) tha(rs+1:rs+nbm(k))'];rs=rs+nbm(k);
         x=x+filter(MF,NF,z(:,k+1));
     end
%
% construct regression matrix
%
nmax=max(max([na+ones(nm,1) nb+nk]'))-1;
M=floor(maxsize/n);nnuu=sum(nbkm)-sum(nkm)+nu;
Rxx=zeros(nma);Ruu=zeros(nnuu);Rxu=zeros(nma,nnuu);Rxy=zeros(nma);
Ruy=zeros(nnuu,nma); F=zeros(n,1);
for k=nmax:M:Ncap
     jj=(k+1:min(Ncap,k+M));
     phix=zeros(length(jj),nma); phiy=phix; phiu=zeros(length(jj),nnuu);
     for kl=1:nma, phiy(:,kl)=-z(jj-kl,1); phix(:,kl)=-x(jj-kl); end
     ss=0;
     for ku=1:nu
            for kl=nkm(ku):nbkm(ku), phiu(:,ss+kl+1-nkm(ku))=z(jj-kl,ku+1);end
             ss=ss+nbkm(ku)-nkm(ku)+1;
     end
     Rxy=Rxy+phix'*phiy;
     Ruy=Ruy+phiu'*phiy;
     Rxu=Rxu+phix'*phiu;
     Ruu=Ruu+phiu'*phiu;
     Rxx=Rxx+phix'*phix;
     if nma>0, F(1:nma)=F(1:nma)+phix'*z(jj,1);end
      F(nma+1:n)=F(nma+1:n)+phiu'*z(jj,1);
end
clear phiu, clear phix, clear phiy,
%
% compute estimate

for j=1:nm
        s=[];rs=0;
        for ku=1:nu
                s=[s,rs+nk(j,ku)-nkm(ku)+1:rs+nb(j,ku)+nk(j,ku)-nkm(ku)];
                rs=rs+nbkm(ku)-nkm(ku)+1;
        end
        RR=[Rxy(1:na(j),1:na(j)) Rxu(1:na(j),s);Ruy(s,1:na(j)) Ruu(s,s)];
        TH=RR\F([1:na(j) nma+s]);
        a=TH(1:na(j))';b=TH(na(j)+1:length(TH))';f=fstab([1 a]);
        [mzv,nzv] = size(zv(:,1));
        yhat=zeros(mzv,nzv);rs=0;
        for ku=1:nu
        bdum=[zeros(1,nk(j,ku)) b(rs+1:rs+nb(j,ku))]; rs=rs+nb(j,ku);
        yhat=yhat+filter(bdum,f,zv(:,ku));
        end
        V(1,j)=(zv(:,1)-yhat)'*(zv(:,1)-yhat)/(Ncap-nmax);
        V(2:nl+1,j)=nn(j,:)';
        if p~=0, V(nl+2,j)=log(cond(RR));,end
end
V(1,nm+1)=Ncap-nmax;