Endogenous Grid Points Method for a model with two endogenous state variables¶

Yuki Yao (Kent)

March 1st 2025

Motivation¶

Many extend Aiyagari-model by adding an extra endogenous state

  • Sweat equity $\kappa$ (Bhandari-McGrattan)
  • Illiquid assets (HANK)
  • Durables $d$ (Guerrieri-Lorenzoni)
  • Human capital $h$ (?)

Optimization with two variables is significantly more challenging

  • Curse of dimensionality
  • No 'bisection' method in more than 1D
  • More time and computational resources

Endogenous Grid Points Method (EGM)¶

EGM is a powerful tool to solve a DP problem

  • Avoid/reduce optimization
  • Achieve higher accuracy
  • But so far, the method applies to a model with one endogenous state ($a^\prime$)

Can we use the EGM for the two dimensional DP problem?

  • can avoid 2D optimization
  • can achieve higher accuracy with fewer points
  • so, gains can be massive in the 2D problem

Caution¶

You need to check if your EGM solution matches with your traditional solution

  • There is no way to proceed without the traditional soluiton

You need to come up with an algorithm specific to your problem

  • painstaking process

Big problem¶

Want to solve a GE model in which

  • Heterogeneity in productivity and states: an extension of Aiyagari economy
  • Occupational choice: households can work for others, or they can run a business
  • Sweat equity/capital: business owners build business relationships and reputation over time

Main computational challenge

  • The problem of business owners has two endogenous state
  • Asset $a$
  • Sweat capital $\kappa$

Small problem¶

Let's zero in on the main challenge

Think of a business owner's problem. The business owner runs a business

  • income flow is $y^{b}\left(z, \kappa \right)$

    • $z$: productivity (stochastic)
    • $\kappa$: owner's sweat capital -- like a client list or business relationship.

  • Building sweat capital $\kappa$ requires effort $e$

  • Also depreciates at a rate $\delta_\kappa$

Dynamic programming with two endogenous states¶

Think of a business owner's problem who stays in business forever

The DP problem is $$ V^{b}\left(a,\kappa; s\right)=\max_{c,e, a^{\prime}\geq0,\kappa^{\prime}}u(c)-g(e)+\beta\sum_{s}\Pi(s^{\prime}|s)V^{b}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) $$ subject to $$ c+a^{\prime} \leq y^{b}\left(z\left(s\right), \kappa\right)+Ra + \bar{T} $$ $$ \kappa^{\prime}=\left(1-\delta_{\kappa}\right)\kappa+e $$

$u(c) = \frac{c^{1-\sigma}}{1-\sigma}$, $g(e) = \frac{e^\eta}{\eta} (\eta > 1)$, $y^{b}\left(z, \kappa \right) = z \kappa^\phi \ (0 < \phi < 1)$

Value Function Iteration with Gridsearch¶

Using the budget constraint and the law of motion, $$ \max_{a^{\prime} \geq 0, \kappa^{\prime}} u( y^{b}\left(\kappa,s\right)+Ra - a^{\prime} )-g(\kappa^{\prime}-\left(1-\delta_{\kappa}\right)\kappa)+\beta\sum_{s}\Pi(s^{\prime}|s)V^{b}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) $$

so we need to max the objective by choosing both $a^\prime$ and $\kappa^\prime$.

How to solve this problem?

For now we use very simple grid search set up grids for $a$ and $\kappa$, $\mathcal{A} = \{a_1, a_2, \cdots, a_{Na}\}$, $\mathcal{K} = \{\kappa_1, \kappa_2, \cdots, \kappa_{N\kappa}\}$

For each $\left(a, \kappa\right) \in \mathcal{A}\times\mathcal{K}$ , find $\left(a^\prime, \kappa^\prime \right) \in \mathcal{A}\times\mathcal{K}$ that maximize the objective

  • Note: You do not necessarily have to do grid search over $\left(a^\prime, \kappa^\prime \right) \in \mathcal{A}\times\mathcal{K}$.

Result from the plain vanilla grid search¶

Plain vanilla grid search is

  • very slow
  • also not accurate since values of the policy function (say $a^\prime$) have to be on the grid

but,...

  • the algorighm is very simple
  • the code often converges without much debugging
  • must be close to the true solution if the grid is very fine

Can we do better (without EGM)?¶

