adding tests for kalman smoother

git-svn-id: https://www.dynare.org/svn/dynare/dynare_v4@2197 ac1d8469-bf42-47a9-8791-bf33cf982152
time-shift
michel 2008-10-26 20:03:48 +00:00
parent 7c6ed357b6
commit 0e125f541e
9 changed files with 572 additions and 0 deletions

View File

@ -0,0 +1,53 @@
var z dw dx dy dc1 dc2;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, uniform_pdf,,, 0.01, 0.1;
end;
varobs dw dx dy z;
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0);
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algo1_mode);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(6,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.dw d.dx d.dy d.z];
disp(max(max(abs(dat(1000:1199,:)-S(:,[2:4 1])))));

View File

@ -0,0 +1,81 @@
var z dw dx dy dc1 dc2;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, uniform_pdf,,, 0.01, 0.1;
end;
varobs dw dx dy z;
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,kalman_algo=2);
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algo1_mode,kalman_algo=2);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(6,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.dw d.dx d.dy d.z];
disp(max(max(abs(dat(1000:1199,:)-S(:,[2:4 1])))));
o1 = load('algo1_results');
obj_endo={'SmoothedVariables'; 'FilteredVariables'; 'UpdatedVariables'};
obj_exo = {'SmoothedShocks';};
nobj_endo = size(obj_endo,1);
nobj_exo = size(obj_exo,1);
for i=1:nobj_endo;
err_endo = zeros(eval(['size(oo_.' obj_endo{i} '.' M_.endo_names(1,:) ',1);']),M_.endo_nbr);
for j=1:M_.endo_nbr;
var1 = eval(['o1.oo_.' obj_endo{i} '.' M_.endo_names(j,:)]);
var2 = eval(['oo_.' obj_endo{i} '.' M_.endo_names(j,:)]);
err_endo(:,j) = var1-var2;
end;
disp(max(max(abs(err_endo))));
end;
err_exo = zeros(200,M_.exo_nbr,nobj_exo);
for i=1:nobj_exo;
err_exo = zeros(size(eval(['oo_.' obj_exo{i} '.' M_.exo_names(1,:)]),1),M_.exo_nbr);
for j=1:M_.exo_nbr;
var1 = eval(['o1.oo_.' obj_exo{i} '.' M_.exo_names(j,:)]);
var2 = eval(['oo_.' obj_exo{i} '.' M_.exo_names(j,:)]);
err_exo(:,j,i) = var1 - var2;
end;
disp(max(max(abs(err_exo))));
end;
disp(max(max(max(abs(err_exo)))));

View File

@ -0,0 +1,56 @@
var z dw dx dy dc1 dc2 w x y;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
w = w(-1) + dw;
x = x(-1) + dx;
y = y(-1) + dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, inv_gamma_pdf,0.01, inf;
end;
varobs w x y;
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,diffuse_filter);
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algo3_mode,diffuse_filter);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2 X.w X.x X.y];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(M_.endo_nbr,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.w d.x d.y];
disp(max(max(abs(dat(1000:1199,:)-S(:,[7:9])))));

View File

