@-DEFINE PARAMETERS, DEFINE SOME USEFUL FUNCTIONS, COMPUTE STEADY STATE-@ library nlsysmt, pgraph; format /rds 6,4; #include nlsysmt.sdf struct nlControl nlc; struct nlOut nlo; nlc = nlControlCreate; nlc.fvtol = 0.000000001; nlc.output = 0; pqgwin many; graphset; _pdate=0; _pnum=2; _plwidth = 6; _plclr = 6; _pltype={6,3,2}; #include hpf; beta = 0.99; delta = 0.02; theta = 0.36; gama = 1.5; rho = 0.95; lambda = 1600; nvar = 3; @-- Order of variables x[1] = capital x[2] = shock x[3] = cons --@ fn rrate(x) = theta*exp(x[2])*x[1]^(theta-1); fn prodn(x) = exp(x[2])*x[1]^theta; fn muc(x) = x[3]^(-gama); @--Define the left and right hand sides of the First Order Conditions--@ fn lhs1(x)=x[1]; fn rhs1(x)=(1-delta)*x[1]+prodn(x)-x[3]; fn lhs2(x)=x[2]; fn rhs2(x)=rho*x[2]; fn lhs3(x)=beta*muc(x)*(1-delta+rrate(x)); fn rhs3(x)=muc(x); @-Define a procedure to check whether FOCs are satisfied-@ proc sse(x); local e; e = zeros(nvar,1); e[1] = lhs1(x)-rhs1(x); e[2] = lhs2(x)-rhs2(x); e[3] = lhs3(x)-rhs3(x); retp(e); endp; @-Make an initial guess for the steady state-@ x0=zeros(nvar,1); x0[1] = 58; x0[2] = 0; x0[3] = 4; @-Compute Steady State-@ {nlc, nlo}=nlsys(nlc, &sse, x0); xbar = nlo.xp; print "steady state capital, shock, consumption"; xbar'; @-COMPUTE DECISION RULES-@ psi=(gradp(&lhs1,xbar)|gradp(&lhs2,xbar)|gradp(&lhs3,xbar)); phi=(gradp(&rhs1,xbar)|gradp(&rhs2,xbar)|gradp(&rhs3,xbar)); A = inv(psi)*phi; { L,Q } = eigv(A); IQ = inv(Q); if (sumc(abs(L).>1)/=1); print "either too few or too many unstable roots"; end; endif; if abs(L[1])>1; ck = -IQ[1,1]/IQ[1,3]; cz = -IQ[1,2]/IQ[1,3]; elseif abs(L[2])>1; ck = -IQ[2,1]/IQ[2,3]; cz = -IQ[2,2]/IQ[2,3]; elseif abs(L[3])>1; ck = -IQ[3,1]/IQ[3,3]; cz = -IQ[3,2]/IQ[3,3]; endif; con_rule = ck~cz; print "coefficients on capital and productivity shock in consumption rule "; con_rule; cap_lom = A[1,1:2] + A[1,3]*con_rule; print "coefficients on capital and productivity shock in capital lom "; cap_lom; z_lom = A[2,1:2] + A[2,3]*con_rule; print "coefficients on capital and productivity shock in prod lom "; z_lom; @-Compute nonstochastic dynamics to steady state-@ horizon = 200; ks = zeros(horizon+1,1); cs = zeros(horizon,1); ps = zeros(horizon,1); ks[1] = 0.1*xbar[1]; i=1; do while i<=horizon; ks[i+1] = cap_lom*( ks[i]-xbar[1]|ps[i]-xbar[2] ) + xbar[1]; cs[i] = con_rule*( ks[i]-xbar[1]|ps[i]-xbar[2] ) + xbar[3]; i=i+1; endo; title("path for capital"); xy(0,ks); title("path for consumption"); xy(0,cs); @--Now want to compute some statistics--@ horizon = 200; stdveps = 0.01; noit = 100; std_Chp = zeros(noit,1); std_Xhp = zeros(noit,1); std_Yhp = zeros(noit,1); cor_Y_X = zeros(noit,1); cor_Y_C = zeros(noit,1); hpmat=hpf(horizon,lambda); dte = 1; do while dte <= noit; ks = zeros(horizon+1,1); cs = zeros(horizon,1); xs = zeros(horizon,1); ps = zeros(horizon+1,1); ys = zeros(horizon,1); inn_seq = stdveps*rndn(1,horizon+1); @ Set initial conditions @ ks[1] = xbar[1]; ps[1] = inn_seq[1]; i=1; do while i<=horizon; ks[i+1] = cap_lom*( ks[i]-xbar[1]|ps[i]-xbar[2] ) + xbar[1]; ps[i+1] = rho*ps[i] + inn_seq[i+1]; cs[i] = con_rule*( ks[i]-xbar[1]|ps[i]-xbar[2] ) + xbar[3]; ys[i] = exp(ps[i])*ks[i]^theta; xs[i] = ks[i+1]-(1-delta)*ks[i]; i=i+1; endo; if (dte==1); title("simulation paths for output, consumption and investment"); _plegctl = {2,4,1,5}; _plegstr = "GDP\0C\0X"; xy(0,ys/meanc(ys)~cs/meanc(cs)~xs/meanc(xs)); endif; @ HP filter logged simulation output @ varhp = hpmat*(ln(cs)~ln(ys)~ln(xs)); std_Chp[dte] = stdc(varhp[.,1]); std_Yhp[dte] = stdc(varhp[.,2]); std_Xhp[dte] = stdc(varhp[.,3]); @ varhp = hpmat*(cs~ys~xs); std_Chp[dte] = stdc(varhp[.,1]./(cs-varhp[.,1])); std_Yhp[dte] = stdc(varhp[.,2]./(ys-varhp[.,2])); std_Xhp[dte] = stdc(varhp[.,3]./(xs-varhp[.,3])); @ corv_C = corrx(varhp[.,1]~varhp[.,2]); cor_Y_C[dte] = corv_C[2,1]; corv_X = corrx(varhp[.,3]~varhp[.,2]); cor_Y_X[dte] = corv_X[2,1]; dte = dte + 1; endo; title("GDP and HP trend in GDP"); _plegstr = "ln(GDP)\0HP trend in ln(GDP)"; xy(0,ln(ys)~ln(ys)-varhp[.,2]); print "ave. % std. dev. of GDP " 100*meanc(std_Yhp); print "ave. % std. dev. of consumption " 100*meanc(std_Chp); print "ave. % std. dev. of investment " 100*meanc(std_Xhp); print "ave. consumption GDP correlation " meanc(cor_Y_C); print "ave. investment GDP correlation " meanc(cor_Y_X); @ --OUTPUT FROM PROGRAM FOLLOWS steady state capital, shock, consumption 48.2992 0.0000 3.0725 coefficients on capital and productivity shock in consumption rule 0.0336 0.9595 coefficients on capital and productivity shock in capital lom 0.9765 3.0791 coefficients on capital and productivity shock in prod lom 0.0000 0.9500 ave. % std. dev. of GDP 1.2754 ave. % std. dev. of consumption 0.4210 ave. % std. dev. of investment 4.1549 ave. consumption GDP correlation 0.9585 ave. investment GDP correlation 0.9918 --@