Smarter grid search

  • $a^\prime$ or $\kappa^\prime$ does not have to be on the grid
  • $\mathcal{A}$ or $\mathcal{K}$ does not need to be uniform (sparse grid)
  • Do not necessarily have to look at all the points

Optimize the objective using a multivariable solver

  • Nelder-Mead, Newton, and their variants
  • and need to check if it is globally optimal (multi-start)

Solve the DP by Newton method not by VFI

  • Note you can combine it with grid search/regular optimization/EGM
  • Much faster than VFI but requires a better initial guess
  • Howard algorithm is a handy way to calculate a Newton step.

Endogenous grid points method¶

Start from $\left(a^\prime, \kappa^\prime \right) \in \mathcal{A}\times\mathcal{K}$, use the Euler equations to

  • back out $c, e, a, \kappa$ that rationalizes $\left(a^\prime, \kappa^\prime \right)$
  • find when the borrowing constraint $a^\prime \geq 0$ binds
  • find $\left(a^\prime, \kappa^\prime \right)$ when $a^\prime \geq 0$ binds
  • and obtain the policy and $V$ on $\left(a, \kappa \right) \in \mathcal{A}\times\mathcal{K}$ by interpolation

Challenge

  • borrowing constraint threshold for $A^\prime (a, \kappa) \geq 0$ depends on not only $a$ but also $\kappa$
  • what is the optimal policy $K^\prime (a, \kappa)$ when the BC binds? (for $A^\prime$, it is zero.)
  • how to interpolate the values in 2D when grid points are on an irregular grid?

If this seems unclear, let's revisit the 1D case for better intuition.

Quick review on the EGM in a generic Aiyagari model¶

Think of a worker's problem.

  • income flow is $w \epsilon \bar{n}$ with $\bar{n} = 1$
  • $w$: equilibrium wage rate
  • $\epsilon$: labor productivity (stochastic)
  • $\bar{n}$: labor supply (for now it is inelastic)

Dynamic programming with one endogenous state¶

Think of a worker's problem who will not run a business

The DP problem is $$ V^w\left(a ; s\right)=\max_{c, a^{\prime}\geq0 } u(c) +\beta\sum_{s^\prime}\Pi(s^{\prime}|s) V^w \left(a^{\prime} ;s^{\prime}\right) $$

subject to $$ c+a^{\prime} \leq w \epsilon(s) \bar{n} + Ra + \bar{T} $$

$$ \kappa^{\prime}= \lambda \kappa $$

with $u(c) = \frac{c^{1-\sigma}}{1-\sigma}$.

To update $V$, we need to calculate the optimal policy $A^\prime \left(a; s\right)$ for all $a \in \mathcal{A}$

and we can update values of $V$ on a grid by applying the max operator.

The usual procedure is to find $A^\prime \left(a; s\right) $ for $a \in \mathcal{A}$ by a solver

Euler equation¶

$$ u^\prime(w \epsilon(s) \bar{n} + Ra - a^{\prime} ) \geq \beta\sum_{s^\prime}\Pi(s^{\prime}|s)V^{w}_{a}\left(a^{\prime} ;s^{\prime}\right)\ \text{w.e. if } a^\prime > 0. $$

Let $\mathcal{A} = \{a_1, a_2, \cdots, a_{Na}\}$,

Rather than putting $a$ on $\mathcal{A}$ and then find $a^\prime$

Put $a^\prime$ on $\mathcal{A}$ and then find $a$

In fact, with a guess for $V^w \left(a, ;\right)$ on the right hand side and $u^\prime(c) = c^{-\sigma}$.

\begin{align*} a & = \frac{1}{R} \left( a^{\prime} - w \epsilon(s) \bar{n} + \left(\beta\sum_{s^\prime}\Pi(s^{\prime}|s)V^{w}_{a}\left(a^{\prime} ;s^{\prime}\right)\right)^{-1/\sigma} \right) \\ & \equiv \bar{A}\left(a\prime; s\right) \end{align*}

if $a^\prime > 0$.

For each $a^\prime_i \in \mathcal{A}$ we obtain a correspoinding $a_i$ as $\left[ \bar{A}\left(a\prime_1; s\right), \cdots, \bar{A}\left(a\prime_{N_a}; s\right) \right]$

The grid points for $a$ are obtained endogenously

Graphical illustration¶

When $a^\prime \geq 0$ is not binding at all, finding $A^\prime \left(a; s\right)$ is straightforward.

