cons_smooth.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "2a9bacb8",
   "metadata": {},
   "source": [
    "# Consumption Smoothing"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "51a257fe",
   "metadata": {},
   "source": [
    "## Overview\n",
    "\n",
    "In this lecture, we’ll study a famous model of the “consumption function” that Milton Friedman [[Friedman, 1956](https://intro.quantecon.org/zreferences.html#id183)] and Robert Hall [[Hall, 1978](https://intro.quantecon.org/zreferences.html#id184)])  proposed to fit some empirical data patterns that the original  Keynesian consumption function  described in this QuantEcon lecture [geometric series](https://intro.quantecon.org/geom_series.html)  missed.\n",
    "\n",
    "We’ll study what is often  called the “consumption-smoothing model.”\n",
    "\n",
    "We’ll use  matrix multiplication and matrix inversion, the same tools that we used in this QuantEcon lecture [present values](https://intro.quantecon.org/pv.html).\n",
    "\n",
    "Formulas presented in  [present value formulas](https://intro.quantecon.org/pv.html) are at the core of the consumption-smoothing model because we shall use them to define a consumer’s “human wealth”.\n",
    "\n",
    "The  key idea that inspired Milton Friedman was that a person’s non-financial income, i.e., his or\n",
    "her wages from working, can be viewed as a dividend stream from ‘‘human capital’’\n",
    "and that standard asset-pricing formulas can be applied to compute\n",
    "‘‘non-financial wealth’’ that capitalizes that  earnings stream.\n",
    "\n",
    ">**Note**\n",
    ">\n",
    ">As we’ll see in this QuantEcon lecture  [equalizing difference model](https://intro.quantecon.org/equalizing_difference.html),\n",
    "Milton Friedman had used this idea  in his PhD thesis at Columbia University,\n",
    "eventually published as [[Kuznets and Friedman, 1939](https://intro.quantecon.org/zreferences.html#id9)] and [[Friedman and Kuznets, 1945](https://intro.quantecon.org/zreferences.html#id10)].\n",
    "\n",
    "It will take a while for a “present value” or asset price explicitly to appear in this lecture, but when it does it will be a key actor."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a96b3476",
   "metadata": {},
   "source": [
    "## Analysis\n",
    "\n",
    "As usual, we’ll start by importing some Python modules."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "efdb64a8",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from collections import namedtuple"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ebedece5",
   "metadata": {},
   "source": [
    "The model describes  a consumer who lives from time $ t=0, 1, \\ldots, T $, receives a stream $ \\{y_t\\}_{t=0}^T $ of non-financial income and chooses a consumption stream $ \\{c_t\\}_{t=0}^T $.\n",
    "\n",
    "We usually think of the non-financial income stream as coming from the person’s earnings from supplying labor.\n",
    "\n",
    "The model  takes a non-financial income stream as an input, regarding it as “exogenous” in the sense that it is  determined outside the model.\n",
    "\n",
    "The consumer faces a gross interest rate of $ R >1 $ that is constant over time, at which she is free to borrow or lend, up to  limits  that we’ll describe below.\n",
    "\n",
    "Let\n",
    "\n",
    "- $ T \\geq 2 $  be a positive integer that constitutes a time-horizon.  \n",
    "- $ y = \\{y_t\\}_{t=0}^T $ be an exogenous  sequence of non-negative non-financial incomes $ y_t $.  \n",
    "- $ a = \\{a_t\\}_{t=0}^{T+1} $ be a sequence of financial wealth.  \n",
    "- $ c = \\{c_t\\}_{t=0}^T $ be a sequence of non-negative consumption rates.  \n",
    "- $ R \\geq 1 $ be a fixed gross one period rate of return on financial assets.  \n",
    "- $ \\beta \\in (0,1) $ be a fixed discount factor.  \n",
    "- $ a_0 $ be a given initial level of financial assets  \n",
    "- $ a_{T+1} \\geq 0 $  be a terminal condition on final assets.  \n",
    "\n",
    "\n",
    "The sequence of financial wealth $ a $ is to be determined by the model.\n",
    "\n",
    "We require it to satisfy  two  **boundary conditions**:\n",
    "\n",
    "- it must  equal an exogenous value  $ a_0 $ at time $ 0 $  \n",
    "- it must equal or exceed an exogenous value  $ a_{T+1} $ at time $ T+1 $.  \n",
    "\n",
    "\n",
    "The **terminal condition** $ a_{T+1} \\geq 0 $ requires that the consumer not leave the model in debt.\n",
    "\n",
    "(We’ll soon see that a utility maximizing consumer won’t want to die leaving positive assets, so she’ll arrange her affairs to make\n",
    "$ a_{T+1} = 0 $.)\n",
    "\n",
    "The consumer faces a sequence of budget constraints that  constrains   sequences $ (y, c, a) $\n",
    "\n",
    "\n",
    "<a id='equation-eq-a-t'></a>\n",
    "$$\n",
    "a_{t+1} = R (a_t+ y_t - c_t), \\quad t =0, 1, \\ldots T \\tag{12.1}\n",
    "$$\n",
    "\n",
    "Equations [(12.1)](#equation-eq-a-t) constitute  $ T+1 $ such budget constraints, one for each $ t=0, 1, \\ldots, T $.\n",
    "\n",
    "Given a sequence $ y $ of non-financial incomes, a large  set of pairs $ (a, c) $ of (financial wealth, consumption) sequences  satisfy the sequence of budget constraints [(12.1)](#equation-eq-a-t).\n",
    "\n",
    "Our model has the following logical flow.\n",
    "\n",
    "- start with an exogenous non-financial income sequence $ y $, an initial financial wealth $ a_0 $, and\n",
    "  a candidate consumption path $ c $.  \n",
    "- use the system of equations [(12.1)](#equation-eq-a-t) for $ t=0, \\ldots, T $ to compute a path $ a $ of financial wealth  \n",
    "- verify that $ a_{T+1} $ satisfies the terminal wealth constraint $ a_{T+1} \\geq 0 $.  \n",
    "  - If it does, declare that the candidate path is **budget feasible**.  \n",
    "  - if the candidate consumption path is not budget feasible, propose a less greedy consumption  path and start over  \n",
    "\n",
    "\n",
    "Below, we’ll describe how to execute these steps using linear algebra – matrix inversion and multiplication.\n",
    "\n",
    "The above procedure seems like a sensible way to find “budget-feasible” consumption paths $ c $, i.e., paths that are consistent\n",
    "with the exogenous non-financial income stream $ y $, the initial financial  asset level $ a_0 $, and the terminal asset level $ a_{T+1} $.\n",
    "\n",
    "In general, there are **many** budget feasible consumption paths $ c $.\n",
    "\n",
    "Among all budget-feasible consumption paths, which one should a consumer want?\n",
    "\n",
    "To answer this question, we shall eventually evaluate alternative budget feasible consumption paths $ c $ using the following utility functional or **welfare criterion**:\n",
    "\n",
    "\n",
    "<a id='equation-welfare'></a>\n",
    "$$\n",
    "W = \\sum_{t=0}^T \\beta^t (g_1 c_t - \\frac{g_2}{2} c_t^2 ) \\tag{12.2}\n",
    "$$\n",
    "\n",
    "where $ g_1 > 0, g_2 > 0 $.\n",
    "\n",
    "When $ \\beta R \\approx 1 $, the fact that the utility function $ g_1 c_t - \\frac{g_2}{2} c_t^2 $ has diminishing marginal utility imparts a preference for consumption that is very smooth.\n",
    "\n",
    "Indeed, we shall see that when $ \\beta R = 1 $ (a condition assumed by Milton Friedman [[Friedman, 1956](https://intro.quantecon.org/zreferences.html#id183)] and Robert Hall [[Hall, 1978](https://intro.quantecon.org/zreferences.html#id184)]),  criterion [(12.2)](#equation-welfare) assigns higher welfare to smoother consumption paths.\n",
    "\n",
    "By **smoother** we mean as close as possible to being constant over time.\n",
    "\n",
    "The preference for smooth consumption paths that is built into the model gives it the  name “consumption-smoothing model”.\n",
    "\n",
    "We’ll postpone verifying our claim that a constant consumption path is optimal when $ \\beta R=1 $\n",
    "by comparing welfare levels that comes from a constant path with ones that involve non-constant paths.\n",
    "\n",
    "Before doing that, let’s dive in and do some calculations that will help us understand how the model works in practice when we provide the consumer with some different streams on non-financial income.\n",
    "\n",
    "Here we use default parameters $ R = 1.05 $, $ g_1 = 1 $, $ g_2 = 1/2 $, and $ T = 65 $.\n",
    "\n",
    "We create a Python **namedtuple** to store these parameters with default values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2855ed1b",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "ConsumptionSmoothing = namedtuple(\"ConsumptionSmoothing\", \n",
    "                        [\"R\", \"g1\", \"g2\", \"β_seq\", \"T\"])\n",
    "\n",
    "def create_consumption_smoothing_model(R=1.05, g1=1, g2=1/2, T=65):\n",
    "    β = 1/R\n",
    "    β_seq = np.array([β**i for i in range(T+1)])\n",
    "    return ConsumptionSmoothing(R, g1, g2, \n",
    "                                β_seq, T)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f411c1b2",
   "metadata": {},
   "source": [
    "## Friedman-Hall consumption-smoothing model\n",
    "\n",
    "A key object is what Milton Friedman called “human” or “non-financial” wealth at time $ 0 $:\n",
    "\n",
    "$$\n",
    "h_0 \\equiv \\sum_{t=0}^T R^{-t} y_t = \\begin{bmatrix} 1 & R^{-1} & \\cdots & R^{-T} \\end{bmatrix}\n",
    "\\begin{bmatrix} y_0 \\cr y_1  \\cr \\vdots \\cr y_T \\end{bmatrix}\n",
    "$$\n",
    "\n",
    "Human or non-financial wealth  at time $ 0 $ is evidently just the present value of the consumer’s non-financial income stream $ y $.\n",
    "\n",
    "Formally it very much resembles the asset price that we computed in this QuantEcon lecture [present values](https://intro.quantecon.org/pv.html).\n",
    "\n",
    "Indeed, this is why Milton Friedman called it “human capital”.\n",
    "\n",
    "By iterating on equation [(12.1)](#equation-eq-a-t) and imposing the terminal condition\n",
    "\n",
    "$$\n",
    "a_{T+1} = 0,\n",
    "$$\n",
    "\n",
    "it is possible to convert a sequence of budget constraints [(12.1)](#equation-eq-a-t) into a single intertemporal constraint\n",
    "\n",
    "\n",
    "<a id='equation-eq-budget-intertemp'></a>\n",
    "$$\n",
    "\\sum_{t=0}^T R^{-t} c_t = a_0 + h_0. \\tag{12.3}\n",
    "$$\n",
    "\n",
    "Equation [(12.3)](#equation-eq-budget-intertemp)  says that the present value of the consumption stream equals the sum of financial and non-financial (or human) wealth.\n",
    "\n",
    "Robert Hall [[Hall, 1978](https://intro.quantecon.org/zreferences.html#id184)] showed that when $ \\beta R = 1 $, a condition Milton Friedman had also  assumed, it is “optimal” for a consumer to smooth consumption by setting\n",
    "\n",
    "$$\n",
    "c_t = c_0 \\quad t =0, 1, \\ldots, T\n",
    "$$\n",
    "\n",
    "(Later we’ll present a “variational argument” that shows that this constant path maximizes\n",
    "criterion [(12.2)](#equation-welfare) when $ \\beta R =1 $.)\n",
    "\n",
    "In this case, we can use the intertemporal budget constraint to write\n",
    "\n",
    "\n",
    "<a id='equation-eq-conssmoothing'></a>\n",
    "$$\n",
    "c_t = c_0  = \\left(\\sum_{t=0}^T R^{-t}\\right)^{-1} (a_0 + h_0), \\quad t= 0, 1, \\ldots, T. \\tag{12.4}\n",
    "$$\n",
    "\n",
    "Equation [(12.4)](#equation-eq-conssmoothing) is the consumption-smoothing model in a nutshell."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "83ff462d",
   "metadata": {},
   "source": [
    "## Mechanics of consumption-smoothing model\n",
    "\n",
    "As promised, we’ll provide step-by-step instructions on how to use linear algebra, readily implemented in Python, to compute all  objects in play in  the consumption-smoothing model.\n",
    "\n",
    "In the calculations below,  we’ll  set default values of  $ R > 1 $, e.g., $ R = 1.05 $, and $ \\beta = R^{-1} $."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3019521f",
   "metadata": {},
   "source": [
    "### Step 1\n",
    "\n",
    "For a $ (T+1) \\times 1 $  vector $ y $, use matrix algebra to compute $ h_0 $\n",
    "\n",
    "$$\n",
    "h_0 = \\sum_{t=0}^T R^{-t} y_t = \\begin{bmatrix} 1 & R^{-1} & \\cdots & R^{-T} \\end{bmatrix}\n",
    "\\begin{bmatrix} y_0 \\cr y_1  \\cr \\vdots \\cr y_T \\end{bmatrix}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bc97a448",
   "metadata": {},
   "source": [
    "### Step 2\n",
    "\n",
    "Compute an  time $ 0 $   consumption $ c_0 $ :\n",
    "\n",
    "$$\n",
    "c_t = c_0 = \\left( \\frac{1 - R^{-1}}{1 - R^{-(T+1)}} \\right) (a_0 + \\sum_{t=0}^T R^{-t} y_t ) , \\quad t = 0, 1, \\ldots, T\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b06aba59",
   "metadata": {},
   "source": [
    "### Step 3\n",
    "\n",
    "Use  the system of equations [(12.1)](#equation-eq-a-t) for $ t=0, \\ldots, T $ to compute a path $ a $ of financial wealth.\n",
    "\n",
    "To do this, we translate that system of difference equations into a single matrix equation as follows:\n",
    "\n",
    "$$\n",
    "\\begin{bmatrix} \n",
    "1 & 0 & 0 & \\cdots & 0 & 0 & 0 \\cr\n",
    "-R & 1 & 0 & \\cdots & 0 & 0 & 0 \\cr\n",
    "0 & -R & 1 & \\cdots & 0 & 0 & 0 \\cr\n",
    "\\vdots  &\\vdots & \\vdots & \\cdots & \\vdots & \\vdots & \\vdots \\cr\n",
    "0 & 0 & 0 & \\cdots & -R & 1 & 0 \\cr\n",
    "0 & 0 & 0 & \\cdots & 0 & -R & 1\n",
    "\\end{bmatrix} \n",
    "\\begin{bmatrix} a_1 \\cr a_2 \\cr a_3 \\cr \\vdots \\cr a_T \\cr a_{T+1} \n",
    "\\end{bmatrix}\n",
    "= R \n",
    "\\begin{bmatrix} y_0 + a_0 - c_0 \\cr y_1 - c_0 \\cr y_2 - c_0 \\cr \\vdots\\cr y_{T-1} - c_0 \\cr y_T - c_0\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "Multiply both sides by the inverse of the matrix on the left side to compute\n",
    "\n",
    "$$\n",
    "\\begin{bmatrix} a_1 \\cr a_2 \\cr a_3 \\cr \\vdots \\cr a_T \\cr a_{T+1} \\end{bmatrix}\n",
    "$$\n",
    "\n",
    "Because we have built into  our calculations that the consumer leaves the model  with exactly zero assets, just barely satisfying the\n",
    "terminal condition that $ a_{T+1} \\geq 0 $, it should turn out   that\n",
    "\n",
    "$$\n",
    "a_{T+1} = 0.\n",
    "$$\n",
    "\n",
    "Let’s verify this with  Python code.\n",
    "\n",
    "First we implement the model with `compute_optimal`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "469085f9",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def compute_optimal(model, a0, y_seq):\n",
    "    R, T = model.R, model.T\n",
    "\n",
    "    # non-financial wealth\n",
    "    h0 = model.β_seq @ y_seq     # since β = 1/R\n",
    "\n",
    "    # c0\n",
    "    c0 = (1 - 1/R) / (1 - (1/R)**(T+1)) * (a0 + h0)\n",
    "    c_seq = c0*np.ones(T+1)\n",
    "\n",
    "    # verify\n",
    "    A = np.diag(-R*np.ones(T), k=-1) + np.eye(T+1)\n",
    "    b = y_seq - c_seq\n",
    "    b[0] = b[0] + a0\n",
    "\n",
    "    a_seq = np.linalg.inv(A) @ b\n",
    "    a_seq = np.concatenate([[a0], a_seq])\n",
    "\n",
    "    return c_seq, a_seq, h0"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e05fa58d",
   "metadata": {},
   "source": [
    "We use an example where the consumer inherits $ a_0<0 $.\n",
    "\n",
    "This  can be interpreted as  student debt with which the consumer begins his or her working life.\n",
    "\n",
    "The non-financial process $ \\{y_t\\}_{t=0}^{T} $ is constant and positive up to $ t=45 $ and then becomes zero afterward.\n",
    "\n",
    "The drop in non-financial income late in life reflects retirement from work."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "31483bfe",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Financial wealth\n",
    "a0 = -2     # such as \"student debt\"\n",
    "\n",
    "# non-financial Income process\n",
    "y_seq = np.concatenate([np.ones(46), np.zeros(20)])\n",
    "\n",
    "cs_model = create_consumption_smoothing_model()\n",
    "c_seq, a_seq, h0 = compute_optimal(cs_model, a0, y_seq)\n",
    "\n",
    "print('check a_T+1=0:', \n",
    "      np.abs(a_seq[-1] - 0) <= 1e-8)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ddcc6ad6",
   "metadata": {},
   "source": [
    "The graphs below  show  paths of non-financial income, consumption, and financial assets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aad3ce3d",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Sequence length\n",
    "T = cs_model.T\n",
    "\n",
    "fig, axes = plt.subplots(1, 2, figsize=(12,5))\n",
    "\n",
    "axes[0].plot(range(T+1), y_seq, label='non-financial income', lw=2)\n",
    "axes[0].plot(range(T+1), c_seq, label='consumption', lw=2)\n",
    "axes[1].plot(range(T+2), a_seq, label='financial wealth', color='green', lw=2)\n",
    "axes[0].set_ylabel(r'$c_t,y_t$')\n",
    "axes[1].set_ylabel(r'$a_t$')\n",
    "\n",
    "for ax in axes:\n",
    "    ax.plot(range(T+2), np.zeros(T+2), '--', lw=1, color='black')\n",
    "    ax.legend()\n",
    "    ax.set_xlabel(r'$t$')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a56af317",
   "metadata": {},
   "source": [
    "Note that $ a_{T+1} = 0 $, as anticipated.\n",
    "\n",
    "We can  evaluate  welfare criterion [(12.2)](#equation-welfare)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9d89a868",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def welfare(model, c_seq):\n",
    "    β_seq, g1, g2 = model.β_seq, model.g1, model.g2\n",
    "\n",
    "    u_seq = g1 * c_seq - g2/2 * c_seq**2\n",
    "    return β_seq @ u_seq\n",
    "\n",
    "print('Welfare:', welfare(cs_model, c_seq))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f800f96d",
   "metadata": {},
   "source": [
    "### Experiments\n",
    "\n",
    "In this section we describe  how a  consumption sequence would optimally respond to different  sequences sequences of non-financial income.\n",
    "\n",
    "First we create  a function `plot_cs` that generates graphs for different instances of the  consumption-smoothing model `cs_model`.\n",
    "\n",
    "This will  help us avoid rewriting code to plot outcomes for different non-financial income sequences."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "287e2edd",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def plot_cs(model,    # consumption-smoothing model      \n",
    "            a0,       # initial financial wealth\n",
    "            y_seq     # non-financial income process\n",
    "           ):\n",
    "    \n",
    "    # Compute optimal consumption\n",
    "    c_seq, a_seq, h0 = compute_optimal(model, a0, y_seq)\n",
    "    \n",
    "    # Sequence length\n",
    "    T = cs_model.T\n",
    "    \n",
    "    fig, axes = plt.subplots(1, 2, figsize=(12,5))\n",
    "    \n",
    "    axes[0].plot(range(T+1), y_seq, label='non-financial income', lw=2)\n",
    "    axes[0].plot(range(T+1), c_seq, label='consumption', lw=2)\n",
    "    axes[1].plot(range(T+2), a_seq, label='financial wealth', color='green', lw=2)\n",
    "    axes[0].set_ylabel(r'$c_t,y_t$')\n",
    "    axes[1].set_ylabel(r'$a_t$')\n",
    "    \n",
    "    for ax in axes:\n",
    "        ax.plot(range(T+2), np.zeros(T+2), '--', lw=1, color='black')\n",
    "        ax.legend()\n",
    "        ax.set_xlabel(r'$t$')\n",
    "    \n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0480bbf0",
   "metadata": {},
   "source": [
    "In the experiments below, please study how consumption and financial asset sequences vary across different sequences for non-financial income."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0dc3253e",
   "metadata": {},
   "source": [
    "#### Experiment 1: one-time gain/loss\n",
    "\n",
    "We first assume a one-time windfall of $ W_0 $ in year 21 of the income sequence $ y $.\n",
    "\n",
    "We’ll make $ W_0 $ big - positive to indicate a one-time windfall, and negative to indicate a one-time “disaster”."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0aa89f35",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Windfall W_0 = 2.5\n",
    "y_seq_pos = np.concatenate([np.ones(21), np.array([2.5]), np.ones(24), np.zeros(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_pos)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "13d6c5e1",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Disaster W_0 = -2.5\n",
    "y_seq_neg = np.concatenate([np.ones(21), np.array([-2.5]), np.ones(24), np.zeros(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_neg)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "61420c84",
   "metadata": {},
   "source": [
    "#### Experiment 2: permanent wage gain/loss\n",
    "\n",
    "Now we assume a permanent  increase in income of $ W $ in year 21 of the $ y $-sequence.\n",
    "\n",
    "Again we can study positive and negative cases"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2bac563d",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Positive permanent income change W = 0.5 when t >= 21\n",
    "y_seq_pos = np.concatenate(\n",
    "    [np.ones(21), 1.5*np.ones(25), np.zeros(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_pos)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2ec3a986",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Negative permanent income change W = -0.5 when t >= 21\n",
    "y_seq_neg = np.concatenate(\n",
    "    [np.ones(21), .5*np.ones(25), np.zeros(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_neg)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2fa73008",
   "metadata": {},
   "source": [
    "#### Experiment 3: a late starter\n",
    "\n",
    "Now we simulate a $ y $ sequence in which a person gets zero for 46 years, and then works and gets 1 for the last 20 years of life (a “late starter”)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8e8e294c",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Late starter\n",
    "y_seq_late = np.concatenate(\n",
    "    [np.ones(46), 2*np.ones(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_late)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "10025102",
   "metadata": {},
   "source": [
    "#### Experiment 4: geometric earner\n",
    "\n",
    "Now we simulate a geometric $ y $ sequence in which a person gets $ y_t = \\lambda^t y_0 $ in first 46 years.\n",
    "\n",
    "We first experiment with $ \\lambda = 1.05 $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9751a241",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Geometric earner parameters where λ = 1.05\n",
    "λ = 1.05\n",
    "y_0 = 1\n",
    "t_max = 46\n",
    "\n",
    "# Generate geometric y sequence\n",
    "geo_seq = λ ** np.arange(t_max) * y_0 \n",
    "y_seq_geo = np.concatenate(\n",
    "            [geo_seq, np.zeros(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_geo)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "14e7dfa3",
   "metadata": {},
   "source": [
    "Now we show the behavior when $ \\lambda = 0.95 $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4fa0e44d",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "λ = 0.95\n",
    "\n",
    "geo_seq = λ ** np.arange(t_max) * y_0 \n",
    "y_seq_geo = np.concatenate(\n",
    "            [geo_seq, np.zeros(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_geo)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "61501c64",
   "metadata": {},
   "source": [
    "What happens when $ \\lambda $ is negative"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "67003411",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "λ = -0.95\n",
    "\n",
    "geo_seq = λ ** np.arange(t_max) * y_0 + 1\n",
    "y_seq_geo = np.concatenate(\n",
    "            [geo_seq, np.ones(20)])\n",
    "\n",
    "plot_cs(cs_model, a0, y_seq_geo)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "776ac943",
   "metadata": {},
   "source": [
    "### Feasible consumption variations\n",
    "\n",
    "We promised to justify  our claim that when $ \\beta R =1 $ as Friedman assumed, a constant consumption play $ c_t = c_0 $ for all $ t $ is optimal.\n",
    "\n",
    "Let’s do that now.\n",
    "\n",
    "The approach we’ll take is  an elementary  example  of the “calculus of variations”.\n",
    "\n",
    "Let’s dive in and see what the key idea is.\n",
    "\n",
    "To explore what types of consumption paths are welfare-improving, we shall create an **admissible consumption path variation sequence** $ \\{v_t\\}_{t=0}^T $\n",
    "that satisfies\n",
    "\n",
    "$$\n",
    "\\sum_{t=0}^T R^{-t} v_t = 0\n",
    "$$\n",
    "\n",
    "This equation says that the **present value** of admissible consumption path variations must be zero.\n",
    "\n",
    "So once again, we encounter a formula for the present value of an “asset”:\n",
    "\n",
    "- we require that the present value of consumption path variations be zero.  \n",
    "\n",
    "\n",
    "Here we’ll restrict ourselves to  a two-parameter class of admissible consumption path variations\n",
    "of the form\n",
    "\n",
    "$$\n",
    "v_t = \\xi_1 \\phi^t - \\xi_0\n",
    "$$\n",
    "\n",
    "We say two and not three-parameter class because $ \\xi_0 $ will be a function of $ (\\phi, \\xi_1; R) $ that guarantees that the variation sequence is feasible.\n",
    "\n",
    "Let’s compute that function.\n",
    "\n",
    "We require\n",
    "\n",
    "$$\n",
    "\\sum_{t=0}^T R^{-t}\\left[ \\xi_1 \\phi^t - \\xi_0 \\right] = 0\n",
    "$$\n",
    "\n",
    "which implies that\n",
    "\n",
    "$$\n",
    "\\xi_1 \\sum_{t=0}^T \\phi_t R^{-t} - \\xi_0 \\sum_{t=0}^T R^{-t} = 0\n",
    "$$\n",
    "\n",
    "which implies that\n",
    "\n",
    "$$\n",
    "\\xi_1 \\frac{1 - (\\phi R^{-1})^{T+1}}{1 - \\phi R^{-1}} - \\xi_0 \\frac{1 - R^{-(T+1)}}{1-R^{-1} } =0\n",
    "$$\n",
    "\n",
    "which implies that\n",
    "\n",
    "$$\n",
    "\\xi_0 = \\xi_0(\\phi, \\xi_1; R) = \\xi_1 \\left(\\frac{1 - R^{-1}}{1 - R^{-(T+1)}}\\right) \\left(\\frac{1 - (\\phi R^{-1})^{T+1}}{1 - \\phi R^{-1}}\\right)\n",
    "$$\n",
    "\n",
    "This is our formula for $ \\xi_0 $.\n",
    "\n",
    "**Key Idea:** if $ c^o $ is a budget-feasible consumption path, then so is $ c^o + v $,\n",
    "where $ v $ is a budget-feasible variation.\n",
    "\n",
    "Given $ R $, we thus have a two parameter class of budget feasible variations $ v $ that we can use\n",
    "to compute alternative consumption paths, then evaluate their welfare.\n",
    "\n",
    "Now let’s compute and plot consumption path variations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7237a5a9",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def compute_variation(model, ξ1, ϕ, a0, y_seq, verbose=1):\n",
    "    R, T, β_seq = model.R, model.T, model.β_seq\n",
    "\n",
    "    ξ0 = ξ1*((1 - 1/R) / (1 - (1/R)**(T+1))) * ((1 - (ϕ/R)**(T+1)) / (1 - ϕ/R))\n",
    "    v_seq = np.array([(ξ1*ϕ**t - ξ0) for t in range(T+1)])\n",
    "    \n",
    "    if verbose == 1:\n",
    "        print('check feasible:', np.isclose(β_seq @ v_seq, 0))     # since β = 1/R\n",
    "\n",
    "    c_opt, _, _ = compute_optimal(model, a0, y_seq)\n",
    "    cvar_seq = c_opt + v_seq\n",
    "\n",
    "    return cvar_seq"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8ca6137f",
   "metadata": {},
   "source": [
    "We visualize variations for $ \\xi_1 \\in \\{.01, .05\\} $ and $ \\phi \\in \\{.95, 1.02\\} $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3ec05cf9",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots()\n",
    "\n",
    "ξ1s = [.01, .05]\n",
    "ϕs= [.95, 1.02]\n",
    "colors = {.01: 'tab:blue', .05: 'tab:green'}\n",
    "\n",
    "params = np.array(np.meshgrid(ξ1s, ϕs)).T.reshape(-1, 2)\n",
    "\n",
    "for i, param in enumerate(params):\n",
    "    ξ1, ϕ = param\n",
    "    print(f'variation {i}: ξ1={ξ1}, ϕ={ϕ}')\n",
    "    cvar_seq = compute_variation(model=cs_model, \n",
    "                                 ξ1=ξ1, ϕ=ϕ, a0=a0, \n",
    "                                 y_seq=y_seq)\n",
    "    print(f'welfare={welfare(cs_model, cvar_seq)}')\n",
    "    print('-'*64)\n",
    "    if i % 2 == 0:\n",
    "        ls = '-.'\n",
    "    else: \n",
    "        ls = '-'  \n",
    "    ax.plot(range(T+1), cvar_seq, ls=ls, \n",
    "            color=colors[ξ1], \n",
    "            label=fr'$\\xi_1 = {ξ1}, \\phi = {ϕ}$')\n",
    "\n",
    "plt.plot(range(T+1), c_seq, \n",
    "         color='orange', label=r'Optimal $\\vec{c}$ ')\n",
    "\n",
    "plt.legend()\n",
    "plt.xlabel(r'$t$')\n",
    "plt.ylabel(r'$c_t$')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "469ac038",
   "metadata": {},
   "source": [
    "We can even use the Python `np.gradient` command to compute derivatives of welfare with respect to our two parameters.\n",
    "\n",
    "(We are actually discovering  the key idea beneath the **calculus of variations**.)\n",
    "\n",
    "First, we define the welfare with respect to $ \\xi_1 $ and $ \\phi $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7b9b3708",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def welfare_rel(ξ1, ϕ):\n",
    "    \"\"\"\n",
    "    Compute welfare of variation sequence \n",
    "    for given ϕ, ξ1 with a consumption-smoothing model\n",
    "    \"\"\"\n",
    "    \n",
    "    cvar_seq = compute_variation(cs_model, ξ1=ξ1, \n",
    "                                 ϕ=ϕ, a0=a0, \n",
    "                                 y_seq=y_seq, \n",
    "                                 verbose=0)\n",
    "    return welfare(cs_model, cvar_seq)\n",
    "\n",
    "# Vectorize the function to allow array input\n",
    "welfare_vec = np.vectorize(welfare_rel)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71565d21",
   "metadata": {},
   "source": [
    "Then we can visualize the relationship between welfare and $ \\xi_1 $ and compute its derivatives"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ca573a8",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "ξ1_arr = np.linspace(-0.5, 0.5, 20)\n",
    "\n",
    "plt.plot(ξ1_arr, welfare_vec(ξ1_arr, 1.02))\n",
    "plt.ylabel('welfare')\n",
    "plt.xlabel(r'$\\xi_1$')\n",
    "plt.show()\n",
    "\n",
    "welfare_grad = welfare_vec(ξ1_arr, 1.02)\n",
    "welfare_grad = np.gradient(welfare_grad)\n",
    "plt.plot(ξ1_arr, welfare_grad)\n",
    "plt.ylabel('derivative of welfare')\n",
    "plt.xlabel(r'$\\xi_1$')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "68c0ac53",
   "metadata": {},
   "source": [
    "The same can be done on $ \\phi $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "af68c244",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "ϕ_arr = np.linspace(-0.5, 0.5, 20)\n",
    "\n",
    "plt.plot(ξ1_arr, welfare_vec(0.05, ϕ_arr))\n",
    "plt.ylabel('welfare')\n",
    "plt.xlabel(r'$\\phi$')\n",
    "plt.show()\n",
    "\n",
    "welfare_grad = welfare_vec(0.05, ϕ_arr)\n",
    "welfare_grad = np.gradient(welfare_grad)\n",
    "plt.plot(ξ1_arr, welfare_grad)\n",
    "plt.ylabel('derivative of welfare')\n",
    "plt.xlabel(r'$\\phi$')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cfe17d59",
   "metadata": {},
   "source": [
    "## Wrapping up the consumption-smoothing model\n",
    "\n",
    "The consumption-smoothing model of Milton Friedman [[Friedman, 1956](https://intro.quantecon.org/zreferences.html#id183)] and Robert Hall [[Hall, 1978](https://intro.quantecon.org/zreferences.html#id184)]) is a cornerstone of modern economics that has important ramifications for the size of the Keynesian  “fiscal policy multiplier” that we  described in\n",
    "QuantEcon lecture [geometric series](https://intro.quantecon.org/geom_series.html).\n",
    "\n",
    "The consumption-smoothingmodel   **lowers** the government expenditure  multiplier relative to  one implied by the original Keynesian consumption function presented in [geometric series](https://intro.quantecon.org/geom_series.html).\n",
    "\n",
    "Friedman’s   work opened the door to an enlightening literature on the aggregate consumption function and associated government expenditure  multipliers that remains  active today."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "48404619",
   "metadata": {},
   "source": [
    "## Appendix: solving difference equations with linear algebra\n",
    "\n",
    "In the preceding sections we have used linear algebra to solve a consumption-smoothing model.\n",
    "\n",
    "The same tools from linear algebra – matrix multiplication and matrix inversion – can be used  to study many other dynamic models.\n",
    "\n",
    "We’ll conclude this lecture by giving a couple of examples.\n",
    "\n",
    "We’ll describe a useful way of representing and “solving” linear difference equations.\n",
    "\n",
    "To generate some $ y $ vectors, we’ll just write down a linear difference equation\n",
    "with appropriate initial conditions and then   use linear algebra to solve it."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2e3ea719",
   "metadata": {},
   "source": [
    "### First-order difference equation\n",
    "\n",
    "We’ll start with a first-order linear difference equation for $ \\{y_t\\}_{t=0}^T $:\n",
    "\n",
    "$$\n",
    "y_{t} = \\lambda y_{t-1}, \\quad t = 1, 2, \\ldots, T\n",
    "$$\n",
    "\n",
    "where  $ y_0 $ is a given  initial condition.\n",
    "\n",
    "We can cast this set of $ T $ equations as a single  matrix equation\n",
    "\n",
    "\n",
    "<a id='equation-eq-first-order-lin-diff'></a>\n",
    "$$\n",
    "\\begin{bmatrix} \n",
    "1 & 0 & 0 & \\cdots & 0 & 0 \\cr\n",
    "-\\lambda & 1 & 0 & \\cdots & 0 & 0 \\cr\n",
    "0 & -\\lambda & 1 & \\cdots & 0 & 0 \\cr\n",
    " \\vdots & \\vdots & \\vdots & \\cdots & \\vdots & \\vdots \\cr\n",
    "0 & 0 & 0 & \\cdots & -\\lambda & 1 \n",
    "\\end{bmatrix} \n",
    "\\begin{bmatrix}\n",
    "y_1 \\cr y_2 \\cr y_3 \\cr \\vdots \\cr y_T \n",
    "\\end{bmatrix}\n",
    "= \n",
    "\\begin{bmatrix} \n",
    "\\lambda y_0 \\cr 0 \\cr 0 \\cr \\vdots \\cr 0 \n",
    "\\end{bmatrix} \\tag{12.5}\n",
    "$$\n",
    "\n",
    "Multiplying both sides of [(12.5)](#equation-eq-first-order-lin-diff)  by the  inverse of the matrix on the left provides the solution\n",
    "\n",
    "\n",
    "<a id='equation-fst-ord-inverse'></a>\n",
    "$$\n",
    "\\begin{bmatrix} \n",
    "y_1 \\cr y_2 \\cr y_3 \\cr \\vdots \\cr y_T \n",
    "\\end{bmatrix} \n",
    "= \n",
    "\\begin{bmatrix} \n",
    "1 & 0 & 0 & \\cdots & 0 & 0 \\cr\n",
    "\\lambda & 1 & 0 & \\cdots & 0 & 0 \\cr\n",
    "\\lambda^2 & \\lambda & 1 & \\cdots & 0 & 0 \\cr\n",
    " \\vdots & \\vdots & \\vdots & \\cdots & \\vdots & \\vdots \\cr\n",
    "\\lambda^{T-1} & \\lambda^{T-2} & \\lambda^{T-3} & \\cdots & \\lambda & 1 \n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix} \n",
    "\\lambda y_0 \\cr 0 \\cr 0 \\cr \\vdots \\cr 0 \n",
    "\\end{bmatrix} \\tag{12.6}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f1d631f9",
   "metadata": {},
   "source": [
    "### Exercise 12.1\n",
    "\n",
    "To get [(12.6)](#equation-fst-ord-inverse), we multiplied both sides of  [(12.5)](#equation-eq-first-order-lin-diff) by  the inverse of the matrix $ A $. Please confirm that\n",
    "\n",
    "$$\n",
    "\\begin{bmatrix} \n",
    "1 & 0 & 0 & \\cdots & 0 & 0 \\cr\n",
    "\\lambda & 1 & 0 & \\cdots & 0 & 0 \\cr\n",
    "\\lambda^2 & \\lambda & 1 & \\cdots & 0 & 0 \\cr\n",
    " \\vdots & \\vdots & \\vdots & \\cdots & \\vdots & \\vdots \\cr\n",
    "\\lambda^{T-1} & \\lambda^{T-2} & \\lambda^{T-3} & \\cdots & \\lambda & 1 \n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "is the inverse of $ A $ and check that $ A A^{-1} = I $"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b420d073",
   "metadata": {},
   "source": [
    "### Second-order difference equation\n",
    "\n",
    "A second-order linear difference equation for $ \\{y_t\\}_{t=0}^T $ is\n",
    "\n",
    "$$\n",
    "y_{t} = \\lambda_1 y_{t-1} + \\lambda_2 y_{t-2}, \\quad t = 1, 2, \\ldots, T\n",
    "$$\n",
    "\n",
    "where now $ y_0 $ and $ y_{-1} $ are two given initial equations determined outside the model.\n",
    "\n",
    "As we did with the first-order difference equation, we can cast this set of $ T $ equations as a single matrix equation\n",
    "\n",
    "$$\n",
    "\\begin{bmatrix} \n",
    "1 & 0 & 0 & \\cdots & 0 & 0 & 0 \\cr\n",
    "-\\lambda_1 & 1 & 0 & \\cdots & 0 & 0 & 0 \\cr\n",
    "-\\lambda_2 & -\\lambda_1 & 1 & \\cdots & 0 & 0 & 0 \\cr\n",
    " \\vdots & \\vdots & \\vdots & \\cdots & \\vdots & \\vdots \\cr\n",
    "0 & 0 & 0 & \\cdots & -\\lambda_2 & -\\lambda_1 & 1 \n",
    "\\end{bmatrix} \n",
    "\\begin{bmatrix} \n",
    "y_1 \\cr y_2 \\cr y_3 \\cr \\vdots \\cr y_T \n",
    "\\end{bmatrix}\n",
    "= \n",
    "\\begin{bmatrix} \n",
    "\\lambda_1 y_0 + \\lambda_2 y_{-1} \\cr \\lambda_2 y_0 \\cr 0 \\cr \\vdots \\cr 0 \n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "Multiplying both sides by  inverse of the matrix on the left again provides the solution."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4e0ed43a",
   "metadata": {},
   "source": [
    "### Exercise 12.2\n",
    "\n",
    "As an exercise, we ask you to represent and solve a **third-order linear difference equation**.\n",
    "How many initial conditions must you specify?"
   ]
  }
 ],
 "metadata": {
  "date": 1738361553.8649197,
  "filename": "cons_smooth.md",
  "kernelspec": {
   "display_name": "Python",
   "language": "python3",
   "name": "python3"
  },
  "title": "Consumption Smoothing"
 },
 "nbformat": 4,
 "nbformat_minor": 5
}