Home > documentation > Smats.m

Smats

PURPOSE ^

SMATS.m

SYNOPSIS ^

function [S0, S1, S2, S3, Sy, Sf, Sz, unique] = Smats(C, GAM0, GAM1, PSI, PPI);

DESCRIPTION ^

 SMATS.m

[S0, S1, S2, S3, Sy, Sf, Sz, unique] = Smats(C, GAM0, GAM1, PSI, PPI);

 System given as
        GAM0*y(t)=GAM1*y(t-1)+C+PSI*e(t)+PPI*eta(t),
 with e(t) an exogenous variable process and eta(t) being endogenously determined
 one-step-ahead expectational errors.  Returned system is
       y(t)=S0 + S1*y(t-1)+S2*e(t)+S3*sunspot(t)+ Sy sum(j=1..inf, Sf^(j-1) Sz e(t+j))  .

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 % SMATS.m
0002 %
0003 %[S0, S1, S2, S3, Sy, Sf, Sz, unique] = Smats(C, GAM0, GAM1, PSI, PPI);
0004 %
0005 % System given as
0006 %        GAM0*y(t)=GAM1*y(t-1)+C+PSI*e(t)+PPI*eta(t),
0007 % with e(t) an exogenous variable process and eta(t) being endogenously determined
0008 % one-step-ahead expectational errors.  Returned system is
0009 %       y(t)=S0 + S1*y(t-1)+S2*e(t)+S3*sunspot(t)+ Sy sum(j=1..inf, Sf^(j-1) Sz e(t+j))  .
0010 
0011 % If div is omitted from argument list, a div>1 is calculated.
0012 % eu(1)=1 for existence, eu(2)=1 for uniqueness.  eu(1)=-1 for
0013 % existence only with not-s.c. z; eu=[-2,-2] for coincident zeros.
0014 
0015 % Date created:
0016 % Modified:  1-DEC-2009
0017 % [add notes here]
0018 %-------------------------------------------------------------------------%
0019 
0020 function [S0, S1, S2, S3, Sy, Sf, Sz, unique] = Smats(C, GAM0, GAM1, PSI, PPI);
0021 
0022 
0023 %% Start Calculations
0024 
0025 eu=[0;0];
0026 
0027 realsmall=1e-6;
0028 
0029 [n k] = size(PPI);
0030 [n l] = size(PSI);
0031 
0032 [A B Q Z VV]=qz(GAM0,GAM1);
0033 
0034 div = 1.01;
0035 nunstab=0;
0036 zxz=0;
0037 
0038 for i=1:n % Determines the dimension of the unstable block of the model. The size of W_1(t)is n-m
0039        if abs(A(i,i)) > 0
0040           divhat=abs(B(i,i))/abs(A(i,i));
0041           if 1+realsmall<divhat && divhat<=div
0042              div=.5*(1+divhat);
0043           end
0044        end
0045    nunstab = nunstab + (abs(B(i,i))>div*abs(A(i,i)));
0046    if abs(A(i,i))<realsmall && abs(B(i,i))<realsmall % Checking coincident zeros
0047       zxz=1;
0048    end
0049 end
0050 
0051 div ;
0052 nunstab;
0053 m = nunstab;
0054 
0055 if ~zxz % If no coincident zeros do QZDIV.m
0056    [A B Q Z]= qzdiv(div,A,B,Q,Z);
0057 end
0058 % vin = 0;
0059 
0060 
0061 %%
0062 gev=[diag(A) diag(B)]; % The ratios diag(B)/diag(A) are the generalized eigenvalues of GAM1 and GAM0
0063 
0064 if zxz
0065    disp('Coincident zeros.  Indeterminacy and/or nonexistence.')
0066    eu=[-2;-2];
0067    % correction added 7/29/2003.  Otherwise the failure to set output
0068    % arguments leads to an error message and no output (including eu).
0069    %disp('Generalized eigenvalues')
0070    %disp(sort(abs(gev(:,2) ./ gev(:,1))));
0071    %disp('Dimension of eta is:')
0072    %disp(k)
0073    S0 = []; S1=[]; S2=[]; S3=[]; unique=[];
0074    return
0075 end
0076 %%
0077 
0078 Q1=Q(1:n-nunstab,:);
0079 Q2=Q(n-nunstab+1:n,:);
0080 Z1 = Z(:,1:n-nunstab)';
0081 Z2 = Z(:,n-nunstab+1:n)';
0082 A22 = A(n-nunstab+1:n,n-nunstab+1:n); % Lambda22 in LS (2003)
0083 B22 = B(n-nunstab+1:n,n-nunstab+1:n); % Omega22 in LS (2003)
0084 etawt=Q2*PPI;
0085 
0086 % Input Matrices for Indeterminate Solution
0087 
0088 A11 = A(1:n-nunstab,1:n-nunstab);
0089 A12 = A(1:n-nunstab,n-nunstab+1:n);
0090 
0091 
0092 % SVD of ETAWT
0093 [U,D,V] = svd(etawt);
0094 md = min(size(D));
0095 bigev = find(diag(D(1:md,1:md))>realsmall); % Partitions the SVD according to the number of non-zero singular values (LS, JEDC - Eq9, p277)
0096 U1 = U(:,bigev); %LS: U.1 ;
0097 V1 = V(:,bigev); %LS: V.1
0098 D11 = D(bigev,bigev); %LS: D11
0099 r = size(D11,1);
0100 V2 = V(:,r+1:end);  
0101 U2 = U(:,r+1:end);  
0102 
0103 eu(1) = length(bigev)>=nunstab;
0104 %%
0105 [Uhat,Dhat,Vhat] = svd(Q1*PPI);
0106 md = min(size(Dhat));
0107 bigev = find(diag(Dhat(1:md,1:md))>realsmall);
0108 Uhat1 = Uhat(:,bigev);
0109 Vhat1 = Vhat(:,bigev);
0110 Dhat11 = Dhat(bigev,bigev);
0111 rr = size(Dhat11,1);
0112 Uhat2 = Uhat(:,rr+1:end); 
0113 Vhat2 = Vhat(:,rr+1:end);
0114 %%
0115 
0116 if eu(1)
0117    disp('Solution Exists')
0118 else
0119    disp('Solution Does Not Exist')
0120    eu
0121    disp('Generalized eigenvalues')
0122    %disp(sort(abs(gev(:,2) ./ gev(:,1))));
0123    %disp('Dimension of eta is:')
0124    %disp(k)
0125    S0 = []; S1=[]; S2=[]; S3=[]; unique=0; Sy=[]; Sf=[]; Sz=[];
0126    return
0127 end
0128 
0129 if  isempty(V2)
0130     unique=1;
0131     disp('The solution is unique');
0132         
0133 else
0134     unique=0;
0135     disp('The solution is indeterminate');
0136 end
0137 
0138 if unique
0139    %disp('solution unique');
0140    eu(2)=1;
0141 else
0142    fprintf(1,'Indeterminacy.  %d loose endog errors.\n',k-r);
0143    disp('Generalized eigenvalues')
0144    disp(gev);
0145    disp(sort(abs(gev(:,2) ./ gev(:,1))));
0146 
0147 end
0148 
0149 %% Check for Uniqueness
0150 if  isempty(Vhat1)
0151     unique=1;
0152 else
0153     unique=norm(Vhat1-V1*V1'*Vhat1)<realsmall*n;
0154 end
0155 
0156 if unique
0157    %disp('Solution Unique');
0158    %disp('Generalized eigenvalues')
0159    %disp(sort(abs(gev(:,2) ./ gev(:,1))));
0160    eu(2)=1;
0161 else
0162    fprintf(1,'Indeterminacy.  %d loose endog errors.\n',size(Vhat1,2)-size(V1,2));
0163    %disp('Generalized eigenvalues')
0164    %disp(sort(abs(gev(:,2) ./ gev(:,1))));
0165    %disp(gev);
0166    %md=abs(diag(a))>realsmall;
0167    %ev=diag(md.*diag(a)+(1-md).*diag(b))\ev;
0168    %disp(ev)
0169 %   return;
0170 end
0171 
0172 %%
0173 PHI = (U1*(D11\V1')*Vhat1*Dhat11*Uhat1')';
0174 
0175 tmat = [eye(n-nunstab) -PHI];
0176 % Calculate matrices which would appear in the solution
0177 G0 = [ tmat * A ; 
0178        zeros(m,n-m) eye(m)];
0179    
0180 G1 = [ tmat * B ; 
0181        zeros(m,n)]; 
0182 
0183 G2 = [zeros(n-m,m);-eye(m)];
0184    
0185 G0I=inv(G0);
0186 
0187 H = Z*G0I;
0188 
0189 %%
0190 
0191 if ~unique
0192     p = 1;  % Set dimension of the sunspot vector
0193     m_1 = 0;%input('Enter 0 for M1 of zeros and 1 for least squares M1: ');
0194     m_2 = 1;%input('Enter 0 for M2 of zeros and 1 for M2 of ones: ');
0195     
0196     if m_1
0197        K1 = Q1*PPI-PHI*Q2*PPI;
0198        K2 = -V1*inv(D11)*U1'*Q2*PSI;
0199        AA = K1*K2;
0200        BB = K1*V2; 
0201        M1 = inv(BB'*BB)*BB'*(-AA);
0202     else
0203        M1 = zeros(k-r,l);
0204     end
0205     if m_2
0206        M2 = ones(k-r,p);
0207     else
0208        M2 = zeros(k-r,p);
0209     end
0210 end
0211 
0212 %% Calculate Solution Matrices
0213 S0 = real(H*[Q1 - PHI * Q2 ; inv(A22-B22)*Q2 ]*C);
0214 S1 = real(H*G1*Z');
0215 Sy = real(H*G2);
0216 Sf = real(B22\A22);
0217 Sz = real(B22\(Q2*PSI));
0218 
0219 if unique  
0220    S2 = H * [ Q1 - PHI * Q2 ; zeros(m,n) ] * PSI;
0221    S2 = real(S2);
0222    
0223    S3 = zeros(0);
0224 else
0225    S2 = H*[Q1-PHI*Q2 ; zeros(m,n)]*PSI + H*[(Q1*PPI-PHI*Q2*PPI)*(-V1*inv(D11)*U1'*Q2*PSI+V2*M1); zeros(m,l)];
0226    S2 = real(S2);
0227    
0228    S3 = H * [ (Q1*PPI - PHI*Q2*PPI)*V2*M2 ; zeros(m,p)];
0229    S3 = real(S3);
0230 end
0231 
0232 
0233 
0234

Generated on Mon 07-Feb-2011 12:06:56 by m2html © 2005