How to plot determinacy/indeterminacy regions?

Hello everyone.

I am trying to replicate the Gali, Lopes-Salido, Valles (2004) paper on “Rule-of-Thumb Consumers…”.

I struggle at least a week to get a figure with determinacy/indeterminacy regions, where two parameters are changing in an interval (0,1).
The mod-file is attached. Although it might be a wrong one, but the problem’s solution shouldn’t depend on it.
One possible solution could be to loop both parameters. But, in the case of looping only one i got this message:

[quote]??? Attempted to access lambdas(0.128708); index must be a positive integer or
logical.

Error in ==> rotc_noshock at 141
M_.params( 5 ) = lambdas(oo_.steady_state(2));

Error in ==> dynare at 132
evalin(‘base’,fname) ;[/quote]

Anyway, I don’t think that would be the right solution path.

So, could anybody help me to crack this?

Thanks a lot.
rotc_noshock.mod (1.68 KB)

I wonder whether my question is far too easy to answer (so that finaly I should find the solution myself) or that nobody knows how to make it. :smiley:

Anyway, I tried to do the sensitivity analysis, but got another error message:

[quote]??? Reference to non-existent field ‘maximum_exo_lag’.

Error in ==> check_model at 22
xlen = M_.maximum_exo_lag+M_.maximum_exo_lead + 1;

Error in ==> stoch_simul at 61
check_model;

Error in ==> rotc_2 at 159
info = stoch_simul(var_list_);

Error in ==> dynare at 132
evalin(‘base’,fname) ;[/quote]

I must confess that I am unexpierienced in programming with matlab, so any technical error messages confuse me to death.
Does anybody know what this message is trying to say to me? :unamused:
The modificated mod-file is attached.
rotc_2.mod (1.66 KB)

Your model does not feature exogenous shocks. How should stoch_simul then work?

Sorry, it is of corse not wise to use stoch_simul. Originally there was no shock introduced in the paper. I just overlooked this line. The particular error from above is solved. Instead I get antother. :angry:

However, is it the right way to achieve the desired graphs? And are there somewhere any examples how to get those?
An example is attached.
Bild 1.pdf (66.9 KB)

Dear wassilo,

I would like to join your problem as I am working on the same paper as you do. Did you already get any further?

I have now the model setup, but I’m also quite puzzled with the papers assumptions on the shocks, in fact, they do not point out which variable is affected by the sunspot shock, so how would you incorporate this? One option would be to assume an AR(1) technology shock, i.e. y = ak(-1)^alphah^(1-alpha) (should be yet log-linearized…)
and a = rho*a(-1)+e; but I dont think that would replicate the paper’s results. So do you have a better idea?

Next, I am still working on some kind of routine to make the plots, seems rather complicated with dynare. Do you have something that works?

I think that in your work you have some erroneous time indices for capital, make sure that k is predetermined and hence k_t -> k(-1) and k_{t+1} -> k

I’ll post my mod file so far, it doesnt work yet as I have no shocks included, please double check, I would appreciate it very much!

Sincerely,

Philipp
glsbasis.mod (4.04 KB)

Were any of you guys able to find a way to plot the determinacy/indeterminacy regions? Any help will be deeply appreciated!

Dear argfm,

I have to admit that I didnt figure out a nice way. There is the globalsensitivity package which is of help, but in this case the model has to be estimated. Otherwise I guess the only way is by calculating the eigenvalue decomposition and tackle the problem from there. So by changing the parameters, the eigenvalues change and you should check for determinacy or not. But this is just a guess, if you have a better solution then Id be glad to hear from you!

Philipp

Take a look at the replication files for Ascari/Sbordone 2014 on my homepage for an example.

I looked at the replication files. Here are two things I am confused about. Thanks!

First, what does options_.qz_criterium control? Seems it should be 1+1e-6 for the code, while 1-1e-6 for dynare_sensitivity.

Second, where is the data from for the variables in varobs? Because if for estimation, it seems to have external data for those variables, but not the case for dynare_sensitivity.

  1. From the manual

[quote] qz_criterium = DOUBLE

Value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving 1^st order problems. Default: 1.000001 (except when estimating with lik_init option equal to 1: the default is 0.999999 in that case; see section Estimation).

[/quote]

In the mod-file, I needed to set it, because I directly call resol, which requires it to be set. The default should be the same of stoch_simul and dynare_sensitivity.
2. There is no data involved here. Rather,

needs to know which the considered “observed variables” are for which some of the sensitivity computations need to be conducted

