Solution Methods for Microeconomic Dynamic Stochastic Optimization Problems

October 31, 2012

Christopher D. Carroll¹

_________________________________________________________________________________

Abstract
These notes describe some tools for solving microeconomic dynamic stochastic optimization problems, and show how to use those tools for efficiently estimating a standard life cycle consumption/saving model using microeconomic data. No attempt is made at a systematic overview of the many possible technical choices; instead, I present a specific set of methods that have proven useful in my own work (and explain why other popular methods, such as value function iteration, are a bad idea). Paired with these notes is Mathematica and Matlab software that solves the problems described in the text.

Keywords: Dynamic Stochastic Optimization, Method of Simulated Moments, Structural Estimation
JEL codes: E21, F41

PDF:

http://econ.jhu.edu/people/ccarroll/SolvingMicroDSOPs.pdf

Slides:

http://econ.jhu.edu/people/ccarroll/SolvingMicroDSOPs-Slides.pdf

Web:

http://econ.jhu.edu/people/ccarroll/SolvingMicroDSOPs/

Archive:

http://econ.jhu.edu/people/ccarroll/SolvingMicroDSOPs.zip

(Contains LaTeX code for this document and software producing figures and results)

¹Carroll: Department of Economics, Johns Hopkins University, Baltimore, MD, http://econ.jhu.edu/people/ccarroll/, ccarroll@jhu.edu, Phone: (410) 516-7602 The notes were originally written for my Advanced Topics in Macroeconomic Theory class at Johns Hopkins University; instructors elsewhere are welcome to use them for teaching purposes. This draft incorporates several improvements related to new results in the paper “Theoretical Foundations of Buffer Stock Saving” (especially tools for approximating the consumption and value functions). Like the last major draft, it also builds on material in “The Method of Endogenous Gridpoints for Solving Dynamic Stochastic Optimization Problems” published in Economics Letters, available at http://econ.jhu.edu/people/ccarroll/EndogenousArchive.zip, and by including sample code for a method of simulated moments estimation of the life cycle model a la Gourinchas and Parker (2002) and Cagetti (2003). Background derivations, notation, and related subjects are treated in my class notes for first year macro, available at http://econ.jhu.edu/people/ccarroll/public/lecturenotes/consumption. I am grateful to several generations of graduate students in helping me to refine these notes, to Marc Chan for help in updating the text and software to be consistent with Carroll (2006), to Kiichi Tokuoka for drafting the section on structural estimation, to Damiano Sandri for exceptionally insightful help in revising and updating the method of simulated moments estimation section, and to Weifeng Wu and Metin Uyanik for revising to be consistent with the ‘method of moderation’ and other improvements. All errors are my own.

1 Introduction

Calculating the optimal amount to save is a remarkably difficult problem. Under empirically realistic assumptions about the nature of risk (and attitudes toward risk), the problem cannot be solved without serious numerical computation. To avoid having to solve this hard problem, economists have shown impressive ingenuity in reformulating the question. Budding graduate students are exposed (often with little motivation) to a host of tricks: Quadratic or Constant Absolute Risk Aversion utility, perfect markets, perfect insurance, perfect foresight, the “timeless” perspective, the restriction of uncertainty to very special kinds,¹ and more.

Explicit or not, the motivation is always to exchange an intractable general problem for a tractable specific alternative. But the burgeoning literature on numerical solutions has shown that the features that yield tractability also profoundly change the solution. The “tricks” are excuses to solve a problem that has defined away the central difficulty: Understanding the proper role of uncertainty in optimal intertemporal choice.

Fortunately, the temptation to use such tricks is waning, thanks to advances in mathematical analysis, increasing computing power, and the growing capabilities of numerical computation software. Together, these tools permit today’s laptop computers to solve the kinds of problems that required supercomputers a decade ago (and, before that, could not be solved at all).

These lecture notes provide a gentle introduction to a particular set of solution tools and show how they can be used to solve some canonical problems in consumption choice and portfolio allocation. Specifically, the notes describe and solve optimization problems for a consumer facing uninsurable idiosyncratic risk to nonfinancial income (e.g., labor or transfer income),² with detailed intuitive discussion of the various mathematical and computational techniques that, together, speed the solution by many orders of magnitude compared to “brute force” methods. The problem is solved with and without liquidity constraints, and the infinite horizon solution is obtained as the limit of the finite horizon solution. After the basic consumption/saving problem with a deterministic interest rate is described and solved, an extension with portfolio choice between a riskless and a risky asset is also solved. Finally, a simple example is presented of how to use these methods (via the statistical ‘method of simulated moments’ or MSM; sometimes called ‘simulated method of moments’ or SMM) to estimate structural parameters like the coefficient of relative risk aversion (a la Gourinchas and Parker (2002) and Cagetti (2003)).

2 The Problem

We are interested in the behavior a consumer whose goal in period $t$ is to maximize discounted utility from consumption over the remainder of a lifetime that ends in period $T$ :

[ T - t ] ∑ n max Et β u (ccct+n) , n=0

(1)

and whose circumstances evolve according to the transition equations³

aaat = mmmt - ccct (2 ) bbb = aaa R (3 ) t+1 t t+1 yyyt+1 = pppt+1θt+1 (4 ) mmmt+1 = bbbt+1 + yyyt+1 (5 )

where

β-pure time discount factor aaat-assets after all actions have been accomplished in period t bbbt+1-‘bank balances ’ (nonhuman wealth ) at the beginning of t + 1 ccct-consumption in period t mmm-‘market resources ’ available for consumption (‘cash -on -hand ’) t pppt+1-‘permanent labor income ’ in period t + 1 Rt+1-interest factor (1 + rt+1 ) from period t to t + 1 yyyt+1-noncapital income in period t + 1.

The exogenous variables evolve as follows:

Rt= R -constant interest factor = 1 + r pppt+1= Γ t+1pppt -permanent labor income dynamics (6 ) 2 2 logθt+n~ N (- σ θ∕2, σθ ) -lognormal transitory shocks ∀ n > 0.

Using the fact about lognormally distributed variables ELogNorm⁴ that if $logΦ ~ N (ϕ, σ2 ) ϕ$ then $log E [Φ ] = ϕ + σ2∕2 ϕ$ , assumption (7) guarantees that $log E [θ] = 0$ which means that $E [θ]$ =1 (the mean value of the transitory shock is 1).

Equation (6) indicates that we are allowing for a predictable average profile of income growth over the lifetime ${ Γ }T 0$ (allowing, for example, for typical career wage paths).⁵

Finally, we assume that the utility function is of the Constant Relative Risk Aversion (CRRA), form, $u (∙) = ∙1- ρ∕(1 - ρ)$ .

As is well known, this problem can be rewritten in recursive (Bellman equation) form

vt(mmmt, pppt) = max u (ccct) + Et [βvt+1 (mmmt+1, pppt+1)] (7 ) ccct

subject to the Dynamic Budget Constraint (DBC) (2)-(5) given above, where $v t$ measures total expected discounted utility from behaving optimally now and henceforth.

3 Normalization

The single most powerful method for speeding the solution of such models is to redefine the problem in a way that reduces the number of state variables (if possible). In the consumption problem here, the obvious idea is to see whether the problem can be rewritten in terms of the ratio of various variables to permanent noncapital (‘labor’) income $pppt$ (henceforth for brevity referred to simply as ‘permanent income.’)

In the last period of life, there is no future, $vT +1 = 0$ , so the optimal plan is to consume everything, implying that

1- ρ mmm T vT (mmmT ,pppT ) = ------. (8 ) 1 - ρ

Now define nonbold variables as the bold variable divided by the level of permanent income in the same period, so that, for example, $m = mmm ∕ppp T T T$ ; and define $vT(mT)=u(mT )$ .⁶ For our CRRA utility function, $u(xy ) = x1 - ρu (y)$ , so equation (8) can be rewritten as

1- ρ 1- ρ v(mmm,ppp ) = ppp1- ρ m-T---= ppp1- ρ Γ 1- ρm-T---= ppp1- ρΓ 1- ρv (m ). TTT T 1 - ρ T- 1 T 1 - ρ T- 1 T T T

Now define a new optimization problem:

v (m ) = max u(c ) + E [β Γ 1- ρv (m )] (9 ) t t ct t t t+1 t+1 t+1 s.t. at = mt - ct mt+1 = (R∕-Γ t+1-)at + θt+1 ◟ ◝ ◜ ◞ ≡Rt+1

The accumulation equation is the normalized version of the transition equation for $mmmt+1$ .⁷ Then it is easy to see that for $t = T - 1$ ,

v (mmm ,ppp ) = ppp1- ρ v (m ) T- 1 T - 1 T- 1 T- 1 T - 1 T- 1

and so on back to all earlier periods. Hence, if we solve the problem (9) which has only a single state variable $m t$ , we can obtain the levels of the value function, consumption, and all other variables of interest simply by multiplying the results by the appropriate function of $ppp t$ , e.g. $ct(mmmt,pppt) = ppptct(mmmt ∕pppt)$ or $1- ρ vt(mmmt, pppt ) = ppp t vt(mt )$ . We have thus reduced the problem from two continuous state variables to one (and thereby enormously simplified its solution).

For some problems it will not be obvious that there is an appropriate ‘normalizing’ variable, but many problems can be normalized if sufficient thought is given. For example, Valencia (2006) shows how a bank’s optimization problem can be normalized by the level of the bank’s productivity.

4 The Usual Theory, and A Bit More Notation

The first order condition for (9) with respect to $ct$ is

′ 1- ρ ′ u (ct) = Et[βRt+1 Γ t+1 v t+1 (mt+1 )] - ρ ′ = Et[βR Γt+1 vt+1(mt+1 )] (10 )

and because the Envelope theorem tells us that

v′t(mt ) = Et [βR Γ - ρ v′t+1 (mt+1 )] (11 ) t+1

we can substitute the LHS of (11) for the RHS of (10) to get

u ′(c ) = v ′(m ) (12 ) t t t

and rolling this equation forward one period yields

′ ′ u (ct+1) = vt+1 (atRt+1 + θt+1) (13 )

while substituting the LHS in equation (10) gives us the Euler equation for consumption

u ′(ct) = Et [βR Γ -t+ ρ1u ′(ct+1)]. (14 )

Now note that in equation (13) neither $mt$ nor $ct$ has any direct effect on $vt+1$ - only the difference between them (i.e. unconsumed market resources or ‘assets’ $at$ ) matters. It is therefore possible (and will turn out to be convenient) to define a function⁸

1- ρ 𝔳t(at ) = Et[β Γt+1 vt+1(Rt+1at + θt+1 )] (15 )

that returns the expected $t + 1$ value associated with ending period $t$ with any given amount of assets. This definition implies that

𝔳 ′(at) = Et[βR Γ - ρv ′ (Rt+1at + θt+1 )] (16 ) t t+1 t+1

or, substituting from equation (13),

′ [ - ρ ′ ] 𝔳t(at) = Et βR Γ t+1u (ct+1 (Rt+1at + θt+1)) . (17 )

Finally, note for future use that the first order condition (10) can now be rewritten as

′ ′ u (ct) = 𝔳 t(mt - ct). (18 )

5 Solving the Next-to-Last Period

The problem in the second-to-last period of life is:

[ 1- ρ ] vT-1(mT-1) = max u(cT - 1) + β ET - 1 ΓT vT ((mT - 1 - cT - 1)RT + θT) , cT-1

and using (1) the fact that $vT = u (c)$ ; (2) the definition of $u (c)$ ; (3) the definition of the expectations operator, and (4) the fact that $Γ T$ is nonstochastic, this becomes

1- ρ ∫ cT- 1 1- ρ ∞ ((mT - 1 - cT - 1)RT + θ )1- ρ vT-1(mT-1) = max -------+ βΓ T -------------------------------dF (θ ) cT-1 1 - ρ 0 1 - ρ

where $F$ is the cumulative distribution function for $θ$ .

In principle, the maximization implicitly defines a function $cT - 1(mT - 1)$ that yields optimal consumption in period $T - 1$ for any given level of resources $mT- 1$ . Unfortunately, however, there is no general analytical solution to this maximization problem, and so for any given $mT - 1$ we must use numerical computational tools to find the $cT - 1$ that maximizes the expression. This is excruciatingly slow because for every potential $cT - 1$ to be considered, the integral must be calculated numerically, and numerical integration is very slow.

5.1 Discretizing the Distribution

