% program2.m
%
% By Dr. Tai-kuang Ho
% Department of Quantitative Finance
% National Tsing Hua University, Taiwan
% ****** Benutzung auf eigene Gefahr! ******
%
% This program solves the small open economy model of Monacelli (2005).
% Monacelli (2005), JMCB
% complete pass-through
% fixed exchange rates
% I use Paul Klein's code (solab.m) to solve the model.
%
% Monacelli (2005) has only domestic and foreign productivity shocks.
% Here I add into his model foreign monetary policy shocks:
% zms=rhozms*zms(-1)+xims


% variables order

%1  qb   qb(+1)=q
%2  sb   sb(+1)=s
%3  rsb  rsb(+1)=rs
%4  z
%5  zs
%6  zms
%7  c
%8  y
%9  r
%10 q
%11 s
%12 pii
%13 piiH
%14 x
%15 mc
%16 ys
%17 piis
%18 rs
%19 xs
%20 mcs

disp('clear');
clear;

% parameters values
betta=0.99;
sigma=1;
eta=1.5;
vphi=3;
thetaH=0.75;
theta =0.75;
rho=0.9;
gama=0.4;

rhor=0.8;
phipii=1.2;
phix=0.5;

rhors=0.8;
phipiis=1.2;
phixs=0.5;

lambdaH=(1-thetaH)*(1-betta*thetaH)/thetaH;
omegas=1+gama*(2-gama)*(sigma*eta-1);

rhozm=0.8;
rhozms=0.8;

% eqs
nvar=20;
nk=6;

A=zeros(nvar,nvar);
B=zeros(nvar,nvar);

A(1,1)=1;
B(1,10)=1;

A(2,2)=1;
B(2,11)=1;

A(3,3)=1;
B(3,18)=1;

A(4,4)=1;
B(4,4)=rho;

A(5,5)=1;
B(5,5)=rho;

A(6,6)=1;
B(6,6)=rhozms;

B(7,7)=-(1-gama);
B(7,8)=1;
B(7,11)=-(2-gama)*gama*eta;
B(7,16)=-gama;

B(8,10)=-1;
B(8,11)=1-gama;

A(9,13)=betta;
B(9,13)=1;
B(9,15)=-lambdaH;

B(10,15)=-1;
B(10,7)=sigma;
B(10,8)=vphi;
B(10,11)=gama;
B(10,4)=-(1+vphi);

B(11,7)=-1;
B(11,16)=1;
B(11,11)=(1/sigma)*(1-gama);

B(12,12)=-1;
B(12,17)=1;
B(12,10)=-1;
B(12,1)=1;

B(13,9)=-1;
B(13,18)=1;

B(14,12)=-1;
B(14,13)=1;
B(14,11)=gama;
B(14,2)=-gama;

B(15,14)=-1;
B(15,8)=1;
B(15,4)=-(omegas*(1+vphi)/(sigma+vphi*omegas));
B(15,16)=-(sigma*(1-omegas)/(sigma+vphi*omegas));

A(16,16)=1;
A(16,17)=1/sigma;
B(16,16)=1;
B(16,18)=1/sigma;

A(17,17)=betta;
B(17,17)=1;
B(17,20)=-(1-theta)*(1-betta*theta)/theta;

B(18,20)=-1;
B(18,16)=(sigma+vphi);
B(18,5)=-(1+vphi);

B(19,18)=-1;
B(19,17)=(1-rhors)*phipiis;
B(19,19)=(1-rhors)*phixs;
B(19,3)=rhors;
B(19,6)=1;

B(20,19)=-1;
B(20,16)=1;
B(20,5)=-(1+vphi)/(sigma+vphi);


[f,p]=solab(A,B,nk);

disp('policy     function'); disp(f);
disp('transition function'); disp(p);


% impulse responses
% use ir.m to compute impulse response function (IRF)
gx=f;
hx=p;

x0=zeros(size(hx,1),1);
x0(end-2)=1; % domestic productivity shock
%x0(end-1)=1; % foreign  productivity shock
%x0(end)=1;   % foreign  monetary shock
IRF=ir(gx,hx,x0,20);

t=(0:size(IRF,1)-1)';
c   =IRF(:,1);
y   =IRF(:,2);
r   =IRF(:,3);
q   =IRF(:,4);
s   =IRF(:,5);
pii =IRF(:,6);
piiH=IRF(:,7);
x   =IRF(:,8);
mc  =IRF(:,9);
ys  =IRF(:,10);
piis=IRF(:,11);
rs  =IRF(:,12);
xs  =IRF(:,13);
mcs =IRF(:,14);
qb  =IRF(:,15);
sb  =IRF(:,16);
rsb =IRF(:,17);
z   =IRF(:,18);
zs  =IRF(:,19);
zms =IRF(:,20);

subplot(3,6,1)
plot(t,c)
title('consumption')

subplot(3,6,2)
plot(t,y)
title('output')

subplot(3,6,3)
plot(t,r)
title('interest rate')

subplot(3,6,4)
plot(t,q)
title('real exchange rate')

subplot(3,6,5)
plot(t,s)
title('terms of trade')

subplot(3,6,6)
plot(t,pii)
title('CPI inflation')

subplot(3,6,7)
plot(t,piiH)
title('producer inflation')

subplot(3,6,8)
plot(t,x)
title('output gap')

subplot(3,6,9)
plot(t,mc)
title('real marginal cost')

subplot(3,6,10)
plot(t,ys)
title('foreign output')

subplot(3,6,11)
plot(t,piis)
title('foreign inflation')

subplot(3,6,12)
plot(t,rs)
title('foreign interest rate')

subplot(3,6,13)
plot(t,xs)
title('foreign output gap')

subplot(3,6,14)
plot(t,mcs)
title('foreign real marginal cost')

subplot(3,6,15)
plot(t,z)
title('domestic productivity shock')

subplot(3,6,16)
plot(t,zs)
title('foreign productivity shocks')

subplot(3,6,17)
plot(t,zms)
title('foreign monetary shocks')

shg


%{
1. I add a lag variable qb, qb(+1)=q.
2. I use the version of interest rate parity that does not involve forward-looking variables.
3. There is no domestic monetary policy, the monetary policy reaction becomes simply r=rs.
%}


