/* =================================================================
This simple program shows you what you have to do in order to solve and
simulate the simple model described in Assignment 2
David Andolfatto
Simon Fraser University
February 2001
=================================================================*/
library nlsys; /* This statement just calls up the nonlinear equations solver */
/* PARAMETER VALUES */
/* There are three key parameters to the model: (1) the curvature parameter on leisre
(eta); (2) the weighting parameter on leisure (psi); and the curvature parameter
on labour in the production technology (theta). In this environment, it can be
demonstrated that the labour supply elasticity in a "steady state" (i.e., z(t)=0)
is equal to:
LSE = eta^(-1)*( 1/lstar -1)
where lstar = fraction of time spent in the nonmarket sector. Typically, we get
measurements like lstar=0.67 and LSE=1.0; plugging these values in, we derive
eta = 2.0.
Regarding the parameter theta, we know that this corresponds to labour's
share of output, so set this equal to 0.70.
Now, from the steady state FOC, we have the expression:
theta*(1-lstar)^(theta-1) = psi*lstar^(-eta)
We can use this expression to solve for the psi that will imply that the average
fraction of time spent in employment is 1/3. */
eta = 2.0;
theta = 0.70;
lstar = 2/3;
psi = ( theta*(1-lstar)^(theta-1) )/( lstar^(-eta) );
/* Now, we'll want to specify a "grid" for our state variable (the technology
parameter). Note that z(t) represents the percent deviation of TFP from
its "trend" (in this case, trend=0). In the data, TFP fluctuates in the order of
three percent (or so) around an HP trend. For simplicity, why don't we
just specify a 5-point grid ranging from plus 2% to minus 2%. Here's how to do
this using GAUSS's seqa procedure (see pg. 1552 of the manual). */
step = 0.01;
zshock = seqa(-0.02,step,5); /* A 5x1 vector of technology parameters */
/* SOLVING FOR THE DECISION RULE */
/* The next step is to solve for the unknown function n(z). We'll do this by
solving for the equilibrium labour input at each point in the grid. Since
the FOC is nonlinear in the labour input, we'll have to use GAUSS's
nonlinear equations solver. */
NZ = rows(zshock);
emp = zeros(NZ,1);
i = 1;
do while i le NZ;
__output=0;
x0 = (1-lstar);
{x,f,g,h} = nlsys(&foc,x0);
/* NLSYS is making a call to a procedure named FOC that I have created below. */
emp[i] = x[1]; /*Store the solution */
i = i+1;
endo;
/* Notice that we have only solved for the function n(z) over a discrete set of points. We can "fill in"
the "missing" parts of the function by linearly interpolating between function evaluations. In this
way, we will have approximated the true solution function with a piecewise linear function. Note
that we can make the approximation as fine as we want simply by specifying a finer and
finer grid. Below, I provide an interpolation procedure called INTERP .*/
/* Suppose, for example, that we want to evaluate our piecewise linear function at a value for z
equal to 0.005 (which is not a grid point). All you have to do now is write:
employment = INTERP(emp,0.005); */
/* SIMULATING THE MODEL */
/* Let us specify the shock process as a first-order Markov process; i.e., z' = rho*z + e. Assume
that the innovation "e" is distributed normally with zero mean and a standard deviation sigma.
Note that for this process, the "long run" standard deviation of the shock "z" is equal to:
std(z) = sigma/sqrt(1-rho^2).
A common parameterization in the RBC literature is to set rho=0.96 and sigma=0.007; in
which case std(z) = 0.025. Notice that our grid goes only from [-0.02,0.02], so that simulating
this process will in good probability take us outside our grid. It is easy enough to just increase
the interval on the grid. But for now, why don't we just specify a smaller value for sigma so that
the realizations for z will (almost always) lie with our current grid. */
rho = 0.96;
sigma = 0.001;
length = 100000; /* Length of simulation (number of periods) */
semp = zeros(length,1); /* Store the simulated employment series in here. */
" ";
"Simulating the Model..."; print;
zz = 0;
t = 1;
do while t le length;
zp = rho*zz + sigma*rndn(1,1);
semp[t] = INTERP(emp,zp);
zz = zp;
t = t+1;
endo;
/* Now let's compute some statistics. */
mnemp = meanc(semp);
sdemp = stdc(ln(semp));
format /rd 6,4;
"Average level of employment: " mnemp;
"Percent standard deviation of employment: " 100*sdemp;
print;
/* Now you can experiment to see how the volatility of employment depends on the model's
parameters, for example, eta. Why don't you explore some of the model's properties? */
/*===============================================================================*/
proc foc(x);
local fn1;
fn1 = exp(zshock[i])*theta*x[1]^(theta-1) - 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( zshock .ge zz );
w1 = (1/step)*( zz - zshock[i-1] );
retp( (1-w1)*rule[i-1] + w1*rule[i] );
endp;
/*================================================================================*/
end;