Our first time-saving step is therefore to construct a discrete approximation to the lognormal distribution that can be used in place of numerical integration. We calculate an $n$ -point approximation as follows.

Define a set of points from $♯0$ to $♯n$ on the $[0, 1]$ interval as the elements of the set $♯={0,1∕n,2∕n, ..., 1}$ .⁹ Call the inverse of the $θ$ distribution $- 1 F$ , and define the points $♯-1i=F-1(♯i)$ . Then the conditional mean of $θ$ in each of the intervals numbered 1 to $n$ is:

∫ - 1 - 1 - 1 ♯i θi ≡ E [θ|♯i- 1 ≤ θ < ♯i ] = ϑ dF (ϑ). (19 ) ♯-i1-1

The method is illustrated in Figure 1. The solid continuous curve represents the “true” CDF $F (θ )$ for a lognormal distribution such that $E[θ]=1$ , $σθ = 0.1$ . The short vertical line segments represent the $n$ equiprobable values of $θ i$ which are used to approximate this distribution.¹⁰

Figure 1: Discrete Approximation to Lognormal Distribution $F$

Recalling our definition of $𝔳t(at )$ , for $t = T - 1$

( ) n 1- ρ 1- ρ 1 ∑ (RT aT - 1 + θi) 𝔳T- 1 (aT- 1) = β ΓT -- --------------------- (20 ) n i=1 1 - ρ

so we can rewrite the maximization problem as

{ 1- ρ } -cT--1- vT-1(mT - 1) = macx + 𝔳T- 1(mT - 1 - cT- 1) . (21 ) T-1 1 - ρ

5.2 The Approximate Consumption and Value Functions

Given a particular value of $m T - 1$ , a numerical maximization routine can now find the $cT- 1$ that maximizes (21) in a reasonable amount of time. The Mathematica program that solves exactly this problem is called 2period.m. (The archive also contains parallel Matlab programs, but these notes will dwell on the specifics of the Mathematica implementation, which is superior in many respects.)

The first thing 2period.m does is to read in the file functions.m which contains definitions of the consumption and value functions; functions.m also defines the function SolveAnotherPeriod which, given the existence in memory of a solution for period $t + 1$ , solves for period $t$ .

The next step is to run the programs setup_params.m, setup_grids.m, setup_shocks.m, respectively. setup_params.m sets values for the parameter values like the coefficient of relative risk aversion. setup_shocks.m calculates the values for the $θ i$ defined above (and puts those values, and the (identical) probability associated with each of them, in the vector variables $θ$ Vals and $θ$ Prob). Finally, setup_grids.m constructs a list of potential values of cash-on-hand and saving, then puts them in the vector variables mVec = aVec = ${0, 1, 2,3, 4}$ respectively. Then 2period.m runs the program setup_lastperiod.m which defines the elements necessary to determine behavior in the last period, in which $c (m ) = m T$ and $vT(m)=u(m )$ .

After all the setup, the only remaining step in 2period.m is to invoke SolveAnotherPeriod, which constructs the solution for period $T-1$ given the presence of the solution for period $T$ (constructed by setup_lastperiod.m).

Because we will always be comparing our solution to the perfect foresight solution, we also construct the variables required to characterize the perfect foresight consumption function in periods prior to $T$ . In particular, we construct the list yExpPDV (which contains the PDV of expected income – ‘expected human wealth’), and yMinPDV which contains the minimum possible discounted value of future income at the beginning of period $T - 1$ (‘minimum human wealth’).¹¹

The perfect foresight consumption function is also constructed (setup_PerfectForesightSolution.m). This program uses the fact that, in Mathematica, functions can be saved as objects using the commands # and &. The # denotes the argument of the function, while the &, placed at the end of the function, tells Mathematica that the function should be saved as an object. In the program, the last period perfect foresight consumption function is saved as an element in the list $cг$ = {(# - 1 + Last[yExpPDV]) Last[ $κ$ Min] &}, where Last[yExpPDV] gives the just-constructed PDV of human wealth at the beginning of $T$ (equal to 1, since current income is included in $h T$ ), and Last[ $κ$ Min] gives the perfect foresight marginal propensity to consume (equal to 1, since it is optimal to spend all resources in the last period). Since # in the code stands in for what was called $m$ in the model, the discounted total wealth is decomposed into discounted non-human wealth # - 1 and discounted human wealth Last[yExpPDV]. The resulting formula then corresponds to $¯c = (m - 1 + h )κ- T T T T$ , which translates to $¯c = m T T$ for $hT=κT=1$ .

The infinite horizon perfect foresight marginal propensity to save

1∕ρ λ = (1 ∕R )(Rβ ) (22 )

is also defined because it will be useful in a number of derivations.¹²

The program then constructs behavior for one iteration back from the last period of life by using the function AddNewPeriodToParamLifeDates. Using the Mathematica command AppendTo, various existing lists (which characterized the solution for period $T$ ) are redefined to include an additional element representing the relevant formulas in the second to last period of life. For example, $κ$ Min now has two elements. The second element, given by 1/(1 + Last[ $λ$ ]/Last[ $κ$ Min]), is the perfect foresight marginal propensity to consume in $t=T-1$ .¹³

Next, the program defines a function $𝔳$ [at_] (in functions_stable.m) which is the exact implementation of (15): It returns the expectation of the value of behaving optimally in period $T$ given any specific amount of assets at the end of $T - 1$ , $a T - 1$ .

The heart of the program is the next expression (in functions.m). This expression loops over the values of the variable mVec, solving the maximization problem (given in equation (21)):

max u[c ] + 𝔳 [mVec [[i ]] -c ] (23 ) c

for each of the $i$ values of mVec (henceforth let’s call these points $mT - 1,i$ ). The maximization routine returns two values: the maximized value, and the value of $c$ which yields that maximized value. When the loop (the Table command) is finished, the variable vAndcList contains two lists, one with the values $vT- 1,i$ and the other with the consumption levels $cT - 1,i$ associated with the $mT - 1,i$ .

5.3 An Interpolated Consumption Function

Now we use the first of the really convenient built-in features of Mathematica. Given a set of points on a function (in this case, the consumption function $c (m ) T-1$ ), Mathematica can create an object called an InterpolatingFunction which when applied to an input $m$ will yield the value of $c$ that corresponds to a linear interpolation of the value of $c$ from the points in the InterpolatingFunction object. We can therefore define an approximation to the consumption function $`c (m ) T- 1 T - 1$ which, when called with an $mT - 1$ that is equal to one of the points in mVec[[i]] returns the associated value of $c T - 1,i$ , and when called with a value of $mT-1$ that is not exactly equal to one of the mVec[[i]], returns the value of $c$ that reflects a linear interpolation between the $c T- 1,i$ associated with the two mVec[[i]] points nearest to $mT - 1$ . Thus if the function is called with $m = 1.75 T - 1$ and the nearest gridpoints are $mj,T-1=1$ and $mk,T - 1 = 2$ then the value of $cT - 1$ returned by the function would be $(0.25c + 0.75c ) j,T - 1 k,T- 1$ . We can define a numerical approximation to the value function $`vT - 1(mT - 1)$ in an exactly analogous way.

Figures 2 and 3 show plots of the $`cT- 1$ and $`vT - 1$ InterpolatingFunctions that are generated by the program 2PeriodInt.m. While the $`c T - 1$ function looks very smooth, the fact that the $`vT - 1$ function is a set of line segments is very evident. This figure provides the beginning of the intuition for why trying to approximate the value function directly is a bad idea (in this context).¹⁴

Figure 2: $cT-1(mT -1)$ (solid) versus $`cT-1(mT -1)$ (dashed)

Figure 3: $vT-1$ (solid) versus $`vT-1(mT -1)$ (dashed)

5.4 Interpolating Expectations

2period.m works well in the sense that it generates a good approximation to the true optimal consumption function. However, there is a clear inefficiency in the program: Since it uses equation (21), for every value of $mT - 1$ the program must calculate the utility consequences of various possible choices of $cT- 1$ as it searches for the best choice. But for any given value of $aT-1$ , there is a good chance that the program may end up calculating the corresponding $𝔳$ many times while maximizing utility from different $mT- 1$ ’s. For example, it is possible that the program will calculate the value of ending the period with $aT - 1 = 0$ dozens of times. It would be much more efficient if the program could make that calculation once and then merely recall the value when it is needed again.

This can be achieved using the same interpolation technique used above to construct a direct numerical approximation to the value function: Define a grid of possible values for saving at time $T - 1$ , $⃗aT - 1$ (aVec in setup_grids.m), designating the specific points $aT- 1,i$ ; for each of these values of $aT - 1,i$ , calculate the vector $⃗ 𝔳T - 1$ as the collection of points $𝔳T-1,i = 𝔳T - 1(aT - 1,i)$ using equation (15); then construct an InterpolatingFunction object $` 𝔳T - 1(aT - 1)$ from the list of points on the function captured in the $⃗aT - 1$ and $⃗𝔳T - 1$ vectors.

Thus, we are now interpolating for the function that reveals the expected value of ending the period with a given amount of assets.¹⁵ The program 2periodIntExp.m solves this problem. Figure 4 compares the true value function to the InterpolatingFunction approximation; the functions are of course identical at the gridpoints chosen for $aT- 1$ and they appear reasonably close except in the region below $m = 1 T- 1$ .

Figure 4: End-Of-Period Value $𝔳T -1(aT-1)$ (solid) versus $` 𝔳T- 1(aT-1)$ (dashed)

Figure 5: $cT-1(mT -1)$ (solid) versus $`cT-1(mT -1)$ (dashed)

Nevertheless, the resulting consumption rule obtained when $`𝔳T- 1(aT - 1)$ is used instead of $𝔳 (a ) T- 1 T - 1$ is surprisingly bad, as shown in figure 5. For example, when $mT - 1$ goes from 2 to 3, $`cT- 1$ goes from about 1 to about 2, yet when $m T- 1$ goes from 3 to 4, $c` T- 1$ goes from about 2 to about 2.05. The function fails even to be strictly concave, which is distressing because Carroll and Kimball (1996) prove that the correct consumption function is strictly concave in a wide class of problems that includes this problem.

5.5 Value Function versus First Order Condition

Loosely speaking, our difficulty is caused by the fact that the consumption choice is governed by the marginal value function, not by the level of the value function (which is what we approximated). To see this, recall that a quadratic utility function exhibits risk aversion because with a stochastic $c$ ,

E [- (c - /c)2] < - (E [c ] - /c)2

(24)

where $c /$ is the ‘bliss point’. However, unlike the CRRA utility function, with quadratic utility the consumption/saving behavior of consumers is unaffected by risk since behavior is determined by the first order condition, which depends on marginal utility, and when utility is quadratic, marginal utility is unaffected by risk:

E [- 2(c - /c)] = - 2 (E [c ] - /c).

(25)

Intuitively, if one’s goal is to accurately capture choices that are governed by marginal value, numerical techniques that approximate the marginal value function will lead to a more accurate approximation to optimal behavior than techniques that approximate the level of the value function.

The first order condition of the maximization problem in period $T - 1$ is:

′ - ρ ′ u(cT- 1) = β ET - 1[ΓT Ru (cT)] (26 ) ( ) ∑ n c- ρ = R β 1- Γ - ρ(R (m - c ) + θ )- ρ . (27 ) T- 1 n T T- 1 T- 1 i i=1

Figure 6: $′ u(c)$ versus $′ ′ `′ `′ 𝔳T-1(3 - c),𝔳T-1(4 - c),𝔳T- 1(3 - c),𝔳T- 1(4 - c)$

The downward-sloping curve in Figure 6 shows the value of $- ρ cT - 1$ for our baseline parameter values for $0 ≤ c ≤ 4 T- 1$ (the horizontal axis). The solid upward-sloping curve shows the value of the RHS of (27) as a function of $c T-1$ under the assumption that $m = 3 T- 1$ . Constructing this figure is rather time-consuming, because for every value of $cT - 1$ plotted we must calculate the RHS of (27). The value of $c T- 1$ for which the RHS and LHS of (27) are equal is the optimal level of consumption given that $m=3 T-1$ , so the intersection of the downward-sloping and the upward-sloping curves gives the optimal value of $cT - 1$ . As we can see, the two curves intersect just below $c = 2 T- 1$ . Similarly, the upward-sloping dashed curve shows the expected value of the RHS of (27) under the assumption that $m = 4 T- 1$ , and the intersection of this curve with $′ u (cT- 1)$ yields the optimal level of consumption if $m=4 T-1$ . These two curves intersect slightly below $c = 2.5 T- 1$ . Thus, increasing $mT- 1$ from 3 to 4 increases optimal consumption by about 0.5.

