diff --git a/matlab/dsge_likelihood.m b/matlab/dsge_likelihood.m index 56c8dc2fa..53d04922c 100644 --- a/matlab/dsge_likelihood.m +++ b/matlab/dsge_likelihood.m @@ -352,7 +352,7 @@ if (kalman_algo == 2) || (kalman_algo == 4) Q = blkdiag(Q,H); R = blkdiag(R,eye(pp)); Pstar = blkdiag(Pstar,H); - Pinf = blckdiag(Pinf,zeros(pp)); + Pinf = blkdiag(Pinf,zeros(pp)); H = zeros(pp,1); mmm = mm+pp; end @@ -395,7 +395,7 @@ switch DynareOptions.lik_init if kalman_algo == 0 kalman_algo = 3; elseif ~((kalman_algo == 3) || (kalman_algo == 4)) - error(['diffuse filter: options_.kalman_algo can only be equal ' ... + error(['The model requires Diffuse filter, but you specified a different Kalman filter. You must set options_.kalman_algo ' ... 'to 0 (default), 3 or 4']) end [Z,T,R,QT,Pstar,Pinf] = schur_statespace_transformation(Z,T,R,Q,DynareOptions.qz_criterium); @@ -436,7 +436,7 @@ switch DynareOptions.lik_init Q = blkdiag(Q,H); R = blkdiag(R,eye(pp)); Pstar = blkdiag(Pstar,H); - Pinf = blckdiag(Pinf,zeros(pp)); + Pinf = blkdiag(Pinf,zeros(pp)); H1 = zeros(pp,1); mmm = mm+pp; end @@ -726,7 +726,7 @@ if (kalman_algo==2) || (kalman_algo==4) Q = blkdiag(Q,H); R = blkdiag(R,eye(pp)); Pstar = blkdiag(Pstar,H); - Pinf = blckdiag(Pinf,zeros(pp)); + Pinf = blkdiag(Pinf,zeros(pp)); H1 = zeros(pp,1); mmm = mm+pp; end diff --git a/matlab/dynare_estimation_1.m b/matlab/dynare_estimation_1.m index 654de281e..1e4eeb34d 100644 --- a/matlab/dynare_estimation_1.m +++ b/matlab/dynare_estimation_1.m @@ -291,9 +291,9 @@ if ~isequal(options_.mode_compute,0) && ~options_.mh_posterior_mode_estimation end parameter_names = bayestopt_.name; if options_.cova_compute || options_.mode_compute==5 || options_.mode_compute==6 - save([M_.fname '_mode.mat'],'xparam1','hh','parameter_names'); + save([M_.fname '_mode.mat'],'xparam1','hh','parameter_names','fval'); else - save([M_.fname '_mode.mat'],'xparam1','parameter_names'); + save([M_.fname '_mode.mat'],'xparam1','parameter_names','fval'); end end diff --git a/matlab/initial_estimation_checks.m b/matlab/initial_estimation_checks.m index f3b9ce06c..642023e0a 100644 --- a/matlab/initial_estimation_checks.m +++ b/matlab/initial_estimation_checks.m @@ -122,6 +122,9 @@ else b=0; fval = 0; end +if DynareOptions.debug + DynareResults.likelihood_at_initial_parameters=fval; +end DynareOptions.analytic_derivation=ana_deriv; if DynareOptions.dsge_var || strcmp(func2str(objective_function),'non_linear_dsge_likelihood') diff --git a/tests/Makefile.am b/tests/Makefile.am index 41cb046e3..e391108cb 100644 --- a/tests/Makefile.am +++ b/tests/Makefile.am @@ -169,6 +169,10 @@ MODFILES = \ kalman_filter_smoother/fs2000_smoother_only.mod \ kalman_filter_smoother/check_variable_dimensions/fs2000.mod \ kalman_filter_smoother/check_variable_dimensions/fs2000_ML.mod \ + kalman/likelihood_from_dynare/fs2000_corr_ME.mod \ + kalman/likelihood_from_dynare/fs2000_corr_ME_missing.mod \ + kalman/likelihood_from_dynare/fs2000_uncorr_ME.mod \ + kalman/likelihood_from_dynare/fs2000_uncorr_ME_missing.mod \ second_order/burnside_1.mod \ second_order/ds1.mod \ second_order/ds2.mod \ @@ -403,6 +407,9 @@ EXTRA_DIST = \ ms-sbvar/archive-files/specification_2v2c.dat \ recursive/data_ca1.m \ kalman_filter_smoother/fsdat_simul.m \ + kalman/likelihood_from_dynare/fsdat_simul_corr_ME_missing.m \ + kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME.m \ + kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME_missing.m \ identification/kim/kim2_steadystate.m \ identification/as2007/as2007_steadystate.m \ estimation/fsdat_simul.m \ diff --git a/tests/kalman/likelihood_from_dynare/fs2000_corr_ME.mod b/tests/kalman/likelihood_from_dynare/fs2000_corr_ME.mod new file mode 100644 index 000000000..304bd7929 --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fs2000_corr_ME.mod @@ -0,0 +1,137 @@ +/* + * This file is based on the cash in advance model described + * Frank Schorfheide (2000): "Loss function-based evaluation of DSGE models", + * Journal of Applied Econometrics, 15(6), 645-670. + * + * The equations are taken from J. Nason and T. Cogley (1994): "Testing the + * implications of long-run neutrality for monetary business cycle models", + * Journal of Applied Econometrics, 9, S37-S70. + * Note that there is an initial minus sign missing in equation (A1), p. S63. + * + * This implementation was written by Michel Juillard. Please note that the + * following copyright notice only applies to this Dynare implementation of the + * model. + */ + +/* + * Copyright (C) 2004-2013 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 . + */ + +var m P c e W R k d n l gy_obs gp_obs y dA; +varexo e_a e_m; + +parameters alp bet gam mst rho psi del theta; + +alp = 0.33; +bet = 0.99; +gam = 0.003; +mst = 1.011; +rho = 0.7; +psi = 0.787; +del = 0.02; +theta=0; + +model; +dA = exp(gam+e_a); +log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m; +-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0; +W = l/n; +-(psi/(1-psi))*(c*P/(1-n))+l/n = 0; +R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W; +1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0; +c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1); +P*c = m; +m-1+d = l; +e = exp(e_a); +y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a)); +gy_obs = dA*y/y(-1); +gp_obs = (P/P(-1))*m(-1)/dA; +end; + +steady_state_model; + dA = exp(gam); + gst = 1/dA; + m = mst; + khst = ( (1-gst*bet*(1-del)) / (alp*gst^alp*bet) )^(1/(alp-1)); + xist = ( ((khst*gst)^alp - (1-gst*(1-del))*khst)/mst )^(-1); + nust = psi*mst^2/( (1-alp)*(1-psi)*bet*gst^alp*khst^alp ); + n = xist/(nust+xist); + P = xist + nust; + k = khst*n; + + l = psi*mst*n/( (1-psi)*(1-n) ); + c = mst/P; + d = l - mst + 1; + y = k^alp*n^(1-alp)*gst^alp; + R = mst/bet; + W = l/n; + ist = y-c; + q = 1 - d; + + e = 1; + + gp_obs = m/dA; + gy_obs = dA; +end; + +varobs gp_obs gy_obs; + +shocks; +var e_a; stderr 0.014; +var e_m; stderr 0.005; +corr gy_obs,gp_obs = 0.5; +end; + +steady; + + +estimated_params; +alp, 0.356; +gam, 0.0085; +del, 0.01; +stderr e_a, 0.035449; +stderr e_m, 0.008862; +corr e_m, e_a, 0; +stderr gp_obs, 1; +stderr gy_obs, 1; +corr gp_obs, gy_obs,0; +end; + +options_.TeX=1; +options_.debug=1; + +%%default +options_.lik_init=1; +estimation(kalman_algo=0,mode_compute=4,order=1,datafile='../../fs2000/fsdat_simul',smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_0=oo_.likelihood_at_initial_parameters; +%%Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=1,mode_file=fs2000_corr_ME_mode,mode_compute=0,order=1,datafile='../../fs2000/fsdat_simul',smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_1=oo_.likelihood_at_initial_parameters; +%%Univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=3,mode_file=fs2000_corr_ME_mode,mode_compute=0,order=1,datafile='../../fs2000/fsdat_simul',smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_3=oo_.likelihood_at_initial_parameters; +%%Diffuse Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=2,mode_file=fs2000_corr_ME_mode,mode_compute=0,datafile='../../fs2000/fsdat_simul',smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_2=oo_.likelihood_at_initial_parameters; +%%Diffuse univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=4,mode_file=fs2000_corr_ME_mode,mode_compute=0,datafile='../../fs2000/fsdat_simul',smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_4=oo_.likelihood_at_initial_parameters; diff --git a/tests/kalman/likelihood_from_dynare/fs2000_corr_ME_missing.mod b/tests/kalman/likelihood_from_dynare/fs2000_corr_ME_missing.mod new file mode 100644 index 000000000..eea44b004 --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fs2000_corr_ME_missing.mod @@ -0,0 +1,137 @@ +/* + * This file is based on the cash in advance model described + * Frank Schorfheide (2000): "Loss function-based evaluation of DSGE models", + * Journal of Applied Econometrics, 15(6), 645-670. + * + * The equations are taken from J. Nason and T. Cogley (1994): "Testing the + * implications of long-run neutrality for monetary business cycle models", + * Journal of Applied Econometrics, 9, S37-S70. + * Note that there is an initial minus sign missing in equation (A1), p. S63. + * + * This implementation was written by Michel Juillard. Please note that the + * following copyright notice only applies to this Dynare implementation of the + * model. + */ + +/* + * Copyright (C) 2004-2013 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 . + */ + +var m P c e W R k d n l gy_obs gp_obs y dA; +varexo e_a e_m; + +parameters alp bet gam mst rho psi del theta; + +alp = 0.33; +bet = 0.99; +gam = 0.003; +mst = 1.011; +rho = 0.7; +psi = 0.787; +del = 0.02; +theta=0; + +model; +dA = exp(gam+e_a); +log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m; +-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0; +W = l/n; +-(psi/(1-psi))*(c*P/(1-n))+l/n = 0; +R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W; +1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0; +c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1); +P*c = m; +m-1+d = l; +e = exp(e_a); +y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a)); +gy_obs = dA*y/y(-1); +gp_obs = (P/P(-1))*m(-1)/dA; +end; + +steady_state_model; + dA = exp(gam); + gst = 1/dA; + m = mst; + khst = ( (1-gst*bet*(1-del)) / (alp*gst^alp*bet) )^(1/(alp-1)); + xist = ( ((khst*gst)^alp - (1-gst*(1-del))*khst)/mst )^(-1); + nust = psi*mst^2/( (1-alp)*(1-psi)*bet*gst^alp*khst^alp ); + n = xist/(nust+xist); + P = xist + nust; + k = khst*n; + + l = psi*mst*n/( (1-psi)*(1-n) ); + c = mst/P; + d = l - mst + 1; + y = k^alp*n^(1-alp)*gst^alp; + R = mst/bet; + W = l/n; + ist = y-c; + q = 1 - d; + + e = 1; + + gp_obs = m/dA; + gy_obs = dA; +end; + +varobs gp_obs gy_obs; + +shocks; +var e_a; stderr 0.014; +var e_m; stderr 0.005; +corr gy_obs,gp_obs = 0.5; +end; + +steady; + + +estimated_params; +alp, 0.356; +gam, 0.0085; +del, 0.01; +stderr e_a, 0.035449; +stderr e_m, 0.008862; +corr e_m, e_a, 0; +stderr gp_obs, 1; +stderr gy_obs, 1; +corr gp_obs, gy_obs,0; +end; + +options_.TeX=1; +options_.debug=1; + +%%default +options_.lik_init=1; +estimation(kalman_algo=0,mode_compute=4,order=1,datafile=fsdat_simul_corr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_0=oo_.likelihood_at_initial_parameters; +%%Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=1,mode_file=fs2000_corr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_corr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_1=oo_.likelihood_at_initial_parameters; +%%Univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=3,mode_file=fs2000_corr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_corr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_3=oo_.likelihood_at_initial_parameters; +%%Diffuse Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=2,mode_file=fs2000_corr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_corr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_2=oo_.likelihood_at_initial_parameters; +%%Diffuse univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=4,mode_file=fs2000_corr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_corr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_4=oo_.likelihood_at_initial_parameters; diff --git a/tests/kalman/likelihood_from_dynare/fs2000_uncorr_ME.mod b/tests/kalman/likelihood_from_dynare/fs2000_uncorr_ME.mod new file mode 100644 index 000000000..f0eb4dea4 --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fs2000_uncorr_ME.mod @@ -0,0 +1,137 @@ +/* + * This file is based on the cash in advance model described + * Frank Schorfheide (2000): "Loss function-based evaluation of DSGE models", + * Journal of Applied Econometrics, 15(6), 645-670. + * + * The equations are taken from J. Nason and T. Cogley (1994): "Testing the + * implications of long-run neutrality for monetary business cycle models", + * Journal of Applied Econometrics, 9, S37-S70. + * Note that there is an initial minus sign missing in equation (A1), p. S63. + * + * This implementation was written by Michel Juillard. Please note that the + * following copyright notice only applies to this Dynare implementation of the + * model. + */ + +/* + * Copyright (C) 2004-2013 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 . + */ + +var m P c e W R k d n l gy_obs gp_obs y dA; +varexo e_a e_m; + +parameters alp bet gam mst rho psi del theta; + +alp = 0.33; +bet = 0.99; +gam = 0.003; +mst = 1.011; +rho = 0.7; +psi = 0.787; +del = 0.02; +theta=0; + +model; +dA = exp(gam+e_a); +log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m; +-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0; +W = l/n; +-(psi/(1-psi))*(c*P/(1-n))+l/n = 0; +R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W; +1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0; +c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1); +P*c = m; +m-1+d = l; +e = exp(e_a); +y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a)); +gy_obs = dA*y/y(-1); +gp_obs = (P/P(-1))*m(-1)/dA; +end; + +steady_state_model; + dA = exp(gam); + gst = 1/dA; + m = mst; + khst = ( (1-gst*bet*(1-del)) / (alp*gst^alp*bet) )^(1/(alp-1)); + xist = ( ((khst*gst)^alp - (1-gst*(1-del))*khst)/mst )^(-1); + nust = psi*mst^2/( (1-alp)*(1-psi)*bet*gst^alp*khst^alp ); + n = xist/(nust+xist); + P = xist + nust; + k = khst*n; + + l = psi*mst*n/( (1-psi)*(1-n) ); + c = mst/P; + d = l - mst + 1; + y = k^alp*n^(1-alp)*gst^alp; + R = mst/bet; + W = l/n; + ist = y-c; + q = 1 - d; + + e = 1; + + gp_obs = m/dA; + gy_obs = dA; +end; + +varobs gp_obs gy_obs; + +shocks; +var e_a; stderr 0.014; +var e_m; stderr 0.005; +end; + +steady; + +//stoch_simul(periods=200, order=1); +//datatomfile('fsdat_simul_uncorr_ME', char('gy_obs', 'gp_obs')); + +estimated_params; +alp, 0.356; +gam, 0.0085; +del, 0.01; +stderr e_a, 0.035449; +stderr e_m, 0.008862; +corr e_m, e_a, 0; +stderr gp_obs, 1; +stderr gy_obs, 1; +//corr gp_obs, gy_obs,0; +end; + +options_.TeX=1; +options_.debug=1; + +%%default +estimation(kalman_algo=0,mode_compute=4,order=1,datafile=fsdat_simul_uncorr_ME,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_0=oo_.likelihood_at_initial_parameters; +%%Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=1,mode_file=fs2000_uncorr_ME_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_1=oo_.likelihood_at_initial_parameters; +%%Univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=3,mode_file=fs2000_uncorr_ME_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_3=oo_.likelihood_at_initial_parameters; +%%Diffuse Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=2,mode_file=fs2000_uncorr_ME_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_2=oo_.likelihood_at_initial_parameters; +%%Diffuse univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=4,mode_file=fs2000_uncorr_ME_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_4=oo_.likelihood_at_initial_parameters; diff --git a/tests/kalman/likelihood_from_dynare/fs2000_uncorr_ME_missing.mod b/tests/kalman/likelihood_from_dynare/fs2000_uncorr_ME_missing.mod new file mode 100644 index 000000000..7ed016b2d --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fs2000_uncorr_ME_missing.mod @@ -0,0 +1,137 @@ +/* + * This file is based on the cash in advance model described + * Frank Schorfheide (2000): "Loss function-based evaluation of DSGE models", + * Journal of Applied Econometrics, 15(6), 645-670. + * + * The equations are taken from J. Nason and T. Cogley (1994): "Testing the + * implications of long-run neutrality for monetary business cycle models", + * Journal of Applied Econometrics, 9, S37-S70. + * Note that there is an initial minus sign missing in equation (A1), p. S63. + * + * This implementation was written by Michel Juillard. Please note that the + * following copyright notice only applies to this Dynare implementation of the + * model. + */ + +/* + * Copyright (C) 2004-2013 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 . + */ + +var m P c e W R k d n l gy_obs gp_obs y dA; +varexo e_a e_m; + +parameters alp bet gam mst rho psi del theta; + +alp = 0.33; +bet = 0.99; +gam = 0.003; +mst = 1.011; +rho = 0.7; +psi = 0.787; +del = 0.02; +theta=0; + +model; +dA = exp(gam+e_a); +log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m; +-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0; +W = l/n; +-(psi/(1-psi))*(c*P/(1-n))+l/n = 0; +R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W; +1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0; +c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1); +P*c = m; +m-1+d = l; +e = exp(e_a); +y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a)); +gy_obs = dA*y/y(-1); +gp_obs = (P/P(-1))*m(-1)/dA; +end; + +steady_state_model; + dA = exp(gam); + gst = 1/dA; + m = mst; + khst = ( (1-gst*bet*(1-del)) / (alp*gst^alp*bet) )^(1/(alp-1)); + xist = ( ((khst*gst)^alp - (1-gst*(1-del))*khst)/mst )^(-1); + nust = psi*mst^2/( (1-alp)*(1-psi)*bet*gst^alp*khst^alp ); + n = xist/(nust+xist); + P = xist + nust; + k = khst*n; + + l = psi*mst*n/( (1-psi)*(1-n) ); + c = mst/P; + d = l - mst + 1; + y = k^alp*n^(1-alp)*gst^alp; + R = mst/bet; + W = l/n; + ist = y-c; + q = 1 - d; + + e = 1; + + gp_obs = m/dA; + gy_obs = dA; +end; + +varobs gp_obs gy_obs; + +shocks; +var e_a; stderr 0.014; +var e_m; stderr 0.005; +end; + +steady; + +//stoch_simul(periods=200, order=1); +//datatomfile('fsdat_simul_uncorr_ME', char('gy_obs', 'gp_obs')); + +estimated_params; +alp, 0.356; +gam, 0.0085; +del, 0.01; +stderr e_a, 0.035449; +stderr e_m, 0.008862; +corr e_m, e_a, 0; +stderr gp_obs, 1; +stderr gy_obs, 1; +//corr gp_obs, gy_obs,0; +end; + +options_.TeX=1; +options_.debug=1; + +%%default +estimation(kalman_algo=0,mode_compute=4,order=1,datafile=fsdat_simul_uncorr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_0=oo_.likelihood_at_initial_parameters; +%%Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=1,mode_file=fs2000_uncorr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_1=oo_.likelihood_at_initial_parameters; +%%Univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=3,mode_file=fs2000_uncorr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_3=oo_.likelihood_at_initial_parameters; +%%Diffuse Multivariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=2,mode_file=fs2000_uncorr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_2=oo_.likelihood_at_initial_parameters; +%%Diffuse univariate Kalman Filter +options_.lik_init=1; +estimation(kalman_algo=4,mode_file=fs2000_uncorr_ME_missing_mode,mode_compute=0,order=1,datafile=fsdat_simul_uncorr_ME_missing,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y gy_obs; +fval_algo_4=oo_.likelihood_at_initial_parameters; diff --git a/tests/kalman/likelihood_from_dynare/fsdat_simul_corr_ME_missing.m b/tests/kalman/likelihood_from_dynare/fsdat_simul_corr_ME_missing.m new file mode 100644 index 000000000..d04535f02 --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fsdat_simul_corr_ME_missing.m @@ -0,0 +1,416 @@ +% Generated data, used by fs2000.mod + +gy_obs =[ + 1.0030045 + 1.0002599 + 0.99104664 + 1.0321162 + 1.0223545 + 1.0043614 + NaN + 1.0092127 + 1.0357197 + 1.0150827 + 1.0051548 + 0.98465775 + 0.99132132 + 0.99904153 + 1.0044641 + 1.0179198 + 1.0113462 + 0.99409421 + 0.99904293 + 1.0448336 + 0.99932433 + 1.0057004 + 0.99619787 + 1.0267504 + 1.0077645 + 1.0058026 + 1.0025891 + 0.9939097 + 0.99604693 + 0.99908569 + 1.0151094 + 0.99348134 + 1.0039124 + 1.0145805 + 0.99800868 + 0.98578138 + 1.0065771 + 0.99843919 + 0.97979062 + 0.98413351 + 0.96468174 + 1.0273857 + 1.0225211 + 0.99958667 + 1.0111157 + 1.0099585 + 0.99480311 + 1.0079265 + 0.98924573 + 1.0070613 + 1.0075706 + 0.9937151 + 1.0224711 + 1.0018891 + 0.99051863 + 1.0042944 + 1.0184055 + 0.99419508 + 0.99756624 + 1.0015983 + 0.9845772 + 1.0004407 + 1.0116237 + 0.9861885 + 1.0073094 + 0.99273355 + 1.0013224 + 0.99777979 + 1.0301686 + 0.96809556 + 0.99917088 + 0.99949253 + 0.96590004 + 1.0083938 + 0.96662298 + 1.0221454 + 1.0069792 + 1.0343996 + 1.0066531 + 1.0072525 + 0.99743563 + 0.99723703 + 1.000372 + 0.99013917 + 1.0095223 + 0.98864268 + 0.98092242 + 0.98886488 + 1.0030341 + 1.01894 + 0.99155059 + 0.99533235 + 0.99734316 + 1.0047356 + 1.0082737 + 0.98425116 + 0.99949212 + 1.0055899 + 1.0065075 + 0.99385069 + 0.98867975 + 0.99804843 + 1.0184038 + 0.99301902 + 1.0177222 + 1.0051924 + 1.0187852 + 1.0098985 + 1.0097172 + 1.0145811 + 0.98721038 + 1.0361722 + 1.0105821 + 0.99469309 + 0.98626785 + 1.013871 + 0.99858924 + 0.99302637 + 1.0042186 + 0.99623745 + 0.98545708 + 1.0225435 + 1.0011861 + 1.0130321 + 0.97861347 + 1.0228193 + 0.99627435 + 1.0272779 + 1.0075172 + 1.0096762 + 1.0129306 + 0.99966549 + 1.0262882 + 1.0026914 + 1.0061475 + 1.009523 + 1.0036127 + 0.99762992 + 0.99092634 + 1.0058469 + 0.99887292 + 1.0060653 + 0.98673557 + 0.98895709 + 0.99111967 + 0.990118 + 0.99788054 + 0.97054709 + 1.0099157 + 1.0107431 + 0.99518695 + 1.0114048 + 0.99376019 + 1.0023369 + 0.98783327 + 1.0051727 + 1.0100462 + 0.98607387 + 1.0000064 + 0.99692442 + 1.012225 + 0.99574078 + 0.98642833 + 0.99008207 + 1.0197359 + 1.0112849 + 0.98711069 + 0.99402748 + 1.0242141 + 1.0135349 + 0.99842505 + 1.0130714 + 0.99887044 + 1.0059058 + 1.0185998 + 1.0073314 + 0.98687706 + 1.0084551 + 0.97698964 + 0.99482714 + 1.0015302 + 1.0105331 + 1.0261767 + 1.0232822 + 1.0084176 + 0.99785167 + 0.99619733 + 1.0055223 + 1.0076326 + 0.99205461 + 1.0030587 + 1.0137012 + 1.0145878 + 1.0190297 + 1.0000681 + 1.0153894 + 1.0140649 + 1.0007236 + 0.97961463 + 1.0125257 + 1.0169503 + NaN + 1.0221185 + +]; + +gp_obs =[ + 1.0079715 + 1.0115853 + 1.0167502 + 1.0068957 + 1.0138189 + 1.0258364 + 1.0243817 + 1.017373 + 1.0020171 + 1.0003742 + 1.0008974 + 1.0104804 + 1.0116393 + 1.0114294 + 0.99932124 + 0.99461459 + NaN + 1.0051446 + 1.020639 + 1.0051964 + 1.0093042 + 1.007068 + 1.01086 + 0.99590086 + 1.0014883 + 1.0117332 + 0.9990095 + 1.0108284 + 1.0103672 + 1.0036722 + 1.0005124 + 1.0190331 + 1.0130978 + 1.007842 + 1.0285436 + 1.0322054 + 1.0213403 + 1.0246486 + 1.0419306 + 1.0258867 + 1.0156316 + 0.99818589 + 0.9894107 + 1.0127584 + 1.0146882 + 1.0136529 + 1.0340107 + 1.0343652 + 1.02971 + 1.0077932 + 1.0198114 + 1.013971 + 1.0061083 + 1.0089573 + 1.0037926 + 1.0082071 + 0.99498155 + 0.99735772 + 0.98765026 + 1.006465 + 1.0196088 + 1.0053233 + 1.0119974 + 1.0188066 + 1.0029302 + 1.0183459 + 1.0034218 + 1.0158799 + 0.98824798 + 1.0274357 + 1.0168832 + 1.0180641 + 1.0294657 + 0.98864091 + 1.0358326 + 0.99889969 + 1.0178322 + 0.99813566 + 1.0073549 + 1.0215985 + 1.0084245 + 1.0080939 + 1.0157021 + 1.0075815 + 1.0032633 + 1.0117871 + 1.0209276 + 1.0077569 + 0.99680958 + 1.0120266 + 1.0017625 + 1.0138811 + 1.0198358 + 1.0059629 + 1.0115416 + 1.0319473 + 1.0167074 + 1.0116111 + 1.0048627 + 1.0217622 + 1.0125221 + 1.0142045 + 0.99792469 + 0.99823971 + 0.99561547 + 0.99850373 + 0.9898464 + 1.0030963 + 1.0051373 + 1.0004213 + 1.0144117 + 0.97185592 + 0.9959518 + 1.0073529 + 1.0051603 + 0.98642572 + 0.99433423 + 1.0112131 + 1.0007695 + 1.0176867 + 1.0134363 + 0.99926191 + 0.99879835 + 0.99878754 + 1.0331374 + 1.0077797 + 1.0127221 + 1.0047393 + 1.0074106 + 0.99784213 + 1.0056495 + 1.0057708 + 0.98817494 + 0.98742176 + 0.99930555 + 1.0000687 + 1.0129754 + 1.009529 + 1.0226731 + 1.0149534 + 1.0164295 + 1.0239469 + 1.0293458 + 1.026199 + 1.0197525 + 1.0126818 + 1.0054473 + 1.0254423 + 1.0069461 + 1.0153135 + 1.0337515 + 1.0178187 + 1.0240469 + 1.0079489 + 1.0186953 + 1.0008628 + 1.0113799 + 1.0140118 + 1.0168007 + 1.011441 + 0.98422774 + 0.98909729 + 1.0157859 + 1.0151586 + 0.99756232 + 0.99497777 + 1.0102841 + 1.0221659 + 0.9937759 + 0.99877193 + 1.0079433 + 0.99667692 + 1.0095959 + 1.0128804 + 1.0156949 + 1.0111951 + 1.0228887 + 1.0122083 + 1.0190197 + 1.0074927 + 1.0268096 + 0.99689352 + 0.98948474 + 1.0024938 + 1.0105543 + 1.014116 + 1.0141217 + 1.0056504 + 1.0101026 + 1.0105069 + 0.99619053 + 1.0059439 + 0.99449473 + 0.99482458 + 1.0037702 + 1.0068087 + 0.99575975 + 1.0030815 + 1.0334014 + 0.99879386 + 0.99625634 + NaN + 0.99233844 + +]; + diff --git a/tests/kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME.m b/tests/kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME.m new file mode 100644 index 000000000..ac2057747 --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME.m @@ -0,0 +1,406 @@ +gy_obs = [ + 1.0089434 + 0.97436837 + 1.0078602 + 0.99728812 + 1.0469033 + 0.98514927 + 1.0130718 + 1.0127905 + 1.0012276 + 1.0207597 + 1.0128382 + 1.0117555 + 1.0093849 + 1.0130848 + 1.0076878 + 1.0151714 + 0.99010787 + 0.9651508 + 1.0075909 + 1.0211352 + 1.0016703 + 1.0067838 + 0.99211778 + 1.0006245 + 1.0164224 + 0.99240467 + 0.98841006 + 1.0021161 + 0.99328172 + 0.99975511 + 0.9894502 + 1.0095403 + 1.0227135 + 0.98474395 + 0.98842474 + 0.99510226 + 1.0003444 + 0.99419424 + 0.98459156 + 1.0006202 + 1.0205175 + 1.0046155 + 0.99261284 + 1.0137533 + 1.0062878 + 0.9882274 + 1.0106694 + 0.99471339 + 1.002965 + 0.99799306 + 1.0134657 + 1.0039836 + 1.0066297 + 1.0084027 + 1.0245597 + 0.9763106 + 1.02033 + 1.0147916 + 1.0219756 + 1.0008057 + 1.0390485 + 1.0210931 + 0.99734292 + 1.0171595 + 1.0130723 + 0.99892094 + 0.98279425 + 1.0066122 + 0.98575479 + 1.0078347 + 1.0036275 + 0.98116612 + 0.99293813 + 0.98842048 + 0.97690963 + 1.0093542 + 1.0027533 + 1.0156453 + 0.99313547 + 1.0004421 + 0.99954572 + 0.98736251 + 1.0238741 + 0.98768174 + 1.000261 + 0.98722229 + 0.98398124 + 1.0072403 + 1.0009303 + 0.99628068 + 0.98538731 + 0.99218841 + 1.009118 + 0.98810044 + 1.0112788 + 1.0004868 + 0.99890858 + 1.0029751 + 1.0219324 + 1.000043 + 1.0058196 + 1.014664 + 1.0044 + 1.0067619 + 1.0008676 + 0.99532428 + 0.99224953 + 0.99444046 + 1.0003366 + 1.0221002 + 0.98855273 + 1.0187089 + 0.98416472 + 1.0006988 + 0.99933767 + 1.0084427 + 0.9910711 + 0.99630044 + 1.0041385 + 0.99578819 + 0.98859148 + 1.0071189 + 1.0057602 + 1.0006798 + 1.0040692 + 0.99357917 + 1.0055212 + 0.99826781 + 1.0294402 + 1.0306182 + 1.0163397 + 0.99544135 + 1.0089258 + 1.0091866 + 1.0031688 + 1.0065311 + 1.0162032 + 1.006835 + 0.98588242 + 1.0031649 + 1.0143694 + 1.0071297 + 1.0151235 + 0.9950707 + 1.0190895 + 1.0036681 + 1.0153039 + 1.0031942 + 0.97146165 + 0.97647363 + 1.0040287 + 1.0074315 + 1.0139851 + 1.0084592 + 1.013138 + 1.0145484 + 1.0008416 + 0.97670707 + 0.98714692 + 1.0106897 + 0.99046031 + 1.0047648 + 1.0064955 + 0.99129614 + 0.97537098 + 0.99551383 + 1.0000844 + 1.0024149 + 1.0294317 + 0.97109038 + 0.98520755 + 1.0043495 + 0.99858097 + 0.99675168 + 1.0112858 + 1.0139032 + 0.99391287 + 1.026372 + 1.003729 + 1.0047125 + 0.99687991 + 1.0218475 + 1.0125423 + 1.0004841 + 0.9992396 + 0.98416027 + 0.99245649 + 1.007447 + 0.99887608 + 0.97908706 + 1.0330613 + 0.99897339 + 0.9891648 + 1.0148316 + 0.99214239 + 0.99630795 + 0.99629966 + 0.99137359 + 0.9809149 + 1.0008659 +]; + +gp_obs = [ + 1.0193403 + 1.0345762 + 1.0011701 + 1.0147224 + 1.008392 + 1.0488327 + 1.0153551 + 1.0099775 + 1.0260561 + 1.0172218 + 1.0014374 + 1.0184572 + 1.0179988 + 1.0060339 + 1.0019537 + 0.99179578 + 1.004346 + 1.0345153 + 1.0004432 + 0.98327074 + 1.0007585 + 1.0034378 + 1.010532 + 1.0121367 + 1.0097161 + 1.0166682 + 1.0089513 + 1.0194821 + 1.0192704 + 1.0220258 + 1.020915 + 1.0176156 + 1.0040707 + 1.0157694 + 1.0357484 + 1.0256259 + 1.0240583 + 1.0095152 + 1.0241605 + 1.0115295 + 1.003636 + 1.0222399 + 1.0250969 + 1.0068969 + 1.0009829 + 1.0166179 + 1.0252018 + 1.0211178 + 0.99867851 + 0.99594002 + 0.9908135 + 0.99762919 + 0.99616309 + 1.0058679 + 0.99323315 + 1.0132879 + 0.98718922 + 0.99739822 + 0.97858593 + 0.99128769 + 0.98624299 + 0.98447966 + 1.0013312 + 0.99189504 + 0.98032699 + 0.99332035 + 1.0129565 + 1.0007785 + 1.0218292 + 1.0030419 + 1.0044453 + 1.0156181 + 1.0040112 + 1.0081137 + 1.0261598 + 1.0053686 + 1.0024674 + 0.99883223 + 1.0224791 + 1.0074723 + 1.0037807 + 1.0348866 + 1.0053664 + 1.0140072 + 1.017359 + 1.0013916 + 1.017887 + 1.008987 + 1.011771 + 1.0201455 + 1.0249464 + 1.0159166 + 1.0162718 + 1.0312397 + 1.0108745 + 1.0132205 + 1.0142484 + 1.0178907 + 1.0065039 + 1.0190304 + 1.0034406 + 1.0053556 + 1.012823 + 1.0009983 + 1.0073148 + 1.0247254 + 1.0140215 + 1.0053603 + 1.006169 + 0.994725 + 1.026685 + 1.0012279 + 1.0160733 + 1.0119851 + 1.0148392 + 0.99760076 + 1.0070377 + 1.0066215 + 0.98130614 + 1.0127043 + 1.0203824 + 1.0067477 + 0.99510728 + 1.0188472 + 1.0100108 + 1.0146874 + 1.0118012 + 1.0111904 + 0.97759194 + 0.99081872 + 0.98425915 + 1.0026496 + 0.98587189 + 0.98648329 + 1.0035766 + 1.0094743 + 0.99460644 + 0.9953724 + 1.0194434 + 1.0065039 + 1.0056522 + 1.0160367 + 1.006524 + 1.0092492 + 0.9864426 + 0.98723638 + 0.9994522 + 1.0026778 + 1.025553 + 1.0030477 + 0.99411719 + 1.0045087 + 0.99375289 + 1.0017609 + 1.0039766 + 0.99976299 + 1.0155671 + 1.0192975 + 1.0135507 + 1.0099869 + 1.0125994 + 1.0050808 + 1.0088531 + 1.0135256 + 1.0322097 + 1.0065808 + 0.99857525 + 1.0008792 + 0.9997691 + 1.02875 + 1.0177818 + 1.0150152 + 1.026416 + 1.0209804 + 1.010633 + 1.009636 + 1.0028257 + 0.9896666 + 1.0094002 + 0.99958414 + 1.0077797 + 0.98933606 + 1.0014885 + 0.99875283 + 1.005051 + 1.016385 + 1.0116282 + 0.99774103 + 1.0101802 + 1.0281101 + 1.0024654 + 1.0174549 + 1.0342738 + 1.0160514 + 1.0238797 + 1.015626 + 1.0075839 + 1.0132055 + 1.0250766 + 1.0175022 +]; + diff --git a/tests/kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME_missing.m b/tests/kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME_missing.m new file mode 100644 index 000000000..9cd0c8da1 --- /dev/null +++ b/tests/kalman/likelihood_from_dynare/fsdat_simul_uncorr_ME_missing.m @@ -0,0 +1,406 @@ +gy_obs = [ + 1.0089434 + 0.97436837 + 1.0078602 + 0.99728812 + 1.0469033 + NaN + 1.0130718 + 1.0127905 + 1.0012276 + 1.0207597 + 1.0128382 + 1.0117555 + 1.0093849 + 1.0130848 + 1.0076878 + 1.0151714 + 0.99010787 + 0.9651508 + 1.0075909 + 1.0211352 + 1.0016703 + 1.0067838 + 0.99211778 + 1.0006245 + 1.0164224 + 0.99240467 + 0.98841006 + 1.0021161 + 0.99328172 + 0.99975511 + 0.9894502 + 1.0095403 + 1.0227135 + 0.98474395 + 0.98842474 + 0.99510226 + 1.0003444 + 0.99419424 + 0.98459156 + 1.0006202 + 1.0205175 + 1.0046155 + 0.99261284 + 1.0137533 + 1.0062878 + 0.9882274 + 1.0106694 + 0.99471339 + 1.002965 + 0.99799306 + 1.0134657 + 1.0039836 + 1.0066297 + 1.0084027 + 1.0245597 + 0.9763106 + 1.02033 + 1.0147916 + 1.0219756 + 1.0008057 + 1.0390485 + 1.0210931 + 0.99734292 + 1.0171595 + 1.0130723 + 0.99892094 + 0.98279425 + 1.0066122 + 0.98575479 + 1.0078347 + 1.0036275 + 0.98116612 + 0.99293813 + 0.98842048 + 0.97690963 + 1.0093542 + 1.0027533 + 1.0156453 + 0.99313547 + 1.0004421 + 0.99954572 + 0.98736251 + 1.0238741 + 0.98768174 + 1.000261 + 0.98722229 + 0.98398124 + 1.0072403 + 1.0009303 + 0.99628068 + 0.98538731 + 0.99218841 + 1.009118 + 0.98810044 + 1.0112788 + 1.0004868 + 0.99890858 + 1.0029751 + 1.0219324 + 1.000043 + 1.0058196 + 1.014664 + 1.0044 + 1.0067619 + 1.0008676 + 0.99532428 + 0.99224953 + 0.99444046 + 1.0003366 + 1.0221002 + 0.98855273 + 1.0187089 + 0.98416472 + 1.0006988 + 0.99933767 + 1.0084427 + 0.9910711 + 0.99630044 + 1.0041385 + 0.99578819 + 0.98859148 + 1.0071189 + 1.0057602 + 1.0006798 + 1.0040692 + 0.99357917 + 1.0055212 + 0.99826781 + 1.0294402 + 1.0306182 + 1.0163397 + 0.99544135 + 1.0089258 + 1.0091866 + 1.0031688 + 1.0065311 + 1.0162032 + 1.006835 + 0.98588242 + 1.0031649 + 1.0143694 + 1.0071297 + 1.0151235 + 0.9950707 + 1.0190895 + 1.0036681 + 1.0153039 + 1.0031942 + 0.97146165 + 0.97647363 + 1.0040287 + 1.0074315 + 1.0139851 + 1.0084592 + 1.013138 + 1.0145484 + 1.0008416 + 0.97670707 + 0.98714692 + 1.0106897 + 0.99046031 + 1.0047648 + 1.0064955 + 0.99129614 + 0.97537098 + 0.99551383 + 1.0000844 + 1.0024149 + 1.0294317 + 0.97109038 + 0.98520755 + 1.0043495 + 0.99858097 + 0.99675168 + 1.0112858 + 1.0139032 + 0.99391287 + 1.026372 + 1.003729 + 1.0047125 + 0.99687991 + 1.0218475 + 1.0125423 + 1.0004841 + 0.9992396 + 0.98416027 + 0.99245649 + 1.007447 + 0.99887608 + 0.97908706 + 1.0330613 + 0.99897339 + 0.9891648 + 1.0148316 + 0.99214239 + 0.99630795 + 0.99629966 + 0.99137359 + NaN + 1.0008659 +]; + +gp_obs = [ + 1.0193403 + 1.0345762 + 1.0011701 + 1.0147224 + 1.008392 + 1.0488327 + 1.0153551 + 1.0099775 + 1.0260561 + 1.0172218 + 1.0014374 + 1.0184572 + 1.0179988 + 1.0060339 + 1.0019537 + 0.99179578 + 1.004346 + 1.0345153 + 1.0004432 + 0.98327074 + 1.0007585 + 1.0034378 + 1.010532 + 1.0121367 + 1.0097161 + 1.0166682 + 1.0089513 + 1.0194821 + 1.0192704 + 1.0220258 + 1.020915 + 1.0176156 + 1.0040707 + 1.0157694 + 1.0357484 + 1.0256259 + 1.0240583 + 1.0095152 + 1.0241605 + 1.0115295 + 1.003636 + 1.0222399 + 1.0250969 + 1.0068969 + 1.0009829 + 1.0166179 + 1.0252018 + 1.0211178 + 0.99867851 + 0.99594002 + 0.9908135 + 0.99762919 + 0.99616309 + 1.0058679 + 0.99323315 + 1.0132879 + 0.98718922 + 0.99739822 + 0.97858593 + 0.99128769 + 0.98624299 + 0.98447966 + 1.0013312 + 0.99189504 + 0.98032699 + 0.99332035 + 1.0129565 + 1.0007785 + 1.0218292 + 1.0030419 + 1.0044453 + 1.0156181 + 1.0040112 + 1.0081137 + 1.0261598 + 1.0053686 + 1.0024674 + 0.99883223 + 1.0224791 + 1.0074723 + 1.0037807 + 1.0348866 + 1.0053664 + 1.0140072 + 1.017359 + 1.0013916 + 1.017887 + 1.008987 + 1.011771 + 1.0201455 + 1.0249464 + 1.0159166 + 1.0162718 + 1.0312397 + 1.0108745 + 1.0132205 + 1.0142484 + 1.0178907 + 1.0065039 + 1.0190304 + 1.0034406 + 1.0053556 + 1.012823 + 1.0009983 + 1.0073148 + 1.0247254 + 1.0140215 + 1.0053603 + 1.006169 + 0.994725 + 1.026685 + 1.0012279 + 1.0160733 + 1.0119851 + 1.0148392 + 0.99760076 + 1.0070377 + 1.0066215 + 0.98130614 + 1.0127043 + 1.0203824 + 1.0067477 + 0.99510728 + 1.0188472 + 1.0100108 + 1.0146874 + 1.0118012 + 1.0111904 + 0.97759194 + 0.99081872 + 0.98425915 + 1.0026496 + 0.98587189 + 0.98648329 + 1.0035766 + 1.0094743 + 0.99460644 + 0.9953724 + 1.0194434 + 1.0065039 + 1.0056522 + 1.0160367 + 1.006524 + 1.0092492 + 0.9864426 + 0.98723638 + 0.9994522 + 1.0026778 + 1.025553 + 1.0030477 + NaN + 1.0045087 + 0.99375289 + 1.0017609 + 1.0039766 + 0.99976299 + 1.0155671 + 1.0192975 + 1.0135507 + 1.0099869 + 1.0125994 + 1.0050808 + 1.0088531 + 1.0135256 + 1.0322097 + 1.0065808 + 0.99857525 + 1.0008792 + 0.9997691 + 1.02875 + 1.0177818 + 1.0150152 + 1.026416 + 1.0209804 + 1.010633 + 1.009636 + 1.0028257 + 0.9896666 + 1.0094002 + 0.99958414 + 1.0077797 + 0.98933606 + 1.0014885 + 0.99875283 + 1.005051 + 1.016385 + 1.0116282 + 0.99774103 + 1.0101802 + 1.0281101 + 1.0024654 + 1.0174549 + 1.0342738 + 1.0160514 + 1.0238797 + 1.015626 + 1.0075839 + 1.0132055 + NaN + 1.0175022 +]; +