aML requires a piecewise-linear duration dependency (baseline hazard). In other words, the user specifies as many nodes (bend points) as desired, and aML estimates the slopes that best fit the data. Unlike Weibull, Gompertz, etc, this barely restricts the functional form of the duration dependency; it is semi-parametric and will adjust to any pattern in the data. Consider the figure below which shows the baseline log-hazard of divorce as a function of marriage duration, i.e., the marriage duration dependency.

Piecewise-Linear Log-Hazard of Divorce

A neat extension, unique to aML, is that it allows for multiple duration dependencies. For example, the hazard of divorce may depend not only on the marriage duration, but also on husband's age, wife's age, calendar time, and possibly more concepts of time. Each may be specified as a piecewise-linear (spline) pattern with its own nodes and slopes. Since all duration dependencies are piecewise-linear, their sum is again piecewise-linear.
Back to top.

Mathematically, there are several equivalent ways to spline a variable. In aML, the definition is such that the coefficients may be directly interpreted as slopes on particular segments of the variable (as opposed to the difference between the current segment's slope and the preceding one, which would result from a "marginal" spline transformation.)
Back to top.

The points are chosen using Gauss-Hermite Quadrature. Suppose you specify a normal distribution with k dimensions and n integration points in each dimension. aML first computes n standardized points which approximate a univariate standard normal distribution. [...] This method is extended to handle k-variate integration by selecting all n^k possibilities of k draws out of n points (with replacement). The appropriate weight for each set is the product of k individual (univariate) weights. These n^k sets of points and weights approximate a k-variate standard normal. In order to approximate other k-variate normals (with zero means), aML premultiplies each set of k draws by the Cholesky decomposition of the covariance matrix of the distribution. The marginal likelihood is computed as the properly weighted accumulation of n^k conditional likelihoods.

             points       weights
     1   -3.324257E+00  2.555784E-03
     2   -1.889176E+00  8.861575E-02
     3   -6.167066E-01  4.088285E-01
     4    6.167066E-01  4.088285E-01
     5    1.889176E+00  8.861575E-02
     6    3.324257E+00  2.555784E-03

The first figure below shows the support points of two univariate normal distributions, approximated by six points each. The size of the circles is indicative of the weights.

Now let's consider a bivariate distribution with covariance matrix . Its Cholesky decomposition is . Premultiplying each pair of univariate support points by this Cholesky matrix yields support points of the correlated bivariate distribution: . The figure below graphs these support points.

One perhaps undesirable implication now becomes clear. Even though the bivariate distribution is symmetric in the sense that residuals may be reversed without consequence, the support points are not symmetric. Residual u₁ is sampled with six unique values; residual u₂ with 36 unique values. This implies that reversal of residuals in model statements will slightly affect the results of estimation. Both are approximations and neither is better than the other. The differences will become smaller as the approximation improves, that is, as the number of support points increases.
Back to top.

If you receive this error message, you should first make sure that the unit of observation is correct. Do you really have huge numbers of outcomes per observation, or is there perhaps something wrong with the IDs that tie outcomes torgether into observations? Check the data summary information (.sum) file, produced by raw2aml, for number of observations and number of branches within observations. Also check the aML output (.out) file for number of observations and outcomes in your model. Is there really something unobserved that is common to all those outcomes in an observation? If not, re-define your ID variable to refer to a smaller unit.

If the data are correctly structured, try running a version of aML that allocates more memory to scratch arrays: bigaml (64MB) or hugeaml (256MB). Run "bigaml" or "hugeaml" in the same way as you would run "aml". These versions support roughly 1,200 to 2,400 outcomes per observation. These larger versions are located in the same directory as your current aml (or aml.exe) executable file.

Please be aware that computations involving very large covariance matrices can be very slow. If nothing appears to happen for a while, try a smaller sample (option obs=...;) just to make sure that the bottleneck is CPU resources.

If even hugeaml results in the above error message, consider ways to reduce the size of your observations. Can you perhaps ignore unobserved heterogeneity at the top level, or can it be replaced by something observed and captured through (dummy) variables? ("Option huber" generates Huber-corrected standard errors, also known as robust standard errors, which offer some protection against clustering of observations.) Or consider keeping a random subset of outcomes; if you cannot identify heterogeneity with, say, 1,000 outcomes per observation, it is probably not there, and you may re-define observations into smaller units.

Finally, beware of aML's default estimates of standard errors. With very many outcomes per observation, your data may end up having very few observations. Standard errors are the square roots of diagonal elements of the parameter estimates' covariance matrix, computed as the inverse of the matrix of second derivatives (Hessian matrix). By default, that Hessian matrix is approximated as minus the sum over observations of the outerproduct of first derivatives (BHHH method). The BHHH approximation is asymptotically correct, but may be poor in finite samples. "Option numerical standard errors" circumvents the BHHH approximation by computing the Hessian as numerical derivatives of (analytical) first derivatives. Use that option to obtain more precise estimates of standard errors. (Perhaps use it for final production runs only; it can be quite CPU-intensive.) Related, "option numerical search" always computes the Hessian numerically, not just to get accurate standard errors upon convergence.
Back to top.

Illustrative Effect on Log-Hazard

define spline Effect; nodes = 6 15; intercept;
hazard model; ...;
   model = ... + durspline(origin=0, ref=Effect) + ...;

One approach is to derive the distribution of the test statistic: . This requires that you know the variances and covariances of , , and . AML writes those out optionally, with the "option variance-covariance matrix". While not difficult, it involves a little algebra. An alternative approach applies a trick: redefine the spline such that the origin is at the point in time where you wish to test for significance, so that the intercept becomes the parameter of interest. In this case:

define spline NewEffect; nodes = -9 0; intercept;
hazard model; ...;
   model = ... + durspline(origin=-15, ref=NewEffect) + ...;

This procedure does not generalize to "kick-in" splines (with "effect=right"), but it does readily generalize to splines that have non- zero origins and to tests at points that do not correspond to nodes. For example, to test whether the effect is significant at time=20:

define spline NewEffect; nodes = -14 -5; intercept;
hazard model; ...;
   model = ... + durspline(origin=-20, ref=NewEffect) + ...;