Now consider the derivative of our function $`𝔳T - 1(aT - 1)$ . Because we have constructed $`𝔳 T- 1$ as a linear interpolation, the slope of $`𝔳T-1(aT-1)$ between any two adjacent points ${aT - 1,i,ai+1,T - 1}$ is constant. The level of the slope immediately below any particular gridpoint is different, of course, from the slope above that gridpoint, a fact which implies that the derivative of $`𝔳 (a ) T- 1 T- 1$ follows a step function.

The solid-line step function in Figure 6 depicts the actual value of $′ `𝔳T-1(3-cT- 1)$ . When we attempt to find optimal values of $cT - 1$ given $m T-1$ using $`𝔳 (a ) T- 1 T- 1$ , the numerical optimization routine will return the $cT - 1$ for which $′ ′ u (cT - 1) = `𝔳T - 1(mT - 1 - cT - 1)$ . Thus, for $m=3 T-1$ the program will return the value of $c T- 1$ for which the downward-sloping $′ u (cT- 1)$ curve intersects with the $′ `𝔳 T- 1(3 - cT - 1)$ ; as the diagram shows, this value is exactly equal to 2. Similarly, if we ask the routine to find the optimal $cT - 1$ for $mT - 1 = 4$ , it finds the point of intersection of $u′(c ) T- 1$ with $`𝔳 ′ (4 - c ) T- 1 T- 1$ ; and as the diagram shows, this intersection is only slightly above 2. Hence, this figure illustrates why the numerical consumption function plotted earlier returned values very close to $cT - 1 = 2$ for both $mT - 1 = 3$ and $m=4 T-1$ .

We would obviously obtain much better estimates of the point of intersection between $′ u (cT- 1)$ and $′ 𝔳 T- 1(mT - 1 - cT - 1)$ if our estimate of $`𝔳′ T-1$ were not a step function. In fact, we already know how to construct linear interpolations to functions, so the obvious next step is to construct a linear interpolating approximation to the expected marginal value of end-of-period assets function $′ 𝔳$ . That is, we calculate

( ) n ′ - ρ 1 ∑ - ρ 𝔳T- 1(aT- 1) = βR Γ T -- (RT aT - 1 + θi) (28 ) n i=1

at the points in aVec yielding $′ ′ { {aT - 1,1,𝔳T - 1,1},{aT - 1,2, 𝔳T- 1,2} ...}$ and construct $`𝔳′T- 1(aT - 1)$ as the linear interpolating function that fits this set of points.

Figure 7: $′ 𝔳T-1(aT-1)$ versus $`′ 𝔳T- 1(aT -1)$

The program file functionsIntExpFOC.m therefore uses the function $𝔳$ a[at_] defined in functions_stable.m as the embodiment of equation (28), and constructs the InterpolatingFunction as described above. The results are shown in Figure 7. The linear interpolating approximation looks roughly as good (or bad) for the marginal value function as it was for the level of the value function. However, Figure 8 shows that the new consumption function (long dashes) is a considerably better approximation of the true consumption function (solid) than was the consumption function obtained by approximating the level of the value function (short dashes).

Figure 8: $cT-1(mT - 1)$ (solid) Versus Two Methods for Constructing $`cT-1(mT - 1)$

5.6 Transformation

However, even the new-and-improved consumption function diverges appallingly from the true solution, especially at lower values of $m$ . That is because the linear interpolation does an increasingly poor job of capturing the nonlinearity of $𝔳 ′T- 1(aT- 1)$ at lower and lower levels of $a$ .

This is where we unveil our next trick. To understand the logic, start by considering the case where $RT = β = Γ T = 1$ and there is no uncertainty (that is, we know for sure that income next period will be $θT = 1$ ). The final Euler equation is then:

- ρ - ρ cT - 1 = cT . (29 )

In the case we are now considering with no uncertainty and no liquidity constraints, the optimizing consumer does not care whether a unit of income is scheduled to be received in the future period $T$ or the current period $T-1$ ; there is perfect certainty that the income will be received, so the consumer treats it as equivalent to a unit of current wealth. Total resources therefore are comprised of two types: current market resources $mT - 1$ and ‘human wealth’ (the PDV of future income) of $𝔥T - 1 = 1$ (where we use the Gothic font to signify that this is the expectation, as of the END of the period, of the income that will be received in future periods; it does not include current income, which has already been incorporated into $mT-1$ ).

The optimal solution is to spend half of total lifetime resources in period $T-1$ and the remainder in period $T$ . Since total resources are known with certainty to be $mT - 1 + 𝔥T - 1 = mT - 1 + 1$ , and since $v′(m)= u′(c ) T-1T-1 T- 1$ this implies that

( ) - ρ v′ (mT - 1) = mT--- 1-+-1 . T- 1 2

(30)

Of course, this is a highly nonlinear function. However, if we raise both sides of (30) to the power $(- 1 ∕ρ)$ the result is a linear function:

′ - 1∕ρ mT - 1 + 1 [vT - 1(mT - 1) ] = -----------. (31 ) 2

This is a specific example of a general phenomenon: A theoretical literature cited in Carroll and Kimball (1996) establishes that under perfect certainty, if the period-by-period marginal utility function is of the form $c- ρ t$ , the marginal value function will be of the form $-ρ (γmt+ζ)$ for some constants ${γ, ζ}$ . This means that if we were solving the perfect foresight problem numerically, we could always calculate a numerically exact (because linear) interpolation. To put this in intuitive terms, the problem we are facing is that the marginal value function is highly nonlinear. But we have a compelling solution to that problem, because the nonlinearity springs largely from the fact that we are raising something to the power $- ρ$ . In effect, we can ‘unwind’ all of the nonlinearity owing to that operation and the remaining nonlinearity will not be nearly so great. Specifically, applying the foregoing insights to the end-of-period value function $𝔳 T- 1$ , we can define

′ - 1∕ρ 𝔠T - 1(aT - 1) ≡ [𝔳 T- 1(aT - 1)] (32 )

which would be linear in the perfect foresight case. Thus, our procedure is to calculate the values of $𝔠T - 1,i$ at each of the $aT - 1,i$ gridpoints, with the idea that we will construct $`𝔠T - 1$ as the interpolating function connecting these points.

5.7 The Self-Imposed ‘Natural’ Borrowing Constraint and the $a T- 1$ Lower Bound

This is the appropriate moment to ask an awkward question that we have neglected until now: How should a function like $`𝔠T - 1$ be evaluated outside the range of points spanned by ${a , ...,a } T- 1,1 T - 1,n$ for which we have calculated the corresponding $𝔠T - 1,i$ gridpoints used to produce our linearly interpolating approximation $`𝔠 T - 1$ (as described in section 5.3)?

The natural answer would seem to be linear extrapolation; for example, we could use

′ `𝔠T-1(aT- 1) = `𝔠T- 1(aT - 1,1) + `𝔠T- 1(aT - 1,1)(aT - 1 - aT- 1,1) (33 )

for values of $aT - 1 < aT - 1,1$ . Unfortunately, this approach will lead us into difficulties. To see why, consider what happens to the true (not approximated) $𝔳T - 1(aT - 1)$ as $aT - 1$ approaches the value $a- = - θR - 1 T - 1 T$ . From (28) we have

( ) ∑n lim𝔳′ (a ) = lim βR Γ - ρ 1- (a R + θ )- ρ(3.4 ) aT-1↓aT-1T - 1 T- 1 aT-1↓aT-1 T n T- 1 T i i=1

But since $θ = θ -- 1$ , exactly at $a = a T- 1 -T - 1$ the first term in the summation would be $(- θ-+ θ1)- ρ = 1 ∕0ρ$ which is infinity. The reason is simple: $-a T - 1$ is the PDV, as of $T - 1$ , of the minimum possible realization of income in period $T$ ( $RT aT - 1 = - θ1$ ). Thus, if the consumer borrows an amount greater than or equal to $θR - 1 -- T$ (that is, if the consumer ends $T - 1$ with $- 1 aT - 1 ≤ - θR T$ ) and then draws the worst possible income shock in period $T$ , he will have to consume zero in period $T$ (or a negative amount), which yields $- ∞$ utility and $∞$ marginal utility (or undefined utility and marginal utility).

These reflections lead us to the conclusion that the consumer faces a ‘self-imposed’ liquidity constraint (which results from the precautionary motive): He will never borrow an amount greater than or equal to $- 1 θRT$ (that is, assets will never reach the lower bound of $a T-1$ ).¹⁶ The constraint is ‘self-imposed’ in the sense that if the utility function were different (say, Constant Absolute Risk Aversion), the consumer would be willing to borrow more than $- 1 θR T$ because a choice of zero or negative consumption in period $T$ would yield some finite amount of utility (though it is very unclear what a proper economic interpretation of negative consumption might be – this is an important reason why CARA utility, like quadratic utility, is increasingly not used for serious quantitative work, though it is still useful for teaching purposes).

This self-imposed constraint cannot be captured well when the $𝔳 ′ T- 1$ function is approximated by a piecewise linear function like $′ `𝔳T- 1$ , because a linear approximation can never reach the correct gridpoint for $′ 𝔳T-1(aT-1)= ∞.$ To see what will happen instead, note first that if we are approximating $𝔳′ T - 1$ the smallest value in aVec must be greater than $aT-1$ (because the expectation for any gridpoint $≤ aT - 1$ is undefined). Then when the approximating $𝔳′ T- 1$ function is evaluated at some value less than the first element in aVec[1], the approximating function will linearly extrapolate the slope that characterized the lowest segment of the piecewise linear approximation (between aVec[1] and aVec[2]), a procedure that will return a positive finite number, even if the requested $aT - 1$ point is below $aT - 1$ . This means that the precautionary saving motive is understated, and by an arbitrarily large amount as the level of assets approaches its true theoretical minimum $a T-1$ .

The foregoing logic demonstrates that the marginal value of saving approaches infinity as $a ↓ a- = - θR - 1 T- 1 T- 1 T$ . But this implies that $′ - 1∕ρ limaT-1↓aT-1𝔠T - 1(aT - 1) = ( 𝔳T- 1(aT - 1)) = 0$ ; that is, as $a$ approaches its minimum possible value, the corresponding amount of $c$ must approach its minimum possible value: zero.

The upshot of this discussion is a realization that all we need to do is to augment each of the $⃗aT - 1$ and $⃗cT - 1$ vectors with an extra point so that the first element in the list used to produce our InterpolatingFunction is ${aT-1,0,cT-1,0} = {aT- - 1,0.}$ .

Figure 9: $𝔠T-1(aT -1)$ versus $` 𝔠T-1(aT-1)$

Figure 9 plots the results (generated by the program 2periodIntExpFOCInv.m). The solid line calculates the exact numerical value of $𝔠 (a ) T- 1 T- 1$ while the dashed line is the linear interpolating approximation $`𝔠T- 1(aT - 1).$ This figure well illustrates the value of the transformation: The true function is close to linear, and so the linear approximation is almost indistinguishable from the true function except at the very lowest values of $aT-1$ .

