NK-EDSGE ECB Inflation

Hello everyone,

I would greatly appreciate help on this issue I have. My NK model is based on Annicchiarico and Di Dio (2015). I want to include phie or emission in the Taylor rule. I have done that already in the code in order to after plot inflation and hopefully with the increase of phie dampen the graph. This does not happen when I plot the graphs. I don’t know what I could change for me to get the desired outcomes. Maybe someone can help!

%% Structural Parameters (quarterly calibration)
beta=0.99;                    % discount factor
alpha=0.33;                   % elasticity of production wrt capital
epsilon=6;                    % elasticity of substitution btw differentiated goods
delta=0.025;                  % depreciation rate
sigma=1;                      % relative risk aversion
phi=1;                        % inverse of Frisch elasticity
g=0.0795;                     % public spending
a=1.2480;                     % exogenous TFP in ss
eta=0.9979;                   % pollution decay
gamma1=0.45;                  % shifter in the emission function
gamma2=0;                     % concavity of emission function
theta1=0.1850;                % shifter in the abatement function
theta2=2.8;                   % convexity in the abatement function
d0=1.3950e-3;                 % constant in damage function
d1=-6.6722e-6;                % 1st order coefficient term in damage function
d2=1.4647e-8;                 % 2nd order coefficient term in damage function
d3=1;                         % damage function shifter
kappaL=19.8413;               % labor disutility shifter
erow=1.3331;                  % row emissions
tau=1;                        % emission tax

%% Steady State
pi=1;                        % inflation targeting 
rr=1/beta;                   % real interest rate
r=pi/beta;                   % nominal interest rate
q=1;                         % marginal value of investment (in terms of lambda)  
rk=1/beta-(1-delta);         % rental rate of capital

options = optimoptions('fsolve','MaxFunEvals',300000,'MaxIter',30000,'TolFun',1e-15);

Z0=[0.7780,800,0.2];                       % Initial guess
Z = fsolve(@(ZZ) find_steady_env(ZZ,eta,erow,gamma1,gamma2,tau,theta1,theta2,d0,d1,d2,d3,a,epsilon,alpha,rk,delta,kappaL,sigma,g,phi),Z0,options);

y=Z(1);                                    % output
x=Z(2);                                    % atmospheric carbon
h=Z(3);                                    % hours

e=x*(1-eta)-erow;                          % emissions
mu=(gamma1*tau/(theta1*theta2)*y^(-gamma2))^(1/(theta2-1)); % abatement
D=d3*(d0+d1*x+d2*x^2);                     % pollution damage
A=(1-D)*a;                                  % TFP 
mc=(epsilon-1)/(epsilon)-(1-gamma2)*gamma1*tau*y^(-gamma2)*(1-mu)-theta1*mu^(theta2); % marginal cost net of abatement and tax
k=alpha*mc*y/rk;                           %  capital
i=delta*k;                                 % investment
w=(1-alpha)*mc*y/h;                        % wage
c=y-i-g-y*theta1*mu^(theta2);              % consumption
nu=e/y;                                    % emission/output ratio

%% Steady state parameters
gss=g;
ass=a;
piss=pi;
rss=r;
yss=y;
xss=x;
hss=h;
tauss=tau;

%% Parameters not affecting the ss
phipi=3;                
phiy=0;               
kappaI=2.48;              
rhoa=0.95;                 
rhog=0.97;                
rhom=0.15;                 
calvo=0.75;               
% adjusment cost coefficient to have the same linear Phillips Curve of the Calvo framework 
kappaP=(epsilon-1)*calvo/(piss^2*(1-calvo)*(1-beta*calvo)); 


% save parameters to use them in Dynare
save par beta alpha delta sigma phi gss ass kappaL piss rss epsilon d0 d1 d2 d3...
    gamma1 gamma2 theta1 theta2 eta erow tauss yss xss hss nu...
     phipi phiy kappaI rhoa rhog rhom kappaP
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%% An Enviromental NK Model %%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% This Dynare code simulates an enviromental NK model as in Annichiarico and Di Dio (2015).

% Author: Valerio Nispi Landi
% First version: 21/08/2020
% This version: 03/04/20201


close all;
warning off

% % Define the enviromental regime
% % tax=1: tax regime. If parameter tau=0 in console_env, no policy regime
% % cap=1: fixed amount of emissions
% % target=1: emission target for unit of output

% Only one of the following can be 1!!!

@#define tax    =  1
@#define cap    =  0
@#define target =  0

%%
%%%%%%%%%%%%%%%%%%%%%%%Endogenous Variables %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
var

c             % consumption 
rk            % rental rate of capital 
rr            % real interest rate
w             % real wage
h             % hours 
y             % output
k             % capital
q             % Tobin Q
i             % investment
lambda        % marginal utility of consumption
r             % nominal interest rate
pi            % inflation 
mc            % net marginal cost
g             % public spending
a             % exogenous TFP
e             % emissions
x             % atmospheric carbon
mu            % abatement
A             % endogenous TFP
tau           % tax/price of emisssion permits

% log variables to have IRFs in percentage deviations from the ss
clog
wlog
hlog
klog
ilog
xlog
elog 
ylog       
;
%%
%%%%%%%%%%%%%%%%%%%%%%%Exogenous Variables%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
varexo 
va           % productivity shock
vg           % public spending shock
vm           % monetary policy shock
;  
    
