Writing the Dynare file

This file writes and runs a Dynare file called DyTruncation.mod to simulate the Ramsey model specified in the notebook Main.ipynb. The calibration can also be found in the same file.

Initialisation

clear
clc
warning('off','all');
warning('off', 'Octave:array-as-logical');

Loading the Julia steady-state and specifying the Dynare execution

We load the matlab file todynare_Truncation.mat created Main.ipynb. We also specify what Dynare does:

Note that in each case, a formatting file is run automatically.

load ../todynare_Truncation.mat

bool_IRF = true;
bool_Trans = !bool_IRF;

Preliminary definitions

We define some variables following the file import.

Nbin  = eco.Nbin;    % total number of bins
PI    = eco.Matab;   % transition matrix 
Sp    = eco.Sp; % size of each bin
K     = eco.A-0*eco.B; % capital
L     = 1; % labor
FKK   = alpha*(alpha-1)*K^(alpha-2)*L^(1-alpha);
FLK   = alpha*(1-alpha)*K^(alpha-1)*L^(-alpha);
Y     = K^alpha*L^(1-alpha); % GDP
I     = delta*K; % investment
P            = ones(Nbin,1);
P(eco.indcc) = 0;

We specify the aggregate shock.

rho_u        = 0.95;
sigma_u      = 0.00312;

Opening the Dynare file

We open the Dynare file DyTruncation.mod that will be written and we specify the precision when writting numerical variables (formatSpec with 12 digits of precision).

fid = fopen('DyTruncation.mod', 'w') ;
formatSpec = '%16.12f';

We also define some variables depending on whether Dynare should run stochastic or perfect-foresight simulations.

if bool_IRF
    str = ['@#define IRF = true\n']; fprintf(fid, str);
    str = ['@#define TRANSITION = false\n']; fprintf(fid, str);
else
    str = ['@#define IRF = false\n']; fprintf(fid, str);
    str = ['@#define TRANSITION = true\n']; fprintf(fid, str);
end
str = ['@#define EXTENDED_PATH_SIMULATION = false\n']; fprintf(fid, str);

Defining endogenous variables

We write the definition of endogenous variables. We use a loop for history-dependent variables.

str = ['var\n'] ; fprintf(fid, str);
str = ['R w FK FL FKK FLK u Z A K GDP I C TT \n'] ; fprintf(fid, str);  % variables
str = ['Rt wt ut Zt Kt GDPt It Ct Tt \n'] ; fprintf(fid, str);

for p = 1:Nbin
    str = ['a',num2str(p),' '] ; fprintf(fid,str);     %  end-of-period
    str = ['at',num2str(p),' '] ; fprintf(fid, str);   % beginning-of-period
    str = ['c',num2str(p),' '] ; fprintf(fid, str);
    str = ['psih',num2str(p),' '] ; fprintf(fid, str);    
    str = ['lambda',num2str(p),' '] ; fprintf(fid, str);    
    str = ['lambdat',num2str(p),' '] ; fprintf(fid, str);    % last period 
    if mod(p,10)==0; % multiple of 10
        str = ['\n'] ;         fprintf(fid, str);
    end;
end;
str = ['; \n\n'] ; fprintf(fid, str);

Defining exogenous variables and parameters

We write the definition of the unique exogenous variable and of model parameters. Their value is taken from the file DyTruncation.mod.

str = ['varexo eps; \n\n'] ; fprintf(fid, str);
str = ['parameters\n'] ; fprintf(fid, str);
str = ['beta alpha theta abar delta gamma rho_u;\n\n'] ; fprintf(fid, str);