Figure 10 similarly shows that when we calculate $``𝔳′ (a ) T- 1 T - 1$ as $- ρ [`𝔠T-1(aT-1)]$ (dashed line) we obtain a much closer approximation to the true function $𝔳 ′ (a ) T- 1 T- 1$ (solid line) than we did in the previous program which did not do the transformation (Figure 7).

Figure 10: $′ 𝔳T-1(aT-1)$ vs. $``′ 𝔳T-1(aT- 1)$ Constructed Using $` 𝔠T-1(aT-1)$

5.8 The Method of Endogenous Gridpoints

Our solution procedure for $cT- 1$ still requires us, for each point in $⃗mT - 1$ (mVect in the code), to use a numerical rootfinding algorithm to search for the value of $cT- 1$ that solves $′ ′ u (cT- 1) = 𝔳T- 1(mT - 1 - cT- 1)$ . Unfortunately, rootfinding is a notoriously slow operation.

Fortunately, our next trick lets us completely skip this computationally burdensome step. The method can be understood by noting that any arbitrary value of $aT - 1,i$ (greater than its lower bound value $aT - 1$ ) will be associated with some marginal valuation as of the end of period $T - 1$ , and the further observation that it is trivial to find the value of $c$ that yields the same marginal valuation, using the first order condition,

′ ′ u (cT- 1,i) = 𝔳T - 1(aT - 1,i) (35 ) cT- 1,i = u′- 1(𝔳 ′ (aT - 1,i)) (36 ) T- 1 = (𝔳′T - 1(aT - 1,i))- 1∕ρ (37 ) ≡ 𝔠T- 1(aT - 1,i) (38 ) ≡ 𝔠T- 1,i. (39 )

But with mutually consistent values of $c T- 1,i$ and $a T- 1,i$ (consistent, in the sense that they are the unique optimal values that correspond to the solution to the problem in a single state), we can obtain the $m T - 1,i$ that corresponds to both of them from

m = c + a . (40 ) T - 1,i T - 1,i T - 1,i

These $mT- 1$ gridpoints are “endogenous” in contrast to the usual solution method of specifying some ex-ante grid of values of $mT - 1$ and then using a rootfinding routine to locate the corresponding optimal $cT-1$ .

Thus, we can generate a set of $mT - 1,i$ and $cT - 1,i$ pairs that can be interpolated between in order to yield $`c(mT - 1 )$ at virtually zero computational cost once we have the $⃗𝔠T- 1$ values in hand!¹⁷ One might worry about whether the ${m, c }$ points obtained in this way will provide a good representation of the consumption function as a whole, but in practice there are good reasons why they work well (basically, this procedure generates a set of gridpoints that is naturally dense right around the parts of the function with the greatest nonlinearity).

Figure 11: $cT- 1(mT -1)$ (solid) versus $`cT-1(mT -1)$ (dashed)

Figure 11 plots the actual consumption function $cT- 1$ and the approximated consumption function $`c T- 1$ derived by the method of endogenous grid points. Compared to the approximate consumption functions illustrated in Figure 8 $`c T- 1$ is quite close to the actual consumption function.

5.9 Improving the $a$ Grid

Thus far, we have arbitrarily used $a$ gridpoints of ${0.,1., 2.,3.,4.}$ (augmented in the last subsection by $aT - 1$ ). But it has been obvious from the figures that the approximated $`𝔠T- 1$ function tends to be farthest from its true value $𝔠T - 1$ at low values of $a$ . Combining this with our insight that $aT-1$ is a lower bound, we are now in position to define a more deliberate method for constructing gridpoints for $aT - 1$ – a method that yields values that are more densely spaced than the uniform grid at low values of $a$ . A pragmatic choice that works well is to find the values such that (1) the last value exceeds the lower bound by the same amount $¯aT- 1$ as our original maximum gridpoint (in our case, 4.); (2) we have the same number of gridpoints as before; and (3) the multi-exponential growth rate (that is, $ee... e$ for some number of exponentiations $n$ ) from each point to the next point is constant (instead of, as previously, imposing constancy of the absolute gap between points).

Figure 12: $𝔠T- 1(aT -1)$ versus $` 𝔠T-1(aT-1)$ , Multi-Exponential aVec

Figure 13: $′ 𝔳T-1(aT-1)$ vs. $``′ 𝔳T-1(aT- 1)$ , Multi-Exponential aVec

The results (generated by the program 2periodIntExpFOCInvEEE.m) are depicted in Figures 12 and 13, which are notably closer to their respective truths than the corresponding figures that used the original grid.

5.10 The Method of Moderation

Unfortunately, this endogenous gridpoints solution is not very well-behaved outside the original range of gridpoints targeted by the solution method. (Though other common solution methods are no better outside their own predefined ranges). Figure 14 demonstrates the point by plotting the amount of precautionary saving implied by a linear extrapolation of our approximated consumption rule (the consumption of the perfect foresight consumer $¯cT- 1$ minus our approximation to optimal consumption under uncertainty, $`cT- 1$ ). Although theory proves that precautionary saving is always positive, the linearly extrapolated numerical approximation eventually predicts negative precautionary saving (at the point in the figure where the extrapolated locus crosses the horizontal axis).

Figure 14: For Large Enough $mT -1$ , Predicted Precautionary Saving is Negative (Oops!)

This error cannot be fixed by extending the upper gridpoint; in the presence of serious uncertainty, the consumption rule will need to be evaluated outside of any prespecified grid (because starting from the top gridpoint, a large enough realization of the uncertain variable will push next period’s realization of assets above that top; a similar argument applies below the bottom gridpoint). While a judicious extrapolation technique can prevent this problem from being fatal (for example by carefully excluding negative precautionary saving), the problem is often dealt with using inelegant methods whose implications for the accuracy of the solution are difficult to gauge.

As a preliminary to our solution, define $𝔥t$ as end-of-period human wealth (the present discounted value of future labor income) for a perfect foresight version of the problem of a ‘risk optimist:’ a consumer who believes with perfect confidence that the shocks will always take the value 1, $θt+n = E[θ ] = 1 ∀ n > 0$ . The solution to a perfect foresight problem of this kind takes the form¹⁸

¯ct(mt ) = (mt + 𝔥t )κt (41 )

for a constant minimal marginal propensity to consume $κt$ given below.

We similarly define $𝔥t$ as ‘minimal human wealth,’ the present discounted value of labor income if the shocks were to take on their worst possible value in every future period $θt+n = θ- ∀ n > 0$ (which we define as corresponding to the beliefs of a ‘pessimist’).

We will call a ‘realist’ the consumer who correctly perceives the true probabilities of the future risks and optimizes accordingly.

A first useful point is that, for the realist, a lower bound for the level of market resources is $mt- = - 𝔥t$ , because if $mt$ equalled this value then there would be a positive finite chance (however small) of receiving $θ = θ- t+n$ in every future period, which would require the consumer to set $ct$ to zero in order to guarantee that the intertemporal budget constraint holds (this is the multiperiod generalization of the discussion in section 5.7 about $a T-1$ ). Since consumption of zero yields negative infinite utility, the solution to realist consumer’s problem is not well defined for values of $m<m tt$ , and the limiting value of the realist’s $c t$ is zero as $mt↓mt$ .

Given this result, it will be convenient to define ‘excess’ market resources as the amount by which actual resources exceed the lower bound, and ‘excess’ human wealth as the amount by which mean expected human wealth exceeds guaranteed minimum human wealth:

=◜-◞◟mt◝ ▴m = m + 𝔥 t t -t ▴ 𝔥t = 𝔥t - 𝔥-. t

We can now transparently define the optimal consumption rules for the two perfect foresight problems, those of the ‘optimist’ and the ‘pessimist.’ The ‘pessimist’ perceives human wealth to be equal to its minimum feasible value $𝔥t$ with certainty, so consumption is given by the perfect foresight solution

c-(m ) = (m + 𝔥 )κ- t t t -t t = ▴mt κt.

The ‘optimist,’ on the other hand, pretends that there is no uncertainty about future income, and therefore consumes

¯c (m ) = (m + 𝔥 - 𝔥 + 𝔥 )κ t t t -t -t t -t = (▴mt + ▴ 𝔥t)κt = ct(mt ) + ▴ 𝔥tκt.

It seems obvious that the spending of the realist will be strictly greater than that of the pessimist and strictly less than that of the optimist. Figure 15 illustrates the proposition for the consumption rule in period $T-1$ .

Figure 15: Moderation Illustrated: $cT-1 < `cT- 1 < ¯cT-1$

Proof is more difficult than might be imagined, but the necessary work is done in Carroll (2011) so we will take the proposition as a fact and proceed by manipulating the inequality:

▴mtκt< ct(mt- + ▴mt ) < (▴mt + ▴ 𝔥t)κt
-▴mtκ> - ct(m-- + ▴mt ) > - (▴mt + ▴ 𝔥t)κ-
t t t
▴𝔥tκt> ( ¯ct(mt- + ▴mt ) - ct(mt- + ▴mt )) > 0
¯ct(mt- + ▴mt ) - ct(mt- + ▴mt )
1> ----------------------------------- > 0
---------------▴-𝔥t κt----------------
◟ ◝ ◜ ◞
≡ ˆϙt

where the fraction in the middle of the last inequality is the ratio of actual precautionary saving (the numerator is the difference between perfect-foresight consumption and optimal consumption in the presence of uncertainty) to the maximum conceivable amount of precautionary saving (the amount that would be undertaken by the pessimist who consumes nothing out of any future income beyond the perfectly certain component).

Defining $μ = log ▴m t t$ (which can range from $- ∞$ to $∞$ ), the object in the middle of the last inequality is

( ¯c (m + eμt) - c (m + eμt)) ˆϙ (μt ) ≡ -t--t-------------t--t--------- , (42 ) t ▴ 𝔥tκt

and we now define

( ) 1----ˆϙt(μt-) ˆχχχt (μt) = log ˆϙ (μ ) (43 ) t t = log (1∕ ˆϙt(μt) - 1) (44 )

which has the virtue that it is linear in the limit as $μt$ approaches $+∞$ .

Given $ˆχχχ$ , the consumption function can be recovered from

------= ˆϙt------ ◜( ◞ ◟ ◝) ------1------- ˆct = ¯ct - 1 + exp (ˆχχχ ) ▴ 𝔥tκt. (45 ) t

Thus, the procedure is to calculate $ˆχχχ t$ at the points $⃗μt$ corresponding to the log of the $▴m⃗t$ points defined above, and then using these to construct an interpolating approximation $`ˆχχχt$ from which we indirectly obtain our approximated consumption rule $` ˆct$ by substituting $` ˆχχχt$ for $ˆχχχ$ in equation (45).

Because this method relies upon the fact that the problem is easy to solve if the decision maker has unreasonable views (either in the optimistic or the pessimistic direction), and because the correct solution is always between these immoderate extremes, we call our solution procedure the ‘method of moderation.’

Results are shown in Figure 16; a reader with very good eyesight might be able to detect the barest hint of a discrepancy between the Truth and the Approximation at the far righthand edge of the figure – a stark contrast with the calamitous divergence evident in Figure 14.

Figure 16: Extrapolated $` ˆcT-1$ Constructed Using the Method of Moderation

5.11 Approximating the Slope Too

Until now, we have calculated the level of consumption at various different gridpoints and used linear interpolation (either directly for $cT- 1$ or indirectly for, say, $ˆχ χχT - 1$ ). But the resulting piecewise linear approximations have the unattractive feature that they are not differentiable at the ‘kink points’ that correspond to the gridpoints where the slope of the function changes discretely.

Carroll (2011) shows that the true consumption function for this problem is ‘smooth:’ It exhibits a well-defined unique marginal propensity to consume at every positive value of $m$ . This suggests that we should calculate, not just the level of consumption, but also the marginal propensity to consume (henceforth $κ$ ) at each gridpoint, and then find an interpolating approximation that smoothly matches both the level and the slope at those points.

This requires us to differentiate (42) and (44), yielding

( ≡κκκt(mt) ) ◜-----◞ ◟-----◝ ˆϙ μ(μt) = (▴ 𝔥t κ-)- 1e μt( κ-- cm (m-- + eμt)) (46 ) t t t t t ( - ˆϙ μ(μ )∕ ˆϙ2 ) ˆχχχ μ(μt) = ----t---t---t-- (47 ) t 1∕ˆϙt (μt) - 1

and (dropping arguments) with some algebra these can be combined to yield

( ) μ ----κt▴mt--▴-𝔥t(-κt --κt-)--- ˆχχχ t = (¯c - c )(¯c - c - κ ▴ 𝔥 ) . (48 ) t t t t -t t

To compute the vector of values of (46) corresponding to the points in $⃗μ t$ , we need the marginal propensities to consume (designated $κ$ ) at each of the gridpoints, $cm t$ (the vector of such values is $⃗κ t$ ). These can be obtained by differentiating the Euler equation (18) (where we define $𝔪t(a ) ≡ 𝔠t(a ) + a$ ):

u′(𝔠 ) = ˆ𝔳a(𝔪 - 𝔠 ) (49 ) t t t t

with respect to $a$ , yielding a marginal propensity to have consumed $a 𝔠$ at each gridpoint:

u′′(𝔠t)𝔠a = ˆ𝔳aa (𝔪t - 𝔠t) (50 ) t t 𝔠at = ˆ𝔳aat (𝔪t - 𝔠t)∕u ′′(𝔠t) (51 )

and the marginal propensity to consume at the beginning of the period is obtained from the marginal propensity to have consumed by noting that

𝔠 = 𝔪 - a 𝔠a + 1 = 𝔪a

which, together with the chain rule $𝔠a = cm 𝔪a$ , yields the MPC from

=𝔪a- m ◜ a◞◟ ◝ a c (𝔠 + 1) = 𝔠 (52 ) cm = 𝔠a∕ (1 + 𝔠a ). (53 )

Designating $`ˆc T - 1$ as the approximated consumption rule obtained using an interpolating polynomial approximation to $ˆχχχ$ that matches both the level and the first derivative at the gridpoints, Figure 17 plots the difference between this latest approximation and the true consumption rule for period $T - 1$ up to the same large value (far beyond the largest gridpoint) used in prior figures. Of course, at the gridpoints the approximation will match the true function; but this figure illustrates that the approximation is quite accurate far beyond the last gridpoint (which is the last point at which the difference touches the horizontal axis). (We plot here the difference between the two functions rather than the level plotted in previous figures, because in levels the approximation error would not be detectable even to the most eagle-eyed reader.)

Figure 17: Difference Between True $cT-1$ and $` ˆcT-1$ Is Minuscule

5.12 Value

Often it is useful to know the value function as well as the consumption rule associated with a problem. Fortunately, many of the tricks used when solving the consumption problem have a direct analogue in approximation of the value function.

Consider the perfect foresight (or ‘optimist’) case in period $T - 1$ :

vT - 1(mT - 1) = u (cT- 1) + βu (cT) ( ) = u (cT- 1) 1 + β( (βT R)1∕ρ)1- ρ ( ) 1∕ρ- 1 = u (cT- 1) 1 + β( βT R) ( ) = u (cT- 1) 1 + (βT R )1∕ρ∕R = u (cT- 1) PDVTt (c)∕cT - 1 ◟------◝◜------◞ ≡ℂTt

where $PDVT(c) t$ is the present discounted value of consumption. A similar function can be constructed recursively for earlier periods, yielding the general expression

v (m ) = u (c )ℂT (54 ) t t t t

which can be transformed as

1∕(1- ρ) Λt ≡ ((1 - ρ)vt) = c (ℂT )1∕(1- ρ) t t

with derivative

¯Λm = (ℂT )1∕(1- ρ)κ-, t t t

and since $T ℂt$ is a constant while the consumption function is linear, $¯Λt$ will also be linear.

We apply the same transformation to the value function for the problem with uncertainty (the realist’s problem) and differentiate

Λ = ((1 - ρ )v (m ))1∕(1- ρ) t t t Λmt = ((1 - ρ )vt(mt ))- 1+1∕(1- ρ) vmt (mt )

and an excellent approximation to the value function can be obtained by calculating the values of $Λ$ at the same gridpoints used by the consumption function approximation, and interpolating among those points.

However, as with the consumption approximation, we can do even better if we realize that the $¯Λ$ function for the optimist’s problem is an upper bound for the $Λ$ function in the presence of uncertainty, and the value function for the pessimist is a lower bound. Analogously to (42), define an upper-case

( ) ¯Λt(mt- + eμt) - Λt (mt-+ eμt) ˆϘt(μt ) = ---------------T--1∕(1--ρ)------ (55 ) ▴ 𝔥tκt(ℂ t )

with derivative (dropping arguments)

μ ˆϘ t = (▴ 𝔥tκt (ℂTt )1∕(1- ρ))- 1eμt (¯Λmt - Λmt ) (56 )

and an upper-case version of the $χχχ$ equation in (44):

( ) 1 - ˆϘt(μt) Xˆt (μt ) = log ------------ (57 ) ˆϘt (μt) ( ) = log 1∕ ˆϘt(μt ) - 1 (58 )

with corresponding derivative

( μ 2) μ - ˆϘ t ∕ˆϘ t ˆX t = ---------- (59 ) 1 ∕ˆϘt - 1

and if we approximate these objects then invert them (as above with the $ˆϙ$ and $ˆχχχ$ functions) we obtain a very high-quality approximation to our inverted value function at the same points for which we have our approximated value function:

= ˆϘ ◜(-------◞◟t-----)◝ ˆΛ = ¯Λ - ------1------- ▴ 𝔥 κ (ℂT )1∕(1- ρ) (60 ) t t 1 + exp (ˆX ) t-t t t

from which we obtain our approximation to the value function and its derivative as

ˆv = u (ˆΛ ) t t ˆvmt = u ′(ˆΛt) ˆΛm mm ′′ m 2 ′ mm ˆv t = u (ˆΛt)(ˆΛ ) + u (ˆΛt)ˆΛ .

Although a linear interpolation that matches the level of $Λ$ at the gridpoints is simple, a Hermite interpolation that matches both the level and the derivative of the $ˆΛt$ function at the gridpoints has the considerable virtue that the $¯vt$ derived from it numerically satisfies the envelope theorem at each of the gridpoints for which the problem has been solved.

If we use the double-derivative calculated above to produce a higher-order Hermite polynomial, our approximation will also match marginal propensity to consume at the gridpoints; this would guarantee that the consumption function generated from the value function would match both the level of consumption and the marginal propensity to consume at the gridpoints; the numerical differences between the newly constructed consumption function and the highly accurate one constructed earlier would be negligible within the grid (though would likely diverge outside it).

5.13 Refinement: A Tighter Upper Bound

Carroll (2011) derives an upper limit $¯κt$ for the MPC as $mt$ approaches its lower bound. Using this fact plus the strict concavity of the consumption function yields the proposition that

ct(mt-+ ▴mt ) < ¯κt▴mt. (61 )

The solution method described above does not guarantee that approximated consumption will respect this constraint between gridpoints, and a failure to respect the constraint can occasionally cause computational problems in solving or simulating the model. Here, we describe a method for constructing an approximation that always satisfies the constraint.

Defining $# mt$ as the ‘cusp’ point where the two upper bounds intersect:

( ) ▴m# + ▴ 𝔥t κt = ¯κt▴m# t t # κt▴ 𝔥t ▴m t = ------------ (1 - κt)κ¯t κ-𝔥t - 𝔥 m# = ---t------t-, t (1 - κt)κ¯t

we want to construct a consumption function for $# mt ∈ (mt, m t ]$ that respects the tighter upper bound:

▴mκ < c (m + ▴m ) < ¯κ ▴m
t-t t -t- t t t
▴mt(¯κt-κt ) > ¯κt▴(mt - ct(mt--+ ▴mt) ) > 0
¯κt▴mt-- ct(mt+-▴mt-)
1 > ▴mt (¯κt- κt) > 0.

Again defining $μt = log ▴mt$ , the object in the middle of the inequality is

¯κt - ct(mt- + eμt)e- μt ˇϙt(μt ) ≡ ------------------------ ¯κt - κt μ ct(m-- + eμt)e - μt - κκκm (m-- + eμt) ˇϙ t (μt ) = -----t-----------------t----t--------. ¯κt - κt

As $mt$ approaches $- mt-$ , $ˇϙt(μt)$ converges to zero, while as $mt$ approaches $+ ∞$ , $ˇϙ (μt) t$ approaches $1$ .

As before, we can derive an approximated consumption function; call it $`ˇct$ . This function will clearly do a better job approximating the consumption function for low values of $mt$ while the previous approximation will perform better for high values of $mt$ .

For middling values of $m$ it is not clear which of these functions will perform better. However, an alternative is available which performs well. Define the highest gridpoint below $# m t$ as $# m¯ˇ t$ and the lowest gridpoint above $m# t$ as $ˆm#- t$ . Then there will be a unique interpolating polynomial that matches the level and slope of the consumption function at these two points. Call this function $˜ct(m )$ .

Using indicator functions that are zero everywhere except for specified intervals,

111 (m ) = 1 if m ≤ ¯ˇm# Lo t 111Mid (m ) = 1 if m¯ˇ#t < m < ˆm#t # 111Hi(m ) = 1 if ˆm-t ≤ m

we can define a well-behaved approximating consumption function

` ` ` `ct = 111Lo ˇct + 111Mid ˜ct + 111Hi ˆct. (62 )

This just says that, for each interval, we use the approximation that is most appropriate. The function is continuous and once-differentiable everywhere, and is therefore well behaved for computational purposes.

To construct the corresponding refined representation of the value function we must first clarify one point: The upper-bound value function that we are constructing will be the one implied by a consumer whose spending behavior is consistent with the refined upper-bound consumption rule.

For $# mt≥m t$ , this consumption rule is the same as before, so the constructed upper-bound value function is also the same. However, for values $m<m# tt$ matters are slightly more complicated.

Start with the fact that at the cusp point,

# # T ¯vt(m t ) = u( ¯ct(m t ))ℂ t # T = u( ▴m t ¯κt)ℂ t .

But for all $mt$ ,

¯vt(m ) = u(¯ct(m )) + ¯𝔳t(m - ¯ct(m )),

and we assume that for the consumer below the cusp point consumption is given by $¯κ▴m t$ so for $m < m# t t$

¯v (m ) = u(κ¯ ▴m ) + 𝔳¯ ((1 - κ¯ )▴m ), t t t t

which is easy to compute because $𝔳t(at) = β¯vt+1 (atR + 1)$ where $¯vt$ is as defined above because a consumer who ends the current period with assets exceeding the lower bound will not expect to be constrained next period. (Recall again that we are merely constructing an object that is guaranteed to be an upper bound for the value that the ‘realist’ consumer will experience.) At the gridpoints defined by the solution of the consumption problem can then construct

¯Λ (m ) = ((1 - ρ)v¯ (m ))1∕(1- ρ) t t

and its derivatives which yields the appropriate vector for constructing $ˇ X$ and $ˇϘ$ . The rest of the procedure is analogous to that performed for the consumption rule and is thus omitted for brevity.

5.14 Extension: A Stochastic Interest Factor

Thus far we have assumed that the interest factor is constant at $R$ . Extending the previous derivations to allow for a perfectly forecastable time-varying interest factor $Rt$ would be trivial. Allowing for a stochastic interest factor is less trivial.

The easiest case is where the interest factor is i.i.d.,

log R ~ N (r + ϕ - σ2 ∕2, σ2) ∀ n > 0 (63 ) t+n r r

where $ϕ$ is the risk premium and the $2 σr ∕2$ adjustment to the mean log return guarantees that an increase in $σ2 r$ constitutes a mean-preserving spread in the level of the return.

This case is reasonably straightforward because Merton (1969) and Samuelson (1969) showed that for a consumer without labor income (or with perfectly forecastable labor income) the consumption function is linear, with a MPC¹⁹

( )1 ∕ρ κ = 1 - β E [R1 - ρ ] (64 ) t t+1

and in this case the previous analysis applies once we substitute this MPC for the one that characterizes the perfect foresight problem without rate-of-return risk.

The more realistic case where the interest factor has some serial correlation is more complex. We consider the simplest case that captures the main features of empirical interest rate dynamics: An AR(1) process. Thus the specification is

r - r = (r - r)γ + ϵ (65 ) t+1 t t+1

where $r$ is the long-run mean log interest factor, $0 < γ < 1$ is the AR(1) serial correlation coefficient, and $ϵt+1$ is the stochastic shock.

The consumer’s problem in this case now has two state variables, $m t$ and $rt$ , and is described by

1- ρ vt(mt,rt) = max u(ct) + Et [βt+1 Γt+1 vt+1(mt+1, rt+1 )] (66 ) ct s.t. at = mt - ct r - r = (r - r)γ + ϵ t+1 t t+1 Rt+1 = exp (rt+1 ) mt+1 = (Rt+1--∕Γ-t+1) at + θt+1. ◟ ◝◜ ◞ ≡Rt+1

We approximate the AR(1) process by a Markov transition matrix using standard techniques. The stochastic interest factor is allowed to take on 11 values centered around the steady-state value $r$ and chosen [how?]. Given this Markov transition matrix, conditional on the Markov AR(1) state the consumption functions for the ‘optimist’ and the ‘pessimist’ will still be linear, with identical MPC’s that are computed numerically. Given these MPC’s, the (conditional) realist’s consumption function can be computed for each Markov state, and the converged consumption rules constitute the solution contingent on the dynamics of the stochastic interest rate process.

In principle, this refinement should be combined with the previous one; further exposition of this combination is omitted here because no new insights spring from the combination of the two techniques.

5.15 Imposing ‘Artificial’ Borrowing Constraints

Optimization problems often come with additional constraints that must be satisfied. Particularly common is an ‘artificial’ liquidity constraint that prevents the consumer’s net worth from falling below some value, often zero.²⁰

With such an additional constraint, the problem is

1- ρ vT-1(mT - 1) = macxT-1 u(cT - 1) + ET - 1[β ΓT vT (mT )] s.t. aT - 1 = mT - 1 - cT - 1 mT = RT aT - 1 + θT aT - 1 ≥ 0.

By definition, the constraint will bind if the unconstrained consumer would choose a level of spending that would violate the constraint. Here, that means that the constraint binds if the $cT - 1$ that satisfies the unconstrained FOC

c - ρ = 𝔳 ′ (m - c ) (67 ) T- 1 T- 1 T- 1 T - 1

is greater than $mT - 1$ . Call $* `cT - 1$ the approximated function returning the level of $c T-1$ that satisfies (67). Then the approximated constrained optimal consumption function will be

`c (m ) = min [m ,`c * (m )]. (68 ) T- 1 T- 1 T - 1 T- 1 T- 1

The introduction of the constraint also introduces a sharp nonlinearity in all of the functions at the point where the constraint begins to bind. As a result, to get solutions that are anywhere close to numerically accurate it is useful to augment the grid of values of the state variable to include the exact value at which the constraint ceases to bind. Fortunately, this is easy to calculate. We know that when the constraint is binding the consumer is saving nothing, which yields marginal value of $′ 𝔳T-1(0)$ . Further, when the constraint is binding, $cT - 1 = mT - 1$ . Thus, the largest value of consumption for which the constraint is binding will be the point for which the marginal utility of consumption is exactly equal to the (expected, discounted) marginal value of saving 0. We know this because the marginal utility of consumption is a downward-sloping function and so if the consumer were to consume $ϵ$ more, the marginal utility of that extra consumption would be below the (discounted, expected) marginal utility of saving, and thus the consumer would engage in positive saving and the constraint would no longer be binding. Thus the level of $mT - 1$ at which the lconstraint stops binding is:²¹

u ′(m ) = 𝔳′ (0) T- 1 T- 1 mT - 1 = (𝔳′T - 1(0)) (- 1∕ρ) = 𝔠T- 1(0). (69 )

Figure 18: Constrained (solid) and Unconstrained (dashed) Consumption

The constrained problem is solved by 2periodIntExpFOCInvPesReaOptCon.m; the resulting consumption rule is shown in Figure 18. For comparison purposes, the approximate consumption rule from Figure 18 is reproduced here as the solid line. The presence of the liquidity constraint requires three changes to the procedures outlined above:

We redefine $𝔥t$ , which now is the PDV of receiving $θt+1 = θ-$ next period and $θt+n = 0 ∀ n > 1$ – that is, the pessimist believes he will receive nothing beyond period $t + 1$
We augment the end-of-period aVec with zero and with a point with a small positive value so that the generated mVec will the binding point $m#$ and a point just above it (so that we can better capture the curvature around that point)
We redefine the optimal consumption rule as in equation (68). This ensures that the liquidity-constrained ‘realist’ will consume more than the redefined ‘pessimist,’ so that we will have $ϙ$ still between $0$ and $1$ and the ‘method of moderation’ will proceed smoothly.

As expected, the liquidity constraint only causes a divergence between the two functions at the point where the optimal unconstrained consumption rule runs into the 45 degree line.

6 Recursion

6.1 Theory

Before we solve for periods earlier than $T - 1$ , we assume for convenience that in each such period a liquidity constraint exists of the kind discussed above, preventing $c$ from exceeding $m$ . This simplifies things a bit because now we can always consider an aVec that starts with zero as its smallest element.

Recall now equations (17) and (18):

′ - ρ ′ 𝔳 t(at) = Et[βR Γt+1u (ct+1 (Rt+1at + θt+1 ))] ′ ′ u (ct) = 𝔳t(mt - ct).

Assuming that the problem has been solved up to period $t + 1$ (and thus assuming that we have an approximated $`c (m ) t+1 t+1$ ), our solution method essentially involves using these two equations in succession to work back progressively from period $T - 1$ to the beginning of life. Stated generally, the method is as follows. (Here, we use the original, rather than the “refined,” method for constructing consumption functions; the generalization of the algorithm below to use the refined method presents no difficulties.)

For the grid of values $at,i$ in aVec $t$ , numerically calculate the values of $𝔠t(at,i)$ and $𝔠′t(at,i)$ ,
$′ - 1∕ρ 𝔠t,i = (𝔳 t(at,i)) , (70) ( [ - ρ - ρ]) - 1∕ρ = β Et RΓ t+1(`ct+1 (Rt+1at,i + θt+1) ) , (71) ′ - 1- 1∕ρ ′′ ϰt,i = - (1 ∕ρ) (𝔳t(at,i) ) 𝔳t(at,i), (72)$
generating vectors of values $⃗𝔠t$ and $⃗ϰt$ (where the letter $ϰ$ is a variant of $κ$ ; we need a variant because $κ$ itself is reserved for the marginal propensity to consume as of the beginning of the period, and here we are calculating the marginal propensity to have consumed).
Construct a corresponding list of values of $ct,i$ and $mt,i$ from $ct,i=𝔠t,i$ and $mt,i = ct,i + at,i$ ; similarly construct a corresponding list of $κt,i$ using equation (53).
Construct a corresponding list of $μt,i$ , the levels and first derivatives of $ϙt,i$ , and the levels and first derivatives of $χt,i$ .
Construct an interpolating approximation $`χ t$ that smoothly matches both the level and the slope at those points.
If we are to approximate the value function, construct a corresponding list of values of $vt,i$ , the levels and first derivatives of $Ϙt,i$ , and the levels and first derivatives of $ˆ Xt,i$ ; and construct an interpolating approximation $Xˆt$ that matches those points.

With $`χt$ in hand, our approximate consumption function is computed directly from the appropriate substitutions in (45) and related equations. With this consumption rule in hand, we can continue the backwards recursion to period $t - 1$ and so on back to the beginning of life.

Note that this loop does not contain steps for constructing $ˆ𝔳′t(at)$ or $ˆv′(m) tt$ . This is because with $ˆ𝔠 (a ) t t$ and $`ˆc (m ) t t$ in hand, we simply define $ˆ𝔳′t(at)= [ˆ𝔠t(at)]- ρ$ and $ˆv ′t(mt ) = u ′(`ˆct(mt ))$ so there is no need to construct interpolating approximations to these functions - they arise ‘free’ (or nearly so) from our constructed $ˆ𝔠t(at )$ and $`ˆc(m) tt$ .

The program multiperiodCon.m²² presents a fairly general and flexible approach to solving problems of this kind. The essential structure of the program is a loop that simply works its way back from an assumed last period of life, using the command AppendTo to record the interpolated $`χt$ functions in the earlier time periods back from the end. For a realistic life cycle problem, it would also be necessary at a minimum to calibrate a nonconstant path of expected income growth over the lifetime that matches the empirical profile; allowing for such a calibration is the reason we have included the ${ Γ }T t$ vector in our computational specification of the problem.

