Added week 8 solutions

.github/workflows/matlab.yml

@@ -71,4 +71,11 @@ jobs:
             cd("progs/matlab")
             VARpDimensionsIllustration
             threeVariableVAROLS
-            threeVariableVARML
+            threeVariableVARML
+      - name: Run week 8 scripts
+        uses: matlab-actions/run-command@v1
+        with:
+          command: |
+            cd("progs/matlab")
+            USOil
+            keatingSR

exercises/svar_irf_solution.tex

@@ -0,0 +1,33 @@
+\begin{enumerate}
+\item There are \(K\) variables and \(K\) structural shocks; hence, there are \(K^2\) IRFs each of length \(H+1\).
+
+\item Let's derive an expression for \(Y_{t+h}\):
+\begin{align*}
+Y_t &= A Y_{t-1} + U_t
+\\
+Y_{t+1} &= A Y_{t} + U_{t+1} = A^2 Y_{t-1} + A U_t + U_{t+1}
+\\
+Y_{t+2} &= A Y_{t+1} + U_{t+2} = A^3 Y_{t-1} + A^2 U_t + A U_{t+1} + U_{t+2}
+\\
+Y_{t+h} &= A^{h+1} Y_{t-1} + \sum_{j=0}^h A^j U_{t+h-j}
+\end{align*}
+Left-multiply by \(J\) (also note that \(J'J=I\)):
+\begin{align*}
+J Y_{t+h} = y_{t+h} = J A^{h+1} Y_{t-1} + \sum_{j=0}^h J A^j \textcolor{red}{J'J} U_{t+h-j} = J A^{h+1} Y_{t-1} + \sum_{j=0}^h J A^j J' u_{t+h-j}
+\end{align*}
+
+\item From the previous exercise:
+\begin{align*}
+y_{t+h} = J A^{h+1} Y_{t-1} + \sum_{j=0}^h J A^j J' \underbrace{u_{t+h-j}}_{B_0^{-1} \varepsilon_{t+h-j}}
+\end{align*}
+Taking the derivative:
+\begin{align*}
+\frac{\partial y_{t+h}}{\partial \varepsilon_{t}'} = \Theta_h = J A^h J' B_0^{-1}
+\end{align*}
+\item Note the following identity: \(p_{t+h}=p_{t-1}+\Delta p_t +\Delta p_{t+1} + \cdots  +\Delta p_{t+h}\).
+That is we simply need to cumulate (\texttt{cumsum}) the impulse responses of the inflation rate to get the impulse response for the price level.
+
+\item A possible implementation:
+\lstinputlisting[style=Matlab-editor,basicstyle=\mlttfamily,title=\lstname]{progs/matlab/irfPlots.m}
+
+\end{enumerate}

exercises/svar_recursive_solution.tex

@@ -0,0 +1,24 @@
+\begin{enumerate}
+\item[1.] The model is identified recursively with the real price of oil ordered first:
+\begin{align*}
+y_t = \begin{pmatrix}\Delta rpoil_t \\ \Delta p_t \\ \Delta gdp_t \end{pmatrix}
+\end{align*}
+such that the real price of oil is \textbf{predetermined} with respect to the U.S. economy.
+In other words, only the structural oil price shock (ordered first) has an immediate effect on the real price of oil,
+the other two structural shocks affect the real price of oil with a delay and not on impact.
+The ordering is thus very important here,
+as our focus is on the effect of an \textbf{unanticipated (exogenous)} increase in the real price of oil.
+Moreover, as we are only interested in one shock,
+the model is only partially identified in that only the oil price shock can be given an economic interpretation.
+
+\item[2./3./4./5.]  \lstinputlisting[style=Matlab-editor,basicstyle=\mlttfamily,title=\lstname]{progs/matlab/USOil.m}
+Here is the helper function for the more general approach:
+\lstinputlisting[style=Matlab-editor,basicstyle=\mlttfamily,title=\lstname]{progs/matlab/USOil_fSR.m}
+\item[6.] In the period under consideration, crude oil is a commodity that is considered vulnerable to negative supply shocks
+due to the political and social volatility of its (Middle East) source.
+In macroeconomics a negative supply shock is an unexpected event that changes the supply of a product,
+resulting in an increase in prices and a decrease in output.
+From the IRFs we see exactly this pattern:
+an unexpected oil price shock creates inflationary pressure on the GDP deflator and a reduction in real GDP in the US\@.
+In this sense, a positive oil price shock (originating e.g.\ from the oil producing countries) indeed acts like a negative domestic supply shock for the U.S. economy.
+\end{enumerate}

exercises/svar_shortrun_solution.tex

@@ -0,0 +1,57 @@
+\begin{enumerate}
+\item The OLS estimate of the reduced-form covariance matrix \(\Sigma_u\) is
+\begin{align*}
+    \widehat{\Sigma}_u = \begin{pmatrix}
+         0.061054  &  -0.015282  &   0.042398  &  0.0037836 \\
+        -0.015282  &    0.52305  &    0.07967  &   0.030602 \\
+         0.042398  &    0.07967  &    0.71689  &    -0.2451 \\
+        0.0037836  &   0.030602  &    -0.2451  &     1.1093 \\
+    \end{pmatrix} 
+\end{align*}
+
+\item In order to identify the structural shocks from \(\Sigma_u = B_0^{-1} \Sigma_\varepsilon B_0^{-1'}\),
+  we require at least \(4(4-1)/2=6\) additional restrictions on \(B_0\) or \(B_0^{-1}\).
+In this exercise, we'll do this on \(B_0\) and not its inverse; that is, we can rewrite the equations in matrix notation:
+\begin{align*}
+    \underbrace{\begin{pmatrix} 1 & 0 & 0 & 0 \\ b_{21,0} & 1 & b_{23,0} & b_{24,0} \\ 0 & 0 & 1 & b_{34,0} \\ b_{41,0} & b_{41,0} & b_{43,0} & 1\end{pmatrix}}_{B_0}
+    \underbrace{\begin{pmatrix}u_t^p\\u_t^{gnp}\\u_t^{i}\\u_t^{m} \end{pmatrix}}_{u_t}
+    =
+    \underbrace{\begin{pmatrix}\varepsilon_t^{AS}\\\varepsilon_t^{IS}\\\varepsilon_t^{MS}\\\varepsilon_t^{MD} \end{pmatrix}}_{\varepsilon_t}
+\end{align*}
+Note that as the diagonal elements of \(B_0\) equal unity, we don't assume that \(E(\varepsilon_t \varepsilon_t')=\Sigma_\varepsilon \) is the identity matrix.
+So this is a different normalization rule as before.
+
+Economically, the above restrictions embody to some extent a baseline IS-LM model:
+\begin{enumerate}
+\item The first equation provides three restrictions by assuming that the price level is predetermined
+except that producers can respond immediately to aggregate supply shocks (e.g.\ an unexpected increase in oil or gas prices).
+This is basically a horizontal aggregate supply (AS) curve; that is why we label the shock with AS\@.
+\item The second equation provides no restrictions as it is assumed that real output responds to all other model variables contemporaneously.
+This equation can be interpreted as an aggregate demand curve, or an IS curve, that is why we label the shock IS\@.
+\item The third equation provides two restrictions by assuming that the interest rate does not react contemporaneously to aggregate measures of output and prices.
+This represents a simple money supply function, according to which the central bank adjusts the rate of interest in relation to the money stock
+and does not immediately observe aggregate output and aggregate prices.
+Of course, this is in contrast to modern monetary policy theory (and practice).
+\item The fourth equation provides one additional restriction by assuming that the first two entries in the last row are identical.
+The underlying idea is that this equation represents a money demand function in which short-run money holdings rise in proportion to NOMINAL GNP\@.
+Moreover, money holdings are allowed to be dependent on the interest rate.
+\end{enumerate}
+In sum, we have \(3+0+2+1=6\) restrictions, that is we have an exactly identified SVAR model,
+identified by short-run exclusion restrictions on \(B_0\).
+As the structure is not recursive, a Cholesky decomposition is not valid and hence we need to rely on a numerical optimizer (or an algorithm for short-run exclusion restrictions).
+
+\item The code might look like this:
+\lstinputlisting[style=Matlab-editor,basicstyle=\mlttfamily,title=\lstname]{progs/matlab/keatingSR.m}
+Here is the helper function to impose the restrictions:
+\lstinputlisting[style=Matlab-editor,basicstyle=\mlttfamily,title=\lstname]{progs/matlab/keatingSR_f.m}
+
+\item An unexpected upward shift of the aggregate supply curve
+\begin{itemize}
+    \item raises the price deflator
+    \item lowers real GNP
+    \item raises the federal funds rate
+    \item lowers the money supply at first, but ultimately raises M1
+\end{itemize}
+This response is not really appropriate for the U.S. economy during the considered sample period,
+making the identifying restrictions rather questionable.
+\end{enumerate}

progs/matlab/USOil.m

@@ -0,0 +1,63 @@
+% -------------------------------------------------------------------------
+% Estimates a SVAR(4) model to identify an oil price shock for the US economy
+% -------------------------------------------------------------------------
+% Willi Mutschler, December 8, 2022
+% willi@mutschler.eu
+% -------------------------------------------------------------------------
+clearvars; clc;close all;
+
+%% data handling
+USOil_data = importdata('../../data/USOil.csv');
+ENDO = USOil_data.data; % note that it already has correct order dlog(poil), dlog(p), dlog(gdp) for recursive identification
+
+%% estimate reduced-form
+nlag = 4;
+opt.const = 1;
+VAR = VARReducedForm(ENDO,nlag,opt);
+
+%% structural identification with Cholesky decomposition
+ENDO = ENDO(:,[1 2 3]); % order dlog(oilprice) first, then dlog(price deflator) and dlog(GDP)
+B0inv_chol = chol(VAR.SigmaOLS([1 2 3],[1 2 3]),'lower');
+% note that the Cholesky decomposition always yields positive diagonal elements, so no normalization needed
+table(B0inv_chol)
+
+%% structural identification with numerical optimization, identification restrictions are put into auxiliary function USOil_fSR
+f = str2func('USOil_fSR');
+StartValueMethod = 1; %0: Use identity matrix, 1: use square root, 2: use cholesky as starting value
+% options for fsolve
+TolX           = 1e-4;  % termination tolerance on the current point
+TolFun         = 1e-9;  % termination tolerance on the function value
+MaxFunEvals    = 50000; % maximum number of function evaluations allowed
+MaxIter        = 1000;  % maximum numberof iterations allowed
+OptimAlgorithm = 'trust-region-dogleg'; % algorithm used in fsolve
+options = optimset('TolX',TolX,'TolFun',TolFun,'MaxFunEvals',MaxFunEvals,'MaxIter',MaxIter,'Algorithm',OptimAlgorithm);
+if StartValueMethod == 0
+    B0inv = eye(size(ENDO,2)); % use identity matrix as starting value
+elseif StartValueMethod == 1
+    B0inv = VAR.SigmaOLS^.5; % use square root of vcov of reduced form as starting value
+elseif StartValueMethod == 2
+    B0inv = chol(VAR.SigmaOLS,'lower'); % use Cholesky decomposition of vcov of reduced form
+end
+f(B0inv,VAR.SigmaOLS)' % test whether function works at initial value (should give you no error)
+
+% call optimization routine fsolve to minimize f
+[B0inv_opt,fval,exitflag,output] = fsolve(f,B0inv,options,VAR.SigmaOLS);
+
+% normalize sign of B0inv such that diagonal elements are positive
+if any(diag(B0inv_opt)<0)
+    idx_negative = diag(B0inv_opt)<0;
+    B0inv_opt(:,(idx_negative==1)) = -1*B0inv_opt(:,(idx_negative==1)); % flip sign
+end
+table(B0inv_opt)
+
+%% compute and plot structural impulse response function
+nSteps = 30;
+cumsumIndicator = [1 0 1]; % the variables in the SVAR are in differences, we want to plot IRFs of the first and last variable by using cumsum
+varNames = ["Real Price of Oil","GDP Deflator Inflation","Real GDP"];
+epsNames = ["Oil Price Shock", "eps2 Shock", "eps3 Shock"];
+IRFpoint_chol = irfPlots(VAR.Acomp,B0inv_chol,nSteps,cumsumIndicator,varNames,epsNames);
+IRFpoint_opt = irfPlots(VAR.Acomp,B0inv_opt,nSteps,cumsumIndicator,varNames,epsNames);
+
+%% both approaches are (numerically) equivalent
+norm(abs(B0inv_chol-B0inv_opt),'Inf') % max abs error
+norm(abs(IRFpoint_chol(:) - IRFpoint_opt(:)),'Inf') % max abs error

progs/matlab/USOil_fSR.m

@@ -0,0 +1,34 @@
+function f = USOil_fSR(B0inv,hatSigmaU)
+% f = USOil_fSR(B0inv,hatSigmaU)
+% -------------------------------------------------------------------------
+% Evaluates the system of nonlinear equations
+%   vech(hatSigmaU) = vech(B0inv*B0inv')
+% subject to specified short-run restrictions
+% -------------------------------------------------------------------------
+% INPUTS
+%   - B0inv     : candidate for short-run impact matrix. [nvars x nvars]
+%   - hatSigmaU : covariance matrix of reduced-form residuals. [nvars x nvars]
+% -------------------------------------------------------------------------
+% OUTPUTS
+%   - f : function value, see below
+% -------------------------------------------------------------------------
+% Willi Mutschler, December 8, 2022
+% willi@mutschler.eu
+% -------------------------------------------------------------------------
+f = [vech(B0inv*B0inv' - hatSigmaU);
+     B0inv(1,2) - 0;
+     B0inv(1,3) - 0;
+     B0inv(2,3) - 0;
+    ];
+
+% % more general way to implement the restrictions
+% % nan means unconstrained, 0 (or any other number) restricts to this number
+% Rshort = [nan   0   0;
+%           nan nan   0;
+%           nan nan nan;
+%           ];
+% selSR = ~isnan(Rshort); % index for short-run restrictions on impact
+% f = [vech(B0inv*B0inv'-hatSigmaU);
+%      B0inv(selSR) - Rshort(selSR);
+%     ];
+end

progs/matlab/keatingSR.m

@@ -0,0 +1,48 @@
+% -------------------------------------------------------------------------
+% Short-run restrictions in the Keating (1992) model.
+% -------------------------------------------------------------------------
+% Willi Mutschler, December 27, 2022
+% willi@mutschler.eu
+% -------------------------------------------------------------------------
+clearvars; clc; close all;
+
+% data handling
+Keating1992 = importdata('../../data/Keating1992.csv');
+ENDO = Keating1992.data;
+[obs_nbr,var_nbr] = size(ENDO);
+
+% estimate reduced-form
+nlag = 4;
+opt.const = 1;
+VAR = VARReducedForm(ENDO,nlag,opt);
+
+% identification restrictions are set up in keatingSR_f_SR
+f = str2func('keatingSR_f');
+
+% options for fsolve
+TolX           = 1e-4;  % termination tolerance on the current point
+TolFun         = 1e-9;  % termination tolerance on the function value
+MaxFunEvals    = 50000; % maximum number of function evaluations allowed
+MaxIter        = 1000;  % maximum numberof iterations allowed
+OptimAlgorithm = 'trust-region-dogleg'; % algorithm used in fsolve
+options = optimset('TolX',TolX,'TolFun',TolFun,'MaxFunEvals',MaxFunEvals,'MaxIter',MaxIter,'Algorithm',OptimAlgorithm);
+
+% initital guess, note that we need candidates for both B_0 and diag(SIG_eps) stored into a single candidate matrix
+B0_diageps= [eye(var_nbr) ones(var_nbr,1)]; % the first 4 columns belong to B_0, the last column are the diagonal elements of SIG_eps
+f(B0_diageps,VAR.SigmaOLS)' % test whether function works at initial value (should give you no error)
+
+% call optimization routine fsolve to minimize f
+[B0_diageps,fval,exitflag,output] = fsolve(f,B0_diageps,options,  VAR.SigmaOLS);
+
+% display results
+B0 = B0_diageps(:,1:var_nbr); % get B0
+SIGeps = B0_diageps(:,var_nbr+1); % these are just the variances
+impact = inv(B0)*diag(sqrt(SIGeps)); % compute impact matrix used to plot IRFs, i.e. inv(B0)*sqrt(SigEps)
+table(impact)
+
+% compute and plot structural impulse response function
+nSteps = 12;
+cumsumIndicator = [1 1 0 1];
+varNames = ["Deflator", "Real GNP", "Federal Funds Rate", "M1"];
+epsNames = ["AS", "IS", "MS", "MD"] + " shock";
+irfPoint = irfPlots(VAR.Acomp,impact,nSteps,cumsumIndicator,varNames,epsNames);

progs/matlab/keatingSR_f.m

@@ -0,0 +1,36 @@
+function f = keatingSR_f(B0_diageps,hatSigmaU)
+% f = keatingSR_f(B0_diageps,hatSigmaU)
+% -------------------------------------------------------------------------
+% Evaluates the system of nonlinear equations vech(hatSigmaU) =
+% vech(B0inv*SIGeps*B0inv') where SIGeps has only values on the diagonal
+% given a candidate for B0 and the diagonal elements of SIGeps
+% subject to the short-run restrictions in the Keating (1992) model.
+% -----------------------------------------------------------------------
+% INPUTS
+%   - B0_diageps   [nvars x (nvars+1)]  candidate matrix for both short-run impact matrix [nvars x nvars] and diagonal elements of SIGeps (nvar x 1) 
+%   - hatSigmaU    [nvars x nvars]      covariance matrix of reduced-form residuals
+% -----------------------------------------------------------------------
+% OUTPUTS
+%   - f : function value, see below
+% -------------------------------------------------------------------------
+% Willi Mutschler, December 14, 2022
+% willi@mutschler.eu
+% -------------------------------------------------------------------------
+
+nvars = size(hatSigmaU,1);  % number of variables
+B0 = B0_diageps(:,1:nvars); % get B0 matrix from candidate matrix
+diageps = diag(B0_diageps(1:nvars,nvars+1)); % get diagonal elements of SIG_eps from candidate matrix and make it a full diagonal matrix
+B0inv = inv(B0); % compute short-run impact matrix
+
+f = [vech(B0inv*diageps*B0inv'-hatSigmaU);
+     B0(1,1) - 1;
+     B0(3,1) - 0;
+     B0(4,1) - B0(4,2);
+     B0(1,2) - 0;
+     B0(2,2) - 1;
+     B0(3,2) - 0;
+     B0(1,3) - 0;
+     B0(3,3) - 1;
+     B0(1,4) - 0;
+     B0(4,4) - 1;
+    ];

week_8.tex

@@ -1,7 +1,7 @@
 % !TEX root = week_8.tex
 \input{exercises/_common_header.tex}
 \Newassociation{solution}{Solution}{week_8_solution}
-\newif\ifDisplaySolutions%\DisplaySolutionstrue
+\newif\ifDisplaySolutions\DisplaySolutionstrue
 
 \begin{document}
 \title{Quantitative Macroeconomics\\~\\Winter 2023/24\\~\\Week 8}