No description has been provided for this image
No description has been provided for this image
No description has been provided for this image

When $a^\prime \geq 0$ is binding, we need to specify the region where $a^\prime \geq 0$ binds

No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image

Misc.¶

Notice that we only have to check the value of $\bar{A} \left(a^\prime = a_1; s\right)$

  • If $\bar{A} \left(a_1; s\right) < a_1 $: not binding
  • Elif $\bar{A} \left(a_1; s\right) \geq a_1 $: it binds for $a_1 \leq a \leq a^0$.

Now, we are back to the 2D problem with sweat equity

Goals¶

The DP problem is $$ V^{b}\left(a,\kappa; s\right)=\max_{c,e, a^{\prime}\geq0,\kappa^{\prime}}u(c)-g(e)+\beta\sum_{s}\Pi(s^{\prime}|s)V^{b}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) $$ subject to $$ c+a^{\prime} \leq y^{b}\left(\kappa,z\left(s\right) \right)+Ra + \bar{T} $$ $$ \kappa^{\prime}=\left(1-\delta_{\kappa}\right)\kappa+e $$

Assume $V\left(a, \kappa; s\right)$ is construcuted by its values on grid $\mathcal{A}\times\mathcal{K}\times S$

VFI requires finding the optimal policy $A^\prime\left(a, \kappa; s\right)$ and $K^\prime\left(a, \kappa; s\right)$ for each $(a, \kappa, s) \in \mathcal{A}\times\mathcal{K}\times S $

Also we want to interpolate $A^\prime\left(a, \kappa; s\right)$ and $K^\prime\left(a, \kappa; s\right)$ using the values on the grid points.

Euler equations¶

Two Euler equations

$$ u^\prime( y^{b}\left(\kappa,s\right)+Ra + \bar{T} - a^{\prime} ) \geq \beta\sum_{s^\prime}\Pi(s^{\prime}|s)V^{b}_{a}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) \text{ w.e. if } a^\prime > 0 $$

$$ g^\prime(\kappa^{\prime}-\left(1-\delta_{\kappa}\right)\kappa) = \beta\sum_{s^\prime}\Pi(s^{\prime}|s)V^{b}_{\kappa}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) $$

Using our specification $u(c) = \frac{c^{1-\sigma}}{1-\sigma}$, $g(e) = \frac{e^\eta}{\eta}$, we can obtain $a$ and $\kappa$.

$$ \begin{align*} a & = \frac{1}{R} \left( a^{\prime} - y^{b}\left(\kappa,s\right) - \bar{T} + \left(\beta\sum_{s^\prime}\Pi(s^{\prime}|s)V^{b}_{a}\left(a^{\prime}, \kappa^{\prime} ;s^{\prime}\right)\right)^{-1/\sigma} \right) \\ & \equiv \bar{A}\left(a\prime, \kappa^{\prime}; s\right) \\ \kappa & = \frac{1}{1-\delta_\kappa} \left( \kappa^\prime - \left( \beta\sum_{s^\prime}\Pi(s^{\prime}|s)V^{b}_{\kappa}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) \right)^{\frac{1}{\eta - 1}} \right) \\ & \equiv \bar{K}\left(a\prime, \kappa^{\prime}; s\right) \end{align*} $$

if $a^\prime > 0$.

For each $\left(a^\prime, \kappa^\prime \right) \in \mathcal{A}\times\mathcal{K}$, we can obtain $\left( \bar{A}\left(a\prime, \kappa^{\prime}; s\right), \bar{K}\left(a\prime, \kappa^{\prime}; s\right) \right)$, a candidate for $\left(a, \kappa \right)$.

but we need to deal with the borrowing constraint $a^\prime \geq 0$.

Implied $\left(a, \kappa \right)$ for $\left(a^\prime, \kappa^\prime \right) \in \mathcal{A} \times \mathcal{K}$¶

First, let's see all the $(a, \kappa)$ points

Magnify it

Which points are feasible and not feasible?

Any point with $\kappa < 0$ is not feasible. $y^b(z, \kappa)$ requires $\kappa \geq 0$

Any point with $a < 0$ is feasible. It is still consistent with $a^\prime \geq 0$

We want to find values for $A^\prime \left(a, \kappa; s\right)$ and $K^\prime \left(a, \kappa; s\right)$ on ★.