6.2 Mathematica Background

Mathematica has several features that are useful in solving the multiperiod problem.

It can treat a user-created function as an object just like a number or a character.
Mathematica uses the ‘list’ as its basic data structure. A Mathematica ‘list’ is a very powerful and flexible data construct. A list of length N in Mathematica can hold essentially anything in each of its N positions - a function, a number, another list, a symbolic expression, or any other object that Mathematica can recognize. The items at position $i$ in a list named ExampleList are retrieved or addressed using the syntax ExampleList[[i]].
The function Apply[FuncName_, DataListName_] takes the function whose name is FuncName (for example, Vt) and the data in DataListName (for example, ${1, 19}$ ) and returns the result that would have been returned by calling the function Vt[1,19].
The function Map[FuncToApply_,DataToApplyItTo_] takes a list of possible arguments to the function FuncToApply and applies that function to each of the elements of that list sequentially. For example, Map[Sin,{1,2,3}] would return a list {Sin[1],Sin[2],Sin[3]}.

6.3 Program Structure

After the usual initializations, the heart of the program works like this.

6.3.1 Iteration

After setting up a variable PeriodsToSolve which defines the total number of periods that the program will solve, the program sets up a “Do[SolveAnotherPeriod,{PeriodsToSolve}]” loop that runs the function SolveAnotherPeriod the number of times corresponding to PeriodsToSolve. Every time SolveAnotherPeriod is run, the interpolated consumption function for one period of life earlier is calculated. The structure of the SolveAnotherPeriod function is as follows:

Add various period-t parameters to their respective lifecycle lists, which is accomplished by calling the AddNewPeriodToParamLifeDates function.
For each $a$ in aVec, construct $𝔠$ as follows:
$( [ - ρ - ρ]) - 1∕ρ 𝔠t(at) = β Et R Γ t+1(`ˆct+1(Rt+1at + θt+1)) (73) ( ) ∑n [ ] - 1∕ρ = β 1- R Γ - ρ(`ˆc (R a + θ ))- ρ . (74) n t+1 t+1 t+1 t i i=1$
Obviously, $𝔠t(a)$ depends on the constructed consumption function from one period later in life. We also construct the corresponding mVec, $κ$ Vec, etc. by calling the AddNewPeriodToSolvedLifeDates function.
For each $m$ in mVec, we can define $▴m$ Vec, find the corresponding optimal consumption vector for a pessimist and an optimist, construct the $ϙ$ and $χ$ vectors, and finally an interpolation function $`χt$ . Similarly we can construct an interpolation function $`Xt$ that approximates the value function. The whole process is done by calling the AddNewPeriodToSolvedLifeDatesPesReaOpt function.
Various period- $t$ functions are derived from $χ`t$ and $`Xt$ (in functions_ConsNVal.m). Note that the liquidity constraint is dealt with by comparing the unconstrained solution $c$ From $χ$ with the 45 degree line.

6.4 Results

As written, the program creates $`χt (μt)$ functions from which the relevant $`ct(mt )$ functions are recovered in any period for any value of $m$ .

Figure 19: Converging $`cT-n (m )$ Functions as $n$ Increases

As an illustration, Figure 19 shows $`cT - n (m )$ for $n = {20, 15, 10, 5,1 }$ . At least one feature of this figure is encouraging: the consumption functions converge as the horizon extends, something that Carroll (2011) shows must be true under certain parametric conditions that are satisfied by the baseline parameter values being used here.

7 Multiple Control Variables

We now consider how to solve problems with multiple control variables. (To reduce notational complexity, in this section we set $Γ t = 1 ∀ t$ .)

7.1 Theory

The new control variable that the consumer can now choose is the portion of the portfolio to invest in risky assets. Designating the gross return on the risky asset as $Rt+1$ , and using $ςt$ to represent the proportion of the portfolio invested in this asset between $t$ and $t + 1$ (restricted here, as often in the literature, to values between 0 and 1, corresponding to an assumption that the consumer cannot be ‘net short’ and cannot issue net equity), the overall return on the consumer’s portfolio between $t$ and $t + 1$ will be:

ℜℜℜt+1 = R (1 - ςt) + Rt+1 ςt (75 ) = R + (Rt+1 - R )ςt (76 )

and the maximization problem is

v (m ) = max u (c ) + β E [v (m )] t t {ct,ςt} t t t+1 t+1 s.t. ℜℜℜt+1 = R + (Rt+1 - R )ςt mt+1 = (mt - ct)ℜℜℜt+1 + θt+1 0 ≤ ςt ≤ 1,

vt(mt) = max u (ct) + Et [βvt+1 ((mt - ct)ℜℜℜt+1 + θt+1 )] {ct,ςt} s.t. 0≤ ςt ≤ 1.

The first order condition with respect to $c t$ is almost identical to that in the single-control problem, equation (10), with the only difference being that the nonstochastic interest factor $R$ is now replaced by $ℜℜℜt+1$ ,

u ′(c ) = β E [ℜℜℜ v ′ (m )], (77 ) t t t+1 t+1 t+1

and the Envelope theorem derivation remains the same, yielding the Euler equation for consumption

′ ′ u (ct) = Et[β ℜℜℜt+1u (ct+1 )]. (78 )

The first order condition with respect to the risky portfolio share is

′ 0 = Et [vt+1(mt+1 )(Rt+1 - R )at] = a E [u′ (c (m ) )(R - R )]. (79 ) t t t+1 t+1 t+1

As before, it will be useful to define $𝔳 t$ as a function that yields the expected $t+1$ value of ending period $t$ in a given state. However, now that there are two control variables, the expectation must be defined as a function of the chosen values of both of those variables, because expected end-of-period value will depend not just on how much the agent saves, but also on how the saved assets are allocated between the risky and riskless assets. Thus we define