@ -0,0 +1,84 @@
var z dw dx dy dc1 dc2 w x y;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
w = w(-1) + dw;
x = x(-1) + dx;
y = y(-1) + dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, inv_gamma_pdf,0.01, inf;
end;
varobs w x y;
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,diffuse_filter);
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algo3_mode,diffuse_filter,kalman_algo=4);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2 X.w X.x X.y];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(M_.endo_nbr,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.w d.x d.y];
disp(max(max(abs(dat(1000:1199,:)-S(:,[7:9])))));
o1 = load('algo3_results');
obj_endo={'SmoothedVariables'; 'FilteredVariables'; 'UpdatedVariables'};
obj_exo = {'SmoothedShocks';};
nobj_endo = size(obj_endo,1);
nobj_exo = size(obj_exo,1);
for i=1:nobj_endo;
err_endo = zeros(eval(['size(oo_.' obj_endo{i} '.' M_.endo_names(1,:) ',1);']),M_.endo_nbr);
for j=1:M_.endo_nbr;
var1 = eval(['o1.oo_.' obj_endo{i} '.' M_.endo_names(j,:)]);
var2 = eval(['oo_.' obj_endo{i} '.' M_.endo_names(j,:)]);
err_endo(:,j) = var1-var2;
end;
disp(max(max(abs(err_endo))));
end;
err_exo = zeros(200,M_.exo_nbr,nobj_exo);
for i=1:nobj_exo;
err_exo = zeros(size(eval(['oo_.' obj_exo{i} '.' M_.exo_names(1,:)]),1),M_.exo_nbr);
for j=1:M_.exo_nbr;
var1 = eval(['o1.oo_.' obj_exo{i} '.' M_.exo_names(j,:)]);
var2 = eval(['oo_.' obj_exo{i} '.' M_.exo_names(j,:)]);
err_exo(:,j,i) = var1 - var2;
end;
disp(max(max(abs(err_exo))));
end;
disp(max(max(max(abs(err_exo)))));

View File

@ -0,0 +1,54 @@
var z dw dx dy dc1 dc2 y;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
y = y(-1) + dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, uniform_pdf,,, 0.01, 0.1;
end;
varobs dw dx y z;
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,diffuse_filter);
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algo3_mode,diffuse_filter);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2 X.y];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(M_.endo_nbr,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.dw d.dx d.y d.z];
disp(max(max(abs(dat(1000:1199,:)-S(:,[2 3 7 1])))));

View File

@ -0,0 +1,57 @@
var z dw dx dy dc1 dc2;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, uniform_pdf,,, 0.01, 0.1;
stderr z, normal_pdf, 0, 0.1;
end;
varobs dw dx dy z;
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0);
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algoH1_mode);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(6,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.dw d.dx d.dy d.z];
X = oo_.SmoothedMeasurementErrors;
ME = [X.dw X.dx X.dy X.z];
disp(max(max(abs(dat(1000:1199,:)-S(:,[2:4 1])-ME))));

View File

@ -0,0 +1,87 @@
var z dw dx dy dc1 dc2;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, uniform_pdf,,, 0.01, 0.1;
stderr z, normal_pdf, 0, 0.1;
end;
varobs dw dx dy z;
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,kalman_algo=2);
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algoH1_mode,kalman_algo=2);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(6,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.dw d.dx d.dy d.z];
X = oo_.SmoothedMeasurementErrors;
ME = [X.dw X.dx X.dy X.z];
disp(max(max(abs(dat(1000:1199,:)-S(:,[2:4 1])-ME))));
o1 = load('algoH1_results');
obj_endo={'SmoothedVariables'; 'FilteredVariables'; 'UpdatedVariables'};
obj_exo = {'SmoothedShocks';};
nobj_endo = size(obj_endo,1);
nobj_exo = size(obj_exo,1);
for i=1:nobj_endo;
err_endo = zeros(eval(['size(oo_.' obj_endo{i} '.' M_.endo_names(1,:) ',1);']),M_.endo_nbr);
for j=1:M_.endo_nbr;
var1 = eval(['o1.oo_.' obj_endo{i} '.' M_.endo_names(j,:)]);
var2 = eval(['oo_.' obj_endo{i} '.' M_.endo_names(j,:)]);
err_endo(:,j) = var1-var2;
end;
disp(max(max(abs(err_endo))));
end;
err_exo = zeros(200,M_.exo_nbr,nobj_exo);
for i=1:nobj_exo;
err_exo = zeros(size(eval(['oo_.' obj_exo{i} '.' M_.exo_names(1,:)]),1),M_.exo_nbr);
for j=1:M_.exo_nbr;
var1 = eval(['o1.oo_.' obj_exo{i} '.' M_.exo_names(j,:)]);
var2 = eval(['oo_.' obj_exo{i} '.' M_.exo_names(j,:)]);
err_exo(:,j,i) = var1 - var2;
end;
disp(max(max(abs(err_exo))));
end;
disp(max(max(max(abs(err_exo)))));
disp(max(max(abs(o1.oo_.SmoothedMeasurementErrors.z - oo_.SmoothedMeasurementErrors.z))));

