Coursework 2026 - Consumption Saving in General Equilibrium

Reproducibility

Set a fixed random seed (for example Random.seed!(2026)) before all stochastic computations (Monte Carlo, shock simulations, ergodic-distribution simulations), and report the seed you used.

Part I : Approximations

Interpolation

Consider the function \(f(x)=\frac{sin(|x|)}{|x|+0.0001}\) defined on \([a,b]\) with \(a=-2\) and \(b=2\).

Consider the grid points \(x_i=a+\frac{i}{N-1}(b-a)\) for \(i=0,\ldots,N-1\).

Define in julia a, b, N and the function f . Compute the vector of grid points s (a vector of size N) then the values y (a vector of size N)that f takes at these grid points.

# your code here

Define a function interp(s,y,x) which approximates the original function f at any \(x\in[-2,2]\) but is defined using only s and y (or if you prefer a,b, N and y). (hint: you can use the library Interpolations.jl)

# your code here

On the same plot, represent (s,y) as a scatter plot, the function f and the interpolated values using interp. For the \(x\)-axis choose \([a-0.1, b+0.1]\). How does interp behave outside of \([a,b]\) (i.e. when it extrapolates)? Can you find a way to make it extrapolate linearly?

# your code here

Discretizing a shock

Consider a random shock \(\epsilon\) that follows a normal law with standard deviation \(\sigma=0.1\). Take function \(f(x) = \frac{1}{1+(x)^2}\).

The goal here is just to compute numerically the expectation

\[\mathbb{E} f(\epsilon)\]

Define the variable \(\sigma\) and the function \(f\) in julia.

# your code here

Monte-Carlo

Take an integer \(K>0\).

Try to evaluate the integral above by running \(K\) random draws i.e. by using the formula \(f^{MC, K} = \frac{1}{K}\sum_{i=1}^K f(\epsilon_i)\) where \(\epsilon_i\) is a random draw. This is the Monte-Carlo approach.

# your code here

Evaluate the standard deviation in \(f^{MC, K}\) by evaluating this function repeatedly. How does the standard deviation depend on \(K\) ?

# your code here

Gauss-Hermite

Discretizing the shock consists in finding real numbers \((x_1, ... x_K)\) and positive real numbers \((w_1, ... w_K)\) such that for a “suitable” function \(f(\epsilon)\) we can write:

\[\mathbb{E} f(\epsilon) \approx \sum_{i=1}^K w_i f(e_i)\]

The Gauss-Hermite quadrature is such a discretization scheme. It is widely documented and there are good Julia libraries for it.

Approximate the expectation above using Gauss-Hermite nodes and weights. What is the error? How does it depend on \(K\) for the function you have chosen? Compare with the Monte-Carlo method.

# your code here

Part II : The Ayiagari Model

In this exercise, we try to solve the Ayiagari model.

The Ayiagari model features a continuum of households indexed by \(i\in[0,1]\), who take wage rate \(\color{red} w\) and gross return on savings \(\color{red} r\) as given.

Each household \(i\) is hit by a random productivity shock (\(\color{green}e^i_t\)), normally distributed, with variance \(\sigma\) and mean \(1\) so that, at the aggregate, the average productivity is \(1\).

Each household supplies \(\color{green}e^i_t\) units of work, and can decide to consume \(c^i_t\) or save out of available income \(m^i_{t}\).

Available income \(m^i_t\) thus follows the law of motion:

\[m^i_{t+1} = {\color{red}w } {\color{green}e^i_{t}} + (m^i_{t} - c^i_{t}) {\color{red}r}\]

There is no borrowing so that end-of-period savings \(a_t^i = m_t^i - c_t^i\) satisfy \(a_t^i\geq 0\) (equivalently \(0<c_t^i\leq m_t^i\)). The household chooses consumption optimally in order to maximize:

\[\mathbb{E}_0 \sum_t \beta^t U(c^i_t)\]

with \(U(x)=\frac{x^{1-\gamma}}{1-\gamma}\) and \(\beta=0.96\), \(\gamma=4.0\).

Denoting end-of-period savings by \(a^i_t = m^i_t - c^i_t\), aggregate savings are turned into productive capital:

\[K = \int_i a^i_t\]

Note that in equilibrium the distribution of savings is invariant, hence the absence of a time-subscript for aggregate capital.

Also, the ergodicity property implies that the invariant distribution across all agents is identical to the distribution of available income levels reached across time by any single agent \(i\) so that we also have:

\[K = \mathbb{E} [ a^i_t ]\]

Capital is rented by competitive firms that all produce with the same Cobb-Douglas technology. Total production is:

\[Y = K^{\alpha} L^{1-\alpha}\]

with \(L=\overline{L}=1\) and \(\alpha=0.3\).