Step 1: Find a threshold for $A^\prime \left(a, \kappa; s\right) \geq 0$¶

As we have seen in the 1D case, need to find thresholds that separetes $\left(a, \kappa \right)$-space

  • In the region with $a^\prime > 0$, $\left(a, \kappa \right)$ is determined by the EE
  • In the region with $a^\prime = 0$, $a^\prime = 0$. What about $\kappa^\prime = $?
  • The thresholds are a set of points $\left(a, \kappa \right)$ that separetes out the two regions.

Fortunately we know the points on the threshold -- all the points from the EE with $a^\prime = 0$

Collect all the $\left(a, \kappa \right) = \left( \bar{A}\left(a\prime = a_1, \kappa^{\prime}; s\right), \bar{K}\left(a\prime = a_1, \kappa^{\prime}; s\right) \right)$ for $\kappa^\prime \in \mathcal{K}$

The threshold is defined by a locus associated with $a^\prime = 0$ with $\kappa ^\prime \in \mathcal{K}$

  • collect all the points with $a^\prime = 0$ with $\kappa ^\prime \in \mathcal{K}$
  • define a locus $A^1 \left(\kappa; s\right)$ by interpolating these points

We can separete $\left(a, \kappa \right)$ - space into two regions

\begin{align*} &\text{If } a \geq A^1(\kappa; s), && \text{not binding} \\ &\text{If } 0 \leq a < A^1(\kappa; s), && \text{binding} \end{align*}

Want to define $A^\prime \left(a, \kappa; s\right)$ and $K^\prime \left(a, \kappa; s\right)$

  • non-binding region $\to$ interpolate the values (we will come back to this point)
  • binding region $\to$ $a^\prime = 0$. What about $\kappa^\prime = $?

Step 2: Find $K^{\prime} \left(a, \kappa; s\right)$ when $a^\prime = 0$¶

We can't extrapolate $\bar{K} \left(a^\prime, \kappa^\prime; s\right)$ to the region with $a^\prime = 0$. It would mean $a^\prime < 0$.

But also notice $a^\prime = 0$ in the binding region.

That means, when the borrowing constraint is binding, the $\kappa$ and $\kappa^\prime$ satisfy $\kappa = \bar{K} \left(a^\prime = 0, \kappa^\prime; s\right)$ regardless of $a$

Those are on the threshold locus with $a = A^1 \left(\kappa; s\right)$

So, for any point $\left(a, \kappa \right)$ in the binding region,

$\kappa^\prime = \bar{K}^\prime \left(a = A^1 \left(\kappa; s\right), \kappa; s\right)$

How to actually calculate $\bar{K}^\prime \left(a = A^1 \left(\kappa; s\right), \kappa; s\right)$ for an arbitrary $\kappa$?

$A^1 \left(\kappa; s\right)$ is constructed by linearly interpolating the points with $a^\prime = 0$ (in $\kappa$)

Just interpolate $\{\left(\kappa, \kappa^\prime \right)\}$ in terms of $\kappa$

Want to evaluate the value of $\kappa^\prime$ at ★

Instead we can take a value at ⬛

As is the case in the 1D EGM, $\kappa^\prime \in \mathcal{K}$ but $\kappa$ is not on the grid

By linear interpolation we can evaluate $\kappa^\prime$ for all $\kappa \in \mathcal{K}$

Step 3: Construct $A^{\prime} \left(a, \kappa; s\right)$ and $K^{\prime} \left(a, \kappa; s\right)$ with interpolation¶

Ultimately we need to construct $A^{\prime} \left(a, \kappa; s\right)$ and $K^{\prime} \left(a, \kappa; s\right)$

  • to evaluate those on $\left(a, \kappa\right) \in \mathcal{A} \times \mathcal{K}$ to update $V$

  • to run a simulation / find a steady state

In the binding region, we already know how to calculate them

Now we need to find a way to calculate them in the non-binding region.

Interpolation on an irregular grid¶

we can approximate $A^{\prime} \left(a, \kappa; s\right)$ by interpolation

  • Conceptually, it is not difficult at all
  • Computationally, it is not as straightforward

Problem: the grid is irrelegular

  • Suppose the grid is regular in that all the $\left(a, \kappa \right)$ is on the mesh of $a$-grid and $\kappa$-grid
  • How to interpolate $A^{\prime} \left(a, \kappa; s\right)$ on the grid?
  • Let's see an example

