@@ -45,7 +45,7 @@ Let's put these algorithm and code optimizations to one side for now.
4545
4646We will use the following imports:
4747
48- ``` {code-cell} ipython
48+ ``` {code-cell} python3
4949import matplotlib.pyplot as plt
5050import numpy as np
5151from scipy.optimize import minimize_scalar, bisect
@@ -142,14 +142,16 @@ The process looks like this:
1421421 . Begin with an array of values $\{ v_0, \ldots, v_I \} $ representing
143143 the values of some initial function $v$ on the grid points $\{ x_0, \ldots, x_I \} $.
1441441 . Build a function $\hat v$ on the state space $\mathbb R_ +$ by
145- linear interpolation, based on these data points.
146- 1 . Obtain and record the value $T \hat v(x_i)$ on each grid point
147- $x_i$ by repeatedly solving the maximization problem in the Bellman
148- equation .
145+ interpolation, based on the interpolation points $ \{ (x_i, v_i) \} $ .
146+ 1 . Insert $\hat v$ into the right hand side of the Bellman equation and
147+ obtain and record the value $T \hat v(x_i)$ on each grid point
148+ $x_i$ .
1491491 . Unless some stopping condition is satisfied, set
150150 $\{ v_0, \ldots, v_I \} = \{ T \hat v(x_0), \ldots, T \hat v(x_I) \} $ and go to step 2.
151151
152- In step 2 we'll use continuous piecewise linear interpolation.
152+ In step 2 we'll use piecewise linear interpolation.
153+
154+
153155
154156### Implementation
155157
@@ -177,26 +179,31 @@ def maximize(g, a, b, args):
177179We'll store the parameters $\beta$ and $\gamma$ and the grid in a
178180` NamedTuple ` called ` Model ` .
179181
182+ We'll also create a helper function called ` create_cake_eating_model ` to store
183+ default parameters and build an instance of ` Model ` .
184+
180185``` {code-cell} python3
181- # Create model data structure
182186class Model(NamedTuple):
183187 β: float
184188 γ: float
185189 x_grid: np.ndarray
186190
187- def create_cake_eating_model(β=0.96, # discount factor
188- γ=1.5, # degree of relative risk aversion
189- x_grid_min=1e-3, # exclude zero for numerical stability
190- x_grid_max=2.5, # size of cake
191- x_grid_size=120):
191+ def create_cake_eating_model(
192+ β: float = 0.96, # discount factor
193+ γ: float = 1.5, # degree of relative risk aversion
194+ x_grid_min: float = 1e-3, # exclude zero for numerical stability
195+ x_grid_max: float = 2.5, # size of cake
196+ x_grid_size: int = 120
197+ ):
192198 """
193199 Creates an instance of the cake eating model.
200+
194201 """
195202 x_grid = np.linspace(x_grid_min, x_grid_max, x_grid_size)
196- return Model(β=β , γ=γ, x_grid= x_grid)
203+ return Model(β, γ, x_grid)
197204```
198205
199- Now we define utility functions that operate on the model:
206+ Here's the CRRA utility function.
200207
201208``` {code-cell} python3
202209def u(c, γ):
@@ -207,33 +214,42 @@ def u(c, γ):
207214 return np.log(c)
208215 else:
209216 return (c ** (1 - γ)) / (1 - γ)
217+ ```
210218
211- def u_prime(c, γ):
212- """
213- First derivative of utility function.
214- """
215- return c ** (-γ)
219+ The next function is the unmaximized right hand side of the Bellman equation.
220+
221+ The array ` v ` is the current guess of $v$, stored as an array on the grid
222+ points.
216223
217- def state_action_value(c, x, v_array, model):
224+ ``` {code-cell} python3
225+ def state_action_value(
226+ c: float, # current consumption
227+ x: float, # the current state (remaining cake)
228+ v: np.ndarray, # current guess of the value function
229+ model: Model # instance of cake eating model
230+ ):
218231 """
219232 Right hand side of the Bellman equation given x and c.
233+
220234 """
235+ # Unpack
221236 β, γ, x_grid = model.β, model.γ, model.x_grid
222- v = lambda x: np.interp(x, x_grid, v_array)
223-
224- return u(c, γ) + β * v(x - c)
237+ # Convert array into function
238+ vf = lambda x: np.interp(x, x_grid, v)
239+ # Return unmaximmized RHS of Bellman equation
240+ return u(c, γ) + β * vf(x - c)
225241```
226242
227243We now define the Bellman operation:
228244
229245``` {code-cell} python3
230- def T(v, model):
246+ def T(
247+ v: np.ndarray, # current guess of the value function
248+ model: Model # instance of cake eating model
249+ ):
231250 """
232251 The Bellman operator. Updates the guess of the value function.
233252
234- * model is an instance of Model
235- * v is an array representing a guess of the value function
236-
237253 """
238254 v_new = np.empty_like(v)
239255
@@ -261,34 +277,42 @@ $x$ grid point.
261277x_grid = model.x_grid
262278v = u(x_grid, model.γ) # Initial guess
263279n = 12 # Number of iterations
264-
265280fig, ax = plt.subplots()
266281
282+ # Initial plot
267283ax.plot(x_grid, v, color=plt.cm.jet(0),
268284 lw=2, alpha=0.6, label='Initial guess')
269285
286+ # Iterate
270287for i in range(n):
271288 v = T(v, model) # Apply the Bellman operator
272289 ax.plot(x_grid, v, color=plt.cm.jet(i / n), lw=2, alpha=0.6)
273290
291+ # One last update and plot
292+ v = T(v, model)
293+ ax.plot(x_grid, v, color=plt.cm.jet(0),
294+ lw=2, alpha=0.6, label='Final guess')
295+
274296ax.legend()
275297ax.set_ylabel('value', fontsize=12)
276298ax.set_xlabel('cake size $x$', fontsize=12)
277299ax.set_title('Value function iterations')
278-
279300plt.show()
280301```
281302
282- To do this more systematically, we introduce a wrapper function called
283- ` compute_value_function ` that iterates until some convergence conditions are
284- satisfied.
303+ To iterate more systematically, we introduce a wrapper function called
304+ ` compute_value_function ` .
305+
306+ It's task is to iterate using $T$ until some convergence conditions are satisfied.
285307
286308``` {code-cell} python3
287- def compute_value_function(model,
288- tol=1e-4,
289- max_iter=1000,
290- verbose=True,
291- print_skip=25):
309+ def compute_value_function(
310+ model: Model,
311+ tol: float = 1e-4,
312+ max_iter: int = 1_000,
313+ verbose: bool = True,
314+ print_skip: int = 25
315+ ):
292316
293317 # Set up loop
294318 v = np.zeros(len(model.x_grid)) # Initial guess
@@ -367,6 +391,9 @@ working with policy functions.
367391We will see that value function iteration can be avoided by iterating on a guess
368392of the policy function instead.
369393
394+ The policy function has less curvature and hence is easier to interpolate than
395+ the value function.
396+
370397These ideas will be explored over the next few lectures.
371398```
372399
@@ -398,14 +425,14 @@ above.
398425Here's the function:
399426
400427``` {code-cell} python3
401- def σ(model, v):
428+ def σ(
429+ v: np.ndarray, # current guess of the value function
430+ model: Model # instance of cake eating model
431+ ):
402432 """
403433 The optimal policy function. Given the value function,
404434 it finds optimal consumption in each state.
405435
406- * model is an instance of Model
407- * v is a value function array
408-
409436 """
410437 c = np.empty_like(v)
411438
@@ -420,7 +447,7 @@ def σ(model, v):
420447Now let's pass the approximate value function and compute optimal consumption:
421448
422449``` {code-cell} python3
423- c = σ(model, v )
450+ c = σ(v, model )
424451```
425452
426453(pol_an)=
0 commit comments