Wage rate and interest rate are then determined by marginal conditions:

\[{\color{red}r}=\alpha \frac{Y}{K}\] \[{\color{red}w}=(1-\alpha) \frac{Y}{L}\]

From the structure of the whole model, it is clear that capital \(K\) determines the whole equilibrium.

Note

Because there is no depreciation, the rental rate determined by the marginal conditions, equals the returns on savings when they are turned into productive capital.

Aggregate production

Create a namedtuple (or a structure) mod to hold all model information.

# your code here

Create a function rates(K, mod)::Tuple{Float64,Float64} that takes aggregate capital and returns the rental rate of capital and the wage rate.

# your code here

Choose an initial level of capital \(K_0\) and corresponding rates \(p_0=(r_0,w_0)\) so that \(\beta < \beta r_0 < 1\). Define the corresponding variables K_0 and p_0.

# your code here

Note

The second part of this inequality ensures that wealth distribution is well defined – agents don’t want to accumulate assets until they have unbounded wealth.

Consumption-Savings Problem

We now solve the individual problem given the vector of rates p. The primary method for this homework is the endogenous grid points method (EGM). If you do not manage to obtain a stable EGM implementation, you can submit a value-function-iteration solution instead.

Given the formulation of the household problem, we look for a decision rule for consumption \(c(m)\) with \(c(m)\in]0,m]\).

Note

Since all agents are identical ex ante, we need to solve only one optimization problem. Hence we drop the \(i\) subscripts.

In theory, the state-space is \(]0,\infty]\) but in practice we approximate the system on \(]0,\overline{A}]\) and adjust \(\overline{A}\) until the solution is well defined.

Note

If \(\color{green}e^i_t\) has a compact support, theoretical considerations show that one can choose such an \(\overline{A}\).

Write the Bellman equation

your text there

Approximate space

To represent functions (\(c()\) and \(V()\)) using a finite number of parameters, we discretize the state-space into a finite grid g, made of \(N\) linearly spaced points between \(\epsilon>0\) and \(\overline{A}\) (the role of \(\epsilon>0\) is to avoid undefined behaviour for \(c=0\)).

We use Gauss-Hermite nodes and weights to approximate the shock \(\epsilon\).

Create a namedtuple dis to represent the discretized model (i.e. the grid and the quadrature).

# your code here

Endogenous Grid Points Method

Let \(m\) denote available income (the state variable) and define end-of-period savings by:

\[a = m-c\]

Then next period’s available income is:

\[m' = {\color{red}r} a + {\color{red}w} e'\]

and the Euler equation for an interior solution becomes:

\[U'(c(m)) = \beta {\color{red}r} \, \mathbb{E}\left[ U'(c(m')) \right]\]

With CRRA utility, this can be rewritten as:

\[c(m) = \left( \beta {\color{red}r} \, \mathbb{E}\left[c(m')^{-\gamma}\right] \right)^{-1/\gamma}\]

Note that given the decision rule tomorrow, the term on the right hand side can be explicitly computed for any given \(a\). For this given \(a\) one can then obtain consumption \(c\) today at the corresponding available income \(m\) using identity \[m = c + a\]. Each iteration thus produces new optimal consumption choices for a grid of available income that is endogenous to the problem, hence the name of the method.

Suggested algorithmic structure for each EGM iteration:

Fix a grid agrid.
For each a in agrid, compute all future income realizations \(m'_{\ell}=r a + w e_{\ell}\).
Interpolate tomorrow’s policy \(c(m')\) at those points.
Compute expected marginal utility and invert Euler to get current \(c\).
Recover endogenous available income \(m = c + a\).
Interpolate the updated policy back onto dis.g.
Enforce borrowing-constraint region explicitly.

Note

The borrowing constraint is \(a\geq 0\). When it binds, the household consumes all available resources, so that \(c(m)=m\).

Write the Euler equation associated with the household problem and explain why, with CRRA utility, it can be inverted to recover current consumption from expected marginal utility tomorrow.

your text here

Create a grid agrid for end-of-period savings. For each a in agrid, define a function that computes the corresponding endogenous available income

\[m = a + \left( \beta {\color{red}r} \sum_{\ell} \omega_{\ell} \, c(m'_{\ell})^{-\gamma} \right)^{-1/\gamma}\]

where \(m'_{\ell} = {\color{red}r} a + {\color{red}w} e_{\ell}\) and \((e_{\ell},\omega_{\ell})\) are the quadrature nodes and weights.

# your code here

Write a function egm_step(cvec, mod, p, dis)::Tuple{Vector, Vector} which takes a consumption policy tomorrow, performs one EGM update, and returns:

the endogenous income grid
the updated consumption values on that grid

Then interpolate this updated policy back onto the original grid dis.g.

In your implementation, report the interpolation/extrapolation choice you use and justify it briefly.

# your code here

Explain how you handle the borrowing constraint for available income levels below the smallest endogenous income grid point. Implement this treatment in your code.

# your code here

Write a function egm(mod, p, dis) which iterates on the consumption rule until convergence. Compare the resulting policy function and execution time with value-function iteration and policy-function iteration.

Use and report a clear convergence criterion, for example:

\[\|c^{(n+1)}-c^{(n)}\|_\infty < \eta\]

with your chosen tolerance \(\eta\) and maximum number of iterations.

# your code here

Optional fallback: Value Function Iteration and Policy Function Iteration

If you do not manage to complete EGM, implement VFI as your fallback method. Policy-function iteration (Howard improvement) is optional but strongly recommended for comparison.

Take an initial guess vvec (preferably concave) for the values on the grid and define the function value_update(m, c, vvec, mod, p, dis)::Float64 which computes the expression inside the max in the Bellman equation for any state m, any acceptable consumption level c and a continuation value, obtained by interpolating the values vvec at any point outside the grid.

# your code here

Use your preferred method for constrained optimization to compute the value \(c\in]0,m]\) which maximizes the function value_update for a given m of your choice.

Note

This is the crucial step. You can test whether it works for different values of m and debug what is going on by making a plot.

# your code here

Write a function bellman_step(vvec, mod, p, dis)::Tuple{Vector, Vector} which takes in a vector representing the value function tomorrow and returns a new value vector and policy vector resulting from the maximization.

# your code here

Write a vfi(mod, p, dis) function, which solves the consumption savings problem by value function iteration. It should return the vector of values and the consumption vector.

# your code here

Check that the solution makes sense: plot the solution, the boundaries, check the extrapolation behaviour…

# your code here

(Optional) Solve the consumption savings problem using policy function iteration (i.e. using Howard improvement steps). Compare execution time with value function iteration.

Note

You might need to define some intermediary functions, like value_step which updates the value at all grid points, for a given vector of policy choices and policy_eval which returns the value vector obtained by iterating value_step many times.

# your code here

Computing the Stable Distribution

Now that we have a consumption rule represented by a vector of consumption values cvec, we would like to compute the corresponding ergodic distribution in order to approximate aggregate capital \(K = \int_i a_i = \int_i (m_i - c(m_i))\).

With Monte Carlo

Write a function transition(m, cvec, mod, p) which returns a new random available income, from an initial m. The consumption choices are defined using cvec and can be interpolated as before.

# your code here

Propose a way to draw \(L=1000\) random points from the ergodic distribution. Plot the result.

# your code here

Compute the average of these values. It corresponds to the capital supply. What is the standard deviation of this method? (you can evaluate the standard deviation by performing the procedure again)

# your code here

With a Markov Chain

Write a function transition(m, c, mod, p, dis)::Vector{Tuple{Float64, Float64}} which returns a vector of next-period available income levels with matching probabilities, obtained from initial income level m, with consumption choice c for the various realizations of the discretized shocks (obtained from the quadrature).

# your code here

Use the transition function to define a stochastic matrix P of size \(N\times N\) such that \(P_{ij}\) represents the probability of reaching grid point \(j\) from grid point i. Use the trembling hand method to deal with income levels that are out of the approximation grid.

Trembling hand

When a given available income \(m\) reached with probability \(w\) falls between two grid points \(m_j\) and \(m_{j+1}\), we consider that it reaches \(m_j\) with probability \(w\frac{m_{j+1}-m}{m_{j+1}-m_j}\) and \(m_{j+1}\) with probability \(w\frac{m-m_j}{m_{j+1}-m_j}\). This trick ensures that the stochastic matrix depends smoothly on the decision rule.

# your code here

Compute the stable distribution \(\mu\) of \(P\). Compute the capital demand.

# your code here

Solve the Model

Using everything you have done so far write a function capital_supply(K, mod, p, dis) which computes the rates corresponding to K, solves the consumption problem, computes the ergodic distribution and the corresponding capital supply.

# your code here

Plot capital supply and capital demand (line \(K=K\)) for different levels of capital. Find the level of capital such that market for capital clears.

# your code here

Implement a root-finding routine for market clearing (for example bisection or secant on excess supply). Report convergence diagnostics (stopping criterion, number of iterations, and final residual).

# your code here

Perform a sensitivity analysis with respect to at least two parameters: \(\beta\) and \(\sigma\). For each case, recompute the equilibrium capital and provide a short economic interpretation.

# your code here

Provide two robustness checks of the equilibrium result:

one with a finer income grid size \(N\)
one with a larger number of shock discretization nodes

Comment on how much your equilibrium capital changes.

# your code here