S´ebastienVillemot Solvingrationalexpectationsmodelsatﬁrstorder:whatDynaredoes DynareWorkingPapersSeries

(1)

Dynare Working Papers Series

http://www.dynare.org/wp/

Solving rational expectations models at first order:

what Dynare does

S´ ebastien Villemot

Working Paper no. 2 April 2011

142, rue du Chevaleret — 75013 Paris — France http://www.cepremap.ens.fr

(2)

Solving rational expectations models at first order:

what Dynare does

^∗

S´ebastien Villemot^†

First version: December 2009 This version: April 2011

Abstract

This paper describes in detail the algorithm implemented in Dynare for computing the first order approximated solution of a nonlinear rational expectations model. The core of the algorithm is a generalized Schur decomposition (also known as the QZ decomposition), as advocated by several authors in the litterature. The contribution of the present paper is to focus on implementation details that make the algorithm more generic and more efficient, especially for large models.

JEL classification: C63; C68

Keywords: Dynare; Numerical methods; Perturbation; Rational expectations

1 Introduction

Perturbation techniques are widely used for solving and estimating rational expectation models and Dynare¹ offers a popular, user-friendly access to these techniques. The purpose of the present paper is to describe in detail the algorithm implemented in Dynare for computing the first order approximated solution of nonlinear rational expectations models.²

This algorithm is based on a generalized Schur decomposition—also known as the QZ decomposition—and is therefore essentially a variation on the methods presented by Klein (2000),Sims (2002) andUhlig(1999).

The contribution of this paper is to present some implementation details that make the algorithm more generic and more efficient for large models. In particular I describe the algorithm for removing the leads and lags of more than one in a nonlinear model. I also describe a way

∗Copyright c 2009, 2011 S´ebastien Villemot. Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.3 or any later version published by the Free Software Foundation; with no Invariant Sections, no Front-Cover Texts, and no Back-Cover Texts.

A copy of the license can be found at: http://www.gnu.org/licenses/fdl.txt

The author acknowledges funding through the Seventh Framework Programme for Research (FP7) of the Eu- ropean Commission’s Socio-economic Sciences and Humanities (SSH) Program under grant agreement SSH-CT- 2009-225149.

†Paris School of Economics and CEPREMAP. E-mail: sebastien.villemot@ens.fr.

1Dynare is a free software platform for handling a wide class of economic models. Seehttp://www.dynare.org andAdjemian et al.(2011) for more details.

2This algorithm is available using the stoch simul command of Dynare. The original implementation of this algorithm was done by Michel Juillard, using MATLAB, and is available in the matlab/dr1.m file which is distributed with Dynare. Another implementation was done by the author, in C++, in the DecisionRules class, in the mex/sources/estimation directory of the source tree. The notations used in the present paper are closer to the C++ implementation than to the MATLAB implementation.

(3)

of reducing the size of the Schur decomposition problem by dealing separately with endogenous variables that appear only at the current date (called static endogenous variables in the following).

It should be noted that Dynare is able to go further than first order and can deliver second and third order approximation of the solution of rational expectations models. These higher order solutions can be computed recursively using the first order solution as a starting point.

For algorithmic details on this issue, the interested reader can refer to Collard and Juillard (2001) or Schmitt-Groh´e and Ur´ıbe(2004).

The rest of this paper is organized as follows. Section 2 presents the class of models to be solved and defines a typology of the endogenous variables. Section3 defines the solution to the model and characterizes its first order approximation. Sections 4 and5 describe the algorithm used to recover this first order approximation.

2 The problem to be solved

2.1 The model

In the following, we consider a dynamic nonlinear rational expectations model composed of several variables (endogenous and exogenous) and several equations. In compact form, the model is written as:

Etf (y_t+1⁺ , y_t, y_t−1⁻ , u_t) = 0 (1) where y_t is the vector of endogenous variables, y_t+1⁺ (resp. y⁻_t−1) is the subset of variables of y_t that appear with a lead (resp. a lag), and ut is the vector of exogenous variables. For obvious identification reasons, the model must contain as many equations as there are endogenous variables; let n be this number.

For the timing of variables, the convention used here is the same as in Dynare: a variable decided at date t should have a time subscript of t. For example, in a typical RBC model, the capital stock used to produce date t output is actually decided at date t − 1, so it will be written as k_t−1 using this convention. Thus accordingly, the law of motion of capital will be kt= (1 − δ)kt−1+ it. Another way of expressing this timing convention is that stock variables should use the “end-of-period” convention.

The vector of exogenous variables u_t (of dimension p) follows a Markov process:

u_t= P (u_t−1, ε_t)

where the ε_t are independent and identically distributed innovations with zero mean and variance-covariance matrix Σ.

Note that the stochastic process adopted here for exogenous variables is more generic than the one allowed in Dynare (which only accepts a white noise for stochastic variables, i.e. ut= εt).

2.2 Typology of variables

All endogenous variables are required to appear at least at one period. However it is not required that all endogenous variables appear at the current period (a weaker condition is actually sufficient, see assumption (5) below).

We define four types of endogenous variables:

Static endogenous variables: those that appear only at the current period. Their number is n^s≤ n, and their indices in the y_t vector are ζ_j^s, j = 1 . . . n^s

(4)

Purely forward endogenous variables: those that appear only at the future period, possibly at the current period, but not at the previous period. Their number is n⁺⁺ ≤ n, and their indices ζ_j⁺⁺, j = 1 . . . n⁺⁺

Purely backward endogenous variables: those that appear only at the previous period, possibly at the current period, but not at the future period. Their number is n⁻⁻ ≤ n, and their indices ζ_j⁻⁻, j = 1 . . . n⁻⁻

Mixed endogenous variables: those that appear both at the future and the previous period, and possibly at the current period. Their number is n^m ≤ n, and their indices ζ_j^m, j = 1 . . . n^m

These four types of variables form a partition of the endogenous variables, and we therefore have:

n^m+ n⁺⁺+ n⁻⁻+ n^s = n We also define:

Forward endogenous variables: the union of mixed and purely forward endogenous variables. Their number is n⁺= n⁺⁺+ n^m, and their indices ζ_j⁺, j = 1 . . . n⁺.

Backward endogenous variables: the union of mixed and purely backward endogenous variables. Their number is n⁻= n⁻⁻+ n^m, and their indices ζ_j⁻, j = 1 . . . n⁻

Dynamic endogenous variables: all the variables except static endogenous variables. Their number is n^d= n − n^s, and their indices ζ_j^d, j = 1 . . . n^d

The seven indices are such that 1 ≤ ζ₁^k< ζ₂^k < . . . < ζ_n^kk ≤ n, where k ∈ {s, +, ++, −, −−, m, d}.

We denote by yt= (y1,t, . . . , yn,t)⁰ the vector of endogenous variables at date t. We denote by y_t^k= (y_ζk

1,t, . . . , y_ζk

nk,t)⁰ a subvector of endogenous variables, where k ∈ {s, +, ++, −, −−, m, d}.

We denote by β⁺_j , j = 1 . . . n^m the indices of mixed endogenous variables inside the ζ_j⁺ sequence, i.e. β⁺_j is such that ζ⁺

β_j⁺ is a mixed endogenous variable. We similarly define β_j⁻ for mixed endogenous variables inside the ζ_j⁻ sequence. We similarly define π_j⁺ (resp. π⁻_j ) for purely forward (resp. purely backward) endogenous variables inside ζ_j⁺ (resp ζ_j⁻).

Finally, the vector of state variables is formed by the union of backward endogenous variables at the previous period and of exogenous variables at the current period, and is therefore of size n⁻+ p.

2.3 Removing extra leads and lags

The form given in equation (1) makes the assumption that endogenous variables appear with at most one lead and one lag, and that exogenous variables appear only at the current period.

This assumption does not imply any loss of generality, since it is easy to transform a nonlinear model with many leads and lags into an equivalent model of the form given in (1), as is detailed below.³

For every variable xt in the original model whose maximum lag is xt−d−k with k > 0 (and d = 1 if x is an endogenous variable or d = 0 if it is an exogenous variable), the transformation is the following:

3The algorithm described in the present section is implemented in the Dynare preprocessor, since version 4.1.

Auxiliary variables are created automatically and will show up at several places in Dynare output; seeAdjemian et al.(2011) for the names of these variables.

(5)

• introduce k new endogenous variables z^j_t, for j ∈ {1, . . . , k};

• add k new equations to the model:

z_t¹ = x_t−d

z_t^j = z_t−1^j−1 for j ∈ {2, . . . , k} ;

• replace all occurrences of x_t−d−j (with j > 0) in the original model by z_t−1^j in the transformed model.

The transformation for variables with a lead is a bit more elaborate because one has to handle the fact that there is an expectation operator in front of all equations. The algorithm is as follows:

• decompose every equation of the original model in the following form:

A +X

i∈I

BiEtCi= 0

where A and the B_i are (possibly nonlinear) expressions containing only current or lagged variables, and the C_i are (possibly nonlinear) expressions which may contain leads; this decomposition is not unique, but one should aim at making the Ci terms as simple as possible;

