%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Basic RBC Model with public debt (no lump-sum tax/transfer)
% rbc_fiscal3_simul.mod
%(Approximation in logs)
% Using chi0na Data, 1978-2013
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%----------------------------------------------------------------
% 
% U(c,n) = log(c) - theta0*n^(1+chi0)/(1+chi0)+h(G)
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%----------------------------------------------------------------
% 0. Housekeeping (close all graphic windows)
%----------------------------------------------------------------


close all;


%----------------------------------------------------------------
% 1. Defining variables and parameters
%----------------------------------------------------------------

var ly lk lc li z lh lw ld lr Rk  lg tk tn  y_obs;
varexo e1 e2 e3 e4 ;

parameters alpha0 beta0 chi0 delta0 theta0 nss tks tns w  k2nss c2nss d2yss rho1 rho2 rho3 rho4   sigma1 sigma2 sigma3 sigma4 gamk gamn;
% Calibration
alpha0 = 0.55;  % capital share in output
delta0 = 0.10;   % depreciation rate
beta0 = 0.95;% subjective discount factor
chi0 = 1/3;
w = 0.2;          % G/Y

tks=0.25;   % capital income tax rate
tns = 0.1;
%load parameterfile
%set_param_value('tns',tns)

nss = 0.3;  % labor supply (hours worked in total time)


k2nss=(alpha0*(1-tks)/(1/beta0-1+delta0))^(1/(1-alpha0));
c2nss = (1-w)*k2nss^alpha0-delta0*k2nss;

theta0  =(1-tns)*(1-alpha0)*k2nss^alpha0/c2nss/nss^(1+chi0);

d2yss = (alpha0*tks+(1-alpha0)*tns-w)/(1/beta0-1);

rho1 = 0.6567;
rho2 = 0.7508;
rho3 = 0.9;
rho4 = 0.9;

sigma1 = 0.0303;
sigma2 = 0.0128;
sigma3 = 0.01;
sigma4 = 0.01;

gamk = 0.015;
gamn = 0.015;


%----------------------------------------------------------------
% 2. Declaring observable variables
%---------------------------------------------------------------

varobs y_obs;


%----------------------------------------------------------------
% 3. Model
%----------------------------------------------------------------

model; 

# Kss = nss*k2nss;
# yss = log(Kss^alpha0*nss^(1-alpha0));
  //Consumption Euler equation 1
  1/exp(lc) = beta0/exp(lc(+1))*(1+(1-tk(+1))*exp(Rk(+1))-delta0);

// Interest rate of public bond 2
exp(lr)=(1-tk)*exp(Rk)-delta0;

// Labor supply 3
 theta0*exp(chi0*lh)= exp(lw)/exp(lc)*(1-tn);

// Labor demand 4
 exp(lw) = exp(z)*(1-alpha0)*exp(lk(-1))^alpha0*exp(lh)^(-alpha0);

//Capital rental rate 5
 exp(Rk) = alpha0*exp(ly)/exp(lk(-1));

//Production function 6
  exp(ly) = exp(z)*(exp(lk(-1))^alpha0)*(exp(lh))^(1-alpha0);

//Resource constraint 7
  exp(lc)+exp(li)+exp(lg) = exp(ly);

//Capital accumulation equation 8
  exp(li) = exp(lk)-(1-delta0)*exp(lk(-1));

// government budget constraint 9

exp(ld)=(1+exp(lr))*exp(ld(-1))+exp(lg)-tk*exp(Rk)*exp(lk)-tn*exp(lw)*exp(lh);

//TFP shock 10
  z = rho1*z(-1)+e1;

// Gov spending shocks 11
  lg = (1-rho2)*(log(w)+yss)+rho2*lg(-1)+e2;

// tax shocks 12-13
 
  tk = (1-rho3)*(tks)+rho3*tk(-1)+ gamk*(exp(ld(-1)-ly)-d2yss)+e3;
  tn = (1-rho4)*(tns)+rho4*tn(-1)+ gamn*(exp(ld(-1)-ly)-d2yss)+e4;

%% observed equation
y_obs = ly-yss;

end;

%-------------------------------------------------------------
% 4. Specifying Steady State
%--------------------------------------------------------------


initval;
  lk = log(nss*k2nss);
  lc = log(nss*c2nss);
  lh = log(nss);  
  li = log(delta0)+lk;
  ly = lk*alpha0+lh*(1-alpha0);
  lw= log(1-alpha0)+ly-log(nss);
  Rk= log(1/beta0-(1-delta0))-log(1-tks);
  lr= log(1/beta0-1);
  ld = log(d2yss)+ly;
  z = 0;
  lg =  log(w)+ly;  
  tk=  tks;
  tn=  tns;
end;

shocks;
var e1 = sigma1^2;
var e2 = sigma2^2;
var e3 = sigma3^2;
var e4 = sigma4^2;
% var e5 = sigma5^2;
end;

steady;
check;
%----------------------------------------------------------------
% 5. Some Results
%----------------------------------------------------------------
%options_.SpectralDensity.trigger = 1;

stoch_simul(hp_filter = 100, order = 1,irf=20)  ly li lc lh lk Rk lw ld z lg  tk tn   ;
save rbc_cn__fiscal3_stoch_simul;