File contents
% This program generates artificial output data, derives the PVM-consistent
% current account and then:
% - estimates the model predicted CA series and compares it to the
% model-consistent simulated series in terms of correlation coefficients
% and variance ratios
% - tests the model both in terms of the regular Wald test of the cross-equation
% restrictions, the Wald test "linearized", and the F-test.
clear all;
lags=8; %number of lags in the VAR
r=0.01; %assumed real interest rate
initial=500; %number of observations to be weeded out to exclude effect of initial conditions
t=120; %sample size
T=initial+t;
nloop=10000; %number of data samples
volra=[]; corr=[]; test=[]; stat=[]; statexact=[]; statlin=[]; testlin=[]; testexact=[]; fcrit=[]; waldcrit=[]; kcoef=[];
% we will now generate nloop samples of T observations for Y using the estimated VAR for each country,
% and btain T observations of model-consistent CA.
% After discarding initial observations, we will use DY and the model-consistent CA to obtain the model predicted series
var = [-0.192310306 0.07640275 0.113610752 0.101462729 0.095366202 0.140838808 0.16289717 -0.017939816 0.063181448 -0.093393605 -0.198562471 0.219022759 -0.156202022 -0.094898027 0.128462974 0.094920512
-0.079870803 -0.011912842 -0.027831189 0.196466741 0.178293504 0.163994162 0.16877331 0.049335969 0.900637022 0.049590797 -0.058251237 -0.066107816 -0.072923869 -0.023516564 0.052925599 0.108542267]; %this is the matrix of VAR coefficients estimated from the country
companionvar = [var(1,:)
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
var(2,:)
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0]; % companion matrix of var
loopcount = 1;
for loop = 1:nloop;
dny=[];x=[]; u1=[]; u2=[]; ca=[]; capred=[];
z2=[]; lhs2=[]; x2=[]; x2p=[]; phi2=[]; res2=[]; x1=[]; phi1=[]; res1=[]; lhs1=[];
V11 =[]; V12=[]; V13=[]; V14=[]; V=[]; k=[];
% we will now create artificial data using the data generating process
dny(1:lags)=1; x(1:lags)=1; ca(1:lags)=0;
sigma=[0.000750832810452 0.000272781813982
0.000272781813982 0.000350984927679];
u=mvnrnd([0 0],sigma,T);
for i=lags+1:T;
dny(i)=var(1,:)*[dny(i-1) dny(i-2) dny(i-3) dny(i-4) dny(i-5) dny(i-6) dny(i-7) dny(i-8) x(i-1) x(i-2) x(i-3) x(i-4) x(i-5) x(i-6) x(i-7) x(i-8)]'+u(i,1); %this generates data for dny according to the estimated VARs
x(i)=var(2,:)*[dny(i-1) dny(i-2) dny(i-3) dny(i-4) dny(i-5) dny(i-6) dny(i-7) dny(i-8) x(i-1) x(i-2) x(i-3) x(i-4) x(i-5) x(i-6) x(i-7) x(i-8)]'+u(i,2);
ca(i) = -[1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]*(1/(1+r))*companionvar*inv(eye(2*lags)-(1/(1+r))*companionvar)*[dny(i) dny(i-1) dny(i-2) dny(i-3) dny(i-4) dny(i-5) dny(i-6) dny(i-7) x(i) x(i-1) x(i-2) x(i-3) x(i-4) x(i-5) x(i-6) x(i-7)]'; % %this is the model consistent current account
end;
z2=[dny(initial+1:T)'-mean(dny(initial+1:T)') ca(initial+1:T)'-mean(ca(initial+1:T)')];
% we are now going to estimate two systems
% 1) regress ca(t)-(1+r)ca(t-1)-dny(t) on the information set at time t-1,
% which gives us the 'exact' test
% 2) a VAR on ca and dny to obtain the predicted series and do the Wald
% test of the cross-equation restrictions of the model
for j=1:lags; x2=[x2 z2(lags+1-j:t-j,1)]; x2p=[x2p z2(lags+1-j:t-j,2)]; end; x2=[x2 x2p];
lhs1 = z2(lags+1:t,2)-(1+r)*z2(lags:t-1,2)-z2(lags+1:t,1);
x1=x2(1:t-lags,:);
phi1 = lhs1'*x1*inv(x1'*x1);
res1 = lhs1-(x1*(phi1'));
sig1 = (res1'*res1)/(size(res1,1)-2*lags);
lhs2 = z2(lags+1:t,:);
phi2 = lhs2'*x2*inv(x2'*x2);
res2 = lhs2-(x2*(phi2'));
companionphi2 = [phi2(1,:)
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
phi2(2,:)
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0]; % companion matrix of phi2
% This will measure the highest eigenvalue of the VAR matrix
eigen=eig(companionphi2);
for i=1:16;
if real(eigen(i))==eigen(i);
eigen(i)=eigen(i);
else;
eigen(i)=NaN;
end;
end;
eigmax(loopcount,1)=max(eigen);
eigmax2=eigmax(~isnan(eigmax));
% we are now constructing V, the variance-covariance matrix of VAR
% coefficients
V11 = (res2(:,1)'*res2(:,1))/(size(res2,1)-2*lags)*inv(x2'*x2);
V12 = (res2(:,1)'*res2(:,2))/(size(res2,1)-2*lags)*inv(x2'*x2);
V21 = (res2(:,2)'*res2(:,1))/(size(res2,1)-2*lags)*inv(x2'*x2);
V22 = (res2(:,2)'*res2(:,2))/(size(res2,1)-2*lags)*inv(x2'*x2);
V = [V11 V12;V21 V22]; diag(V);
k = -[1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]*(1/(1+r))*companionphi2*inv(eye(2*lags)-(1/(1+r))*companionphi2);
kcoef(loopcount,:)=k;
capred(1:lags) = 0;
for i=lags+1:t;
capred(i) = k*[z2(i,1) z2(i-1,1) z2(i-2,1) z2(i-3,1) z2(i-4,1) z2(i-5,1) z2(i-6,1) z2(i-7,1) z2(i,2) z2(i-1,2) z2(i-2,2) z2(i-3,2) z2(i-4,2) z2(i-5,2) z2(i-6,2) z2(i-7,2)]';
end;
% this constructs the Jacobian of derivatives of K with respect to the VAR
% parameters
for i = 1:4*lags;
a=[]; varvec=[]; kp=[];
a = reshape(transpose(phi2),1,4*lags);
a(i) = a(i)+0.000001;
varvec = [a(1:2*lags)
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
a(2*lags+1:4*lags)
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0];
kp = -[1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]*(1/(1+r))*varvec*inv(eye(2*lags)-(1/(1+r))*varvec);
J(1:2*lags,i) = (transpose(kp)-transpose(k))/0.000001;
end;
% We will now construct the 'exact' test...
statexact(loopcount,1) = ((lhs1'*lhs1-res1'*res1)/(2*lags))/((res1'*res1)/(size(res1,1)-2*lags)); %this is the F-statistic
fcrit(loopcount,1)=fcdf(statexact(loopcount,1),2*lags,size(res1,1)-2*lags);
% and the wald test
stat(loopcount,1)=(k-[0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0])*inv((J*V*J'))*(k-[0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0])';
waldcrit(loopcount,1)=chi2cdf(stat(loopcount,1),2*lags);
% and the linearized wald test
R=[-eye(2*lags) eye(2*lags)];
statlin(loopcount,1)=(R*[phi2(1,:)';phi2(2,:)']-[0 0 0 0 0 0 0 0 1+r 0 0 0 0 0 0 0]')'*inv(R*V*R')*(R*[phi2(1,:)';phi2(2,:)']-[0 0 0 0 0 0 0 0 1+r 0 0 0 0 0 0 0]');
waldlin(loopcount,1)=chi2cdf(statlin(loopcount,1),2*lags);
if waldcrit(loopcount,1) <= 0.95;
test(loopcount,1) = 1;
else;
test(loopcount,1) = 0;
end;
if waldlin(loopcount,1) <= 0.95;
testlin(loopcount,1) = 1;
else;
testlin(loopcount,1) = 0;
end;
if fcrit(loopcount,1) <= 0.95;
testexact(loopcount,1) = 1;
else;
testexact(loopcount,1) = 0;
end;
% this computes the variance ratio and correlation coefficient between predicted and actual
z=corrcoef(capred(lags+1:t),z2(lags+1:t,2));
corr(loopcount,1)=z(1,2);
volra(loopcount,1)=(std(capred(lags+1:t))/std(z2(lags+1:t,2)))^2;
loopcount=loopcount+1;
end; % end of loops
Click here to get the file