View File

@ -0,0 +1,62 @@
var z dw dx dy dc1 dc2 w x y;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
model(linear);
dw = rho_w*dw(-1)+a1*(dc1(-1))+e_w;
dx = rho_x*dx(-1)+a2*(dc1(-1))+e_x;
dy = rho_y*dy(-1)+a3*(dc2(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
dc1 = dc1(-1)+dx-b*dy-c*dw;
dc2 = dc2(-1)+dx-b*dy;
w = w(-1) + dw;
x = x(-1) + dx;
y = y(-1) + dy;
end;
estimated_params;
rho_w, normal_pdf, 0.5,0.2;
rho_x, normal_pdf, 0.5,0.2;
rho_y, normal_pdf, 0.5,0.2;
rho_z, normal_pdf, 0.8,0.2;
a1, normal_pdf, 0.1,0.2;
a2, normal_pdf, -0.1,0.2;
a3, normal_pdf, 0.1,0.2;
b , normal_pdf, 1,0.2;
c , normal_pdf, 1,0.2;
stderr e_w, uniform_pdf,,, 0.01, 0.1;
stderr e_x, uniform_pdf,,, 0.01, 0.1;
stderr e_y, uniform_pdf,,, 0.01, 0.1;
stderr e_z, inv_gamma_pdf,0.01, inf;
stderr w, inv_gamma_pdf, 0.01,inf;
end;
varobs w x y;
estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,diffuse_filter);
//estimation(datafile=data,first_obs=1000,nobs=200,mh_replic=0,mode_compute=0,mode_file=algo3_mode,diffuse_filter);
stoch_simul(irf=0);
//checking smoother consistency
X = oo_.SmoothedVariables;
S = [X.z X.dw X.dx X.dy X.dc1 X.dc2 X.w X.x X.y];
X = oo_.SmoothedShocks;
E = [X.e_w X.e_x X.e_y X.e_z];
A = oo_.dr.ghx;
B = oo_.dr.ghu;
err = zeros(M_.endo_nbr,200);
for t=2:200;
err(:,t) = S(t,:)'-A*S(t-1,:)'-B*E(t,:)';
end;
disp(max(max(abs(err))));
d=load('data');
dat = [d.w d.x d.y];
X = oo_.SmoothedMeasurementErrors;
ME = [X.w X.x X.y];
disp(max(max(abs(dat(1000:1199,:)-S(:,[7:9])-ME))));

View File

@ -0,0 +1,38 @@
var w x y z dw dx dy;
varexo e_w e_x e_y e_z;
parameters rho_w rho_x rho_y rho_z a1 a2 a3 b c;
rho_w = 0.5;
rho_x = 0.5;
rho_y = 0.5;
rho_z = 0.8;
a1 = 0.1;
a2 = -0.1;
a3 = 0.1;
b = 1;
c = 1;
model(linear);
dw = rho_w*dw(-1)+a1*(x(-1)-b*y(-1)-c*w(-1))+e_w;
dx = rho_x*dx(-1)+a2*(x(-1)-b*y(-1)-c*w(-1))+e_x;
dy = rho_y*dy(-1)+a3*(x(-1)-b*y(-1))+e_y;
z = rho_z*z(-1)+dw-dx+e_z;
w = w(-1)+dw;
x = x(-1)+dx;
y = y(-1)+dy;
end;
shocks;
var e_w; stderr 0.05;
var e_x; stderr 0.05;
var e_y; stderr 0.05;
var e_z; stderr 0.05;
end;
stoch_simul(periods=2000,irf=0,simul_seed=7);
plot([w x y z]);
save data w x y z dw dx dy;