{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Make graphs appear within notebook\n",
    "%matplotlib inline\n",
    "\n",
    "# Import numerical computations library\n",
    "import numpy as np\n",
    "\n",
    "# Import plotting library\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Import dataframe management library\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Restore old behavior of rounding default axis ranges\n",
    "from matplotlib import rcParams\n",
    "\n",
    "rcParams['axes.autolimit_mode'] = 'round_numbers'\n",
    "rcParams['axes.xmargin'] = 0\n",
    "rcParams['axes.ymargin'] = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def k_star(s=0.2, δ=0.05, n=0.01, g=0.02, α=1/3):\n",
    "    return (s/(δ+n+g+n*g))**(1/(1-α))\n",
    "\n",
    "def f(k, α=1/3):\n",
    "    return k**α"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def c_star(s=0.2, δ=0.05, n=0.01, g=0.02, α=1/3):\n",
    "    return (1-s)*f(k_star(s, δ, n, g, α), α)\n",
    "\n",
    "def c_k(k, s=0.2, α=1/3):\n",
    "    return (1-s)*f(k, α)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "δ=0.05\n",
    "n=0.01\n",
    "g=0.02\n",
    "α=1/3\n",
    "\n",
    "ss = np.linspace(0, 1, 1000)\n",
    "\n",
    "plt.plot(100*ss, c_star(ss, δ, n, g, α))\n",
    "plt.vlines(100*α, 0, c_star(α, δ, n, g, α), lw=1, linestyle='--')\n",
    "\n",
    "plt.title('Steady state consumption as function of saving rate')\n",
    "plt.ylabel('Steady state consumption')\n",
    "plt.xlabel('Saving rate (%)')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def k_next(k, s=0.2, δ=0.05, n=0.01, g=0.02, α=1/3):\n",
    "    return (s*f(k, α)+(1-δ)*k)/((1+n)*(1+g))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "T = 100\n",
    "k_t = np.zeros(T)\n",
    "c_t = np.zeros(T)\n",
    "\n",
    "print('Initial saving rate below golden rule saving rate')\n",
    "\n",
    "s_old = 0.3\n",
    "s_GR  = α\n",
    "\n",
    "k_t[:10] = k_star(s=s_old)\n",
    "c_t[:10] = c_k(k_t[:10], s=s_old)\n",
    "for t in range(10, T):\n",
    "    k_t[t] = k_next(k_t[t-1], s=s_GR)\n",
    "    c_t[t] = c_k(k_t[t], s=s_GR)\n",
    "\n",
    "plt.plot(k_t)\n",
    "plt.hlines(k_star(s=s_GR), 0, T, lw=1, linestyle='--')\n",
    "plt.title('Capital after an increase in saving rate')\n",
    "plt.show()\n",
    "\n",
    "plt.plot(c_t)\n",
    "plt.hlines(c_star(s=s_GR), 0, T, lw=1, linestyle='--')\n",
    "plt.title('Consumption after an increase in saving rate')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "T = 100\n",
    "k_t = np.zeros(T)\n",
    "c_t = np.zeros(T)\n",
    "\n",
    "print('Initial saving rate above golden rule saving rate')\n",
    "\n",
    "s_old = 0.35\n",
    "s_GR  = α\n",
    "\n",
    "k_t[:10] = k_star(s=s_old)\n",
    "c_t[:10] = c_k(k_t[:10], s=s_old)\n",
    "for t in range(10, T):\n",
    "    k_t[t] = k_next(k_t[t-1], s=s_GR)\n",
    "    c_t[t] = c_k(k_t[t], s=s_GR)\n",
    "\n",
    "plt.plot(k_t)\n",
    "plt.hlines(k_star(s=s_GR), 0, T, lw=1, linestyle='--')\n",
    "plt.title('Capital after a decrease in saving rate')\n",
    "plt.show()\n",
    "\n",
    "plt.plot(c_t)\n",
    "plt.hlines(c_star(s=s_GR), 0, T, lw=1, linestyle='--')\n",
    "plt.title('Consumption after a decrease in saving rate')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Read dataset\n",
    "pwt = pd.read_stata('PWT/pwt90.dta')\n",
    "\n",
    "# Modify the dataset so that it is easier to work with\n",
    "\n",
    "# Store country names and codes for later use\n",
    "countries = pwt['country']\n",
    "countries = countries.drop_duplicates()\n",
    "\n",
    "countrycodes = pwt['countrycode']\n",
    "countrycodes = countrycodes.drop_duplicates()\n",
    "\n",
    "# Set MultiIndex\n",
    "pwt.set_index(['countrycode', 'year'], inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "(100*pwt.loc['TWN']['csh_i']).plot(lw=2, label='Taiwan')\n",
    "(100*pwt.loc['KOR']['csh_i']).plot(lw=2, label='South Korea')\n",
    "(100*pwt.loc['SGP']['csh_i']).plot(lw=2, label='Singapore')\n",
    "(100*pwt.loc['HKG']['csh_i']).plot(lw=2, label='Hong Kong')\n",
    "\n",
    "plt.xlim(1950, 2014)\n",
    "plt.ylim(0, 70)\n",
    "\n",
    "plt.title('Investment rates in \"growth miracle\" countries')\n",
    "plt.xlabel('Year')\n",
    "plt.ylabel('Investment share of GDP (%)')\n",
    "\n",
    "plt.legend(loc='lower right', ncol=2, frameon=False)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
