{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "lambda=\n", "[1.09007 1.09007 0.5]\n" ] } ], "source": [ "# Solving and simulating a baseline New Keynesian model adapted from Martin Ellison's Matlab code: \n", "# http://users.ox.ac.uk/~exet2581/recursive/recursive.html\n", "# tested in Julia 0.6.1.1\n", "# this code is part of chapter 3, \"Solving and Simulating DSGE Models\" from the book: \"Introduction to Quantitative Macroeconomics using Julia\"\n", "# Academic Press - Elsevier\n", "# for comments, email at: petre(dot)caraiani(at)gmail(dot)com\n", "\n", "#Calibrated parameter values\n", "beta=0.99;\n", "sigma=1.0;\n", "chi=1.55;\n", "eta=0;\n", "theta=2.064;\n", "omega=0.5;\n", "alpha=3;\n", "delta=1.5;\n", "rho=0.5;\n", "\n", "kappa=(1-omega)*(1-beta*omega)/(alpha*omega);\n", "\n", "#Number of predetermined and jump variables\n", "n=1;#predetermined variables\n", "m=2;#jump variables\n", "nu=1;\n", "\n", "cutoff=1.0;\n", "\n", "#Define state space matrices\n", "\n", "A0=zeros(3,3);\n", "A0[1,1]=1;\n", "A0[2,2]=1;\n", "A0[2,3]=sigma^-1;\n", "A0[3,3]=beta;\n", "\n", "A1=zeros(3,3);\n", "A1[1,1]=rho;\n", "A1[2,1]=sigma^-1;\n", "A1[2,2]=1;\n", "A1[2,3]=sigma^-1*delta;\n", "A1[3,2]=-kappa;\n", "A1[3,3]=1;\n", "\n", "B0=zeros(3,1);\n", "B0[1,1]=1;\n", "\n", "#Calculate alternative state space matrices\n", "\n", "A=inv(A0)*A1;\n", "B=inv(A0)*B0;\n", "\n", "lambda=eigvals(A);\n", "\n", "println(\"lambda=\");\n", "println(transpose(real(lambda)));\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "BK condition\n", "Number of jump variables:2\n", "Number of unstable roots:2\n", "BK satistifed\n" ] } ], "source": [ "#Blanchard-Kahn Conditions\n", "\n", "bk_n = 0;\n", "\n", "for i = 1:length(lambda)\n", " if abs.(lambda[i]) > 1.0\n", " bk_n = bk_n+1; \n", " end \n", "end;\n", "\n", "\n", "println(\"BK condition\");\n", "println(\"Number of jump variables:\",m);\n", "println(\"Number of unstable roots:\",bk_n);\n", "\n", "if bk_n==m;\n", " println(\"BK satistifed\")\n", "elseif bk_n>n\n", " println(\"Too many unstable roots\")\n", "else\n", " println(\"Too few unstable roots\")\n", "end;\n", "\n" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "#Sort eigenvectors and eigenvalues\n", "v,w = eig(A)\n", "r=[abs.(v) w']\n", "\n", "for i = 1:(m+n)\n", " for j = i+1:(m+n)\n", " if real(r[i,1])> real(r[j,1])\n", " tmp = r[i,:];\n", " r[i,:]=r[j,:];\n", " r[j,:]=tmp;\n", " elseif real(r[i,1])== real(r[j,1])\n", " if imag(r[i,1])> imag(r[j,1])\n", " tmp = r[i,:];\n", " r[i,:]=r[j,:];\n", " r[j,:]=tmp; \n", " end; \n", " end;\n", " end; \n", "end;\n", "\n", "#r=sortrows(r, by=x->real(x[1]));" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "lam= r[:,1];\n", "M=r[:,2:4]';\n", "\n", "P=inv(M);\n", "lam=diagm(lam);\n", "\n", "#solutions\n", "LAMBDA1=lam[1:n,1:n];\n", "LAMBDA2=lam[(n+1):n+m,(n+1):n+m];\n", "\n", "P11 = P[1:n,1:n];\n", "P12 = P[1:n,(n+1):n+m];\n", "P21 = P[(n+1):n+m,1:n];\n", "P22 = P[(n+1):n+m,(n+1):n+m];\n", "\n", "R=P*B;\n", "\n", "#decision rules\n", "C11=real(inv(P11-P12*inv(P22)*P21)*LAMBDA1*(P11-P12*inv(P22)*P21));\n", "C12=real(inv(P11-P12*inv(P22)*P21)*R[1]);\n", "C21=real(-inv(P22)*P21);\n" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "#IRF\n", "horiz=24;\n", "e=zeros(horiz,1);\n", "x=zeros(horiz,2);\n", "eps=zeros(horiz,1);\n", "eps[1,1]=1;\n", "\n", "for i = 1:horiz-1\n", "e[i+1,:]=C11*e[i,1]+C12*eps[i,1]\n", "\n", "x[i+1,:]=C21*e[i,1]\n", "end;\n", "\n", "\n", "#Pkg.add(\"Plots\")\n", "using Plots\n", "pyplot()\n", "#plotly() # Choose the Plotly.jl backend for web interactivity\n", "plot(x[2:24,2],linewidth=1,label=\"Inflation\",title=\"Impulse Response Function to a monetary policy shock\")\n", "plot!(x[2:24,1],linewidth=1,label=\"Output Gap\")" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "#Simulation and moments\n", "#variable names\n", "name=['x',`pi`];\n", "periods=500;\n", "drop=250;\n", "rep=100;\n", "nlag=4;\n", "nvars=size(name,1);\n", "#matrices store results\n", "stdst=zeros(rep,nvars);\n", "corst=zeros(rep,nvars*nvars);\n", "auto=zeros(rep,nvars*nlag);" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Simulation results\n", "Average standard deviation of output gap:0.04060552607066903\n", "Average standard deviation of inflation:0.013535175356889672\n" ] } ], "source": [ "for i = 1:rep\n", "eps=sqrt(0.01)*rand(1,periods); \n", "x=zeros(periods,2);\n", "e=zeros(periods,1);\n", "\n", "#compute variables \n", "e[1,1]=0; \n", "x[1,:]=0; \n", "for j = 1:periods-1\n", "e[j+1,:]=C11*e[j]+C12*eps[j]\n", "x[j+1,:]=C21*e[j,1]\n", " \n", "end\n", "#drop burnins\n", "stdst[i,1]=std(x[:,1]) \n", "stdst[i,2]=std(x[:,2]) \n", "end\n", "\n", "#Display results\n", "println(\"Simulation results\");\n", "println(\"Average standard deviation of output gap:\",mean(stdst[:,1]));\n", "println(\"Average standard deviation of inflation:\",mean(stdst[:,2]));\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "anaconda-cloud": {}, "kernelspec": { "display_name": "Julia 0.6.1", "language": "julia", "name": "julia-0.6" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "0.6.1" } }, "nbformat": 4, "nbformat_minor": 2 }