• for every C_i where there is a lead of 2 or more on an endogenous variable, or a lead on an exogenous variable:

– let k be the minimal number of periods so that C_i^(−k) has at most one lead on endogenous variables and no lead on exogenous variables (where C_i^(−k) stands for the transformation of C_i where all variables have been lagged by k periods);

– introduce k new endogenous variables z^j_t, for j ∈ {1, . . . , k};

– add k new equations to the model:

(

z_t¹= EtC_i^(−k)

z_t^j = Etz_t+1^j−1 for j ∈ {2, . . . , k} ;

– replace all occurrences of EtC_i in the original model by Etz_t+1^k in the transformed model.

It is straightforward to see that this transformed model is under the form given in (1). And by the law of iterated expectations, it is equivalent to the original one.

3 The solution and its first order approximation

We first define the deterministic steady state of the model as the vector (¯y, ¯u, ¯ε) satisfying:

¯ ε = 0

¯

u = P (¯u, ¯ε) f (¯y⁺, ¯y, ¯y⁻, ¯u) = 0

(6)

Finding the deterministic steady state involves the resolution of a multivariate nonlinear system.⁴ Then, finding the rational expectation solution of the model means finding the policy functions (also known as decision rules), which give current endogenous variables as a function of state variables:

y_t= g(y_t−1⁻ , u_t)

Note that, by definition of the deterministic steady state, we have ¯y = g(¯y⁻, ¯u).

The function g is characterized by the following functional equation:

Etfg⁺(g⁻(y⁻_t−1, ut), ut+1), g(y_t−1⁻ , ut), y_t−1⁻ , ut = 0 (2) where g⁺ (resp. g⁻) is the restriction of g to forward (resp. backward) endogenous variables.

In the general case, this functional equation cannot be solved exactly, and one has to resort to numerical techniques to get an approximated solution. The remainder of this paper describes the first order perturbation technique implemented in Dynare. Let:

f_y⁺ = ∂f

∂y_t+1⁺ , f_y⁰ = ∂f

∂y_t, f_y⁻ = ∂f

∂y_t−1⁻ , fu= ∂f

∂u_t g_y = ∂g

∂y_t−1⁻ , g_u = ∂g

∂ut

P_u= ∂P

∂u_t−1, P_ε= ∂P

∂ε_t where the derivatives are taken at ¯y, ¯u and ¯ε.

The first order approximation of the policy function is therefore:

ˆ

g(y_t−1⁻ , ut) = ¯y + gyyˆ⁻_t−1+ guuˆt

where ˆy_t−1⁻ = y⁻_t−1− ¯y⁻, ˆut= ut− ¯u, and gy and gu are unknowns at this stage.

A first order approximation of (2) around ¯y and ¯u gives:

f (¯y⁺, ¯y, ¯y⁻, ¯u) + f_y⁺[g⁺_y(g_y⁻yˆ_t−1⁻ + g⁻_uuˆ_t) + g⁺_uEt[P_uuˆ_t+ P_εε_t+1]]

+ f_y⁰(gyyˆ_t−1⁻ + guuˆt) + f_y⁻yˆ⁻_t−1+ fuuˆt= 0 where g⁺_y, g_y⁻, g_u⁻, g⁺_u are the derivatives of the restrictions of g with obvious notation. Com- puting the expectancy term, taking into account the property of the deterministic steady state, and reorganizing the terms, we obtain:

(f_y⁺g_y⁺g_y⁻+ f_y⁰gy+ f_y⁻)ˆy_t−1⁻ + (f_y⁺g⁺_yg⁻_u + f_y⁺g_u⁺Pu+ fygu+ fu)ˆut= 0 (3) In the next sections, we exploit this equation in order to recover the unknown coefficients g_u and g_y.

4 Recovering g

_y

Taking into account the term multiplying ˆy⁻_t−1, equation (3) imposes:

f_y⁺g⁺_yg⁻_y + f_y⁰gy+ f_y⁻= 0

4Dynare offers efficient ways of performing this task, but this is out of the scope of the present paper.

(7)

This amounts to:

f_y+yˆ_t+1⁺ + f_y0yˆ_t+ f_y⁻yˆ⁻_t−1= 0 (4) Let S be the n × n^ssubmatrix of f_y0 where only the columns for static endogenous variables are kept, i.e. Si,j = f_y⁰_,i,ζ^s

j. A QR decomposition⁵ gives S = QR where Q is an n × n orthogonal matrix and R an n × n^s upper triangular matrix.

For the model to be identified, we assume that:

rank(R) = n^s. (5)

Thus, equation (4) can be rewritten as:

A⁺yˆ⁺_t+1+ A⁰yˆ_t+ A⁻yˆ⁻_t−1= 0 (6) where A⁺= Q⁰f_y⁺, A⁰ = Q⁰f_y0and A⁻= Q⁰f_y⁻. By construction, columns of static endogenous variables in A⁰ are zero in their lower part: ∀i > n^s, ∀j ≤ n^s, A⁰_i,ζs

j = 0.

4.1 Non-static endogenous variables Taking only the n^d lower rows of system (6), we get:

A˜⁺yˆ_t+1⁺ + Ã⁰⁺yˆ_t⁺+ Ã⁰⁻yˆ⁻_t + Ã⁻yˆ⁻_t−1= 0 (7) where Ã⁺ (resp. Ã⁻) contains the last n^d rows of A⁺ (resp. A⁻). Matrices Ã⁰⁺ and Ã⁰⁻ can be defined in two ways, depending on where we deal with mixed endogenous variables:

• ˜A⁰⁺ is a submatrix of A⁰ where only the last n^d rows and the columns for forward endogenous variables are kept ( ˜A⁰⁺_i,j = A⁰_n_s_+i,ζ+

j

), and ˜A⁰⁻ is such that purely backward columns are taken from A⁰ ( ˜A⁰⁻

i,π⁻_j = A⁰_n_s_+i,ζ−−

j

), and the rest is zero;

• ˜A⁰⁻ is a submatrix of A⁰ where only the last n^d rows and the columns for backward endogenous variables are kept ( ˜A⁰⁻_i,j = A⁰

n^s+i,ζ_j⁻), and ˜A⁰⁺ is such that purely forward columns are taken from A⁰ ( ˜A⁰⁺

i,π⁺_j = A⁰

n^s+i,ζ_j⁺⁺), and the rest is zero.

Note that in equation (7), static endogenous variables no longer appear.

The structural state space representation of (7) is:

A˜⁰⁻ A˜⁺ I⁻ 0

| {z }

D

ˆy_t⁻ ˆ y_t+1⁺

=− ˜A⁻ − ˜A⁰⁺

0 I⁺

| {z }

E

ˆy_t−1⁻ ˆ y_t⁺

where I⁻ is an n^m× n⁻ selection matrix for mixed endogenous variables, such that I⁻

i,β⁻_i = 1, and zero otherwise. Similarly, I⁺is an n^m× n⁺matrix, such that I⁺

i,β_i⁺ = 1, and zero otherwise.

Therefore, D and E are square matrices of size n⁺⁺+ n⁻⁻+ 2n^m. Using the fact that ˆy_t+1⁺ = g⁺_yyˆ⁻_t , this can be rewritten as:

DI_n⁻ g⁺_y

ˆ

y⁻_t = EI_n⁻ g_y⁺

ˆ

y_t−1⁻ (8)

5SeeGolub and Van Loan(1996, section 5.2).

(8)

where I_n⁻ is the identity matrix of size n⁻.

A generalized Schur decomposition (also known as the QZ decomposition) of the pencil (D, E) is performed:⁶

D = QT Z

E = QSZ

where T is upper triangular, S quasi upper triangular, and Q and Z are orthogonal matrices.

The decomposition is done is such a way that stable generalized eigenvalues (modulus less than 1) are in the upper left corner of T and S.

Matrices T and S are block decomposed so that the upper left block of both matrices is square and contains generalized eigenvalues of modulus less than 1, and the lower right block is square and contains generalized eigenvalues of modulus strictly greater than 1.

Equation (8) can be rewritten as:

T11 T12

0 T₂₂

Z11 Z12

Z₂₁ Z₂₂

I_n⁻ g_y⁺

ˆ

y_t⁻=S11 S12

0 S₂₂

Z11 Z12

Z₂₁ Z₂₂

I_n⁻ g⁺_y

ˆ

y⁻_t−1 (9) where T11and S11 are square and contain stable generalized eigenvalues, while T22 and S22 are square and contain explosive generalized eigenvalues.

To exclude explosive trajectories, we impose:

Z11 Z12

Z₂₁ Z₂₂

I_n⁻ g_y⁺

=X 0

(10) which implies:

g_y⁺= −(Z₂₂)⁻¹Z₂₁

Note that the squareness of Z₂₂is theBlanchard and Kahn(1980) order condition (i.e. the requirement to have as many explosive eigenvalues as forward or mixed endogenous variables), and the non-singularity of Z22 is theBlanchard and Kahn (1980) rank condition.

Using equation (10) and the fact that ˆy⁻_t = g_y⁻yˆ⁻_t−1, equation (9) implies:

T₁₁ T12

0 T22

X 0

g⁻_y =S₁₁ S12

0 S22

X 0

Then, using the fact that solving equation (10) for X gives X = (Z₁₁⁰ )⁻¹, the upper part of this system gives the solution for g_y⁻:

g⁻_y = X⁻¹T₁₁⁻¹S11X = Z₁₁⁰ T₁₁⁻¹S11(Z₁₁⁰ )⁻¹

Note that mixed variables appear in both g⁺ and g⁻: the corresponding lines will be equal across the two matrices by construction.

4.2 Static endogenous variables

The n^s upper lines of equation (6) can be written as:

A˘⁺yˆ⁺_t+1+ ˘A^0dyˆ_t^d+ ˘A^0syˆ_t^s+ ˘A⁻yˆ⁻_t−1= 0 (11) where ˘A⁺ (resp. A˘⁻) contains the first n^s rows of A⁺ (resp. A⁻). Matrix ˘A^0s (resp. A˘^0d) contains the first n^s rows and only the static (resp. non-static) columns of A⁰. Recall that ˘A^0s

6SeeGolub and Van Loan(1996, section 7.7) for theoretical and practical details on this decomposition.

(9)

is a square upper triangular matrix by construction, and it is invertible because of assumption (5).

Equation (11) can be rewritten as:

A˘⁺g⁺_yg⁻_yyˆ_t−1⁻ + ˘A^0dg_y^dyˆ_t−1⁻ + ˘A^0syˆ_t^s+ ˘A⁻yˆ_t−1⁻ = 0

where g_y^d, the restriction of g_y to non-static endogenous variables, is obtained by combining g⁺_y and g⁻_y. We therefore have:

g_y^s= −h ˘A^0si−1

˘A⁺g_y⁺g_y⁻+ ˘A^0dg^d_y+ ˘A⁻

5 Recovering g

_u

Equation (3) restricted to ˆut imposes:

f_y+g_y⁺g_u⁻+ f_y+g_u⁺P_u+ f_yg_u+ f_u = 0, and be rewritten as:

(f_y+g_y⁺J⁻+ f_y)g_u+ f_y+J⁺g_uP_u+ f_u = 0

where J⁻ (resp J⁺) is an n⁻× n matrix (resp. n⁺× n matrix) selecting only the backward (resp. forward) endogenous variables. In the particular case solved by Dynare, where Pu = 0, the solution to this equation is:

gu = −(f_y⁺g_y⁺J⁻+ fy)⁻¹fu

In the general case, this equation is a specialized Sylvester equation, which can be solved using the algorithm proposed byKamenik (2003)⁷.

References

Adjemian, S., H. Bastani, M. Juillard, F. Mihoubi, G. Perendia, M. Ratto, and S. Villemot (2011): “Dynare: Reference Manual, Version 4,” Dynare Working Papers 1, CEPREMAP.

Blanchard, O. J. and C. M. Kahn (1980): “The Solution of Linear Difference Models under Rational Expectations,” Econometrica, 48, 1305–11.

Collard, F. and M. Juillard (2001): “A Higher-Order Taylor Expansion Approach to Simulation of Stochastic Forward-Looking Models with an Application to a Nonlinear Phillips Curve Model,” Computational Economics, 17, 125–39.

Golub, G. H. and C. F. Van Loan (1996): Matrix Computations, The John Hopkins Uni- versity Press, third ed.

Kamenik, O. (2003): “Solution of Specialized Sylvester Equation,” Tech. rep., Manuscript.

Klein, P. (2000): “Using the generalized Schur form to solve a multivariate linear rational expectations model,” Journal of Economic Dynamics and Control, 24, 1405–1423.

7This paper is distributed with Dynare, in the sylvester.pdf file under the documentation directory.

(10)

Schmitt-Groh´e, S. and M. Ur´ıbe (2004): “Solving dynamic general equilibrium models using a second-order approximation to the policy function,” Journal of Economic Dynamics and Control, 28, 755–775.

Sims, C. A. (2002): “Solving Linear Rational Expectations Models,” Computational Eco- nomics, 20, 1–20.

Uhlig, H. (1999): “A toolkit for analysing nonlinear dynamic stochastic models easily,”

in Computational Methods for the Study of Dynamic Economics, ed. by R. Marimon and A. Scott, Oxford University Press, 30–61.