Theoretical Question: compact representation of third-order perturbation solution

Stochastic simulations
1

Hi,
the Dynare framework assumes that all variables have at max only on lead and one lag, so our model framework is (as in Villemot (2011)):

E_t f(y_{t+1}^{+}, y_{t}, y_{t-1}^{-},u_t)

where y_{t} are the n endogenous variables, y_{t-1}^{-} the variables that appear with a lag, and y_{t+1}^{+} the variables that appear with a lead (as in the M_.lead_lag_incidence matrix). I am wondering whether there is a neat way to represent the third-order approximation ghxxx without the use of tensors. Here is what I mean:

First-order Approximation

It can be shown that at first-order the perturbation approach for g_y=ghx satisfies the following equation:

0 = f_{y^+} g_y^{+} g_y^{-} + f_{y^0} g_y + f_{y^-}

where the notation follows Villemot (2011). This can be written more compactly as

0 = f_y \cdot z

where f_y is g1, i.e the first derivative of the dynamic model, re-ordered to the DR-order, and

z =\begin{pmatrix} I\\g_y\\g_y^+ g_y^{-} \end{pmatrix}

Second-order Approximation

At second order it can be shown, see e.g. the korder documentation by Kamenik and Julliard (or the code in dyn_second_order_solver.m) that the perturbation approach for g_{yy}=ghxx satisfies the following equation

0 = f_{y^0} g_{yy} + f_{y^+} g_{yy}^+ (g_y^- \otimes g_y^-) + f_{y^+} g_y^+ g_{yy}^- + f_{y^-y^-} (I \otimes I) + f_{y^-y^0} (g_y \otimes I) + f_{y^-y^+} (g_y^+ g_y^- \otimes I) + f_{y^0y^-}(I \otimes g_y) + f_{y^0y^0}(g_y \otimes g_y) + f_{y^0y^+} (g_y^+ g_y^- \otimes g_y) + f_{y^+y^-} (I \otimes g_y^+ g_y^-) + f_{y^+y^0} (g_y \otimes g_y^+ g_y^-) + f_{y^+y^+} (g_y^+ g_y^- \otimes g_y^+ g_y^-)

or more compactly

0 = f_{yy} \cdot \left(z \otimes z \right) + f_y \cdot zz

where f_{yy} is g2, i.e the first derivative of the dynamic model, re-ordered to the DR-order, and

zz=\begin{pmatrix} 0\\g_{yy}\\g_{yy}^+ (g_y^{-} \otimes g_y^{-}) + g_y^+ g_{yy}^{-} \end{pmatrix}

Third-order Approximation (HERE IS MY QUESTION)

Is there a similar compact representation? I came up with the following:

0 = f_{yyy} \cdot \left(z \otimes z \otimes z \right) + f_{yy} \cdot \left( zz \otimes z + z \otimes zz \right) + f_y \cdot zzz

where f_{yyy} is g3, i.e the first derivative of the dynamic model, re-ordered to the DR-order, and

zzz=\begin{pmatrix} 0\\g_{yyy}\\g_{yyy}^+ (g_y^{-} \otimes g_y^{-} \otimes g_y^{-}) + g_{yy}^+ \cdot ( g_{yy}^{-} \otimes g_{y}^{-} + g_{y}^{-} \otimes g_{yy}^{-} ) + g_{y}^{+}g_{yyy}^{-} \end{pmatrix}

But when testing this with a small script file (see below) this was not correct. In the script I also have the long expression with the individual terms, but this is also not working. What (or which terms or permutations) am I missing?

Illustration CodePerturbationScript.m (9.1 KB)

I have written a small script that illustrates the above equations for any mod file after running stoch_simul(order=3). But unfortunately for the third-order I am doing something wrong either conceptually or in the script…

Any help would be appreciated!
Willi

2

Hi Willi,

the matrix representation exists only for order 1 and 2. Consider the composition of function
y=h(x) and z=g(y)=g(h(x))

First order:

z_{x_a} = \sum_i g_i h^i_a

Second order

z_{x_ax_b} = \sum_i \sum_j g_{ij}h^i_ah^j_b + \sum_i g_i h^i_{ab}

Third order

z_{x_ax_bx_c} = \sum_i\sum_j\sum_kg_{ijk}h^i_ah^j_bh^k_c + \sum_i \sum_j g_{ij}\left(h^i_{ab}h^j_c +h^i_{ac}h^j_b+h^i_ah^j_{bc}\right) + \sum_i g_i h^i_{abc}

There is no way to represent the double sum term with matrix multiplications and Kronecker products.
If you consider a Taylor expansion,multiplying the derivatives by x or \hat x, you can simplify using symmetry. This breaks down if you consider a function with several vectors as arguments

3

Hey @MichelJuillard ,

thank you for the clarification, but I think I have to disagree, there is one :slight_smile:

Andrew Binning gave me the tip to look into Oren Levintal’s “Fifth-Order Perturbation Solution to DSGE Models,” Journal of Economic Dynamics and Control , 80 (2017), and its online appendix. There it is shown that using Matlab’s ipermute function, one can compute a permutation sum matrix \Omega_1, that provides the representation (basically going from tensor to matrix). Thus, defining

zzz = \begin{pmatrix}0\\ g_{yyy}\\ g_{yyy}^+ (g_y^- \otimes g_y^- \otimes g_y^-) + g_{yy}^+ (g_y^- \otimes g_{yy}^-) \Omega_1 + g_y^+ g_{yyy}^- \end{pmatrix}

the compact representation becomes:

0=f_{yyy} \cdot (z \otimes z \otimes z) + f_{yy} \cdot (z \otimes zz) \Omega_1 + f_y \cdot zzz

where \Omega_1 is computed by the following matlab code (provided by Oren Levinal’s replication codes)

M=reshape(1:length(ic)^3,1,length(ic),length(ic),length(ic));
Ix=speye(length(ic)^3);
OMEGA1=Ix(:,ipermute(M,[1,3,4,2]))+Ix(:,ipermute(M,[1,2,4,3]))+Ix(:,ipermute(M,[1,2,3,4])); %matrix that sums permuted second-order tensors (reshaped to matrix)

Just to be complete, here is my updated and corrected script:
PerturbationScript.m (3.5 KB)

4

This is right: you can introduce an additional permutation matrix.