str = ['alpha','   = ',num2str(alpha,formatSpec),';\n']; fprintf(fid, str);
str = ['beta','    = ',num2str(beta,formatSpec),';\n']; fprintf(fid, str);
str = ['theta','   = ',num2str(theta,formatSpec),';\n']; fprintf(fid, str);
str = ['abar','    = ',num2str(min(eco.aep),formatSpec),';\n']; fprintf(fid, str);
str = ['delta','   = ',num2str(delta,formatSpec),';\n']; fprintf(fid, str);
str = ['gamma','   = ',num2str(gamma,formatSpec),';\n']; fprintf(fid, str);
str = ['rho_u','   = ',num2str(rho_u),';\n']; fprintf(fid, str);

Writing equations: the Dynare model part

We write the model equations. We label them using equa, which simplifies debugging. We use a threshold tol to discard very small truncated histories.

equa = 1;
tol  = 10^-18;  % threshold for transition, to avoid to consider very small transitions

str  = ['%%%%%%%%%%%%%%%%%%%%%%%%']; fprintf(fid, str);
str  = ['\n\n'] ; fprintf(fid, str);

str  = ['model;\n\n']; fprintf(fid, str);

We start with truncated-history equations and loop over all truncated histories.

for h = 1:Nbin
    nsi = eco.ytype(h);

Budget constraints

    str = ['c',num2str(h),' = ',num2str((eco.ytype(h)),formatSpec), '*w + R*','at',num2str(h),'-a', num2str(h),'- TT;']; fprintf(fid, str);
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;

Beginning-of-period wealth ãyN

    str = ['at',num2str(h),' = 10^-18 ']; fprintf(fid, str);
    for hi = 1:Nbin
        if (PI(h,hi))>tol   % hi to h
          str = ['+ ',num2str(Sp(hi)*PI(h,hi)/Sp(h),formatSpec),'*a',num2str(hi),'(-1)']; fprintf(fid, str);
        end
    end
    str = [';']; fprintf(fid, str);str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;

Euler equations

    if P(h)==1 % not constrained
        str = [num2str(eco.xsyn1(h),formatSpec),'*(c',num2str(h),')^(-gamma) =', ...
            num2str(eco.Res(h),formatSpec),' + beta*(R(+1))*(']; fprintf(fid, str);
        for hp = 1:Nbin
            if (PI(hp,h)*eco.xsyn1(hp))>tol
                str = ['+ ',num2str((PI(hp,h))*eco.xsyn1(hp),formatSpec),'*(c',num2str(hp),'(+1))^-gamma']; fprintf(fid, str);
            end
        end
        str = [');']; fprintf(fid, str);str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
    else
        str = ['a',num2str(h),' = ', num2str(eco.aep(h)/Sp(h),formatSpec),';\n']; fprintf(fid, str);
        str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
    end

Definition of the social valuation of liquidity ψ and planner’s FOC with respect to savings

    if P(h)==1
        % Unconstrainted truncated histories
        str = ['psih',num2str(h),'=',num2str(eco.xsyn1(h),formatSpec) ,'*(c',num2str(h),')^(-gamma) - (lambda',num2str(h),...
         '- R*lambdat',num2str(h),')*(',num2str(-gamma*eco.xsyn2(h),formatSpec),')*((c',num2str(h),')^(-gamma-1));'];
        fprintf(fid, str);str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
        str = ['psih',num2str(h),'  = beta*(R(+1))*(']; fprintf(fid, str);
        for hp = 1:Nbin
            if (PI(hp,h))>tol   
                str = ['+ ',num2str((PI(hp,h)),formatSpec),'*psih',num2str(hp),'(+1)']; fprintf(fid, str);
            end
        end 
        str = [') + ']; fprintf(fid, str);
        str = ['beta*(10^-18']; fprintf(fid, str);
        for hp = 1:Nbin
            str = ['+ (psih',num2str(hp),'(+1))*',num2str(Sp(hp),formatSpec),'*((at',num2str(hp),...
                '(+1))*FKK(+1) + FLK(+1)*' ,num2str((eco.ytype(hp)),formatSpec),')',]; 
            fprintf(fid, str);
        end 
        str = [') + ']; fprintf(fid, str);
        str = ['beta*FKK(+1)*(10^-18']; fprintf(fid, str);
        for hp = 1:Nbin
            str = ['+',num2str(Sp(hp),formatSpec),'*(lambdat',num2str(hp),'(+1))*' , ...
                num2str(eco.xsyn1(hp),formatSpec) ,'*(c',num2str(hp),'(+1))^-gamma']; fprintf(fid, str);
        end 
        str = [');']; fprintf(fid, str);str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;   
    else
        # Credit-constrained credit-constrained 
        str = ['psih',num2str(h),'=',num2str(eco.xsyn1(h),formatSpec) ,'*(c',num2str(h),')^(-gamma) + ',...
         'R*lambdat',num2str(h),'*(',num2str(-gamma*eco.xsyn2(h),formatSpec),')*((c',num2str(h),')^(-gamma-1));'];
        fprintf(fid, str);str = [' %% equation ', num2str(equa),'\n']; 
        fprintf(fid, str); equa = equa+1;
        str = ['lambda',num2str(h),' = ',num2str(Pl.lambda(h),formatSpec),';\n']; fprintf(fid, str);
        str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
    end