%%
%%%%%%%%%%%%%%%%%%%%%%%Parameters%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

parameters

beta         % discount factor
alpha        % elasticity of production wrt capital
epsilon      % elasticity of substitution btw differentiated goods
delta        % depreciation rate
sigma        % relative risk aversion
phi          % inverse of Frisch elasticity
eta          % pollution decay
gamma1       % shifter in the emission function
gamma2       % concavity of emission function
theta1       % shifter in the abatement function
theta2       % convexity in the abatement function
d0           % constant in damage function
d1           % 1st order coefficient term in damage function
d2           % 2nd order coefficient term in damage function
d3           % damage function shifter
kappaL       % labor disutility shifter
erow         % row emissions
tauss        % emission tax
nu           % emission target

gss
ass 
piss
rss 
yss 
xss 
hss 

kappaP        % price adjustment costs
kappaI        % investment adjustment cost (as in CEE). If 0, q is constant 
phipi         % mp response to inflation
phiy          % mp response to output
rhoa          % tfp persistence 
rhog          % public spending persistence
rhom          % monetary policy inertia
;

load par;  % load mat file created in console
for jj=1:length(M_.param_names)
set_param_value(M_.param_names{jj},eval(M_.param_names{jj})); 
end;


%%
%%%%%%%%%%%%%%%%%%%%%%%Non-Linear Model%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

model;
% Households
lambda=c^-sigma;  
1=beta*lambda(1)/lambda*r/pi(1);
1=beta*lambda(1)/lambda*(rk(1)+(1-delta)*q(1))/q;
kappaL*h^(phi)=w*lambda;        
k=(1-delta)*k(-1)+(1-kappaI/2*(i/i(-1)-1)^2)*i;
1=q*(1-kappaI/2*(i/i(-1)-1)^2-kappaI*(i/i(-1)-1)*i/i(-1))+kappaI*beta*lambda(1)/lambda*q(1)*(i(1)/i-1)*(i(1)/i)^2;    
% Firms
y=A*k(-1)^alpha*h^(1-alpha);
(1-alpha)*mc*y=w*h;
alpha*mc*y=rk*k(-1);
mu=(tau*gamma1/(theta1*theta2)*y^(-gamma2))^(1/(theta2-1));
(pi-piss)*pi=beta*(lambda(1)/lambda*y(1)/y*pi(1)*(pi(1)-piss))+epsilon/kappaP*(mc+(1-gamma2)*gamma1*tau*(1-mu)*y^(-gamma2)+theta1*mu^(theta2)-(epsilon-1)/epsilon);
% Market clearing
y=c+i+g+(kappaP/2*(pi-piss)^2)*y+theta1*mu^(theta2)*y;
% Monetary Policy
r/(rss)=((pi/piss)^(phipi)*(e/steady_state(e))^(phie))^(1-rhom)*(r(-1)/rss)^(rhom)*exp(vm);
% Enviromental policy
@#if tax
tau=tauss;
@#endif
@#if cap
e=xss*(1-eta)-erow;
@#endif
@#if target
e=nu*y;
@#endif
% Environment
x=eta*x(-1)+e+erow;
e=(1-mu)*gamma1*y^(1-gamma2);
A=(1-d3*(d0+d1*x+d2*x^2))*a;              
% Shocks
log(a)=(1-rhoa)*log(ass)+rhoa*log(a(-1))+va;  
log(g)=(1-rhog)*log(gss)+rhog*log(g(-1))+vg;  
% Auxiliary variables
rr=r/pi(1);
clog=log(c); wlog=log(w); hlog=log(h); klog=log(k); ilog=log(i); 
elog=log(e); xlog=log(x); ylog=log(y);
end;

%% Steady State

steady_state_model;
tau=tauss;
pi=piss;                       
rr=1/beta;                   
r=pi/beta;                   
q=1;   
a=ass;
g=gss;
rk=1/beta-(1-delta);        
y=yss;
h=hss;
x=xss;
e=x*(1-eta)-erow;                          
mu=(gamma1*tau/(theta1*theta2)*y^(-gamma2))^(1/(theta2-1)); 
A=(1-d3*(d0+d1*x+d2*x^2))*a;              
mc=(epsilon-1)/(epsilon)-(1-gamma2)*gamma1*tau*y^(-gamma2)*(1-mu)-theta1*mu^(theta2);
k=alpha*mc*y/rk;                           
i=delta*k;                                 
w=(1-alpha)*mc*y/h;                        
c=y-i-g-y*theta1*mu^(theta2);               
lambda=c^(-sigma);


clog=log(c);
wlog=log(w);
hlog=log(h);
klog=log(k); 
ilog=log(i); 
elog=log(e);
xlog=log(x); 
ylog=log(y);
 
end;

steady;
check;

%% Shocks
shocks;
%var va; stderr 0.01;     
%var vg; stderr 0.01;    
var vm; stderr 0.0025; 
end;


% %% IRFs
 stoch_simul(irf=20,order=1,nograph);
 cap_oo_irfs=oo_.irfs;
 save cap_irf cap_oo_irfs

Can you please elaborate what exactly the issue is and where one can see it. You also did not provide the find_steady_env.m file.

You still did not provide the file called in

Note that the Taylor rule does not describe the equilibrium behavior but is one equation describing the equilibrium. There can be inintuitive nonlinear effects. I would recommend to start with a simpler model and then proceed step by step.