Thanks!

  1. I tried both 1.000001 and 0.999999. Appears stoch_simul works for both, but dynare _sensitivity only accepts 0.999999.

  2. How can I know which are the considered variables for a specific sensitivity analysis?

[quote=“jpfeifer”]1. From the manual

[quote] qz_criterium = DOUBLE

Value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving 1^st order problems. Default: 1.000001 (except when estimating with lik_init option equal to 1: the default is 0.999999 in that case; see section Estimation).

[/quote]

In the mod-file, I needed to set it, because I directly call resol, which requires it to be set. The default should be the same of stoch_simul and dynare_sensitivity.
2. There is no data involved here. Rather,

needs to know which the considered “observed variables” are for which some of the sensitivity computations need to be conducted[/quote]

For this particular exercise, you just need to specify some variable as observed. It does not matter which one, because none of the capabilities relying on varobs are actually triggered.

Dear Mr. Pfeifer,

I have copied your chapter 3 mod file from Gali (reduced it to 5 eq) and added a debt equation to get the indeterminacy Graphs à la Leeper (1991). However, I get the following error message:

In an assignment A(:slight_smile: = B, the number of elements in A and B must be the same.

Error in dyn_first_order_solver (line 250)
info(2) = temp’*temp;

Error in stochastic_solvers (line 267)
[dr,info] = dyn_first_order_solver(jacobia_,M_,dr,options_,task);

Error in resol (line 144)
[dr,info] = stochastic_solvers(dr,check_flag,M,options,oo);

Error in NKM_ch3_3eq_debt (line 205)
[dr,info]=resol(0,M_,options_,oo_);

Error in dynare (line 235)
evalin(‘base’,fname) ;

Which I do not understand why, since if I run the model with some pre-specified value, it works. Apparently, the set_param_value also works. I am a bit confused…

I have attached the .mod-file and qould appreciate any comments!

Best!

NKM_ch3_3eq_debt.mod (5.0 KB)

You need to add

if isempty(options_.qz_criterium)
    options_.qz_criterium = 1+1e-6;
end

before the loop

Dear Mr Pfeifer,

I have a question following these messages, following your comment and based on the code above, what is the command to plot the determinacy/indeterminacy graphs ?
The code provided by Palumenis just run the interaction part but does not plot the graphs.
Thanks in advance.
Best,
Hugo

That depends. At DSGE_mod/Ascari_Sbordone_2014.mod at master · JohannesPfeifer/DSGE_mod · GitHub
I use

    figure('Name','Figure 11: The Determinacy Region and Trend Inflation')
    contourf(phi_pi_mat,phi_y_mat,Z_plot_total,5);
    colormap(hot);

For instance if I use that code :

phi_pih_vec=linspace(0.5,2,100);
Psi_vec=linspace(0,0.5,100);
[phi_pih_mat,Psi_mat]=meshgrid(phi_pih_vec,Psi_vec);
 
info_mat=NaN(size(phi_pih_mat));
if isempty(options_.qz_criterium)
    options_.qz_criterium = 1+1e-6;
end
for phi_pih_iter=1:length(phi_pih_vec)
phi_pih_iter
    for Psi_iter=1:length(Psi_vec)
        set_param_value('phi_pih',phi_pih_mat(phi_pih_iter,Psi_iter));
        set_param_value('Psi',Psi_mat(phi_pih_iter,Psi_iter));
        [dr,info]=resol(0,M_,options_,oo_);
        info_mat(phi_pih_iter,Psi_iter)=info(1);
    end
end

Then the above code should apply.

Could you please show me how you would implement that in my code ? I do not have the Z_plot_total command in my code.
Thank you very much.
Best regards.
Hugo

Simply adjust the naming. My full code was

        info_mat=NaN(size(phi_pi_mat));
        for phi_pi_iter=1:length(phi_pi_vec)
            for phi_y_iter=1:length(phi_y_vec)
                set_param_value('phi_pi',phi_pi_mat(phi_pi_iter,phi_y_iter));
                set_param_value('phi_y',phi_y_mat(phi_pi_iter,phi_y_iter));
                [dr,info]=resol(0,M_,options_,oo_);
                info_mat(phi_pi_iter,phi_y_iter)=info(1);
            end
        end
        Z_plot=zeros(size(info_mat));
        Z_plot(info_mat==0)=1;
        figure
        contourf(phi_pi_mat,phi_y_mat,Z_plot,1)
        Z_plot_total(info_mat==0)=trend_inflation_iter;