Previous-period Lagrange multiplier λ̃yN

    str = ['lambdat',num2str(h),' =']; fprintf(fid, str);
    for hi = 1:Nbin
        if PI(h,hi)>tol
          str = ['+ ',num2str((PI(h,hi))*Sp(hi)/Sp(h),formatSpec),'*lambda',num2str(hi),'(-1)']; 
          fprintf(fid, str);
        end
    end
    str = [';']; fprintf(fid, str);str = [' %% equation ', num2str(equa),'\n']; 
    fprintf(fid, str); equa = equa+1;
end

Planner’s FOC with respect to lump-sum tax T

str = ['theta*(TT)^(theta-1) =     ']; fprintf(fid, str); 
for hi = 1:Nbin
      if (Sp(hi))>tol
          str = [' +',num2str(Sp(hi),formatSpec) ,'*psih',num2str(hi) ]; fprintf(fid, str);
      end
     if mod(hi,10)==0;str = ['\n'] ;fprintf(fid, str); end;
end             
str = [';']; fprintf(fid, str);   str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;

Prices and aggregate quantities

str = ['w = (1-alpha)*(Z^(1/(1-alpha)))*((R- 1+ delta)/alpha)^(alpha/(alpha-1));']; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;   
str = ['(R-1+delta)/(alpha*Z) = (K(-1)/1)^(alpha-1);']; fprintf(fid, str); 
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['C + I + TT = GDP;'] ; fprintf(fid, str);
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1; 
str = ['Z = 1 + u;'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['u = rho_u*u(-1) + eps;'] ; fprintf(fid, str);    
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['A = K;'] ; fprintf(fid, str);
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
    str = ['\n']; fprintf(fid, str);

Derivatives of F

We specify a number of derivatives of F: FL, FK, FKL, and FKK.

str = ['FL = (1 - alpha)*Z*(K(-1)/1)^alpha;'] ; fprintf(fid, str);   
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1; 
str = ['FK = alpha*Z*(K(-1)/1)^(alpha-1);'] ; fprintf(fid, str);   
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1; 
str = ['FLK = (1-alpha)*alpha*Z*(K(-1))^(alpha-1);'] ; fprintf(fid, str);   
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1; 
str = ['FKK = alpha*(alpha-1)*Z*(K(-1))^(alpha-2);'] ; fprintf(fid, str);   
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1; 
    str = ['\n']; fprintf(fid, str);

Aggregate consumption C

We sum over oall truncated-history consumption levels.

str = ['C ='] ; fprintf(fid, str);
for h = 1:Nbin
    str = ['+',num2str(Sp(h),formatSpec),'*c',num2str(h)] ; fprintf(fid, str);
    if mod(h,5)==0; 
        str = ['\n'] ;         
        fprintf(fid, str);
    end;
end;
str = [';']; fprintf(fid, str);
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;

Investment I

str = ['I = K - (1 - delta)*K(-1);']; fprintf(fid, str); 
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;

Aggregate savings A

str = ['A = '] ; fprintf(fid, str);
for h = 1:Nbin
        str = ['+', num2str(Sp(h),formatSpec),'*a',num2str(h)] ; fprintf(fid, str);
    if mod(p,10)==10; tr = ['\n']; fprintf(fid, str); end;
end;
str = [';']; fprintf(fid, str); 
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;

Variales in deviations

str = ['ut = 100*u;'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['Kt = 100*(K/',num2str(K,formatSpec),'-1);'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['Ct = 100*(C/',num2str(Ctot,formatSpec),'-1);'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['wt = 100*(w/',num2str(eco.w,formatSpec),'-1);'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['Rt = 100*(R - ',num2str(eco.R,formatSpec),');'] ; fprintf(fid, str); 
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['Tt = 100*(TT/',num2str(eco.TT,formatSpec),'-1);'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['Zt = ut;'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['GDPt = 100*(GDP/',num2str(Y,formatSpec),'-1);'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['It = 100*(I/',num2str(I,formatSpec),'-1);'] ; fprintf(fid, str);  
    str = [' %% equation ', num2str(equa),'\n']; fprintf(fid, str); equa = equa+1;
str = ['\n'] ;         fprintf(fid, str);

Closing the model part

str = ['end;\n\n'] ; fprintf(fid, str);
str = ['%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%'] ; fprintf(fid, str);
str = ['\n\n'] ; fprintf(fid, str);

Solving the model: the Dynare steady state part

%str = ['steady_state_model;\n\n'] ; fprintf(fid, str);
str = ['initval;\n\n'] ; fprintf(fid, str);fflush(fid);

for h = 1:Nbin
    str = ['a',num2str(h),' = ',num2str(eco.aep(h)/Sp(h),formatSpec),';\n'] ; fprintf(fid, str);
    str = ['at',num2str(h),' = ',num2str(eco.abp(h)/Sp(h),formatSpec),';\n']; fprintf(fid, str);
    str = ['c',num2str(h),' = ',num2str(eco.cp(h)/Sp(h),formatSpec),';\n']; fprintf(fid, str);
    str = ['lambda',num2str(h),' = ',num2str(1*Pl.lambda(h),formatSpec),';\n']; fprintf(fid, str);    
    str = ['lambdat',num2str(h),' = ',num2str(1*Pl.lambdat(h),formatSpec),';\n']; fprintf(fid, str);    
    str = ['psih',num2str(h),' = ',num2str(Pl.psih(h),formatSpec),';\n']; fprintf(fid, str);     
end;

str = ['R = ', num2str(eco.R,formatSpec),';\n'] ; fprintf(fid, str);
str = ['w = ', num2str(eco.w,formatSpec),';\n'] ; fprintf(fid, str);
str = ['FK = ', num2str((eco.R-1+delta),formatSpec),';\n'] ; fprintf(fid, str);
str = ['FL = ', num2str(1*eco.w,formatSpec),';\n'] ; fprintf(fid, str);

str = ['FKK = ', num2str(FKK,formatSpec),';\n'] ; fprintf(fid, str);
str = ['FLK = ', num2str(FLK,formatSpec),';\n'] ; fprintf(fid, str);

str = ['K = ', num2str(K,formatSpec),';\n'] ; fprintf(fid, str);
str = ['Z = ', num2str(1),';\n'] ; fprintf(fid, str);
str = ['A = ', num2str(eco.A,formatSpec),';\n'] ; fprintf(fid, str);
str = ['u = ', num2str(0),';\n'] ; fprintf(fid, str);
str = ['GDP = ', num2str(Y,formatSpec),';\n'] ; fprintf(fid, str);

str = ['I = ', num2str(I,formatSpec),';\n'] ; fprintf(fid, str);
str = ['C = ', num2str(Ctot,formatSpec),';\n'] ; fprintf(fid, str);
str = ['TT = ', num2str(eco.TT,formatSpec),';\n'] ; fprintf(fid, str);

str = ['ut = 0;\n'] ; fprintf(fid, str);
str = ['Kt = 0;\n'] ; fprintf(fid, str);
str = ['Ct  = 0;\n'] ; fprintf(fid, str);
str = ['wt = 0;\n'] ; fprintf(fid, str);
str = ['Rt = 0;\n'] ; fprintf(fid, str);
str = ['Zt = 0;\n'] ; fprintf(fid, str);
str = ['GDPt = 0;\n'] ; fprintf(fid, str);
str = ['It = 0;\n'] ; fprintf(fid, str);
str = ['Tt = 0;\n'] ; fprintf(fid, str);

str = ['end;\n\n'] ; fprintf(fid, str); % Closing the steady-state / initval blok

The simulation part

IRFs

We simulate IRFs if bool_IRF is set to true.

str = ['@#if IRF \n\n']; fprintf(fid, str);

str = ['resid;\n\n'] ; fprintf(fid, str);
%str = ['steady(nocheck);\n\n'] ; fprintf(fid, str);
%str = ['check;\n\n'] ; fprintf(fid, str);

%str = ['options_.solve_tolf=10^-6;\n\n'] ; fprintf(fid, str);
str = ['steady(maxit = 100,solve_algo=3);\n\n'] ; fprintf(fid, str);
str = ['steady;\n\n'] ; fprintf(fid, str);

str = ['shocks;\n var eps; stderr ',num2str(sigma_u,formatSpec),';\n end;\n'] ; fprintf(fid, str);

str = ['options_.TeX=1;\n\n'] ; fprintf(fid, str);
str = ['stoch_simul (order=1,irf=200, periods=10000) Zt GDPt Ct Kt It Rt wt Tt Z GDP C K I R w TT;\n\n'] ; fprintf(fid, str);

str = ['@#endif\n\n']; fprintf(fid, str);

Transitions

We simulate a transition if bool_IRF is set to false.

str = ['@#if TRANSITION\n\n']; fprintf(fid, str);

str = ['endval;\n\n'] ; fprintf(fid, str);
for h = 1:Nbin   
    str = ['lambda',num2str(h),' = ',num2str(Pl.lambda(h),formatSpec),';\n'] ; fprintf(fid, str);         
end;
str = ['end;\n\n'] ; fprintf(fid, str);

str = ['resid;\n\n'] ; fprintf(fid, str);
%str = ['steady(nocheck);\n\n'] ; fprintf(fid, str);
%str = ['check;\n\n'] ; fprintf(fid, str);

%str = ['options_.solve_tolf=10^-6;\n\n'] ; fprintf(fid, str);
str = ['steady(maxit = 100,solve_algo=3);\n\n'] ; fprintf(fid, str);
str = ['steady;\n\n'] ; fprintf(fid, str);

str = ['perfect_foresight_setup (periods=400);\n\n']; fprintf(fid, str);
str = ['perfect_foresight_solver(robust_lin_solve);\n\n']; fprintf(fid, str);

str = ['rplot TT ;\n'] ; fprintf(fid, str);
str = ['plot (C./GDP) ;\n'] ; fprintf(fid, str);
str = ['plot (TT./GDP) ;\n'] ; fprintf(fid, str);

str = ['@#endif\n\n']; fprintf(fid, str);

Closing the file

Writing the mod file is completed and we close the mod file.

fclose(fid);

Running Dynare

We run the Dynare file we have just created.

dynare DyTruncation noclearall

Post-process files

We run the post-process files to save data and report some statistics. The two files are: post-dynare-stoch-simul.m and post-dynare-perfect-foresight.m.

if bool_IRF
    run post-dynare-stoch-simul.m
else
    run post-dynare-perfect-foresight.m
end