var_forecast: use example1 in forecast, add code to use estimation via rfvar3
parent
d46d107837
commit
dc7fca7ece
|
@ -216,108 +216,3 @@ dimy = size(ydum);
|
||||||
ydum = reshape(permute(ydum,[1 3 2]),dimy(1)*dimy(3),nv);
|
ydum = reshape(permute(ydum,[1 3 2]),dimy(1)*dimy(3),nv);
|
||||||
xdum = reshape(permute(xdum,[1 3 2]),dimy(1)*dimy(3),nx);
|
xdum = reshape(permute(xdum,[1 3 2]),dimy(1)*dimy(3),nx);
|
||||||
breaks = breaks(1:(end-1));
|
breaks = breaks(1:(end-1));
|
||||||
|
|
||||||
|
|
||||||
function var=rfvar3(ydata,lags,xdata,breaks,lambda,mu)
|
|
||||||
%function var=rfvar3(ydata,lags,xdata,breaks,lambda,mu)
|
|
||||||
% This algorithm goes for accuracy without worrying about memory requirements.
|
|
||||||
% ydata: dependent variable data matrix
|
|
||||||
% xdata: exogenous variable data matrix
|
|
||||||
% lags: number of lags
|
|
||||||
% breaks: rows in ydata and xdata after which there is a break. This allows for
|
|
||||||
% discontinuities in the data (e.g. war years) and for the possibility of
|
|
||||||
% adding dummy observations to implement a prior. This must be a column vector.
|
|
||||||
% Note that a single dummy observation becomes lags+1 rows of the data matrix,
|
|
||||||
% with a break separating it from the rest of the data. The function treats the
|
|
||||||
% first lags observations at the top and after each "break" in ydata and xdata as
|
|
||||||
% initial conditions.
|
|
||||||
% lambda: weight on "co-persistence" prior dummy observations. This expresses
|
|
||||||
% belief that when data on *all* y's are stable at their initial levels, they will
|
|
||||||
% tend to persist at that level. lambda=5 is a reasonable first try. With lambda<0,
|
|
||||||
% constant term is not included in the dummy observation, so that stationary models
|
|
||||||
% with means equal to initial ybar do not fit the prior mean. With lambda>0, the prior
|
|
||||||
% implies that large constants are unlikely if unit roots are present.
|
|
||||||
% mu: weight on "own persistence" prior dummy observation. Expresses belief
|
|
||||||
% that when y_i has been stable at its initial level, it will tend to persist
|
|
||||||
% at that level, regardless of the values of other variables. There is
|
|
||||||
% one of these for each variable. A reasonable first guess is mu=2.
|
|
||||||
% The program assumes that the first lags rows of ydata and xdata are real data, not dummies.
|
|
||||||
% Dummy observations should go at the end, if any. If pre-sample x's are not available,
|
|
||||||
% repeating the initial xdata(lags+1,:) row or copying xdata(lags+1:2*lags,:) into
|
|
||||||
% xdata(1:lags,:) are reasonable subsititutes. These values are used in forming the
|
|
||||||
% persistence priors.
|
|
||||||
|
|
||||||
% Original file downloaded from:
|
|
||||||
% http://sims.princeton.edu/yftp/VARtools/matlab/rfvar3.m
|
|
||||||
|
|
||||||
[T,nvar] = size(ydata);
|
|
||||||
nox = isempty(xdata);
|
|
||||||
if ~nox
|
|
||||||
[T2,nx] = size(xdata);
|
|
||||||
else
|
|
||||||
T2 = T;
|
|
||||||
nx = 0;
|
|
||||||
xdata = zeros(T2,0);
|
|
||||||
end
|
|
||||||
% note that x must be same length as y, even though first part of x will not be used.
|
|
||||||
% This is so that the lags parameter can be changed without reshaping the xdata matrix.
|
|
||||||
if T2 ~= T, error('Mismatch of x and y data lengths'),end
|
|
||||||
if nargin < 4
|
|
||||||
nbreaks = 0;
|
|
||||||
breaks = [];
|
|
||||||
else
|
|
||||||
nbreaks = length(breaks);
|
|
||||||
end
|
|
||||||
breaks = [0;breaks;T];
|
|
||||||
smpl = [];
|
|
||||||
for nb = 1:nbreaks+1
|
|
||||||
smpl = [smpl;[breaks(nb)+lags+1:breaks(nb+1)]'];
|
|
||||||
end
|
|
||||||
Tsmpl = size(smpl,1);
|
|
||||||
X = zeros(Tsmpl,nvar,lags);
|
|
||||||
for is = 1:length(smpl)
|
|
||||||
X(is,:,:) = ydata(smpl(is)-(1:lags),:)';
|
|
||||||
end
|
|
||||||
X = [X(:,:) xdata(smpl,:)];
|
|
||||||
y = ydata(smpl,:);
|
|
||||||
% Everything now set up with input data for y=Xb+e
|
|
||||||
|
|
||||||
% Add persistence dummies
|
|
||||||
if lambda ~= 0 || mu > 0
|
|
||||||
ybar = mean(ydata(1:lags,:),1);
|
|
||||||
if ~nox
|
|
||||||
xbar = mean(xdata(1:lags,:),1);
|
|
||||||
else
|
|
||||||
xbar = [];
|
|
||||||
end
|
|
||||||
if lambda ~= 0
|
|
||||||
if lambda>0
|
|
||||||
xdum = lambda*[repmat(ybar,1,lags) xbar];
|
|
||||||
else
|
|
||||||
lambda = -lambda;
|
|
||||||
xdum = lambda*[repmat(ybar,1,lags) zeros(size(xbar))];
|
|
||||||
end
|
|
||||||
ydum = zeros(1,nvar);
|
|
||||||
ydum(1,:) = lambda*ybar;
|
|
||||||
y = [y;ydum];
|
|
||||||
X = [X;xdum];
|
|
||||||
end
|
|
||||||
if mu>0
|
|
||||||
xdum = [repmat(diag(ybar),1,lags) zeros(nvar,nx)]*mu;
|
|
||||||
ydum = mu*diag(ybar);
|
|
||||||
X = [X;xdum];
|
|
||||||
y = [y;ydum];
|
|
||||||
end
|
|
||||||
end
|
|
||||||
|
|
||||||
% Compute OLS regression and residuals
|
|
||||||
[vl,d,vr] = svd(X,0);
|
|
||||||
di = 1./diag(d);
|
|
||||||
B = (vr.*repmat(di',nvar*lags+nx,1))*vl'*y;
|
|
||||||
u = y-X*B;
|
|
||||||
xxi = vr.*repmat(di',nvar*lags+nx,1);
|
|
||||||
xxi = xxi*xxi';
|
|
||||||
|
|
||||||
var.B = B;
|
|
||||||
var.u = u;
|
|
||||||
var.xxi = xxi;
|
|
||||||
|
|
|
@ -0,0 +1,122 @@
|
||||||
|
function var=rfvar3(ydata,lags,xdata,breaks,lambda,mu)
|
||||||
|
%function var=rfvar3(ydata,lags,xdata,breaks,lambda,mu)
|
||||||
|
% This algorithm goes for accuracy without worrying about memory requirements.
|
||||||
|
% ydata: dependent variable data matrix
|
||||||
|
% xdata: exogenous variable data matrix
|
||||||
|
% lags: number of lags
|
||||||
|
% breaks: rows in ydata and xdata after which there is a break. This allows for
|
||||||
|
% discontinuities in the data (e.g. war years) and for the possibility of
|
||||||
|
% adding dummy observations to implement a prior. This must be a column vector.
|
||||||
|
% Note that a single dummy observation becomes lags+1 rows of the data matrix,
|
||||||
|
% with a break separating it from the rest of the data. The function treats the
|
||||||
|
% first lags observations at the top and after each "break" in ydata and xdata as
|
||||||
|
% initial conditions.
|
||||||
|
% lambda: weight on "co-persistence" prior dummy observations. This expresses
|
||||||
|
% belief that when data on *all* y's are stable at their initial levels, they will
|
||||||
|
% tend to persist at that level. lambda=5 is a reasonable first try. With lambda<0,
|
||||||
|
% constant term is not included in the dummy observation, so that stationary models
|
||||||
|
% with means equal to initial ybar do not fit the prior mean. With lambda>0, the prior
|
||||||
|
% implies that large constants are unlikely if unit roots are present.
|
||||||
|
% mu: weight on "own persistence" prior dummy observation. Expresses belief
|
||||||
|
% that when y_i has been stable at its initial level, it will tend to persist
|
||||||
|
% at that level, regardless of the values of other variables. There is
|
||||||
|
% one of these for each variable. A reasonable first guess is mu=2.
|
||||||
|
% The program assumes that the first lags rows of ydata and xdata are real data, not dummies.
|
||||||
|
% Dummy observations should go at the end, if any. If pre-sample x's are not available,
|
||||||
|
% repeating the initial xdata(lags+1,:) row or copying xdata(lags+1:2*lags,:) into
|
||||||
|
% xdata(1:lags,:) are reasonable subsititutes. These values are used in forming the
|
||||||
|
% persistence priors.
|
||||||
|
|
||||||
|
% Original file downloaded from:
|
||||||
|
% http://sims.princeton.edu/yftp/VARtools/matlab/rfvar3.m
|
||||||
|
|
||||||
|
% Copyright (C) 2003-2007 Christopher Sims
|
||||||
|
% Copyright (C) 2007-2012 Dynare Team
|
||||||
|
%
|
||||||
|
% This file is part of Dynare.
|
||||||
|
%
|
||||||
|
% Dynare is free software: you can redistribute it and/or modify
|
||||||
|
% it under the terms of the GNU General Public License as published by
|
||||||
|
% the Free Software Foundation, either version 3 of the License, or
|
||||||
|
% (at your option) any later version.
|
||||||
|
%
|
||||||
|
% Dynare is distributed in the hope that it will be useful,
|
||||||
|
% but WITHOUT ANY WARRANTY; without even the implied warranty of
|
||||||
|
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
||||||
|
% GNU General Public License for more details.
|
||||||
|
%
|
||||||
|
% You should have received a copy of the GNU General Public License
|
||||||
|
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
|
||||||
|
|
||||||
|
[T,nvar] = size(ydata);
|
||||||
|
nox = isempty(xdata);
|
||||||
|
if ~nox
|
||||||
|
[T2,nx] = size(xdata);
|
||||||
|
else
|
||||||
|
T2 = T;
|
||||||
|
nx = 0;
|
||||||
|
xdata = zeros(T2,0);
|
||||||
|
end
|
||||||
|
% note that x must be same length as y, even though first part of x will not be used.
|
||||||
|
% This is so that the lags parameter can be changed without reshaping the xdata matrix.
|
||||||
|
if T2 ~= T, error('Mismatch of x and y data lengths'),end
|
||||||
|
if nargin < 4
|
||||||
|
nbreaks = 0;
|
||||||
|
breaks = [];
|
||||||
|
else
|
||||||
|
nbreaks = length(breaks);
|
||||||
|
end
|
||||||
|
breaks = [0;breaks;T];
|
||||||
|
smpl = [];
|
||||||
|
for nb = 1:nbreaks+1
|
||||||
|
smpl = [smpl;[breaks(nb)+lags+1:breaks(nb+1)]'];
|
||||||
|
end
|
||||||
|
Tsmpl = size(smpl,1);
|
||||||
|
X = zeros(Tsmpl,nvar,lags);
|
||||||
|
for is = 1:length(smpl)
|
||||||
|
X(is,:,:) = ydata(smpl(is)-(1:lags),:)';
|
||||||
|
end
|
||||||
|
X = [X(:,:) xdata(smpl,:)];
|
||||||
|
y = ydata(smpl,:);
|
||||||
|
% Everything now set up with input data for y=Xb+e
|
||||||
|
|
||||||
|
% Add persistence dummies
|
||||||
|
if lambda ~= 0 || mu > 0
|
||||||
|
ybar = mean(ydata(1:lags,:),1);
|
||||||
|
if ~nox
|
||||||
|
xbar = mean(xdata(1:lags,:),1);
|
||||||
|
else
|
||||||
|
xbar = [];
|
||||||
|
end
|
||||||
|
if lambda ~= 0
|
||||||
|
if lambda>0
|
||||||
|
xdum = lambda*[repmat(ybar,1,lags) xbar];
|
||||||
|
else
|
||||||
|
lambda = -lambda;
|
||||||
|
xdum = lambda*[repmat(ybar,1,lags) zeros(size(xbar))];
|
||||||
|
end
|
||||||
|
ydum = zeros(1,nvar);
|
||||||
|
ydum(1,:) = lambda*ybar;
|
||||||
|
y = [y;ydum];
|
||||||
|
X = [X;xdum];
|
||||||
|
end
|
||||||
|
if mu>0
|
||||||
|
xdum = [repmat(diag(ybar),1,lags) zeros(nvar,nx)]*mu;
|
||||||
|
ydum = mu*diag(ybar);
|
||||||
|
X = [X;xdum];
|
||||||
|
y = [y;ydum];
|
||||||
|
end
|
||||||
|
end
|
||||||
|
|
||||||
|
% Compute OLS regression and residuals
|
||||||
|
[vl,d,vr] = svd(X,0);
|
||||||
|
di = 1./diag(d);
|
||||||
|
B = (vr.*repmat(di',nvar*lags+nx,1))*vl'*y;
|
||||||
|
u = y-X*B;
|
||||||
|
xxi = vr.*repmat(di',nvar*lags+nx,1);
|
||||||
|
xxi = xxi*xxi';
|
||||||
|
|
||||||
|
var.B = B;
|
||||||
|
var.u = u;
|
||||||
|
var.xxi = xxi;
|
||||||
|
end
|
|
@ -0,0 +1,51 @@
|
||||||
|
// Example 1 from Collard's guide to Dynare
|
||||||
|
var y, c, k, a, h, b;
|
||||||
|
varexo e, u;
|
||||||
|
|
||||||
|
verbatim;
|
||||||
|
% I want these comments included in
|
||||||
|
% example1.m 1999q1 1999y
|
||||||
|
%
|
||||||
|
var = 1;
|
||||||
|
end;
|
||||||
|
|
||||||
|
parameters beta, rho, alpha, delta, theta, psi, tau;
|
||||||
|
|
||||||
|
alpha = 0.36;
|
||||||
|
rho = 0.95;
|
||||||
|
tau = 0.025;
|
||||||
|
beta = 0.99;
|
||||||
|
delta = 0.025;
|
||||||
|
psi = 0;
|
||||||
|
theta = 2.95;
|
||||||
|
|
||||||
|
phi = 0.1;
|
||||||
|
|
||||||
|
model;
|
||||||
|
c*theta*h^(1+psi)=(1-alpha)*y;
|
||||||
|
k = beta*(((exp(b)*c)/(exp(b(+1))*c(+1)))
|
||||||
|
*(exp(b(+1))*alpha*y(+1)+(1-delta)*k));
|
||||||
|
y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
|
||||||
|
k = exp(b)*(y-c)+(1-delta)*k(-1);
|
||||||
|
a = rho*a(-1)+tau*b(-1) + e;
|
||||||
|
b = tau*a(-1)+rho*b(-1) + u;
|
||||||
|
end;
|
||||||
|
|
||||||
|
initval;
|
||||||
|
y = 1.08068253095672;
|
||||||
|
c = 0.80359242014163;
|
||||||
|
h = 0.29175631001732;
|
||||||
|
k = 11.08360443260358;
|
||||||
|
a = 0;
|
||||||
|
b = 0;
|
||||||
|
e = 0;
|
||||||
|
u = 0;
|
||||||
|
end;
|
||||||
|
|
||||||
|
shocks;
|
||||||
|
var e; stderr 0.009;
|
||||||
|
var u; stderr 0.009;
|
||||||
|
var e, u = phi*0.009*0.009;
|
||||||
|
end;
|
||||||
|
|
||||||
|
stoch_simul(order=1, periods=200);
|
|
@ -0,0 +1,53 @@
|
||||||
|
// Example 1 from Collard's guide to Dynare
|
||||||
|
var y, c, k, a, h, b;
|
||||||
|
varexo e, u;
|
||||||
|
|
||||||
|
verbatim;
|
||||||
|
% I want these comments included in
|
||||||
|
% example1.m 1999q1 1999y
|
||||||
|
%
|
||||||
|
var = 1;
|
||||||
|
end;
|
||||||
|
|
||||||
|
parameters beta, rho, alpha, delta, theta, psi, tau;
|
||||||
|
|
||||||
|
alpha = 0.36;
|
||||||
|
rho = 0.95;
|
||||||
|
tau = 0.025;
|
||||||
|
beta = 0.99;
|
||||||
|
delta = 0.025;
|
||||||
|
psi = 0;
|
||||||
|
theta = 2.95;
|
||||||
|
|
||||||
|
phi = 0.1;
|
||||||
|
|
||||||
|
var_model(model_name=my_var_est, order=1) y c;
|
||||||
|
|
||||||
|
model;
|
||||||
|
c*theta*h^(1+psi)=(1-alpha)*y;
|
||||||
|
k = beta*(((exp(b)*c)/(exp(b(+1))*c(1)))
|
||||||
|
*(exp(b(+1))*alpha*var_expectation(y, model_name=my_var_est)+(1-delta)*k));
|
||||||
|
y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
|
||||||
|
k = exp(b)*(y-c)+(1-delta)*k(-1);
|
||||||
|
a = rho*a(-1)+tau*b(-1) + e;
|
||||||
|
b = tau*a(-1)+rho*b(-1) + u;
|
||||||
|
end;
|
||||||
|
|
||||||
|
initval;
|
||||||
|
y = 1.08068253095672;
|
||||||
|
c = 0.80359242014163;
|
||||||
|
h = 0.29175631001732;
|
||||||
|
k = 11.08360443260358;
|
||||||
|
a = 0;
|
||||||
|
b = 0;
|
||||||
|
e = 0;
|
||||||
|
u = 0;
|
||||||
|
end;
|
||||||
|
|
||||||
|
shocks;
|
||||||
|
var e; stderr 0.009;
|
||||||
|
var u; stderr 0.009;
|
||||||
|
var e, u = phi*0.009*0.009;
|
||||||
|
end;
|
||||||
|
|
||||||
|
stoch_simul(order=1, periods=200);
|
|
@ -1,30 +1,45 @@
|
||||||
clearvars
|
clearvars
|
||||||
clearvars -global
|
clearvars -global
|
||||||
|
close all
|
||||||
|
|
||||||
!rm nkm_saved_data.mat
|
!rm nkm_saved_data.mat
|
||||||
!rm my_var_est.mat
|
!rm my_var_est.mat
|
||||||
!rm nkm_var_saved_data.mat
|
!rm nkm_var_saved_data.mat
|
||||||
|
|
||||||
%% Call Dynare without VAR
|
%% Call Dynare without VAR
|
||||||
dynare nkm.mod
|
dynare example1.mod
|
||||||
save('nkm_saved_data.mat')
|
save('nkm_saved_data.mat')
|
||||||
|
|
||||||
%% Estimate VAR using simulated data
|
%% Estimate VAR using simulated data
|
||||||
disp('VAR Estimation');
|
disp('VAR Estimation');
|
||||||
Y = oo_.endo_simul(2:3, 2:end);
|
|
||||||
Z = [ones(1, size(Y,2)); oo_.endo_simul(2:3, 1:end-1)];
|
|
||||||
|
|
||||||
% OLS
|
|
||||||
B = Y*transpose(Z)/(Z*transpose(Z));
|
|
||||||
|
|
||||||
%% save values
|
% MLS estimate of mu and B (autoregressive_matrices)
|
||||||
|
% Y = mu + B*Z
|
||||||
|
% from New Introduction to Multiple Time Series Analysis
|
||||||
|
Y = oo_.endo_simul(1:2, 2:end);
|
||||||
|
Z = [ ...
|
||||||
|
ones(1, size(Y,2)); ...
|
||||||
|
oo_.endo_simul(1:2, 1:end-1); ...
|
||||||
|
];
|
||||||
|
%B = Y*Z'*inv(Z*Z');
|
||||||
|
B = Y*Z'/(Z*Z');
|
||||||
mu = B(:, 1);
|
mu = B(:, 1);
|
||||||
autoregressive_matrices{1} = B(:, 2:end);
|
autoregressive_matrices{1} = B(:, 2:end);
|
||||||
|
|
||||||
|
% Sims
|
||||||
|
% (provides same result as above)
|
||||||
|
% var = rfvar3(Y',1,zeros(size(Y')),0,5,2)
|
||||||
|
% mu = var.B(:, 1);
|
||||||
|
% autoregressive_matrices{1} = var.B(:, 2:end);
|
||||||
|
|
||||||
|
%% save values
|
||||||
save('my_var_est.mat', 'mu', 'autoregressive_matrices');
|
save('my_var_est.mat', 'mu', 'autoregressive_matrices');
|
||||||
|
|
||||||
%% Call Dynare with VAR
|
%% Call Dynare with VAR
|
||||||
clearvars
|
clearvars
|
||||||
clearvars -global
|
clearvars -global
|
||||||
dynare nkm_var.mod
|
dynare example1_var.mod
|
||||||
save('nkm_var_saved_data.mat')
|
save('nkm_var_saved_data.mat')
|
||||||
|
|
||||||
%% compare values
|
%% compare values
|
||||||
|
@ -35,11 +50,19 @@ zerotol = 1e-12;
|
||||||
nv = load('nkm_saved_data.mat');
|
nv = load('nkm_saved_data.mat');
|
||||||
wv = load('nkm_var_saved_data.mat');
|
wv = load('nkm_var_saved_data.mat');
|
||||||
|
|
||||||
assert(max(max(abs(nv.y - wv.y))) < zerotol);
|
ridx = 3;
|
||||||
assert(max(max(abs(nv.pi - wv.pi))) < zerotol);
|
cidx = 2;
|
||||||
assert(max(max(abs(nv.i - wv.i))) < zerotol);
|
exo_names = nv.M_.exo_names;
|
||||||
|
endo_names = nv.M_.endo_names;
|
||||||
|
|
||||||
fn = fieldnames(nv.oo_.irfs);
|
for i = 1:length(exo_names)
|
||||||
for i=1:length(fn)
|
figure('Name', ['Shock to ' exo_names(i)]);
|
||||||
assert(max(max(abs(nv.oo_.irfs.(fn{i}) - (wv.oo_.irfs.(fn{i}))))) < zerotol);
|
for j = 1:length(endo_names)
|
||||||
|
subplot(ridx, cidx, j);
|
||||||
|
hold on
|
||||||
|
title(endo_names(j));
|
||||||
|
plot(nv.oo_.irfs.([endo_names(j) '_' exo_names(i)]));
|
||||||
|
plot(wv.oo_.irfs.([endo_names(j) '_' exo_names(i)]), '--');
|
||||||
|
hold off
|
||||||
|
end
|
||||||
end
|
end
|
Loading…
Reference in New Issue