𝔳 (a ,ς ) = E [βv (m )] t t t t t+1 t+1

which has derivatives

a m 𝔳 t = Et[β ℜℜℜt+1v t+1 (mt+1 )] 𝔳 ς = E [β (R - R)vm (m )]a t t t+1 t+1 t+1 t

implying that the first order conditions (78) and (79) and can be rewritten

u′(c ) = 𝔳a(m - c ,ς ) (80 ) t t t t t 0 = 𝔳ςt(at, ςt). (81 )

7.2 Application

Our first step is to specify the stochastic process for $Rt+1$ . We follow the common practice of assuming that returns are lognormally distributed, $2 2 log R ~ N (ϕ + r - σ ϕ∕2, σϕ )$ where $ϕ$ is the equity premium over the returns $r$ available on the riskless asset.²³

As with labor income uncertainty, it is necessary to discretize the rate-of-return risk in order to have a problem that is soluble in a reasonable amount of time. We follow the same procedure as for labor income uncertainty, generating a set of $m$ equiprobable shocks to the rate of return; in a slight abuse of notation, we will designate the portfolio-weighted return (contingent on the chosen portfolio share in equity, and potentially contingent on any other aspect of the consumer’s problem) simply as $ℜℜℜi,j$ (where dependence on $i$ is allowed to permit the possibility of nonzero correlation between the return on the risky asset and the shock to labor income (for example, in recessions the stock market falls and labor income also declines).

The direct expressions for the derivatives of $𝔳 t$ are

( ) n m a 1 ∑ ∑ - ρ 𝔳t(at,ςt) = β ---- ℜℜℜi,j (ct+1(ℜℜℜi,j at + θi)) (82 ) mn i=1 j=1 ( ) n m ς 1 ∑ ∑ - ρ 𝔳t(at,ςt) = β ---- (Ri,j - R )(ct+1 (ℜℜℜi,jat + θi)) .(83 ) mn i=1 j=1

Writing these equations out explicitly makes a problem very apparent: For every different combination of ${at, ςt}$ that the routine wishes to consider, it must perform two double-summations of $m × n$ terms. Once again, there is an inefficiency if it must perform these same calculations many times for the same or nearby values of ${at, ςt}$ , and again the solution is to construct an approximation to the derivatives of the $𝔳$ function.

Details of the construction of the interpolating approximation are given below; assume for the moment that we have the approximations $ˆ𝔳at$ and $ˆ𝔳 ςt$ in hand and we want to proceed. As noted above, nonlinear equation solvers (including those built into Mathematica) can find the solution to a set of simultaneous equations. Thus we could ask Mathematica to solve

c-t ρ = ˆ𝔳at(mt - ct,ςt) (84 ) ς 0 = ˆ𝔳t(mt - ct,ςt) (85 )

simultaneously for $c$ and $ς$ at the set of potential $mt$ values defined in mVec. However, multidimensional constrained maximization problems are difficult and sometimes quite slow to solve. There is a better way. Define the problem

𝔴 (a ) = max 𝔳 (a , ς) (86 ) t t ςt t t t s.t. 0 ≤ ςt ≤ 1 (87 )

where the bar accent on $𝔳$ indicates that this is the $𝔳$ that has been optimized with respect to all of the arguments other than the one still present ( $a t$ ). We solve this problem for the set of gridpoints in aVec and use the results to construct the interpolating function $`𝔴a(at) t$ .²⁴ With this function in hand, we can use the first order condition from the single-control problem

c- ρ = `𝔴a (mt - ct) t t

to solve for the optimal level of consumption as a function of $mt$ . Thus we have transformed the multidimensional optimization problem into a sequence of two simple optimization problems for which solutions are much easier and more reliable.

Note the parallel between this trick and the fundamental insight of dynamic programming: Dynamic programming techniques transform a multi-period (or infinite-period) optimization problem into a sequence of two-period optimization problems which are individually much easier to solve; we have done the same thing here, but with multiple dimensions of controls rather than multiple periods.

7.3 Implementation

The program which solves the constrained problem with multiple control variables is multicontrolCon.m.

Some of the functions defined in multicontrolCon.m correspond to the derivatives of $𝔳t (at,ςt)$ .

The first function definition that does not resemble anything in multiperiod.m is $ς$ Raw[at_]. This function, for its input value of $at$ , calculates the value of the portfolio share $ςt$ which satisfies the first order condition (85), tests whether the optimal portfolio share would violate the constraints, and if so resets the portfolio share to the constrained optimum. The function returns the optimal value of the portfolio share itself, $* ςt$ , from which the functions $¯ a 𝔳 t(at)$ and $ˆςt(at)$ will be constructed.

As $ˆςt(at)$ can be constructed by $ς$ Raw[at_], $¯a 𝔳t(at )$ is constructed by another newly defined function $𝔳$ aOpt[at_], where the naming convention is obviously that ‘Opt’ stands for ‘Optimized.’ With $¯a 𝔳t(at)$ in hand (as well as the appropriately redefined $¯𝔳t(at)$ and $¯𝔳ata(at )$ ) the analysis is essentially identical to that for the standard multiperiod problem with a single control variable.

The structure of the program in detail is as follows. First, perform the usual initializations. Then initialize $ς$ Vec and the other variables specific to the multiple control problem.²⁵ In particular, there are now three kinds of functions: those with both $a t$ and $ςt$ as arguments, those with just $at$ , and those with $mt$ .

Once the setup is complete, the heart of the program is the following.

Construct $ς 𝔳t(at,ςt)$ using the usual calculation over the tensor defined by the combinations of the elements of aVec and $ς$ Vec.
For any level of saving at, the function $ς$ Raw[at_] performs a rootfinding operation²⁶
$0 = 𝔳 ς(at,ςt) (88) t s.t. 0 ≤ ςt ≤ 1 (89)$
and generates the corresponding optimal portfolio share $ς*t$ .
Construct the function $𝔴$ a[at_]
$𝔴a (a ) ≡ 𝔳a(a ,ς*(a )) (90) t t t t t t$
where $* ςt (at)$ is computed by $ς$ Raw[at_].
Using $a 𝔴 t(at) ≡$ $𝔴$ a[at_] and the redefined $𝔴t (at)$ and $𝔴aa(a) tt$ (in place of $𝔳a (a ) ≡ t t$ $𝔳$ a[at_] in multiperiod.m), follow the same procedures as in multiperiod.m to generate $`c(m) t$ .

7.4 Results

Figure 20 plots the first-period consumption function generated by the program; qualitatively it does not look much different from the consumption functions generated by the program without portfolio choice. Figure 21 plots the optimal portfolio share as a function of the level of assets. This figure exhibits several interesting features. First, even with a coefficient of relative risk aversion of 6, an equity premium of only 4 percent, and an annual standard deviation in equity returns of 15 percent, the optimal choice is for the agent to invest a proportion 1 (100 percent) of the portfolio in stocks (instead of the safe bank account with riskless return $R$ ) is at values of $at$ less than about 2. Second, the proportion of the portfolio kept in stocks is declining in the level of wealth - i.e., the poor should hold all of their meager assets in stocks, while the rich should be cautious, holding more of their wealth in safe bank deposits and less in stocks. This seemingly bizarre (and highly counterfactual) prediction reflects the nature of the risks the consumer faces. Those consumers who are poor in measured financial wealth are likely to derive a high proportion of future consumption from their labor income. Since by assumption labor income risk is uncorrelated with rate-of-return risk, the covariance between their future consumption and future stock returns is relatively low. By contrast, persons with relatively large wealth will be paying for a large proportion of future consumption out of that wealth, and hence if they invest too much of it in stocks their consumption will have a high covariance with stock returns. Consequently, they reduce that correlation by holding some of their wealth in the riskless form.

Figure 20: $c(m1 )$ With Portfolio Choice

Figure 21: Portfolio Share in Risky Assets in First Period $ς(a)$

8 The Infinite Horizon

All of the solution methods presented so far have involved period-by-period iteration from an assumed last period of life, as is appropriate for life cycle problems. However, if the parameter values for the problem satisfy certain conditions (detailed in Carroll (2011)), the consumption rules (and the rest of the problem) will converge to a fixed rule as the horizon (remaining lifetime) gets large, as illustrated in Figure 19. Furthermore, Deaton (1991), Carroll (1992; 1997) and others have argued that the ‘buffer-stock’ saving behavior that emerges under some further restrictions on parameter values is a good approximation of the behavior of typical consumers over much of the lifetime. Methods for finding the converged functions are therefore of interest, and are dealt with in this section.

Of course, the simplest such method is to solve the problem as specified above for a large number of periods. This is feasible, but there are much faster methods.

8.1 Convergence

In solving an infinite-horizon problem, it is necessary to have some metric that determines when to stop because a solution that is ‘good enough’ has been found.

A natural metric is defined by the unique ‘target’ level of wealth that Carroll (2011) proves will exist in problems of this kind: The $mˇ$ such that

Et [mt+1 ∕mt ] = 1 if mt = mˇ

(91)

where the $∨$ accent is meant to signify that this is the value that other $m$ ’s ‘point to.’

Given a consumption rule $c(m )$ it is straightforward to find the corresponding $mˇ$ . So for our problem, a solution is declared to have converged if the following criterion is met: $| ˇmt+1 - ˇmt | < ϵ$ , where $ϵ$ is a very small number and measures our degree of convergence tolerance.

Similar criteria can obviously be specified for other problems. However, it is always wise to plot successive function differences and to experiment a bit with convergence criteria to verify that the function has converged for all practical purposes.

8.2 Coarse then Fine $θ$ Vec

The speed of solution is roughly proportionate²⁷ to the number of points used in approximating the distribution of shocks. At least 3 gridpoints should probably be used as an initial minimum, and my experience is that increasing the number of gridpoints beyond 7 generally yields only very small changes in the solution. The program multiperiodCon_infhor.m begins with three gridpoints, and then solves for successively finer $θ$ Vec.

9 Structural Estimation

This section describes how to use the methods developed above to structurally estimate a life-cycle consumption model, following closely the work of Cagetti (2003).²⁸ The key idea of structural estimation is to look for the parameter values (for the time preference rate, relative risk aversion, or other parameters) which lead to the best possible match between simulated and empirical moments. (The code for the structural estimation is in the self-contained subfolder StructuralEstimation in the Matlab and Mathematica directories.)

9.1 Life Cycle Model

The decision problem for the household at age $t$ is:

{ [ ] } ∑ T ( ) max u (ccct) + Et ℶs- t Πsi=t+1 ˆβi/𝔇/i u (cccs) s=t+1

(92)

subject to the constraints

aaa = mmm - ccc s s s mmms+1 = Raaas + Ys+1 Ys+1 = ppps+1 θs+1 ppps+1 = Γ s+1ppps Ψs+1

where

/𝔇:probability alive (not dead ) until age s given alive at age s - 1 /s ˆβs:time-varying discount factor between age s - 1 and s Ψs:mean -one shock to permanent income ℶ:time-invariant discount factor

and all the other variables are defined as in section 2.

Households start life at age $s = 25$ and live with probability 1 until retirement ( $s = 65$ ). Thereafter the survival probability shrinks every year and agents are dead by $s = 91$ as assumed by Cagetti. Note that in addition to a typical time-invariant discount factor $ℶ$ , there is a time-varying discount factor $βˆs$ in (92) which captures the effect of time-varying demographic variables (e.g. changes in family size).

Transitory and permanent shocks are distributed as follows:

( | 0 with probability ℘ > 0 { Ξs= θs∕℘ with probability (1 - ℘), where log θs ∽ N (- σ2θ∕2, σ(2θ9)3) |( 2 2 logψs∽N ( - σψ ∕2,σ ψ) (94)

where $℘$ is the probability of unemployment (and unemployment shocks are turned off after retirement).

The parameter values for the shocks are taken from Carroll (1992), $℘ = 0.5∕100$ , $σ θ = 0.1$ , and $σψ=0.1$ .²⁹ The income growth profile $Γ s$ is from Carroll (1997) and the values of $//𝔇s$ and $βˆs$ are obtained from Cagetti (2003) (Figure 22).³⁰ The interest rate is assumed to equal $1.03$ . The model parameters are included in Table 1.

Figure 22: Time Varying Parameters

Table 1: Parameter Values

-----------------------
σθ0.1 Carroll (1992)
σψ0.1 Carroll (1992)
℘0.005 Carroll (1992)
Γsfigure 22 Carroll (1997)
ˆ
βs,//𝔇sfigure 22 Cagetti (2003)
R1.03-----Cagetti-(2003)--

The parameters $ℶ$ and $ρ$ are structurally estimated following the procedure described below.

9.2 Estimation

When economists say that they are performing “structural estimation” of a model like this, they mean that they have devised a formal procedure for searching for values for the parameters $ℶ$ and $ρ$ at which some measure of the model’s outcome (like “median wealth by age”) is as close as possible to an empirical measure of the same thing. Here, we choose to match the median of the wealth to permanent income ratio across 7 age groups, from age $26 - 30$ up to $56-60$ .³¹ The choice of matching the medians rather the means is motivated by the fact that the wealth distribution is much more concentrated at the top than the model is capable of explaining using a single set of parameter values. This means that in practice one must pick some portion of the population who one wants to match well; since the model has little hope of capturing the behavior of Bill Gates, but might conceivably match the behavior of Homer Simpson, we choose to match medians rather than means.

As explained in section 3, it is convenient to work with the normalized version the model which can be written as:

{ } v(m)= max u (c ) + ℶ//𝔇 βˆ E [(ψ Γ )1- ρv (m )] tt ct t t+1 t+1 t t+1 t+1 t+1 t+1 s.t. at= mt - ct ( ) m= a ----R----- + θ t+1 t ψt+1 Γ t+1 t+1 ◟----◝ ◜-----◞ ≡Rt+1

with the first order condition:

u′(ct)= ℶ//𝔇t+1 ˆβt+1R Et [u ′ (ψt+1 Γ t+1ct+1 (atRt+1 + θt+1))] . (95 )

The first step is to solve for the consumption functions at each age using the routines included in the setup_ConsFn.m file. We need to discretize the shock distribution and solve for the policy functions by backward induction using equation (95) following the procedure in sections 5 and 6 (ConstructcFuncLife). The latter routine is slightly complicated by the fact that we are considering a life-cycle model and therefore the growth rate of permanent income, the probability of death, the time-varying discount factor and the distribution of shocks will be different across the years. We thus must ensure that at each backward iteration the right parameter values are used.

Once we have the age varying consumption functions, we can proceed to generate the simulated data and compute the simulated medians using the routines defined in the setup_Sim.m file. We first have to draw the shocks for each agent and period. This involves discretizing the shock distribution for as many points as the number of agents we want to simulate (ConstructShockDistribution). We then randomly permute this shock vector as many times as we need to simulate the model for, thus obtaining a time varying shock for each agent (ConstructSimShocks). This is much more time efficient than drawing at each time from the shock distribution a shock for each agent, and also ensures a stable distribution of shocks across the simulation periods even for a small number of agents. (Similarly, in order to speed up the process, at each backward iteration we compute the consumption function and other variables as a vector at once.) Then, following Cagetti (2003), we initialize the wealth-to-income ratio of agents at age $25$ by randomly assigning the equal probability values to $0.17$ , $0.50$ and $0.83$ and run the simulation (Simulate). In particular we consider a population of agents at age 25 and follow their consumption and wealth accumulation dynamics as they reach the age of $60$ , using the appropriate age-specific consumption functions and the age-varying parameters. The simulated medians are obtained by taking the medians of the wealth to income ratio of the $7$ age groups.

Given these simulated medians, we can estimate the model by calculating empirical medians and measure the model’s success by calculating the difference between the empirical median and the actual median. Specifically, defining $ξ$ as the set of parameters to be estimated (in the current case $ξ={ρ,ℶ }$ ), we could search for the parameter values which solve

∑7 τ τ minξ |ς - s (ξ)| (96 ) τ=1

where $ςτ$ and $sτ$ are respectively the empirical and simulated medians of the wealth to permanent income ratio for age group $τ$ .

A drawback of proceeding in this way is that it treats the empirically estimated medians as though they reflected perfect measurements of the truth. Imagine, however, that one of the age groups happened to have (in the consumer survey) four times as many data observations as another age group; then we would expect the median to be more precisely estimated for the age group with more observations; yet (96) assigns equal importance to a deviation between the model and the data for all age groups.

We can get around this problem (and a variety of others) by instead minimizing a slightly more complex object:

∑N min ωi |ςτi - sτ(ξ) | ξ i

(97)

where $ωi$ is the weight of household $i$ in the entire population,³² and $τ ςi$ is the empirical wealth-to-permanent-income ratio of household $i$ whose head belongs to age group $τ$ . $ωi$ is needed because unequal weight is assigned to each observation in the Survey of Consumer Finances (SCF). The absolute value is used since the formula is based on the fact that the median is the value that minimizes the sum of the absolute deviations from itself.

The actual data are taken from several waves of the SCF and the medians and means for each age category are plotted in figure 23. More details on the SCF data are included in appendix A.

Figure 23: Wealth to Permanent Income Ratios from SCF (means (dashed) and medians (solid))

The key function to perform structural estimation is defined in the setup_Estimation.m file as follows:

GapEmpiricalSimulatedMedians [ρ,ℶ ]:= [ ConstructcFuncLife [ρ, ℶ ]; Simulate; ∑ N ωi |ςiτ- sτ( ξ)| i ];

For a given pair of the parameters to be estimated, the GapEmpiricalSimulatedMedians routine therefore:

solves for the consumption functions by calling ConstructcFuncLife
simulates the data and computes the simulated medians by calling Simulate
returns the value of equation (97)

We delegate the task of finding the coefficients that minimize the GapEmpiricalSimulatedMedians function to the Mathematica built-in numerical minimizer FindMinimum. This task can be quite time demanding and rather problematic if the GapEmpiricalSimulatedMedians function has very flat regions or sharp features. It is thus wise to verify the accuracy of the solution, for example by experimenting with a variety of alternative starting values for the parameter search.

Finally the standard errors are computed by bootstrap using the routines in the setup_Bootstrap.m file.³³ This involves:

drawing new shocks for the simulation
drawing a random sample (with replacement) of actual data from the SCF
obtaining new estimates for $ρ$ and $ℶ$

We repeat the above procedure several times (Bootstrap) and take the standard deviation for each of the estimated parameters across the various bootstrap iterations.

The file StructEstimation.m produces our $ρ$ and $ℶ$ estimates with standard errors using 10,000 simulated agents.³⁴ Results are reported in Table 2.³⁵ Figure 24 shows the contour plot of the GapEmpiricalSimulatedMedians function and the parameter estimates. The contour plot shows equally spaced isoquants of the GapEmpiricalSimulatedMedians function, i.e. the pairs of $ρ$ and $ℶ$ which lead to the same deviations between simulated and empirical medians (equivalent values of equation (97)). We can thus interestingly see that there is a large rather flat region, or more formally speaking there exists a broad set of parameter pairs which leads to similar simulated wealth to income ratios. Intuitively, the flatter and larger is this region, the harder it is for the structural estimation procedure to precisely identify the parameters.

Table 2: Estimation Results

Figure 24: Contour Plot (larger values are shown lighter) with ${ρ,ℶ }$ Estimates (red dot).

10 Conclusion

There are many alternative choices that can be made for solving microeconomic dynamic stochastic optimization problems. The set of techniques, and associated programs, described in these notes represents an approach that I have found to be powerful, flexible, and efficient, but other problems may require other techniques. For a much broader treatment of many of the issues considered here, see Judd (1998).

Appendices

A Further Details on SCF Data

Data used in the estimation is constructed using the SCF 1992, 1995, 1998, 2001 and 2004 waves. The definition of wealth is net worth including housing wealth, but excluding pensions and social securities. The data set contains only households whose heads are aged 26-60 and excludes singles, following Cagetti (2003).³⁶ Furthermore, the data set contains only households whose heads are college graduates. The total sample size is 4,774.

In the waves between 1995 and 2004 of the SCF, levels of normal income are reported. The question in the questionnaire is "About what would your income have been if it had been a normal year?" We consider the level of normal income as corresponding to the model’s theoretical object $P$ , permanent noncapital income. Levels of normal income are not reported in the 1992 wave. Instead, in this wave there is a variable which reports whether the level of income is normal or not. Regarding the 1992 wave, only observations which report that the level of income is normal are used, and the levels of income of remaining observations in the 1992 wave are interpreted as the levels of permanent income.

Normal income levels in the SCF are before-tax figures. These before-tax permanent income figures must be rescaled so that the median of the rescaled permanent income of each age group matches the median of each age group’s income which is assumed in the simulation. This rescaled permanent income is interpreted as after-tax permanent income. This rescaling is crucial since in the estimation empirical profiles are matched with simulated ones which are generated using after-tax permanent income (remember the income process assumed in the main text). Wealth / permanent income ratio is computed by dividing the level of wealth by the level of (after-tax) permanent income, and this ratio is used for the estimation.³⁷

References

ATTANASIO, O.P., J. BANKS, C. MEGHIR, AND G. WEBER (1999): “Humps and Bumps in Lifetime Consumption,” Journal of Business & Economic Statistics, 17(1), 22–35.

CAGETTI, MARCO (2003): “Wealth Accumulation Over the Life Cycle and Precautionary Savings,” Journal of Business and Economic Statistics, 21(3), 339–353.

CARROLL, CHRISTOPHER D. (1992): “The Buffer-Stock Theory of Saving: Some Macroeconomic Evidence,” Brookings Papers on Economic Activity, 1992(2), 61–156, http://econ.jhu.edu/people/ccarroll/BufferStockBPEA.pdf.

(1997): “Buffer Stock Saving and the Life Cycle/Permanent Income Hypothesis,” Quarterly Journal of Economics, CXII(1), 1–56, http://econ.jhu.edu/people/ccarroll/BSLCPIH.zip.

(2006): “The Method of Endogenous Gridpoints for Solving Dynamic Stochastic Optimization Problems,” Economics Letters, pp. 312–320, http://econ.jhu.edu/people/ccarroll/EndogenousGridpoints.pdf.

(2011): “Theoretical Foundations of Buffer Stock Saving,” Manuscript, Department of Economics, Johns Hopkins University, http://econ.jhu.edu/people/ccarroll/papers/BufferStockTheory.

(Current): “Math Facts Useful for Graduate Macroeconomics,” Online Lecture Notes.

CARROLL, CHRISTOPHER D., AND MILES S. KIMBALL (1996): “On the Concavity of the Consumption Function,” Econometrica, 64(4), 981–992, http://econ.jhu.edu/people/ccarroll/concavity.pdf.

CARROLL, CHRISTOPHER D., AND ANDREW A. SAMWICK (1997): “The Nature of Precautionary Wealth,” Journal of Monetary Economics, 40(1), 41–71.

DEATON, ANGUS S. (1991): “Saving and Liquidity Constraints,” Econometrica, 59, 1221–1248.

DEN HAAN, WOUTER J, AND ALBERT MARCET (1990): “Solving the Stochastic Growth Model by Parameterizing Expectations,” Journal of Business and Economic Statistics, 8(1), 31–34, Available at http://ideas.repec.org/a/bes/jnlbes/v8y1990i1p31-34.html.

GOURINCHAS, PIERRE-OLIVIER, AND JONATHAN PARKER (2002): “Consumption Over the Life Cycle,” Econometrica, 70(1), 47–89.

HOROWITZ, JOEL L. (2001): “The Bootstrap,” in Handbook of Econometrics, ed. by James J. Heckman, and Edward Leamer, vol. 5. Elsevier/North Holland.

JUDD, KENNETH L. (1998): Numerical Methods in Economics. The MIT Press, Cambridge, Massachusetts.

KOPECKY, KAREN A., AND RICHARD M.H. SUEN (2010): “Finite State Markov-Chain Approximations To Highly Persistent Processes,” Review of Economic Dynamics, 13(3), 701–714, http://www.karenkopecky.net/RouwenhorstPaper.pdf.

MERTON, ROBERT C. (1969): “Lifetime Portfolio Selection under Uncertainty: The Continuous Time Case,” Review of Economics and Statistics, 50, 247–257.

PALUMBO, MICHAEL G (1999): “Uncertain Medical Expenses and Precautionary Saving Near the End of the Life Cycle,” Review of Economic Studies, 66(2), 395–421, Available at http://ideas.repec.org/a/bla/restud/v66y1999i2p395-421.html.

SAMUELSON, PAUL A. (1969): “Lifetime Portfolio Selection by Dynamic Stochastic Programming,” Review of Economics and Statistics, 51, 239–46.

VALENCIA, FABIAN (2006): “Banks’ Financial Structure and Business Cycles,” Ph.D. thesis, Johns Hopkins University.