{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Solución del Modelo Neoclásico de Crecimiento usando el Shooting Algorithm\n", "\n", "Mauricio M. Tejada\n", "\n", "Departamento de Economía, Universidad Alberto Hurtado\n", "\n", "---" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## El Problema de Optimización " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "El problema de optimización del planificador central es el siguiente:\n", "\n", "\\begin{eqnarray*}\n", "\\max U_{0} & = & \\sum_{t=0}^{\\infty}\\beta^{t} g_L^t \\left[\\frac{(\\hat{c}_{t}g_A^t)^{1-\\sigma}-1}{1-\\sigma}\\right]\\\\\n", "s.a & & \\hat{c}_{t}+\\hat{i}_{t}=\\hat{k}_{t}^{\\alpha}\\\\\n", " & & \\hat{k}_{t+1}g_Ag_L=\\hat{i}_{t}+(1-\\delta)\\hat{k}_{t}\\\\\n", " & & \\hat{k}_{0}\\,dado.\n", "\\end{eqnarray*}\n", "\n", "con $L_0=A_0=1$, $g_L=1+n$ y $g_A=1+x$. Alternativamente:\n", "\n", "\\begin{eqnarray*}\n", "\\max U_{0} & = & \\sum_{t=0}^{\\infty}\\beta^{t}g_L^t\\left[\\frac{(\\hat{c}_{t}g_A^t)^{1-\\sigma}-1}{1-\\sigma}\\right]\\\\\n", "s.a & & \\hat{c}_{t}+\\hat{k}_{t+1}g_A g_L=\\hat{k}_{t}^{\\alpha}+(1-\\delta)\\hat{k}_{t}\\\\\n", " & & k_{0}\\,dado.\n", "\\end{eqnarray*}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "De las CPO tenemos el siguiente sistema de ecuaciones en diferencias (no lineal):\n", "\n", "\n", "\\begin{eqnarray*}\n", "\\hat{c}_{t}^{-\\sigma} & = & \\beta g_A^{-\\sigma}\\hat{c}_{t+1}^{-\\sigma}\\left[\\alpha \\hat{k}_{t+1}^{\\alpha-1}+1-\\delta\\right]\\\\\n", "\\hat{c}_{t}+\\hat{k}_{t+1} g_Ag_L& = & \\hat{k}_{t}^{\\alpha}+1-\\delta\\hat{k}_{t}\n", "\\end{eqnarray*}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Estado Estacionario " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "El estado estacionario resuelve:\n", "\\begin{eqnarray*}\n", "1 & = & \\beta g_A^{-\\sigma}\\left[\\alpha \\hat{k}^{*\\alpha-1}+1-\\delta\\right]\\\\\n", "\\hat{c}^{*} & = & \\hat{k}^{*\\alpha}+(1-\\delta-g_A g_L) \\hat{k}^{*}\n", "\\end{eqnarray*}\n", "\n", "La solución puede ser hallada algebraicamente:\n", "\n", "\\begin{eqnarray*}\n", "k^{*} & = & \\left[\\frac{1-\\beta g_A^{-\\sigma}(1-\\delta)}{\\beta g_A^{-\\sigma}\\alpha }\\right]^{\\frac{1}{\\alpha-1}}\\\\\n", "c^{*} & = & A\\left[\\frac{1-\\beta g_A^{-\\sigma}(1-\\delta)}{\\beta g_A^{-\\sigma} \\alpha }\\right]^{\\frac{\\alpha}{\\alpha-1}}+(1-\\delta-g_A g_L)\\left[\\frac{1-\\beta g_A^{-\\sigma}(1-\\delta)}{\\beta g_A^{-\\sigma}\\alpha }\\right]^{\\frac{1}{\\alpha-1}}\n", "\\end{eqnarray*}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## El Algorítmo de Solución" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "1. Dado $\\hat{k}_0$, conjeturar $\\hat{c}_0$\n", "2. Calcular $\\hat{k}_1,\\hat{k}_2,...,\\hat{k}_t$ y $\\hat{c}_1,\\hat{c}_2,...,\\hat{c}_t$ usando las CPO y para un $t$ grande.\n", " - En este paso se requiere resolver el sistema (no lineal) para hallar $\\hat{k}_{t+1}$ y $\\hat{c}_{t+1}$ dados $\\hat{k}_{t}$ y $\\hat{c}_{t}$.\n", "3. ¿$\\hat{k}_t$ converge a $\\hat{k}^*$? Si esto se cumple tenemos la solución, caso contrario es necesario actualizar $\\hat{c}_0$:\n", " - Si $\\hat{k}_t \\rightarrow \\infty$, existe sobre-ahorro y por tanto es necesario incrementar $\\hat{c}_0$.\n", " - Si $\\hat{k}_t \\rightarrow 0$, existe sub-ahorro y por tanto es necesario incrementar $\\hat{c}_0$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Solución del Modelo" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# Paquetes necesarios\n", "using NLsolve, Plots" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Estado Estacionario (Senda de Crecimiento Balanceado)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# Parámetros\n", "α = 0.3\n", "β = 0.98\n", "σ = 1.5\n", "δ = 0.2\n", "gA = 1 + 0.02\n", "gL = 1 + 0.01\n", "\n", "params = (α, β, σ, δ, gA, gL);" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "sistema_ee (generic function with 2 methods)" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Hallamos el estado estacionario del modelo\n", "\n", "function sistema_ee(x,p)\n", " α, β, σ, δ, gA, gL = p\n", " \n", " k = x[1]\n", " c = x[2]\n", " \n", " fee = zeros(2)\n", " \n", " fee[1] = 1 - β*(gA^(-σ))*(α*(k^(α-1))+1-δ)\n", " fee[2] = (k^α) + (1-δ-gA*gL)*k - c \n", " \n", " return fee\n", "end" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Results of Nonlinear Solver Algorithm\n", " * Algorithm: Trust-region with dogleg and autoscaling\n", " * Starting Point: [1.5, 1.0]\n", " * Zero: [1.2888784172449008, 0.7824048915854173]\n", " * Inf-norm of residuals: 0.000000\n", " * Iterations: 4\n", " * Convergence: true\n", " * |x - x'| < 0.0e+00: false\n", " * |f(x)| < 1.0e-08: true\n", " * Function Calls (f): 5\n", " * Jacobian Calls (df/dx): 5" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "solucion_ee = nlsolve(x -> sistema_ee(x,params), [1.5, 1.0], inplace = false)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Solución Numérica\n", "k^*=1.2888784172449008\n", "c^*=0.7824048915854173\n" ] } ], "source": [ "kstar = solucion_ee.zero[1]\n", "cstar = solucion_ee.zero[2]\n", "\n", "println(\"Solución Numérica\")\n", "println(\"k^*=$kstar\") \n", "println(\"c^*=$cstar\") " ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Solución Algebraica\n", "k^*=1.2888784172449266\n", "c^*=0.7824048915854149\n" ] } ], "source": [ "# Verificamos la solución hallada con NLsolve:\n", "kstar2 = ((1-β*gA^(-σ)*(1-δ))/(β*gA^(-σ)*α))^(1/(α-1))\n", "cstar2 = kstar2^α + (1-δ-gA*gL)*kstar2\n", "\n", "println(\"Solución Algebraica\")\n", "println(\"k^*=$kstar2\") \n", "println(\"c^*=$cstar2\") " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Dinámica de Transición" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "k0=0.6444392086224504\n" ] } ], "source": [ "# Estado Inicial\n", "k0 = 0.5*kstar\n", "println(\"k0=$k0\") \n", "\n", "# Conjetura inicial para le consumo\n", "cL = 0.000001*cstar\n", "cH = (kstar^α) + (1-δ)*kstar\n", "\n", "# Parámetros del método de solución\n", "T = 100\n", "tol = 1e-4\n", "step = 0.1\n", "niter = 200;\n", "flag = 0;" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Utilizamos el método de bisección para actualizar la conjetura del consumo inicial en el algoritmo de solución." ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Solución encontrada!\n", "El consumo inicial es: 0.5992856861680562\n" ] } ], "source": [ "k = zeros(T)\n", "c = zeros(T)\n", "\n", "for i in 1:niter\n", "\n", " c0 = (cL+cH)/2\n", " c[1] = c0\n", " #println([cL, cH])\n", " k[1] = k0\n", " \n", " for t in 1:T-1\n", " \n", " k[t+1] = (k[t]^α + (1-δ)*k[t] - c[t])/(gA*gL)\n", " \n", " if k[t+1]<0\n", " flag = 1 # Consumo inicial muy alto\n", " break\n", " elseif k[t+1]>(kstar^α) + (1-δ)*kstar\n", " flag = 2 # Consumo inicial muy bajo\n", " break\n", " else\n", " c[t+1] = ((c[t]^σ)*β*(gA^(-σ))*(α*(k[t+1]^(α-1)) + 1 - δ))^(1/σ)\n", " end\n", " \n", " if t == T-1\n", " if abs(k[T]-kstar)/kstar < tol && abs(c[T]-cstar)/cstar < tol\n", " flag = 10\n", " elseif k[T]-kstar < -tol\n", " flag = 1 # Consumo inicial muy alto\n", " elseif k[T]-kstar > tol\n", " flag = 2 # Consumo inicial muy bajo\n", " end\n", " end\n", " \n", " end \n", " \n", " if flag == 1\n", " cH = copy(c0) # Eliminamos la parte alta del intervalo\n", " elseif flag == 2\n", " cL = copy(c0) # Eliminamos la parte bajo del intervalo\n", " elseif flag == 10\n", " println(\"Solución encontrada!\")\n", " println(\"El consumo inicial es: $c0\")\n", " break\n", " end\n", " \n", "end\n" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "TG = 30\n", "plot(c[1:TG], title = \"Trayectoria del Consumo\", legend = false,\n", " xlabel = \"Periodos\", ylabel = \"c[t]\", linewidth = 2)" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot(k[1:TG], title = \"Trayectoria del Capital\", legend = false,\n", " xlabel = \"Periodos\", ylabel = \"k[t]\", linewidth = 2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Crecimiento del PIB Per Cápita" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Analizamos el proceso de transición del crecimiento del PIB per cápita hacia su senda de crecimiento balanceado." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "per = 1:length(k)\n", "y = k.^α\n", "ystar = kstar.^α\n", "ypc = y.*gA.^per\n", "ypc_star = ystar.*gA.^per\n", "\n", "plot(per[1:TG],[ypc[1:TG] ypc_star[1:TG]], title = \"PIB Per Cápita\",\n", " xlabel = \"Periodos\", ylabel = \"PIBpc[t]\", label = [\"PIB Per Capita\" \"SCB\"],\n", " linewidth = 2)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Δypc = ((ypc[2:end]./ypc[1:end-1]).-1)*100\n", "Δypc_star = ((ypc_star[2:end]./ypc_star[1:end-1]).-1)*100\n", "\n", "plot(per[1:TG],[Δypc[1:TG] Δypc_star[1:TG]], title = \"Crecimiento del PIB Per Cápita\",\n", " xlabel = \"Periodos\", ylabel = \"Porcentaje\", label = [\"PIB Per Capita\" \"SCB\"],\n", " linewidth = 2, ylim = (0,10))" ] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.6.0", "language": "julia", "name": "julia-1.6" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.6.0" } }, "nbformat": 4, "nbformat_minor": 4 }