12
$\begingroup$

I am trying to implement the Nelder-Mead algorithm for optimising a function. The wikipedia page about Nelder-Mead is surprisingly clear about the entire algorithm, except for its stopping criterion. There it sadly says:

Check for convergence[clarification needed].

I tried and tested a couple of criteria myself:

  • Stop if $f(x_{N+1}) - f(x_1) < \epsilon$ where $\epsilon$ is small and where $x_i$ is the $i$-th vertex of the simplex, ordered from low ($f(x_1)$) to high ($f(x_{N+1})$) function values. In other words, when the maximum value of the simplex is almost equal to the minimum value. I found this didn't work properly, since this gives no guarantees about what the function does inside the simplex. Example, consider the function: $$f(x) = x^2$$This is of course trivial to optimise, but let's say we do this with NM, and let our two simplex points be $x_1 = -1$ and $x_2=1$. The algorithm would converge here without finding its optimum.

  • The second option involves evaluating the centroid of the simplex: stop if $|f(x_1) - f(x_c)| < \epsilon$. This assumes that if the lowest point of the simplex and the centroid have such similar values, the simplex is sufficiently small to call convergence.

Is this a proper way to check for convergence? Or is there an established way to check for this? I couldn't find any sources on this, since most search-hits focus on the complexity of the algorithm.

$\endgroup$
4
  • $\begingroup$ 1. It's not clear to me why you're comparing what happens at $x_{N+1}$ with $x_1$; surely you'd want to compare it with what happens at $x_N$. 2. convergence checks are a particularly tricky area in a lot of optimization; you need that the function is not changing much, but if the arguments are changing rapidly (even if the function is barely changing) you may not have converged, so people often use criteria that look at both. There's also the issue of whether you use a relative or an absolute criterion (neither is enough - e.g. a relative test when you're very close to 0 won't get triggered) $\endgroup$
    – Glen_b
    Oct 2 '17 at 23:57
  • 3
    $\begingroup$ The best stopping criterion for Nelder Mead is before you start. $\endgroup$ Oct 3 '17 at 0:25
  • $\begingroup$ Just to avoid confusion w.r.t notation in @Glen_b's comment ... I believe the subscripts here refer to the vertices of the simplex, not the iteration of the algorithm. So that the first convergence criterion proposed in this question, compares the lowest and highest function values of vertices in the $N$-dimensional parameter space ... it's not explicitly stated in the question, but the description of the algorithm on the linked wikipedia page (and in original paper) order the $N+1$ vertices from lowest function value to highest. $\endgroup$
    – Nate Pope
    Oct 3 '17 at 5:59
  • $\begingroup$ @NatePope That was my intention yes, I'll add clarification to the question.\ $\endgroup$
    – JAD
    Oct 3 '17 at 6:21
6
+50
$\begingroup$

The account of this "downhill simplex algorithm" in the original versions of Numerical Recipes is particularly lucid and helpful. I will therefore quote relevant parts of it. Here is the background:

In one-dimensional minimization, it was possible to bracket a minimum... . Alas! There is no analogous procedure in multidimensional space. ... The best we can do is give our algorithm a starting guess; that is, an $N$-vector of independent variables as the first point to try. The algorithm is then supposed to make its own way downhill through the unimaginable complexity of an $N$-dimensional topography until it encounters an (at least local) minimum.

The downhill simplex method must be started not just with a single point, but with $N+1$ points, defining an initial simplex. [You can take these points to be an initial starting point $P_0$ along with] $$P_i = P_0 + \lambda e_i\tag{10.4.1}$$ where the $e_i$'s are $N$ unit vectors and where $\lambda$ is a constant which is your guess of the problem's characteristic length scale. ...

Most steps just [move] the point of the simplex where the function is largest ("highest point") through the opposite face of the simplex to a lower point. ...

Now for the issue at hand, terminating the algorithm. Note the generality of the account: the authors provide intuitive and useful advice for terminating any multidimensional optimizer and then show specifically how it applies to this particular algorithm. The first paragraph answers the question before us:

Termination criteria can be delicate ... . We typically can identify one "cycle" or "step" of our multidimensional algorithm. It is then possible to terminate when the vector distance moved in that step is fractionally smaller in magnitude than some tolerance TOL. Alternatively, we could require that the decrease in the function value in the terminating step be fractionally smaller than some tolerance FTOL. ...

Note well that either of the above criteria might be fooled by a single anomalous step that, for one reason or another, failed to get anywhere. Therefore, it is frequently a good idea to restart a multidimensional minimization routine at a point where it claims to have found a minimum. For this restart, you should reinitialize any ancillary input quantities. In the downhill simplex method, for example, you should reinitialize $N$ of the $N+1$ vertices of the simplex again by equation $(10.4.1)$, with $P_0$ being one of the vertices of the claimed minimum.

Restarts should never be very expensive; your algorithm did, after all, converge to the restart point once, and now you are starting the algorithm already there.

[Pages 290-292.]

