{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Deterministic Optimal Economic Growth Model\n",
"\n",
"**Randall Romero Aguilar, PhD**\n",
"\n",
"This demo is based on the original Matlab demo accompanying the Computational Economics and Finance 2001 textbook by Mario Miranda and Paul Fackler.\n",
"\n",
"Original (Matlab) CompEcon file: **demdp06.m**\n",
"\n",
"Running this file requires the Python version of CompEcon. This can be installed with pip by running\n",
"\n",
" !pip install compecon --upgrade\n",
"\n",
"Last updated: 2022-Oct-23\n",
"
"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": false
},
"source": [
"## About\n",
"\n",
"Welfare maximizing social planner must decide how much society should consume and invest. Model is of special interest because it has a known closed-form solution.\n",
"\n",
"- States\n",
" - s stock of wealth\n",
"- Actions\n",
" - k capital investment\n",
"- Parameters\n",
" - β capital production elasticity\n",
" - δ discount factor"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import pandas as pd\n",
"import matplotlib.pyplot as plt\n",
"from compecon import BasisChebyshev, DPmodel, DPoptions, qnwnorm"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Model parameters\n",
"\n",
"Assuming that the marginal productivity of capital is $\\beta=0.5$ and the discount factor is $\\delta=0.9$"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"β, δ = 0.5, 0.9"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Analytic results\n",
"\n",
"The steady-state values for this model are"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"sstar = (β * δ) ** (β / (1 - β)) # steady-state wealth\n",
"kstar = β * δ * sstar # steady-state capital investment\n",
"vstar = np.log(sstar - kstar) / (1 - δ) # steady-state value\n",
"pstar = 1 / (sstar * (1 - β * δ)) # steady-state shadow price\n",
"b = 1 / (1 - δ * β)\n",
"\n",
"print('\\n\\nSteady-State')\n",
"for var, value in zip(['Wealth','Investment','Value','Shadow price'], [sstar,kstar,vstar,pstar]):\n",
" print(f'\\t{var:12s} = {value:8.4f}')"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Wealth 0.450000\n",
"Investment 0.202500\n",
"Value -13.963447\n",
"Shadow Price 4.040404\n",
"dtype: float64"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"sstar = (β * δ) ** (β / (1 - β)) # steady-state wealth\n",
"kstar = β * δ * sstar # steady-state capital investment\n",
"vstar = np.log(sstar - kstar) / (1 - δ) # steady-state value\n",
"pstar = 1 / (sstar * (1 - β * δ)) # steady-state shadow price\n",
"b = 1 / (1 - δ * β)\n",
"\n",
"ss = pd.Series([sstar,kstar,vstar,pstar], index=['Wealth','Investment','Value','Shadow Price'])\n",
"ss"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The true value function is\n",
"\\begin{equation}\n",
"V(s) = v^* + \\frac{1}{1-\\delta\\beta}\\left(\\log(s) -\\log(s^*)\\right)\n",
"\\end{equation}"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"def vtrue(wealth): # analytic value function\n",
" return vstar + b * (np.log(wealth) - np.log(sstar))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The true policy function is\n",
"\\begin{equation}\n",
"k(s) = \\delta\\beta s\n",
"\\end{equation}"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"def ktrue(wealth): #analytic policy function\n",
" return δ*β*wealth"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Numeric results"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### State space\n",
"The state variable is s=\"Wealth\", which we restrict to $0\\in[0.2, 1.0]$. \n",
"\n",
"Here, we represent it with a Chebyshev basis, with $n=15$ nodes."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"n, smin, smax = 15, 0.2, 1.0\n",
"basis = BasisChebyshev(n, smin, smax, labels=['Wealth'])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Action space\n",
"The choice variable k=\"Investment\" must be nonnegative."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"def bounds(s, i=None, j=None):\n",
" return np.zeros_like(s), s[:]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Reward function\n",
"The reward function is the utility of consumption=$s-k$."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"def reward(s, k, i=None, j=None):\n",
" sk = s - k\n",
" u = np.log(sk)\n",
" ux= - sk ** -1\n",
" uxx = - sk ** -2\n",
" return u, ux, uxx"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### State transition function\n",
"Next period, wealth will be equal to production from available initial capital $k$, that is $s' = k^\\beta$"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"def transition(s, k, i=None, j=None, in_=None, e=None):\n",
" g = k ** β\n",
" gx = β * k **(β - 1)\n",
" gxx = (β - 1) * β * k ** (β - 2)\n",
" return g, gx, gxx"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Model structure\n",
"The value of wealth $s$ satisfies the Bellman equation\n",
"\\begin{equation*}\n",
"V(s) = \\max_k\\left\\{\\log(s-k) + \\delta V(k^\\beta) \\right\\}\n",
"\\end{equation*}\n",
"\n",
"To solve and simulate this model,use the CompEcon class `DPmodel`"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [],
"source": [
"growth_model = DPmodel(basis, reward, transition, bounds,\n",
" x=['Investment'],\n",
" discount=δ)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Solving the model\n",
"\n",
"Solving the growth model by collocation, using *Newton* algorithm and a maximum of 20 iterations"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Solving infinite-horizon model collocation equation by Newton's method\n",
"iter change time \n",
"------------------------------\n",
" 0 4.7e-07 0.0040\n",
" 1 1.6e-12 0.0080\n",
"Elapsed Time = 0.01 Seconds\n"
]
}
],
"source": [
"options = dict(show=True,\n",
" algorithm='newton',\n",
" maxit=20)\n",
"\n",
"snodes = growth_model.Value.nodes\n",
"S = growth_model.solve(vtrue(snodes), ktrue(snodes), **options)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"`DPmodel.solve` returns a pandas `DataFrame` with the following data:"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" Wealth | \n",
" value | \n",
" resid | \n",
" Investment | \n",
"
\n",
" \n",
" | Wealth | \n",
" | \n",
" | \n",
" | \n",
" | \n",
"
\n",
" \n",
" \n",
" \n",
" | 0.200000 | \n",
" 0.200000 | \n",
" -15.437865 | \n",
" 2.766030e-07 | \n",
" 0.090000 | \n",
"
\n",
" \n",
" | 0.205369 | \n",
" 0.205369 | \n",
" -15.389699 | \n",
" -2.099164e-07 | \n",
" 0.092416 | \n",
"
\n",
" \n",
" | 0.210738 | \n",
" 0.210738 | \n",
" -15.342776 | \n",
" -2.488242e-07 | \n",
" 0.094832 | \n",
"
\n",
" \n",
" | 0.216107 | \n",
" 0.216107 | \n",
" -15.297033 | \n",
" -1.102768e-07 | \n",
" 0.097248 | \n",
"
\n",
" \n",
" | 0.221477 | \n",
" 0.221477 | \n",
" -15.252412 | \n",
" 5.644805e-08 | \n",
" 0.099664 | \n",
"
\n",
" \n",
"
\n",
"
\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" Wealth | \n",
" value | \n",
" resid | \n",
" Investment | \n",
" shadow price | \n",
" error | \n",
" D.Wealth | \n",
"
\n",
" \n",
" | Wealth | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
"
\n",
" \n",
" \n",
" \n",
" | 0.200000 | \n",
" 0.200000 | \n",
" -15.437865 | \n",
" 2.766030e-07 | \n",
" 0.090000 | \n",
" 9.090758 | \n",
" 4.241867e-07 | \n",
" 0.100000 | \n",
"
\n",
" \n",
" | 0.205369 | \n",
" 0.205369 | \n",
" -15.389699 | \n",
" -2.099164e-07 | \n",
" 0.092416 | \n",
" 8.853206 | \n",
" -3.099038e-08 | \n",
" 0.098631 | \n",
"
\n",
" \n",
" | 0.210738 | \n",
" 0.210738 | \n",
" -15.342776 | \n",
" -2.488242e-07 | \n",
" 0.094832 | \n",
" 8.627699 | \n",
" -4.010325e-08 | \n",
" 0.097210 | \n",
"
\n",
" \n",
" | 0.216107 | \n",
" 0.216107 | \n",
" -15.297033 | \n",
" -1.102768e-07 | \n",
" 0.097248 | \n",
" 8.413362 | \n",
" 1.253504e-07 | \n",
" 0.095739 | \n",
"
\n",
" \n",
" | 0.221477 | \n",
" 0.221477 | \n",
" -15.252412 | \n",
" 5.644805e-08 | \n",
" 0.099664 | \n",
" 8.209399 | \n",
" 3.151235e-07 | \n",
" 0.094220 | \n",
"
\n",
" \n",
"
\n",
"
\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" Wealth | \n",
" Investment | \n",
" value | \n",
" value_Wealth | \n",
" Wealth_next | \n",
" shadow price | \n",
" error | \n",
" D.Wealth | \n",
"
\n",
" \n",
" | Wealth | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
"
\n",
" \n",
" \n",
" \n",
" | 0.202191 | \n",
" 0.202191 | \n",
" -0.020528 | \n",
" -15.014820 | \n",
" 4.444953 | \n",
" 0.202192 | \n",
" 4.444953 | \n",
" 0.403234 | \n",
" 2.977733e-07 | \n",
"
\n",
" \n",
" | 0.219577 | \n",
" 0.219577 | \n",
" -0.004880 | \n",
" -14.937786 | \n",
" 4.416570 | \n",
" 0.219578 | \n",
" 4.416570 | \n",
" 0.330284 | \n",
" 2.769166e-07 | \n",
"
\n",
" \n",
" | 0.253590 | \n",
" 0.253590 | \n",
" 0.025731 | \n",
" -14.788512 | \n",
" 4.361045 | \n",
" 0.253590 | \n",
" 4.361045 | \n",
" 0.217716 | \n",
" 2.361148e-07 | \n",
"
\n",
" \n",
" | 0.302742 | \n",
" 0.302742 | \n",
" 0.069968 | \n",
" -14.576129 | \n",
" 4.280804 | \n",
" 0.302742 | \n",
" 4.280804 | \n",
" 0.107984 | \n",
" 1.771510e-07 | \n",
"
\n",
" \n",
" | 0.364886 | \n",
" 0.364886 | \n",
" 0.125897 | \n",
" -14.313256 | \n",
" 4.179353 | \n",
" 0.364886 | \n",
" 4.179353 | \n",
" 0.031397 | \n",
" 1.026024e-07 | \n",
"
\n",
" \n",
"
\n",
"