2D linear interpolation with a mesh grid¶

We just need to find two adjacent grid points of $a$-grid, and another two adjacent grid points of $\kappa$-grid

No different than interpolating a function in 1D

Challenge with interpolation on an irregular grid¶

Now the grid is irregular. You now understand that it may not be always straightforward to interpolate a function

Suppose you want to use a quadrangular that encompasses a point to interpolate a function

  • First, split the space into quadrangulars (not necessarily unique)
  • Then you need to find the quadrangular
  • Both steps require algorithms

You may not want to use a quadrangular but a triangle.

  • Notice a point with $\kappa < 0$ is not feasible (but $a<0$ is feasible)
  • Sometimes there is no quadrangular a point resides in.

Popular approach in the literature¶

One algorithm to split the space into triangles is called Delaunay triangulation

  • there is a package to implement the triangulation
  • the package also searches for a triangle that contains a point
  • Ludwig and Schön (2014) implemented the algorithem to implement the 2D EGM

Other method papers on 2D-EGM

  • Hintermaier and Koeniger (2010): Hybrid method ($a$ and $\kappa^\prime$ are on the exo. grid. Solve for $a^\prime$)
  • Ludwig and Schön (2014): use Delaunay triangulation
  • White (2015): use the irregular quadrangulars with an efficient locating algorithm

Here I propose a different approach, using 1D linear interpolation

From Endogenous Grid to Exogenous Grid¶

From the Euler equation, we have obtained $\left(a, \kappa \right)$ for all $\left(a^\prime , \kappa^\prime \right) \in \mathcal{A} \times \mathcal{K}$

To update $V$, need to evaluate $A^{\prime} \left(a, \kappa; s\right)$ and $K^{\prime} \left(a, \kappa; s\right)$ on $\left(a, \kappa \right) \in \mathcal{A} \times \mathcal{K}$

By some 2D interpolation scheme, we could evaluate $\left(a^\prime , \kappa^\prime \right)$ for all $\left(a, \kappa \right) \in \mathcal{A} \times \mathcal{K}$ in one step

  • as we have seen, this is actually not that simple

My proposal: move variables onto the exo. grid one by one.

  • First, move $\kappa$ onto $\mathcal{K}$, keeping $a^\prime$ on $\mathcal{A}$

  • Next, move $a$ onto $\mathcal{A}$, keeping $\kappa$ on $\mathcal{K}$

  • Voilà!

Each locus represents all the points with $a^\prime = a_n \in \mathcal{A}$

Denote them $a = A^{n} \left(\kappa; s\right)$ for $n \in \left[1, \cdots, N_a\right]$

but each grid point on a locus is associated with $\kappa^\prime = \kappa_n \in \mathcal{K}$

we can interpolate the grid points to find the value of the function for $\kappa \in \mathcal{K}$

  • $A^{\prime} \left(A^{n} \left(\kappa; s\right), \kappa; s\right)$: it does not depend on $\kappa$
  • $K^{\prime} \left(A^{n} \left(\kappa; s\right), \kappa; s\right)$: interpolate the value in terms of $\kappa$

Now add points with $\left(a^\prime, \kappa\right) \in \mathcal{A}\times\mathcal{K}$

Now the grid points $\kappa \in \mathcal{K}$ and $a^\prime \in \mathcal{A}$

On the horizontal line $\kappa = \kappa_2$, we can interpolate the values in terms of $a$

Taking stock¶

Repeating the process for all $\kappa \in \mathcal{K}$, we obtain $A^{\prime} \left(a, \kappa; s\right)$ and $K^{\prime} \left(a, \kappa; s\right)$ on the gridpoint in the non-binding region.

Now we can evaluate the right hand side of

$$ V^{b}\left(a,\kappa; s\right)=\max_{c,e, a^{\prime}\geq0,\kappa^{\prime}}u(c)-g(e)+\beta\sum_{s}\Pi(s^{\prime}|s)V^{b}\left(a^{\prime},\kappa^{\prime};s^{\prime}\right) $$ on the grid, and update $V^b$.

We can also evaluate $A^{\prime} \left(a, \kappa; s\right)$ and $K^{\prime} \left(a, \kappa; s\right)$ using the points on the mesh grid.

  • used to run a simulation/construct a transition matrix
  • interpolation on a mesh 2D grid is straightforward

Now let's see the code in action.