Add unit tests for HP filter, bandpass filter, and spectral density

time-shift
Johannes Pfeifer 2015-10-12 20:33:09 +02:00
parent 40de494568
commit 2c63ca8843
4 changed files with 353 additions and 1 deletions

View File

@ -7,7 +7,10 @@ MODFILES = \
estimation/MH_recover/fs2000_recover.mod \
estimation/t_proposal/fs2000_student.mod \
estimation/TaRB/fs2000_tarb.mod \
moments/example1_var_decomp \
moments/example1_var_decomp.mod \
moments/example1_hp_test.mod \
moments/example1_bp_test.mod \
moments/test_AR1_spectral_density.mod \
gsa/ls2003.mod \
gsa/ls2003a.mod \
gsa/cod_ML_morris/cod_ML_morris.mod \

View File

@ -0,0 +1,116 @@
/*
* Example 1 from F. Collard (2001): "Stochastic simulations with DYNARE:
* A practical guide" (see "guide.pdf" in the documentation directory).
*/
/*
* Copyright (C) 2001-2010 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/>.
*/
var y, c, k, a, h, b;
varexo e, u;
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;
steady(solve_algo=4,maxit=1000);
stoch_simul(order=1,nofunctions,irf=0,bandpass_filter=[6 32],hp_ngrid=8192);
oo_filtered_all_shocks_theoretical=oo_;
stoch_simul(order=1,nofunctions,periods=1000000);
oo_filtered_all_shocks_simulated=oo_;
shocks;
var e; stderr 0;
var u; stderr 0.009;
var e, u = phi*0.009*0;
end;
stoch_simul(order=1,nofunctions,irf=0,periods=0);
oo_filtered_one_shock_theoretical=oo_;
stoch_simul(order=1,nofunctions,periods=5000000);
oo_filtered_one_shock_simulated=oo_;
if max(abs((1-diag(oo_filtered_one_shock_simulated.var)./(diag(oo_filtered_all_shocks_simulated.var)))*100-oo_filtered_all_shocks_theoretical.variance_decomposition(:,1)))>2
error('Variance Decomposition wrong')
end
if max(max(abs(oo_filtered_all_shocks_simulated.var-oo_filtered_all_shocks_theoretical.var)))>2e-4;
error('Covariance wrong')
end
if max(max(abs(oo_filtered_one_shock_simulated.var-oo_filtered_one_shock_theoretical.var)))>1e-4;
error('Covariance wrong')
end
for ii=1:options_.ar
autocorr_model_all_shocks_simulated(:,ii)=diag(oo_filtered_all_shocks_simulated.autocorr{ii});
autocorr_model_all_shocks_theoretical(:,ii)=diag(oo_filtered_all_shocks_theoretical.autocorr{ii});
autocorr_model_one_shock_simulated(:,ii)=diag(oo_filtered_one_shock_simulated.autocorr{ii});
autocorr_model_one_shock_theoretical(:,ii)=diag(oo_filtered_one_shock_theoretical.autocorr{ii});
end
if max(max(abs(autocorr_model_all_shocks_simulated-autocorr_model_all_shocks_theoretical)))>2e-2;
error('Correlation wrong')
end
if max(max(abs(autocorr_model_one_shock_simulated-autocorr_model_one_shock_theoretical)))>2e-2;
error('Correlation wrong')
end

View File

