From 96d716343d9f9c27d2266d6014f82a1a1daa4d12 Mon Sep 17 00:00:00 2001 From: Houtan Bastani Date: Wed, 15 Nov 2017 12:36:12 +0100 Subject: [PATCH] add pooled_fgls (with test) --- matlab/ols/pooled_fgls.m | 83 ++++++++++ matlab/ols/pooled_ols.m | 86 ++++++++-- .../panel_var/panel_var_diff_NB_2country.mod | 147 +++++++++++++++++ .../panel_var_diff_NB_simulation_test.mod | 151 ++++++++++++++++++ tests/ECB/pooled_fgls/run_simulation_test.m | 47 ++++++ 5 files changed, 498 insertions(+), 16 deletions(-) create mode 100644 matlab/ols/pooled_fgls.m create mode 100755 tests/ECB/panel_var/panel_var_diff_NB_2country.mod create mode 100644 tests/ECB/pooled_fgls/panel_var_diff_NB_simulation_test.mod create mode 100644 tests/ECB/pooled_fgls/run_simulation_test.m diff --git a/matlab/ols/pooled_fgls.m b/matlab/ols/pooled_fgls.m new file mode 100644 index 000000000..dd31faadf --- /dev/null +++ b/matlab/ols/pooled_fgls.m @@ -0,0 +1,83 @@ +function pooled_fgls(ds, param_common, param_regex) +% function pooled_fgls(ds, param_common, param_regex) +% Run Pooled FGLS +% Apply parameter values found to corresponding parameter values in the +% other blocks of the model +% +% INPUTS +% ds [dseries] data to use in estimation +% param_common [cellstr] List of values to insert into param_regex, +% e.g. country codes {'FR', 'DE', 'IT'} +% param_regex [cellstr] Where '*' should be replaced by the first +% value in param_common +% +% OUTPUTS +% none +% +% SPECIAL REQUIREMENTS +% dynare must be run with the option: json=parse + +% Copyright (C) 2017 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 . + +global M_ oo_ + +pooled_ols(ds, param_common, param_regex, true, 'pooled_fgls'); + +neq = length(fieldnames(oo_.pooled_fgls.resid)); +nobs = length(oo_.pooled_fgls.sample_range); +for i = 1:neq + ui = oo_.pooled_fgls.resid.(oo_.pooled_fgls.residnames{i}); + for j = i:neq + uj = oo_.pooled_fgls.resid.(oo_.pooled_fgls.residnames{j}); + M_.Sigma(i, j) = (ui'*uj)/nobs; + M_.Sigma(j, i) = M_.Sigma(i, j); + end +end + +kLeye = kron(chol(M_.Sigma), eye(nobs)); +[q, r] = qr(kLeye*oo_.pooled_fgls.X, 0); +oo_.pooled_fgls.beta = r\(q'*kLeye*oo_.pooled_fgls.Y); + +param_names_trim = cellstr(M_.param_names); +regexcountries = ['(' strjoin(param_common(1:end),'|') ')']; +pbeta = oo_.pooled_fgls.pbeta; +assigned_idxs = false(size(pbeta)); +for i = 1:length(param_regex) + beta_idx = strcmp(pbeta, strrep(param_regex{i}, '*', oo_.pooled_fgls.country_name)); + assigned_idxs = assigned_idxs | beta_idx; + value = oo_.pooled_fgls.beta(beta_idx); + assert(~isempty(value)); + M_.params(~cellfun(@isempty, regexp(param_names_trim, ... + strrep(param_regex{i}, '*', regexcountries)))) = value; +end +idxs = find(assigned_idxs == 0); +values = oo_.pooled_fgls.beta(idxs); +names = pbeta(idxs); +assert(length(values) == length(names)); +for i = 1:length(idxs) + M_.params(strcmp(param_names_trim, names{i})) = values(i); +end + +oo_.pooled_fgls = rmfield(oo_.pooled_fgls, 'X'); +oo_.pooled_fgls = rmfield(oo_.pooled_fgls, 'Y'); +oo_.pooled_fgls = rmfield(oo_.pooled_fgls, 'varcovar'); +oo_.pooled_fgls = rmfield(oo_.pooled_fgls, 'residnames'); +oo_.pooled_fgls = rmfield(oo_.pooled_fgls, 'pbeta'); +oo_.pooled_fgls = rmfield(oo_.pooled_fgls, 'country_name'); + +end diff --git a/matlab/ols/pooled_ols.m b/matlab/ols/pooled_ols.m index eccf1704a..10f3be656 100644 --- a/matlab/ols/pooled_ols.m +++ b/matlab/ols/pooled_ols.m @@ -1,15 +1,19 @@ -function pooled_ols(ds, param_common, param_regex) -% function pooled_ols(ds, param_common, param_regex) +function pooled_ols(ds, param_common, param_regex, overlapping_dates, save_structure_name) +% pooled_ols(ds, param_common, param_regex, overlapping_dates, save_structure_name) % Run Pooled OLS % Apply parameter values found to corresponding parameter values in the % other blocks of the model % % INPUTS -% ds [dseries] data to use in estimation -% param_common [cellstr] List of values to insert into param_regex, -% e.g. country codes {'FR', 'DE', 'IT'} -% param_regex [cellstr] Where '*' should be replaced by the first -% value in param_common +% ds [dseries] data to use in estimation +% param_common [cellstr] List of values to insert into param_regex, +% e.g. country codes {'FR', 'DE', 'IT'} +% param_regex [cellstr] Where '*' should be replaced by the first +% value in param_common +% overlapping_dates [bool] if the dates across the equations should +% overlap +% save_structure_name [string] Name of structure in oo_ to save results in +% (pooled_ols by default) % % OUTPUTS % none @@ -38,6 +42,7 @@ global M_ oo_ % Check input arguments assert(~isempty(ds) && isdseries(ds), 'The first argument must be a dseries'); + if isempty(param_common) && isempty(param_regex) disp('Performing OLS instead of Pooled OLS...') dyn_ols(ds); @@ -46,12 +51,24 @@ end assert(~isempty(param_common) && iscellstr(param_common), 'The second argument must be a cellstr'); assert(~isempty(param_regex) && iscellstr(param_regex), 'The third argument must be a cellstr'); +if nargin < 4 + overlapping_dates = false; +else + assert(islogical(overlapping_dates) && length(overlapping_dates) == 1, 'The fourth argument must be a bool'); +end + +if nargin < 5 + save_structure_name = 'pooled_ols'; +else + assert(ischar(save_structure_name), 'The fifth argument must be a string'); +end + +%% Read JSON jsonfile = [M_.fname '_original.json']; if exist(jsonfile, 'file') ~= 2 error('Could not find %s! Please use the json=parse option (See the Dynare invocation section in the reference manual).', jsonfile); end -%% Read JSON jsonmodel = loadjson(jsonfile); jsonmodel = jsonmodel.model; [lhs, rhs, lineno] = getEquationsByTags(jsonmodel); @@ -78,6 +95,8 @@ pbeta = {}; Y = []; X = []; startidxs = zeros(length(lhs), 1); +startdates = cell(length(lhs), 1); +enddates = cell(length(lhs), 1); residnames = cell(length(lhs), 1); for i = 1:length(lhs) rhs_ = strsplit(rhs{i}, {'+','-','*','/','^','log(','ln(','log10(','exp(','(',')','diff('}); @@ -161,14 +180,49 @@ for i = 1:length(lhs) lp = min(ydata.lastobservedperiod, xjdata.lastobservedperiod); startidxs(i) = length(Y) + 1; + startdates{i} = fp; + enddates{i} = lp; Y(startidxs(i):startidxs(i)+lp-fp, 1) = ydata(fp:lp).data; X(startidxs(i):startidxs(i)+lp-fp, pidxs) = xjdata(fp:lp).data; end +if overlapping_dates + maxfp = max([startdates{:}]); + minlp = min([enddates{:}]); + nobs = minlp - maxfp; + newY = zeros(nobs*length(lhs), 1); + newX = zeros(nobs*length(lhs), columns(X)); + newstartidxs = zeros(size(startidxs)); + newstartidxs(1) = 1; + for i = 1:length(lhs) + if i == length(lhs) + yds = dseries(Y(startidxs(i):end), startdates{i}); + xds = dseries(X(startidxs(i):end, :), startdates{i}); + else + yds = dseries(Y(startidxs(i):startidxs(i+1)-1), startdates{i}); + xds = dseries(X(startidxs(i):startidxs(i+1)-1, :), startdates{i}); + end + newY(newstartidxs(i):newstartidxs(i) + nobs, 1) = yds(maxfp:minlp).data; + newX(newstartidxs(i):newstartidxs(i) + nobs, :) = xds(maxfp:minlp, :).data; + if i ~= length(lhs) + newstartidxs(i+1) = newstartidxs(i) + nobs + 1; + end + end + Y = newY; + X = newX; + startidxs = newstartidxs; + oo_.(save_structure_name).sample_range = maxfp:minlp; + oo_.(save_structure_name).residnames = residnames; + oo_.(save_structure_name).Y = Y; + oo_.(save_structure_name).X = X; + oo_.(save_structure_name).pbeta = pbeta; + oo_.(save_structure_name).country_name = country_name; +end + %% Estimation % Estimated Parameters [q, r] = qr(X, 0); -oo_.pooled_ols.beta = r\(q'*Y); +oo_.(save_structure_name).beta = r\(q'*Y); % Assign parameter values back to parameters using param_regex & param_common param_names_trim = cellfun(@strtrim, num2cell(M_.param_names(:,:),2), 'Uniform', 0); @@ -177,28 +231,28 @@ assigned_idxs = false(size(pbeta)); for i = 1:length(param_regex) beta_idx = strcmp(pbeta, strrep(param_regex{i}, '*', country_name)); assigned_idxs = assigned_idxs | beta_idx; - value = oo_.pooled_ols.beta(beta_idx); + value = oo_.(save_structure_name).beta(beta_idx); assert(~isempty(value)); M_.params(~cellfun(@isempty, regexp(param_names_trim, ... strrep(param_regex{i}, '*', regexcountries)))) = value; end idxs = find(assigned_idxs == 0); -values = oo_.pooled_ols.beta(idxs); +values = oo_.(save_structure_name).beta(idxs); names = pbeta(idxs); assert(length(values) == length(names)); for i = 1:length(idxs) M_.params(strcmp(param_names_trim, names{i})) = values(i); end -residuals = Y - X * oo_.pooled_ols.beta; +residuals = Y - X * oo_.(save_structure_name).beta; for i = 1:length(lhs) if i == length(lhs) - oo_.pooled_ols.resid.(residnames{i}) = residuals(startidxs(i):end); + oo_.(save_structure_name).resid.(residnames{i}) = residuals(startidxs(i):end); else - oo_.pooled_ols.resid.(residnames{i}) = residuals(startidxs(i):startidxs(i+1)-1); + oo_.(save_structure_name).resid.(residnames{i}) = residuals(startidxs(i):startidxs(i+1)-1); end - oo_.pooled_ols.varcovar.(['eq' num2str(i)]) = oo_.pooled_ols.resid.(residnames{i})*oo_.pooled_ols.resid.(residnames{i})'; + oo_.(save_structure_name).varcovar.(['eq' num2str(i)]) = oo_.(save_structure_name).resid.(residnames{i})*oo_.(save_structure_name).resid.(residnames{i})'; idx = find(strcmp(residnames{i}, M_exo_names_trim)); - M_.Sigma_e(idx, idx) = var(oo_.pooled_ols.resid.(residnames{i})); + M_.Sigma_e(idx, idx) = var(oo_.(save_structure_name).resid.(residnames{i})); end end diff --git a/tests/ECB/panel_var/panel_var_diff_NB_2country.mod b/tests/ECB/panel_var/panel_var_diff_NB_2country.mod new file mode 100755 index 000000000..610826743 --- /dev/null +++ b/tests/ECB/panel_var/panel_var_diff_NB_2country.mod @@ -0,0 +1,147 @@ +// --+ options: json=compute +-- + +/* REMARK +** ------ +** +** You need to have the first line on top of the mod file. The options defined on this line are passed +** to the dynare command (you can add other options, separated by spaces or commas). The option defined +** here is mandatory for the decomposition. It forces Dynare to output another representation of the +** model in JSON file (additionaly to the matlab files) which is used here to manipulate the equations. +*/ + +var +U2_Q_YED +U2_G_YER +U2_STN +U2_ESTN +U2_EHIC +DE_Q_YED +DE_G_YER +DE_EHIC +; + +varexo +res_U2_Q_YED +res_U2_G_YER +res_U2_STN +res_U2_ESTN +res_U2_EHIC +res_DE_Q_YED +res_DE_G_YER +res_DE_EHIC +; + +parameters +u2_q_yed_ecm_u2_q_yed_L1 +u2_q_yed_ecm_u2_stn_L1 +u2_q_yed_u2_g_yer_L1 +u2_q_yed_u2_stn_L1 +u2_g_yer_ecm_u2_q_yed_L1 +u2_g_yer_ecm_u2_stn_L1 +u2_g_yer_u2_q_yed_L1 +u2_g_yer_u2_g_yer_L1 +u2_g_yer_u2_stn_L1 +u2_stn_ecm_u2_q_yed_L1 +u2_stn_ecm_u2_stn_L1 +u2_stn_u2_q_yed_L1 +u2_stn_u2_g_yer_L1 +u2_estn_u2_estn_L1 +u2_ehic_u2_ehic_L1 + +de_q_yed_ecm_de_q_yed_L1 +de_q_yed_ecm_u2_stn_L1 +de_q_yed_de_g_yer_L1 +de_q_yed_u2_stn_L1 +de_g_yer_ecm_de_q_yed_L1 +de_g_yer_ecm_u2_stn_L1 +de_g_yer_de_q_yed_L1 +de_g_yer_de_g_yer_L1 +de_g_yer_u2_stn_L1 +de_ehic_de_ehic_L1 +; + +u2_q_yed_ecm_u2_q_yed_L1 = -0.82237516589315 ; +u2_q_yed_ecm_u2_stn_L1 = -0.323715338568976 ; +u2_q_yed_u2_g_yer_L1 = 0.0401361895021084 ; +u2_q_yed_u2_stn_L1 = 0.058397703958446 ; +u2_g_yer_ecm_u2_q_yed_L1 = 0.0189896046977421 ; +u2_g_yer_ecm_u2_stn_L1 = -0.109597659887432 ; +u2_g_yer_u2_q_yed_L1 = 0.0037667967632025 ; +u2_g_yer_u2_g_yer_L1 = 0.480506381923644 ; +u2_g_yer_u2_stn_L1 = -0.0722359286123494 ; +u2_stn_ecm_u2_q_yed_L1 = -0.0438500662608356 ; +u2_stn_ecm_u2_stn_L1 = -0.153283917138772 ; +u2_stn_u2_q_yed_L1 = 0.0328744983772825 ; +u2_stn_u2_g_yer_L1 = 0.292121949736756 ; +u2_estn_u2_estn_L1 = 1 ; +u2_ehic_u2_ehic_L1 = 1 ; + +de_q_yed_ecm_de_q_yed_L1 = -0.822375165893149 ; +de_q_yed_ecm_u2_stn_L1 = -0.323715338568977 ; +de_q_yed_de_g_yer_L1 = 0.0401361895021082 ; +de_q_yed_u2_stn_L1 = 0.0583977039584461 ; +de_g_yer_ecm_de_q_yed_L1 = 0.0189896046977422 ; +de_g_yer_ecm_u2_stn_L1 = -0.109597659887433 ; +de_g_yer_de_q_yed_L1 = 0.00376679676320256; +de_g_yer_de_g_yer_L1 = 0.480506381923643 ; +de_g_yer_u2_stn_L1 = -0.0722359286123494 ; +de_ehic_de_ehic_L1 = 1 ; + + +model; + +diff(U2_Q_YED) = u2_q_yed_ecm_u2_q_yed_L1 * (U2_Q_YED(-1) - U2_EHIC(-1)) + + u2_q_yed_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + u2_q_yed_u2_g_yer_L1 * diff(U2_G_YER(-1)) + + u2_q_yed_u2_stn_L1 * diff(U2_STN(-1)) + + res_U2_Q_YED ; + +diff(U2_G_YER) = u2_g_yer_ecm_u2_q_yed_L1 * (U2_Q_YED(-1) - U2_EHIC(-1)) + + u2_g_yer_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + u2_g_yer_u2_q_yed_L1 * diff(U2_Q_YED(-1)) + + u2_g_yer_u2_g_yer_L1 * diff(U2_G_YER(-1)) + + u2_g_yer_u2_stn_L1 * diff(U2_STN(-1)) + + res_U2_G_YER ; + +diff(U2_STN) = u2_stn_ecm_u2_q_yed_L1 * (U2_Q_YED(-1) - U2_EHIC(-1)) + + u2_stn_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + u2_stn_u2_q_yed_L1 * diff(U2_Q_YED(-1)) + + u2_stn_u2_g_yer_L1 * diff(U2_G_YER(-1)) + + res_U2_STN ; + +U2_ESTN = u2_estn_u2_estn_L1 * U2_ESTN(-1) + + res_U2_ESTN ; + +U2_EHIC = u2_ehic_u2_ehic_L1 * U2_EHIC(-1) + + res_U2_EHIC ; + +diff(DE_Q_YED) = de_q_yed_ecm_de_q_yed_L1 * (DE_Q_YED(-1) - DE_EHIC(-1)) + + de_q_yed_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + de_q_yed_de_g_yer_L1 * diff(DE_G_YER(-1)) + + de_q_yed_u2_stn_L1 * diff(U2_STN(-1)) + + res_DE_Q_YED ; + +diff(DE_G_YER) = de_g_yer_ecm_de_q_yed_L1 * (DE_Q_YED(-1) - DE_EHIC(-1)) + + de_g_yer_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + de_g_yer_de_q_yed_L1 * diff(DE_Q_YED(-1)) + + de_g_yer_de_g_yer_L1 * diff(DE_G_YER(-1)) + + de_g_yer_u2_stn_L1 * diff(U2_STN(-1)) + + res_DE_G_YER ; + +DE_EHIC = de_ehic_de_ehic_L1 * DE_EHIC(-1) + + res_DE_EHIC ; + +end; + +pooled_fgls(dseries('data_var.mat'), ... + {'de','u2'}, ... + {'*_q_yed_ecm_*_q_yed_L1', ... + '*_q_yed_ecm_u2_stn_L1', ... + '*_q_yed_*_g_yer_L1', ... + '*_q_yed_u2_stn_L1', ... + '*_g_yer_ecm_*_q_yed_L1', ... + '*_g_yer_ecm_u2_stn_L1', ... + '*_g_yer_*_q_yed_L1', ... + '*_g_yer_*_g_yer_L1', ... + '*_g_yer_u2_stn_L1', ... + '*_ehic_*_ehic_L1'}); diff --git a/tests/ECB/pooled_fgls/panel_var_diff_NB_simulation_test.mod b/tests/ECB/pooled_fgls/panel_var_diff_NB_simulation_test.mod new file mode 100644 index 000000000..db41db854 --- /dev/null +++ b/tests/ECB/pooled_fgls/panel_var_diff_NB_simulation_test.mod @@ -0,0 +1,151 @@ +// --+ options: json=compute +-- + +/* REMARK +** ------ +** +** You need to have the first line on top of the mod file. The options defined on this line are passed +** to the dynare command (you can add other options, separated by spaces or commas). The option defined +** here is mandatory for the decomposition. It forces Dynare to output another representation of the +** model in JSON file (additionaly to the matlab files) which is used here to manipulate the equations. +*/ + +var +U2_Q_YED +U2_G_YER +U2_STN +U2_ESTN +U2_EHIC +DE_Q_YED +DE_G_YER +DE_EHIC + +; + +varexo +res_U2_Q_YED +res_U2_G_YER +res_U2_STN +res_U2_ESTN +res_U2_EHIC +res_DE_Q_YED +res_DE_G_YER +res_DE_EHIC +; + +parameters +u2_q_yed_ecm_u2_q_yed_L1 +u2_q_yed_ecm_u2_stn_L1 +u2_q_yed_u2_g_yer_L1 +u2_q_yed_u2_stn_L1 +u2_g_yer_ecm_u2_q_yed_L1 +u2_g_yer_ecm_u2_stn_L1 +u2_g_yer_u2_q_yed_L1 +u2_g_yer_u2_g_yer_L1 +u2_g_yer_u2_stn_L1 +u2_stn_ecm_u2_q_yed_L1 +u2_stn_ecm_u2_stn_L1 +u2_stn_u2_q_yed_L1 +u2_stn_u2_g_yer_L1 +u2_estn_u2_estn_L1 +u2_ehic_u2_ehic_L1 + +de_q_yed_ecm_de_q_yed_L1 +de_q_yed_ecm_u2_stn_L1 +de_q_yed_de_g_yer_L1 +de_q_yed_u2_stn_L1 +de_g_yer_ecm_de_q_yed_L1 +de_g_yer_ecm_u2_stn_L1 +de_g_yer_de_q_yed_L1 +de_g_yer_de_g_yer_L1 +de_g_yer_u2_stn_L1 +de_ehic_de_ehic_L1 + + +; + +u2_q_yed_ecm_u2_q_yed_L1 = -0.82237516589315 ; +u2_q_yed_ecm_u2_stn_L1 = -0.323715338568976 ; +u2_q_yed_u2_g_yer_L1 = 0.0401361895021084 ; +u2_q_yed_u2_stn_L1 = 0.058397703958446 ; +u2_g_yer_ecm_u2_q_yed_L1 = 0.0189896046977421 ; +u2_g_yer_ecm_u2_stn_L1 = -0.109597659887432 ; +u2_g_yer_u2_q_yed_L1 = 0.0037667967632025 ; +u2_g_yer_u2_g_yer_L1 = 0.480506381923644 ; +u2_g_yer_u2_stn_L1 = -0.0722359286123494 ; +u2_stn_ecm_u2_q_yed_L1 = -0.0438500662608356 ; +u2_stn_ecm_u2_stn_L1 = -0.153283917138772 ; +u2_stn_u2_q_yed_L1 = 0.0328744983772825 ; +u2_stn_u2_g_yer_L1 = 0.292121949736756 ; +u2_estn_u2_estn_L1 = 1 ; +u2_ehic_u2_ehic_L1 = 1 ; + +de_q_yed_ecm_de_q_yed_L1 = -0.822375165893149 ; +de_q_yed_ecm_u2_stn_L1 = -0.323715338568977 ; +de_q_yed_de_g_yer_L1 = 0.0401361895021082 ; +de_q_yed_u2_stn_L1 = 0.0583977039584461 ; +de_g_yer_ecm_de_q_yed_L1 = 0.0189896046977422 ; +de_g_yer_ecm_u2_stn_L1 = -0.109597659887433 ; +de_g_yer_de_q_yed_L1 = 0.00376679676320256; +de_g_yer_de_g_yer_L1 = 0.480506381923643 ; +de_g_yer_u2_stn_L1 = -0.0722359286123494 ; +de_ehic_de_ehic_L1 = 1 ; + + +model(linear); + +diff(U2_Q_YED) = u2_q_yed_ecm_u2_q_yed_L1 * (U2_Q_YED(-1) - U2_EHIC(-1)) + + u2_q_yed_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + u2_q_yed_u2_g_yer_L1 * diff(U2_G_YER(-1)) + + u2_q_yed_u2_stn_L1 * diff(U2_STN(-1)) + + res_U2_Q_YED ; + +diff(U2_G_YER) = u2_g_yer_ecm_u2_q_yed_L1 * (U2_Q_YED(-1) - U2_EHIC(-1)) + + u2_g_yer_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + u2_g_yer_u2_q_yed_L1 * diff(U2_Q_YED(-1)) + + u2_g_yer_u2_g_yer_L1 * diff(U2_G_YER(-1)) + + u2_g_yer_u2_stn_L1 * diff(U2_STN(-1)) + + res_U2_G_YER ; + +diff(U2_STN) = u2_stn_ecm_u2_q_yed_L1 * (U2_Q_YED(-1) - U2_EHIC(-1)) + + u2_stn_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + u2_stn_u2_q_yed_L1 * diff(U2_Q_YED(-1)) + + u2_stn_u2_g_yer_L1 * diff(U2_G_YER(-1)) + + res_U2_STN ; + +U2_ESTN = u2_estn_u2_estn_L1 * U2_ESTN(-1) + + res_U2_ESTN ; + +U2_EHIC = u2_ehic_u2_ehic_L1 * U2_EHIC(-1) + + res_U2_EHIC ; + +diff(DE_Q_YED) = de_q_yed_ecm_de_q_yed_L1 * (DE_Q_YED(-1) - DE_EHIC(-1)) + + de_q_yed_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + de_q_yed_de_g_yer_L1 * diff(DE_G_YER(-1)) + + de_q_yed_u2_stn_L1 * diff(U2_STN(-1)) + + res_DE_Q_YED ; + +diff(DE_G_YER) = de_g_yer_ecm_de_q_yed_L1 * (DE_Q_YED(-1) - DE_EHIC(-1)) + + de_g_yer_ecm_u2_stn_L1 * (U2_STN(-1) - U2_ESTN(-1)) + + de_g_yer_de_q_yed_L1 * diff(DE_Q_YED(-1)) + + de_g_yer_de_g_yer_L1 * diff(DE_G_YER(-1)) + + de_g_yer_u2_stn_L1 * diff(U2_STN(-1)) + + res_DE_G_YER ; + +DE_EHIC = de_ehic_de_ehic_L1 * DE_EHIC(-1) + + res_DE_EHIC ; + + + +end; + +shocks; +var res_U2_Q_YED = 0.005; +var res_U2_G_YER = 0.005; +var res_U2_STN = 0.005; +var res_U2_ESTN = 0.005; +var res_U2_EHIC = 0.005; +var res_DE_Q_YED = 0.005; +var res_DE_G_YER = 0.005; +var res_DE_EHIC = 0.005; +end; + diff --git a/tests/ECB/pooled_fgls/run_simulation_test.m b/tests/ECB/pooled_fgls/run_simulation_test.m new file mode 100644 index 000000000..f2a94c933 --- /dev/null +++ b/tests/ECB/pooled_fgls/run_simulation_test.m @@ -0,0 +1,47 @@ +close all + +dynare panel_var_diff_NB_simulation_test.mod; + +NSIMS = 1000; + +calibrated_values = M_.params; +Sigma_e = M_.Sigma_e; + +options_.bnlms.set_dynare_seed_to_default = false; + +M_endo_names_trim = cellstr(M_.endo_names); +nparampool = length(M_.params); +BETA = zeros(NSIMS, nparampool); +for i=1:NSIMS + i + firstobs = rand(3, length(M_endo_names_trim)); + M_.params = calibrated_values; + M_.Sigma_e = Sigma_e; + simdata = simul_backward_model(dseries(firstobs, dates('1995Q1'), M_endo_names_trim), 10000); + simdata = simdata(simdata.dates(5001:6000)); + pooled_fgls(simdata, ... + {'de','u2'}, ... + {'*_q_yed_ecm_*_q_yed_L1', ... + '*_q_yed_ecm_u2_stn_L1', ... + '*_q_yed_*_g_yer_L1', ... + '*_q_yed_u2_stn_L1', ... + '*_g_yer_ecm_*_q_yed_L1', ... + '*_g_yer_ecm_u2_stn_L1', ... + '*_g_yer_*_q_yed_L1', ... + '*_g_yer_*_g_yer_L1', ... + '*_g_yer_u2_stn_L1', ... + '*_ehic_*_ehic_L1'}); + BETA(i, :) = M_.params'; + oldsim = simdata; +end + +mean(BETA)' - calibrated_values + +for i=1:nparampool + figure + hold on + title(strrep(M_.param_names(i,:), '_', '\_')); + histogram(BETA(:,i),50); + line([calibrated_values(i) calibrated_values(i)], [0 NSIMS/10], 'LineWidth', 2, 'Color', 'r'); + hold off +end \ No newline at end of file