Perturbation AIM (updated
4/3)
Perturbation AIM is an algorithm that solves a linear or nonlinear
dynamic, discrete-time system of rational expectations equations to
arbitrarily high order, in the sense of a Taylor series
approximation. The algorithm requires that your system of equations
has at least one steady-state point around which to approximate a
solution, that a solution exists in at least a neighborhood of that
steady-state point, and that the solution and original system of
equations be analytic (or at least sufficiently smooth) in at least a
neighborhood of that steady-state point. The solution does not
necessarily need to be first-order stable around the steady state:
unit roots are handled automatically, and even explosive first-order
roots are acceptable but require the modeler to specify which
explosive roots are to be allowed (this requires the user to make
minor modifications to the code).
An important feature of the algorithm is that, so long as the solution
to your problem is analytic, Perturbation AIM provides
a global solution in the
sense that it is guaranteed to converge to the true solution over the
whole domain of convergence of the Taylor series, and this domain is
often very large. Moreover, the convergence
is uniform on compact
subsets, which implies that there exists a finite order of approximation n which achieves any desired level of accuracy over
any given (possibly large) compact set. Note the contrast here with
fixed-order approximation methods such as log-linearization or
second-order approximation, which only provide local solutions.
For an introduction to the problem and details of the solution
algorithm, read the following paper:
| “Higher-Order Perturbation
Solutions to Dynamic, Discrete-Time
Rational Expectations Models” (Eric Swanson, Gary Anderson, and Andrew
Levin). |
(abstract)
(PDF)
|
We have implemented the Perturbation AIM algorithm in
Mathematica. While second-order approximations can be done
fairly easily in any computer language, generalization to higher
orders is dramatically simpler in Mathematica, and Mathematica’s
large library of built-in routines for symbolic differentiation and
solving linear systems of equations has the following advantages:
- Allows for much simpler and intuitive code (less than 200
lines), making the implementation of the algorithm more transparent
and more likely to be free of human coding error (bugs).
- Allows the algorithm to manipulate symbolic expressions and thus
proceed in the same way and produce the same output as if the user
were performing it by hand, making the implementation of the algorithm
and the generated output more intuitive.
- Eliminates the need for user-written differentiation and linear
solution routines, which are likely to be more amateur than built-in
Mathematica routines and hence more likely to suffer from bugs and
numerical inaccuracies.
- Allows the use of symbolic differentiation, improving numerical
accuracy, particularly at higher orders.
- Allows the computation of all coefficients to arbitrarily high
precision, eliminating inaccuracies from machine-precision
computations that cumulate with each recursive step of the
algorithm.
In addition to the internal advantages of PerturbationAIM described
above, there are a number of advantages from the economic modeler's
point of view:
- No need to categorize variables as “predetermined”
or “non-predetermined.” Knowing that all variables
dated t−1 or earlier
are observed and all variables dated t are to be solved is sufficient for the computer to
determine which linear combinations of variables are
“predetermined” and which “non-predetermined.”
Anderson and Moore (1985), Sims (2000), and Collard and Juillard
(2001) also make this observation. This feature of the code becomes
increasingly convenient as the number of variables in the model
increases.
- No need for the modeler to substitute out identities. The
computer is completely
capable of handling these substitutions.
- Ability to solve models with arbitrarily long lags and leads
(including shocks dated later than t
or t+1), with no need for the
modeler to transform the model by hand. Note also that the usual trick
for linear models of defining an auxiliary variable zt ≡ Etxt+1,
and then considering zt+1,
fails for nonlinear models because F(Etxt+1) is unequal to EtF(xt+1) in general. Perturbation AIM
handles multiple leads correctly in the nonlinear as well as the linear
case.
- Ability to solve models with non-normally distributed shocks.
Perturbation AIM allows the modeler to substitute in any set of moments
for the stochastic shocks of the model, allowing consideration of
models with essentially any shock distribution (for which the moments
are finite).
New Version of the Code (Perturbation AIM
2.7) (updated 4/3)
The current version of the code is Perturbation AIM 2.7. Relative to
version 2.6, perturbationAIM 2.7 contains an internal optimization
that speeds up the code considerably; these modifications are small in
a coding sense, but the gain in speed is large enough that it merits
being released to the public as a new version.
Relative to versions 2.5 and older, Perturbation AIM 2.7 and 2.6
contain some numerical improvements (better precision and numerical
stability) and some user interface enhancements (catch a few more user
input errors, report estimated precision), but otherwise the syntax
and function output of versions 2.7 and 2.6 are identical to version
2.5.
Relative to version 2.4, Perturbation AIM 2.5, 2.6, and 2.7 feature:
- A new AIMPrecision environment variable, which allows the user to
easily specify the internal numerical precision that Mathematica should
use to solve the model. Although Perturbation AIM has always had
the
capability of using aribtrary-precision arthmetic, the new environment
variable and other internal enhancements make this feature of the code
much more user-friendly.
- Greater numerical precision and robustness. For example, when
steady-state or first-order coefficients are found to be essentially
zero (less than the AIMZeroTol environment variable), then Perturbation
AIM now imposes exact zeros for those coefficients and
re-solves with the hard zero constraint imposed. This increases the
numerical precision of the solution.
- The ability to call Perturbation AIM with model equations and
parameter values as separate arguments. This makes parameter estimation
with Perturbation AIM much faster, since now the derivatives of the
model only need to be computed once and then the various parameter
vectors can be substituted into those derivatives repeatedly.
- Improved compatibility with Mathematica 6.0. The previous version
of Perturbation AIM generated occasional warnings regarding deprecated
function calls. Those function calls have been updated in the current
version of the code.
Some examples that utilize these new features of the code are provided
in the “Using the Code” section below.
Obtaining the Code
Click on one of the following links to download Perturbation AIM (if
your browser loads it in as a text file, then just save the text file
to your hard drive):
For serious research work, you will want to use the most recent
version of the code (version 2.7). The previous versions (especially
2.4 and the “simple” version 2.4) are provided for those
who wish to understand how the algorithm works. The
“simple” version has about 200 lines of code (plus 100
lines of comments) and has been designed to be more transparent than
the full version, hence easier to understand. The “full”
version 2.4 is functionally identical to the “simple”
version 2.4 but is a little more complicated and less transparent due
to the implementation of some additional internal tricks that speed up
computation and reduce the likelihood of numerical inaccuracies,
particularly for larger models and higher orders of
approximation. Versions 2.5 and 2.6 of the code have even more
internal tricks, optimizations, and enhancements that improve the
numerical robustness and flexibility of the code.
No matter which version of the code you use, you
are strongly encouraged to
use its arbitrary-precision features when computing your final
results, as experience has shown that numerical problems are
surpisingly likely to occur, in any higher-order code written in any
computer language. Iteratively solving huge systems of linear
equations is inherently going to lose precision at every iterative
step, not to mention potential condition number problems,
near-singularities, and numerical inaccuracies that are increasingly
likely to arise as one considers more and more curved
models. Solutions computed using the arbitrary-precision features of
the code are much more
likely to be free of these kinds of numerical problems.
All versions of the code above require Mathematica version
5—while it is not difficult to modify the routines to work in
earlier versions of Mathematica, those earlier version
are much slower at solving
systems of linear equations. Perturbation AIM 2.5 and 2.6 were
designed to work with Mathematica 6 and may not work with earlier
versions, although they should still work with Mathematica 7. Note
that Mathematica 7 introduces parallel computing procedures that I
will try to incorporate into the perturbation AIM code at some point
in the near future.
Using the Code
To run the code, you must first start up Mathematica, and then type
from a command line (including the “<<”)
<<perturbationAIMsimple.m
If you are using Windows rather than Unix, you may first have to tell
Mathematica which directory or folder you are working in. For example,
on my laptop, the command:
SetDirectory["C:\Documents and
Settings\Eric T. Swanson\My Documents\papers\pert"]
changes the current working directory to the one in quotes.
Alternatively, you can use Mathematica’s pulldown menus to input
the perturbationAIM.m file. This will define all of the functions you
need to solve rational expectations models to arbitrarily high
order. Of course, you need to define the equations and
parameters of your model. The following example files (particularly
the first one) should serve as useful templates and include many
comments, tips, and hints:
| etsSG.m |
standard stochastic growth model.
|
testpert.m
|
suite of toy model examples
showing some of the capabilities of the code.
|
rssmodel.m
|
baseline model from Rudebusch
and Swanson (2008 JME), shows use of AIMPrecision environment variable.
(updated 2/4)
|
rssmodeliterate.m
|
same model as above, but
efficiently calls Perturbation AIM repeatedly with a range of
parameter values. (updated 2/4)
|
etsSGarbp.m
|
stochastic growth model with
solution coefficients computed to 50 significant digits (this is for
version 2.4 of the code; in version 2.5, just use the model etsSG.m
above and set AIMPrecision = 70 and AIMZeroTol = 10-20 to
get the same result.)
|
To run the example file etsSG.m, simply type:
<<etsSG.m
at a Mathematica command prompt (or use the pulldown menus). In order
to play with the output, substitute in numerical values for the
moments of your shocks, and so on, you will probably want to become
familiar with basic Mathematica commands and
syntax. The Mathematica
Documentation Center provides complete and searchable
documentation for all of Mathematica,
including The
Mathematica Book, which can also be purchased in hardcopy from the
usual places.
Some Caveats
Using numerical computation methods is always a bit of an art, and the
Perturbation AIM code is no exception. With standard,
machine-precision computation, numerical difficulties can easily
arise, and in fact our experience has been that they occur with a
surprising degree of frequency—for example, even in the simple
stochastic growth model above, a few of the machine-precision
coefficients in the third-order expansion are incorrect, as you can
verify by checking them against their arbitrary-precision
counterparts. Moreover, our experience has been that including
representative-agent welfare as one of the equations of the model (a
fairly common application of second-order methods) has a tendency to
increase numerical problems, partly because steady-state welfare is
typically 100 times larger than the steady-state values of the other
variables in the model.
The simplest, most reliable way to avoid numerical inaccuracies in
your answers is to use the aribtrary-precision feature of the code,
and to experiment a little with the level of precision (using the
AIMPrecision environment variable) to verify that the solution does
not change meaningfully. The user should also experiment with varying
the AIMZeroTol environment variable, which is the smallest numerical
value that the code treats as being economically meaningful, as
opposed to being simply an artifact of finite-precision
arithmetic. The default value of AIMZeroTol is 10-10. Note
that, as with all machine-precision or finite-precision computations,
a smaller value for this number does not necessarily lead to more
accurate results, and can even lead to less accurate results.
Thus, users of Perturbation AIM are strongly
encouraged to compute their final results using the
arbitrary-precision feature of the code. While this comes at some cost
in terms of computation time, the cost is generally small compared to
the bottleneck in the program, which experience has shown to be the
computation of symbolic derivatives.