@ -0,0 +1,157 @@
/*
* Example 1 from F. Collard (2001): "Stochastic simulations with DYNARE:
* A practical guide" (see "guide.pdf" in the documentation directory).
*/
/*
* Copyright (C) 2001-2010 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/>.
*/
var y, c, k, a, h, b;
varexo e, u;
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;
steady(solve_algo=4,maxit=1000);
options_.hp_ngrid=2048*4;
options_.bandpass.indicator=0;
options_.bandpass.passband=[6 32];
stoch_simul(order=1,nofunctions,hp_filter=1600,irf=0);
total_var_filtered=diag(oo_.var);
oo_filtered_all_shocks=oo_;
stoch_simul(order=1,nofunctions,hp_filter=0,periods=5000000,nomoments);
options_.nomoments=0;
oo_unfiltered_all_shocks=oo_;
[junk, y_filtered]=sample_hp_filter(y,1600);
[junk, c_filtered]=sample_hp_filter(c,1600);
[junk, k_filtered]=sample_hp_filter(k,1600);
[junk, a_filtered]=sample_hp_filter(a,1600);
[junk, h_filtered]=sample_hp_filter(h,1600);
[junk, b_filtered]=sample_hp_filter(b,1600);
verbatim;
total_std_all_shocks_filtered_sim=std([y_filtered c_filtered k_filtered a_filtered h_filtered b_filtered])
cov_filtered_all_shocks=cov([y_filtered c_filtered k_filtered a_filtered h_filtered b_filtered])
acf(1,:)=autocorr([y_filtered ],5)';
acf(2,:)=autocorr([c_filtered ],5)';
acf(3,:)=autocorr([k_filtered ],5)';
acf(4,:)=autocorr([a_filtered ],5)';
acf(5,:)=autocorr([h_filtered ],5)';
acf(6,:)=autocorr([b_filtered ],5)';
autocorr_filtered_all_shocks=acf(:,2:end);
end;
shocks;
var e; stderr 0;
var u; stderr 0.009;
var e, u = phi*0.009*0;
end;
stoch_simul(order=1,nofunctions,hp_filter=1600,irf=0,periods=0);
total_var_filtered_one_shock=diag(oo_.var);
oo_filtered_one_shock=oo_;
stoch_simul(order=1,nofunctions,hp_filter=0,periods=5000000,nomoments);
oo_unfiltered_one_shock=oo_;
[junk, y_filtered]=sample_hp_filter(y,1600);
[junk, c_filtered]=sample_hp_filter(c,1600);
[junk, k_filtered]=sample_hp_filter(k,1600);
[junk, a_filtered]=sample_hp_filter(a,1600);
[junk, h_filtered]=sample_hp_filter(h,1600);
[junk, b_filtered]=sample_hp_filter(b,1600);
verbatim;
total_std_one_shock_filtered_sim=std([y_filtered c_filtered k_filtered a_filtered h_filtered b_filtered])
cov_filtered_one_shock=cov([y_filtered c_filtered k_filtered a_filtered h_filtered b_filtered])
acf(1,:)=autocorr([y_filtered ],5)';
acf(2,:)=autocorr([c_filtered ],5)';
acf(3,:)=autocorr([k_filtered ],5)';
acf(4,:)=autocorr([a_filtered ],5)';
acf(5,:)=autocorr([h_filtered ],5)';
acf(6,:)=autocorr([b_filtered ],5)';
autocorr_filtered_one_shock=acf(:,2:end);
end;
if max(abs((1-(total_std_one_shock_filtered_sim.^2)./(total_std_all_shocks_filtered_sim.^2))*100-oo_filtered_all_shocks.variance_decomposition(:,1)'))>2
error('Variance Decomposition wrong')
end
if max(max(abs(oo_filtered_all_shocks.var-cov_filtered_all_shocks)))>1e-4;
error('Covariance wrong')
end
if max(max(abs(oo_filtered_one_shock.var-cov_filtered_one_shock)))>5e-5;
error('Covariance wrong')
end
for ii=1:options_.ar
autocorr_model_all_shocks(:,ii)=diag(oo_filtered_all_shocks.autocorr{ii});
autocorr_model_one_shock(:,ii)=diag(oo_filtered_one_shock.autocorr{ii});
end
if max(max(abs(autocorr_model_all_shocks-autocorr_filtered_all_shocks)))>1e-2;
error('Covariance wrong')
end
if max(max(abs(autocorr_model_one_shock-autocorr_filtered_one_shock)))>1e-2;
error('Covariance wrong')
end

View File

@ -0,0 +1,76 @@
var white_noise ar1;
varexo e;
parameters phi;
phi=0.9;
model;
white_noise=e;
ar1=phi*ar1(-1)+e;
end;
shocks;
var e = 1;
end;
options_.SpectralDensity.trigger=1;
options_.bandpass.indicator=0;
options_.hp_ngrid=2048;
stoch_simul(order=1,nofunctions,hp_filter=0,irf=0,periods=1000000);
white_noise_sample=white_noise;
theoretical_spectrum_white_noise=1^2/(2*pi); %Hamilton (1994), 6.1.9
if max(abs(oo_.SpectralDensity.density(strmatch('white_noise',M_.endo_names,'exact'),:)-theoretical_spectrum_white_noise))>1e-10
error('Spectral Density is wrong')
end
theoretical_spectrum_AR1=1/(2*pi)*(1^2./(1+phi^2-2*phi*cos(oo_.SpectralDensity.freqs))); %Hamilton (1994), 6.1.13
if max(abs(oo_.SpectralDensity.density(strmatch('ar1',M_.endo_names,'exact'),:)-theoretical_spectrum_AR1'))>1e-10
error('Spectral Density is wrong')
end
stoch_simul(order=1,nofunctions,hp_filter=1600,irf=0,periods=0);
lambda=options_.hp_filter;
Kalman_gain=(4*lambda*(1 - cos(oo_.SpectralDensity.freqs)).^2 ./ (1 + 4*lambda*(1 - cos(oo_.SpectralDensity.freqs)).^2));
theoretical_spectrum_white_noise_hp_filtered=1^2/(2*pi)*Kalman_gain.^2; %Hamilton (1994), 6.1.9
if max(abs(oo_.SpectralDensity.density(strmatch('white_noise',M_.endo_names,'exact'),:)-theoretical_spectrum_white_noise_hp_filtered'))>1e-10
error('Spectral Density is wrong')
end
theoretical_spectrum_AR1_hp_filtered=1/(2*pi)*(1^2./(1+phi^2-2*phi*cos(oo_.SpectralDensity.freqs))).*Kalman_gain.^2; %Hamilton (1994), 6.1.13
if max(abs(oo_.SpectralDensity.density(strmatch('ar1',M_.endo_names,'exact'),:)-theoretical_spectrum_AR1_hp_filtered'))>1e-10
error('Spectral Density is wrong')
end
options_.hp_filter=0;
stoch_simul(order=1,nofunctions,bandpass_filter=[6 32],irf=0);
theoretical_spectrum_white_noise=repmat(theoretical_spectrum_white_noise,1,options_.hp_ngrid);
passband=oo_.SpectralDensity.freqs>=2*pi/options_.bandpass.passband(2) & oo_.SpectralDensity.freqs<=2*pi/options_.bandpass.passband(1);
if max(abs(oo_.SpectralDensity.density(strmatch('white_noise',M_.endo_names,'exact'),passband)-theoretical_spectrum_white_noise(passband)))>1e-10
error('Spectral Density is wrong')
end
if max(abs(oo_.SpectralDensity.density(strmatch('white_noise',M_.endo_names,'exact'),~passband)-0))>1e-10
error('Spectral Density is wrong')
end
if max(abs(oo_.SpectralDensity.density(strmatch('ar1',M_.endo_names,'exact'),passband)-theoretical_spectrum_AR1(passband)'))>1e-10
error('Spectral Density is wrong')
end
if max(abs(oo_.SpectralDensity.density(strmatch('ar1',M_.endo_names,'exact'),~passband)-0))>1e-10
error('Spectral Density is wrong')
end
% [pow,f]=psd(a_sample,1024,1,[],512);
% figure
% plot(f,pow/(2*pi))
%
% % figure
% % [pow,f]=psd(a_sample,1000,1,[],500);
% % plot(f(3:end)*2*pi,pow(3:end)/(2*pi));