Numerical Maximum Likelihood with ARMA Constraints

Why Numerical Optimization Needs Care

An ARMA likelihood is a function of several parameters whose effects interact: changing an autoregressive coefficient changes the innovations, and changing a moving-average coefficient changes them again. A numerical optimizer can search this function, but an unconstrained search may propose coefficients that violate causality or invertibility. Even when every candidate is admissible, a poor starting point or an inadequately checked stopping rule can produce an unreliable result.

The goal here is a practical and mathematically controlled workflow. We define the likelihood being optimized, map unconstrained numerical coordinates to valid ARMA parameters, choose starting values, and check whether an apparent optimum is credible. The worked calculations use a conditional Gaussian likelihood with a specified initial-value approximation; the same optimization principles apply when the likelihood is evaluated by the innovations algorithm.

Model and Objective Function

Use the lag-polynomial convention \(y_t-\mu=\phi(B)^{-1}\theta(B)\varepsilon_t\), where \(\phi(B)=1-\phi_1B-\cdots-\phi_pB^p\) and \(\theta(B)=1+\theta_1B+\cdots+\theta_qB^q\). The innovations satisfy \(\varepsilon_t\sim N(0,\sigma^2)\). The usual likelihood interpretation assumes Gaussian innovations; a conditional Gaussian likelihood is also a useful objective when Gaussianity is treated as a working assumption. Causality requires the roots of \(\phi(z)\) to lie outside the unit circle, and invertibility requires the roots of \(\theta(z)\) to lie outside it.

For a conditional likelihood, fix an initialization rule for the unobserved initial values and compute innovations recursively. Write the resulting residuals as \(e_1,\ldots,e_m\), with \(m\) the number of contributions included in the chosen conditional objective. Conditional on the initialization, the Gaussian log-likelihood is

where \(\boldsymbol\beta\) collects the mean and ARMA coefficients. For fixed coefficients, let \(S(\boldsymbol\beta)=\sum_{t=1}^{m}e_t(\boldsymbol\beta)^2\). If \(S(\boldsymbol\beta)>0\), maximizing over positive \(\sigma^2\) gives \(\widehat{\sigma}^2=S(\boldsymbol\beta)/m\). Substitution yields a profiled likelihood, so maximizing it over the ARMA coefficients is equivalent to minimizing \(S(\boldsymbol\beta)\). If \(S(\boldsymbol\beta)=0\), this formula does not apply: the likelihood increases without bound as \(\sigma^2\downarrow0\), and there is no positive-variance maximum-likelihood estimate at those coefficients.

Transforming Parameters to Preserve Stability

For AR(1), causality requires \(|\phi_1|<1\). An unconstrained coordinate \(\eta\in\mathbb R\) can be mapped into this interval by \(\phi_1=\tanh(\eta)\). The map is smooth and one-to-one, with derivative \(d\phi_1/d\eta=1-\tanh^2(\eta)\). It guarantees an admissible coefficient for every finite \(\eta\).

For higher orders, independently applying \(\tanh\) to each ordinary AR coefficient is not sufficient: the stationary region is jointly constrained. One standard approach uses reflection coefficients \(\kappa_j\in(-1,1)\). In one common recursion, begin with the order-one coefficient \(\phi_{1,1}=\kappa_1\), then construct coefficients at order \(r\) by \(\phi_{r,r}=\kappa_r\) and \(\phi_{r,j}=\phi_{r-1,j}-\kappa_r\phi_{r-1,r-j}\), for \(j=1,\ldots,r-1\). This parameterization maps the reflection coefficients into a causal AR polynomial. Set each \(\kappa_j=\tanh(\eta_j)\) to give the optimizer unconstrained coordinates. For the MA polynomial under the convention above, apply the same stable-polynomial construction to coefficients \(a_j=-\theta_j\), then set \(\theta_j=-a_j\); this enforces invertibility.

A variance parameter can be represented as \(\sigma^2=\exp(\lambda)\), which ensures positivity. If the variance is profiled out analytically, however, there is no need for the optimizer to search over \(\lambda\). These transformations re-express the likelihood optimization; they do not add a Jacobian term to a maximum-likelihood objective.

Starting Values, Gradients, and Stopping Checks

Starting values should be plausible, not merely convenient. Moment estimates can provide an initial ARMA approximation, and a simpler nested model can supply additional candidate starts. For a particular data set, use more than one start when feasible: a small grid of plausible values, perturbations around a moment estimate, or a combination of these. Transform every start into the optimizer's unconstrained coordinates. If the mean is estimated, initialize it with a sensible estimate of the series level, and keep its treatment consistent across evaluations.

Gradient-based algorithms use derivatives of the objective to propose improvements. Derivatives can be computed analytically, by automatic differentiation, or by carefully checked finite differences. Scaling matters: a poorly scaled objective or parameter can lead to small numerical steps even when a meaningful improvement remains. After optimization, do not rely solely on a software message saying “converged.” Check that the objective is finite, the transformed coefficients map to admissible values, and the gradient is small at the reported solution. When using transformed parameters, also inspect the gradient with respect to the original coefficients where practical. Near a transformation boundary, the derivative of \(\tanh(\eta)\) is close to zero, so a small transformed-coordinate gradient alone may be misleading.

For several starts, compare the final objective values and parameter values. Similar best objective values with materially different parameters can indicate a flat region or competing solutions; a substantially worse endpoint suggests that a run found a poorer optimum. Restarting from the best candidate, tightening tolerances, or checking derivatives can help distinguish numerical trouble from a genuine feature of the likelihood. The goal is not to claim global optimality from a single run, but to gather evidence that the reported solution is stable and best among the candidates examined.

A Short Justification for Profiling the Variance

For fixed coefficients with \(S>0\), differentiate the log-likelihood with respect to \(\sigma^2\):

Setting this derivative to zero gives \(m\sigma^2=S\), hence \(\widehat{\sigma}^2=S/m\). The second derivative at that value is negative, so it is a maximum. If \(S=0\), the log-likelihood is \(-\frac{m}{2}\log(2\pi\sigma^2)\), which tends to \(+\infty\) as \(\sigma^2\downarrow0\); there is then no finite maximizer with positive variance. This is why the positive-sum-of-squares condition is essential when profiling.

Common Pitfalls and Practical Interpretation

Numerical optimization is part of statistical modeling, not a black box that makes the modeling decisions automatically. In particular, the likelihood must correspond to a declared initialization or exact-likelihood calculation, and all runs being compared must use the same definition.

• Bounding each ordinary coefficient separately. Having every coefficient between \(-1\) and \(1\) does not generally make a higher-order AR or MA polynomial stable. Use a parameterization designed to enforce the joint root constraints, and verify the resulting roots or equivalent stability conditions.
• Trusting a single start or a status message. One successful-looking run does not rule out a better candidate or a numerical failure. Compare starts and objective values, inspect gradients, and verify that the result is admissible.
• Using the wrong sign in the residual recursion. Under the stated lag-polynomial convention, AR(1) residuals are \(e_t=y_t-\phi_1y_{t-1}\). Substituting a negative coefficient requires careful sign handling; compute residuals directly from the model before evaluating the likelihood.
• Applying the variance formula when the residual sum of squares is zero. The formula \(S/m\) is the positive-variance MLE only for \(S>0\). At zero residual sum of squares the likelihood is unbounded as the variance approaches zero.