/* ================================================================= This simple program solves and simulates a simple model that features exogenous fiscal spending shocks (financed with a lump-sum tax). David Andolfatto Simon Fraser University January 2004 =================================================================*/ library nlsys; /* ==================================================================== Preferences: U(c,l) = ln(c)+psi*(1-eta)^(-1)*( l^(1-eta) - 1) ) Technology: y = An^theta Fiscal policy: ln(g') = (1-rho)*ln(gbar) + rho*ln(g) + e, where e is distributed N(0,sigma^2). gbar here correponds to the average level of real, per capita government purchases. You can estimate the parameters (rho, gbar, sigma) by running a simple regression (on a stationary time series for government purchases). ======================================================================*/ /* Model Calibration */ gbar = 1.0; rho = 0.5; sigma = 0.01; stdg = sigma/(1-rho^2)^(0.5); /* This is the (ergodic) percent standard deviation in g. */ eta = 1; /*calibrated to labor supply elasticity */ theta = 0.65; /* calibrated to labor's share */ nstar = 0.30; /* calibrated to data on hours */ gyratio = 0.20; /* ratio of government purchases to GDP -- calibrated from data) */ /* From the theoretical restriction characterizing the equilibrium level of employment, we have: theta*nstar^(theta-1) = (1-gbar)*nstar^theta*psi*(1-nstar)^(-eta) We can use this restriction to solve for the psi that calibrates the model to the average level of employment; ie. */ psi = theta*nstar^(theta-1) / ( (1-gyratio)*nstar^theta*(1-nstar)^(-eta) ); /* Finally, we can compute the TFP parameter */ A = gbar/( gyratio*nstar^theta ); /* Solving for the equilibrium level of employment */ /* First, specify the grid over g */ lbound = -3*stdg; lbound = exp(lbound); ubound = 3*stdg; ubound = exp(ubound); NG = 11; /* Number of grid points */ step = (ubound - lbound)/(NG -1); gshock = seqa(lbound,step,NG); emp = zeros(NG,1); i = 1; do while i le NG; __output=0; x0 = nstar; {x,f,g,h} = nlsys(&foc,x0); emp[i] = x[1]; i = i+1; endo; /* Print the solution to the screen */ print; "Solution: "; gshock~emp; print; "Press RETURN to continue."; cons; /* Simulating the Model */ length = 100; /* This number should be the length of the sample period that you used to calibrate your model */ semp = zeros(length,1); sgdp = semp; swage = semp; sgov = semp; reps = 100; /* Number of replications in Monte Carlo simulation */ sdgdp = zeros(reps,1); sdemp = sdgdp; sdwage = sdgdp; sdgov = sdgdp; corgdpemp = sdgdp; corgdpwage = sdgdp; corgdpgov = sdgdp; " "; "Simulating the Model..."; print; iter = 1; do while iter le reps; gg = gbar; t = 1; do while t le length; gp = gbar^(1-rho)*gg^rho*exp(sigma*rndn(1,1)); sgov[t] = gp; semp[t] = ln( INTERP(emp,gp) ); sgdp[t] = ln(A) + theta*semp[t]; swage[t] = ln(theta) + ln(A) + (theta-1)*semp[t]; gg = gp; t = t+1; endo; sdgdp[iter] = stdc(sgdp); sdemp[iter] = stdc(semp); sdwage[iter] = stdc(swage); sdgov[iter] = stdc( sgov ); temp = corrx( sgdp~semp ); corgdpemp[iter] = temp[1,2]; temp = corrx( sgdp~swage ); corgdpwage[iter] = temp[1,2]; temp = corrx( sgdp~sgov ); corgdpgov[iter] = temp[1,2]; iter = iter+1; endo; /* Now let's compute some statistics. */ msdgdp = meanc(sdgdp); msdemp = meanc(sdemp); msdwage = meanc(sdwage); msdgov = meanc(sdgov); mcorgdpemp = meanc( corgdpemp ); mcorgdpwage = meanc( corgdpwage ); mcorgdpgov = meanc( corgdpgov ); format /rd 6,4; print; "Covariance Properties of Simulated Time-Series"; print; " Percent Standard Deviations (standard errors in brackets):"; print; " GDP: " 100*msdgdp "( " 100*stdc(sdgdp) ")"; " Employment: " 100*msdemp "( " 100*stdc(sdemp) ")"; " Real wage: " 100*msdwage "( " 100*stdc(sdwage) ")"; " Govt: " 100*msdgov "( " 100*stdc(sdgov) ")"; print; " Correlations (standard errors in brackets): "; print; " Cor(GDP, Employment): " mcorgdpemp "( " stdc(corgdpemp) ")"; " Cor(GDP, Real Wage): " mcorgdpwage "( " stdc(corgdpwage) ")"; " Cor(GDP, Govt): " mcorgdpgov "( " stdc(corgdpgov) ")"; print; print; print; end; /*===============================================================================*/ proc foc(x); local fn1; fn1 = theta*A*x[1]^(theta-1) - (A*x[1]^theta - gshock[i])*psi*(1-x[1])^(-eta); retp( fn1 ); endp; /*================================================================================*/ /* Interpolation Procedure: Here, I want a function that takes as input some arbitrary shock value (not necessarily on the grid) and return the value implied by the decision rule (an NZx1 vector of points). */ proc INTERP(rule,zz); local i, w1; i = maxindc( gshock .ge zz ); i = maxc(i|2); w1 = (1/step)*( zz - gshock[i-1] ); retp( (1-w1)*rule[i-1] + w1*rule[i] ); endp; /*================================================================================*/ end;