The code accompanying this text in Numerical Recipes clarifies the meaning of "fractionally smaller": the difference between values $x$ and $y$ (either values of the argument or values of the function) is "fractionally smaller" than a threshold $T\gt 0$ when

$$\frac{|x| - |y|}{f(x,y)} = 2\frac{|x|-|y|}{|x| + |y|} \lt T\tag{1}$$

with $f(x,y) = (|x|+|y|)/2$.

The left hand side of $(1)$ is sometimes known as the "relative absolute difference." In some fields it is expressed as a percent, where it is called the "relative percent error." See the Wikipedia article on Relative change and difference for more options and terminology.

Reference

William H. Press et al., Numerical Recipes: The Art of Scientific Computing. Cambridge University Press (1986). Visit http://numerical.recipes/ for the latest editions.

$\endgroup$
2
  • 1
    $\begingroup$ Thanks for the insight about restarting. I thought this was just running the algorithm from different starting points, but there actually seems to be more to that. $\endgroup$
    – JAD
    Oct 8 '17 at 18:37
  • $\begingroup$ I hadn't thought before about the possible meanings of "restarting." In the present context, I might have used a term like "polishing" for the "restart," but maybe that's not quite right either. The kind of "restart" advocated for the simplex method may be rather special to it. $\endgroup$
    – whuber
    Oct 8 '17 at 18:43
9
$\begingroup$

Not a complete answer, but too long for a comment and may put you on the right track.

This topic is briefly treated on page 171 of "Compact Numerical Methods for Computers" 2nd Ed., by John C. Nash. And happens to be the reference cited for the Nelder-Mead routine implemented in R's optim() function. Quoting the relevant part:

The thorniest question concerning minimisation algorithms must, therefore, be addressed: when has the minimum been found? Nelder and Mead suggest the 'standard error' of the function values: $$\mathrm{test} = \left[ \left( \sum_{i=1}^{n+1}[S(b_i) - \bar{S}]^2 \right) / n \right]^{1/2}$$ Where $$\bar{S} = \sum_{i=1}^{n+1} S(b_i)/(n+1).$$

I'll interrupt to clarify that $S(.)$ is the function being minimized, the $b$ are the $n+1$ points that define the $n$-dimensional simplex; the point with the highest function value is $b_H$ and the point with the lowest function value is $b_L$. Nash continues:

The procedure is taken to have converged when the test value falls below a preassigned tolerance. In the statistical applications which interested Nelder and Mead, this approach is reasonable. However, the author has found that this criterion can cause premature termination of the procedure on problems with fairly flat areas on the function surface. In a statistical context one might want to stop if such a such a region were encountered, but presuming the minimum is sought, it seems logical to use the simpler test for equality between $S(b_L)$ and $S(b_H)$, that is, a test for equal height of all points in the simplex.

A quick look at the source of optim() indicates that it uses the difference between the highest and lowest function values (of the points defining the simplex) to determine convergence: if (VH <= VL + convtol || VL <= abstol) break; Where VH is the high value and VL the low value. This comes with the caveat that I took a very quick look at the source, and am probably missing something.

Now, your option (1) appears to be the second approach advocated by Nash. He also discusses the problem you encountered:

Finally, it is still possible to converge at a point which is not the minimum. If, for instance, the $(n+1)$ points of the simplex are all in one plane (which is a line in two dimensions), the simplex can only move in $(n-1)$ directions in the $n$-dimensional space and may not be able to proceed towards the minimum. O'Neill (1971), in a FORTRAN implementation of the Nelder-Mead ideas, tests the function value at either side of the supposed minimum along each of the parameter axes. If any function value is found lower than the current supposed minimum, then the procedure is restarted.

The original references that Nash refers to here are:

Nelder JA, Mead R. 1965. A simplex method for function minimization. The Computer Journal 7: 308-313.

O'Neill R. 1971. Algorithm AS 47: Function minimization using a simplex procedure. Applied Statistics 20: 338-345.

$\endgroup$
3
$\begingroup$

A straw man: track only the lowest corner $$f_{\text{min}}(t) \equiv \text{min}_{\text{all corners}} \ f(X_i, t) $$ with the "patience" stopping rule:

# stop when you run out of patience, no improvement for say 10 iterations in a row:
if t > tbest + patience:
    message = "iter %d: f %g >= fbest %g" ...
    return message, fbest, xbest

Monitoring all $n+1$ corners is definitely useful in checking for reasonable coordinate scaling, constraints, N-M initial simplex. Whether tracking all corners can improve the combination of

  1. the problem: rough terrain, perhaps with bad scaling or silly constraints
  2. the algorithm, balancing exploring and moving (and software complexity)
  3. the stopping rule proper

remains to be seen — real test cases welcome.

(A real Stopiter class has many stop conditions, in addition to patience; simplest is wall clock time.)

See also:
NLopt: many algorithms including Nelder-Mead, easy to compare. See especially the notes on Comparing algorithms
Downhill: patience stopping with min_improvement

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.