B.D. McCullough 1307common method of achieving this level of accuracy is to use FORTRAN with a
multiple precision pre-processor (e.g., Bailey, 1993). An alternative is to employ
the software package Mathematica, which can combine symbolic calculation with
extended precision computation to produce benchmark-quality results. For exam-
ple, on the NIST StRD ANOVA tests, 32-bit word double precision can do no better
than four or five digits of accuracy. Mathematica can, on these problems, return
a full 15 digits of accuracy (McCullough, 2000c). For this analysis, we therefore
employ Mathematica.
Here we consider only Box and Jenkins (1976) (henceforth BJ) methods based on
least squares for two reasons. First, these are by far the most widely-used methods;
only two of the packages considered by NAM offered exact maximum likelihood.
Second, there is a definitive reference for the procedures – BJ. There exist many
alternative methods for computing maximum likelihood methods, and computing
a maximum likelihood benchmark for ARIMA estimation constitutes a project in
itself.
This is not an unimportant topic. ARMA estimation and forecasting is a mainstay
of applied time series analysis. Note that forecasts have not been benchmarked: this
needs to be done, too. Recently, Yalta and Jenal (2009) attempted to double-check,
by hand, the forecasts coming from an ARMA procedure. They could not reproduce
the program’s results. On further investigation, they found that not only were the
forecasts incorrect, even the ARMA coefficients were incorrect. (For a benchmark,
they used the approach of having several packages give the same answer to the
same problem – the package in question did not agree with all the other packages.)
Sub-section 28.6.2 describes definitions and notations. Sub-section 28.6.3
presents CLS results; analytic derivatives are employed for the ARMA(1,1) case, and
comparison with numerical derivatives sheds light on the choice of the differencing
interval for the computation of numerical derivatives. Then, using extended preci-
sion computation, a benchmark for CLS is presented. Sub-section 28.6.4 presents
a benchmark for ULS.
28.6.1 Definitions and notation
The ARMA(p,q) model can be written as:
at=w ̃t−φ 1 w ̃t− 1 −φ 2 w ̃t− 2 −...−φpw ̃t−p+θ 1 at− 1 +θ 2 at− 2 +...θqat−q, (28.4)wherewt=∇dztandw ̃t=wt−μwithE[wt]=μ. In general, whend>0itis
assumed thatμ=0. However, we takeμ =0 as an additional parameter to be
estimated. The form in which the equation is written affects the intercept term
in an ARMA model. For example, BJ write( 1 −#(B))yt=c+( 1 −(B))at, while
an alternate formulation, adopted here, takes( 1 −#(B))(yt−μ)=( 1 −(B))at.
Comparing the two formulations, it is obvious thatμ=c/( 1 −
∑p
j= 1 φp). In the
sequel, only the casep=1,q=1 is treated.
The variation on the Marquardt algorithm proposed by BJ, and which is imple-
mented in several packages, is extremely simple, and focuses on minimizing a sum