Numerical Optimisation

What numerical optimisation does

It is used when an optimum cannot be obtained conveniently in closed form, including model calibration, likelihood estimation, treatment allocation and value-maximising design. This page explains gradients, curvature, line searches, derivative-free and global methods, stopping rules, constraints and the checks needed before an optimiser's result can be trusted.

An optimiser searches the mathematical problem supplied to it. Convergence of the algorithm does not prove that the objective represents the right health-economic question, that the optimum is global or that the input evidence is credible.

The objective defines the target of the search

An optimisation problem begins with decision or parameter vector (x) and objective function (f(x)). A minimisation formulation is conventional because a maximisation problem can be converted by changing the sign. Bounds or other constraints restrict the points the algorithm may consider.

An unconstrained problem is:

$$\min_{x\in\mathbb R^n} f(x)$$

A bound-constrained problem is:

$$\min_x f(x) \quad\text{subject to}\quad \ell_i\le x_i\le u_i$$

Common health-economic objectives include:

The objective's scale and units should be interpretable. Arbitrary weights can change the optimum even when every component looks reasonable in isolation.

A local optimum may not be the global optimum

A local minimum has no better point in a sufficiently small neighbourhood, while a global minimum has no better feasible point anywhere in the search domain. Convex problems allow strong global conclusions; non-convex problems may contain several local minima, saddle points and flat regions. One successful run can therefore identify only one basin of attraction.

Point (x^*) is a local minimum if there is a neighbourhood (\mathcal N) such that:

$$f(x^*)\le f(x) \quad\text{for all }x\in\mathcal N$$

It is a global minimum if:

$$f(x^*)\le f(x) \quad\text{for every feasible }x$$

Likelihoods for mixture models, calibrated transmission models and policy problems with thresholds can be strongly non-convex. Multiple starting points or global-search methods are needed when the decision cannot tolerate an unrecognised local solution.

Gradients identify directions of local change

The gradient collects the first derivatives of the objective with respect to every variable. It points in the direction of steepest local increase under the Euclidean norm, so its negative points toward steepest local decrease. Gradient-based methods use this information to propose the next iterate.

$$\nabla f(x)= \begin{bmatrix} \partial f/\partial x_1 \\ \vdots \\ \partial f/\partial x_n \end{bmatrix}$$

A basic gradient-descent update is:

$$x_{k+1}=x_k-\alpha_k\nabla f(x_k)$$

where (\alpha_k>0) is the step size. A step that is too small can make progress slow, while a step that is too large can overshoot or increase the objective.

A worked gradient-descent example

Consider the one-dimensional objective:

$$f(x)=(x-3)^2+2$$

Its derivative is:

$$f'(x)=2(x-3)$$

Starting at (x_0=0) with constant step size (\alpha=0.25):

$$x_1=0-0.25[-6]=1.5$$ $$x_2=1.5-0.25[-3]=2.25$$ $$x_3=2.25-0.25[-1.5]=2.625$$

The iterates approach the unique minimum (x^=3), where (f(x^)=2). The example is well-conditioned and convex; real health-economic objectives may be noisier, constrained or multi-modal.

Line searches choose a step along a direction

Many algorithms first choose a search direction (p_k) and then select a step length (\alpha_k). A line search seeks sufficient improvement without requiring exact minimisation along the line. Conditions such as Armijo or Wolfe balance decrease with useful progress.

The update is:

$$x_{k+1}=x_k+\alpha_k p_k$$

A descent direction satisfies:

$$\nabla f(x_k)^{\mathsf T}p_k<0$$

The Armijo sufficient-decrease condition is:

$$f(x_k+\alpha_kp_k) \le f(x_k)+c_1\alpha_k \nabla f(x_k)^{\mathsf T}p_k$$

for (0<c_1<1). The line-search parameters affect efficiency and robustness but should not change the substantive optimum in a well-solved problem.

Curvature can accelerate convergence

Second-order methods use curvature information from the Hessian matrix. Near a well-behaved minimum, Newton's method can converge rapidly because it rescales the gradient according to local curvature. The Hessian must be invertible and appropriately positive definite for the ordinary minimisation step to be a descent direction.

The Hessian is:

$$H(x)=\nabla^2f(x)$$

Newton's step solves:

$$H(x_k)p_k=-\nabla f(x_k)$$

and updates:

$$x_{k+1}=x_k+p_k$$

Explicitly calculating (H^{-1}) is generally unnecessary and less stable than solving the linear system. A singular or indefinite Hessian can require damping, regularisation or a trust-region method.

Quasi-Newton methods approximate curvature

Quasi-Newton methods build an approximation to the Hessian or its inverse from changes in gradients. They often achieve faster convergence than gradient descent without the full cost of exact second derivatives. BFGS and limited-memory BFGS are widely used examples.

Define:

$$s_k=x_{k+1}-x_k$$

and:

$$y_k=\nabla f(x_{k+1})-\nabla f(x_k)$$

The BFGS update maintains a curvature approximation when the curvature condition:

$$s_k^{\mathsf T}y_k>0$$

is satisfied. Limited-memory methods store only recent vector pairs, making them suitable for high-dimensional problems.

Trust regions restrict the local model

Trust-region methods approximate the objective near the current point and restrict the step to a region where that approximation is considered reliable. The region expands when predicted and actual improvement agree and contracts when they do not. This can be more robust than an unrestricted Newton step in difficult curvature.

A quadratic local model is:

$$m_k(p)= f(x_k)+ \nabla f(x_k)^{\mathsf T}p+ \frac{1}{2}p^{\mathsf T}B_kp$$

subject to:

$$\|p\|\le\Delta_k$$

The agreement ratio is:

$$\rho_k= \frac{f(x_k)-f(x_k+p_k)} {m_k(0)-m_k(p_k)}$$

A poor ratio signals that the local approximation was unreliable over the proposed step.

Derivatives can be analytic, automatic or numerical

Optimisation quality often depends on derivative accuracy. Analytic derivatives can be efficient but error prone to derive and code. Automatic differentiation applies chain rules through the implemented computation, while finite differences estimate derivatives from objective evaluations.

A central finite-difference approximation is:

$$\frac{\partial f}{\partial x_j} \approx \frac{f(x+he_j)-f(x-he_j)}{2h}$$

Choosing (h) creates a trade-off. A large step produces truncation error, while a very small step produces cancellation and floating-point error.

Derivative validation can compare analytic or automatic derivatives with finite differences at several interior points. Agreement should be assessed relative to scale rather than through exact decimal equality.

Derivative-free methods use objective evaluations

Derivative-free algorithms are useful when derivatives are unavailable, unreliable or the objective is produced by a simulation. Examples include Nelder-Mead, pattern search, Powell methods and model-based trust-region algorithms. They do not remove the need for a smooth or informative objective.

Derivative-free methods can require many evaluations and scale poorly with dimension. Stochastic simulation noise can make apparent improvements indistinguishable from random variation, so common random numbers, replication control or noise-aware algorithms may be needed.

Nelder-Mead maintains a simplex of (n+1) points in (n) dimensions and transforms it through reflection, expansion, contraction and shrinkage. It is useful in some low-dimensional problems but does not provide a universal convergence guarantee for arbitrary objectives.

Global methods explore multiple regions

Global optimisation methods attempt to reduce dependence on one initial point and one local basin. They include multi-start local search, branch-and-bound, simulated annealing, evolutionary algorithms and deterministic global methods. No general method can guarantee an efficient global solution for every nonlinear problem.

A practical multi-start approach is:

1. Define the feasible search domain using defensible parameter bounds.
2. Generate diverse starting points through a grid, random design or space-filling design.
3. Run the local optimiser from every starting point.
4. Cluster converged solutions by parameter values and objective value.
5. Investigate competing optima and retain uncertainty about non-identification.

Repeated convergence to the same point increases confidence but is not a mathematical proof of global optimality. The starting design must adequately cover the relevant domain.

Constraints require specialised handling

Bound, equality and inequality constraints restrict the search. Algorithms may use projected gradients, active sets, sequential quadratic programming, interior-point methods, augmented Lagrangians or penalties. The returned point must be evaluated against the original constraints, not only the transformed objective.

A general constrained problem is:

$$\min_x f(x)$$

subject to:

$$g_j(x)\le0, \qquad h_k(x)=0$$

At a regular solution, KKT residuals provide first-order optimality diagnostics. Integer or binary variables require discrete or mixed-integer methods rather than continuous optimisation followed by rounding.

Scaling affects steps and stopping rules

Variables can differ by many orders of magnitude, such as probabilities near 0.001 and costs near 100,000. Poor scaling can make gradient components, finite-difference steps and constraint tolerances incomparable. Reparameterisation can improve conditioning and convergence.

Common transformations include:

• Log-transform positive parameters: (z=\log x).
• Logit-transform probabilities: (z=\log[p/(1-p)]).
• Centre and scale unconstrained variables.
• Divide objective components by interpretable scales.
• Scale constraints so comparable violations have comparable numerical magnitude.

Results must be transformed back correctly, including uncertainty and derivatives. A transformed optimisation should represent the same substantive problem rather than silently changing prior or penalty assumptions.

Stopping criteria should reflect several diagnostics

An optimiser can stop because the gradient is small, the step is small, the objective changes little, constraints are satisfied or an iteration limit is reached. No single criterion is sufficient in every problem. A small objective change can occur on a flat plateau far from a meaningful optimum.

Typical diagnostics include:

$$\|\nabla f(x_k)\|\le\epsilon_g$$ $$\|x_{k+1}-x_k\|\le \epsilon_x(1+\|x_k\|)$$ $$|f(x_{k+1})-f(x_k)| \le \epsilon_f(1+|f(x_k)|)$$

Constraint violations should also fall below problem-specific tolerances. A stop caused by maximum iterations, failed line search or numerical error is not ordinary convergence.

Calibration is an optimisation problem with identification risks

Health-economic models are often calibrated by minimising the discrepancy between model outputs and external targets. Several parameter combinations may fit the targets similarly, creating flat ridges or multiple minima. A small objective value does not prove that parameters are uniquely identified.

A weighted least-squares calibration objective is:

$$f(\theta)= \sum_{k=1}^{K}w_k \left[m_k(\theta)-y_k\right]^2$$

where (m_k(\theta)) is model output, (y_k) is target and (w_k) is a prespecified weight. Weights may reflect target variance, importance or scale, but each interpretation creates a different objective.

Calibration validation should examine parameter profiles, Hessian or Jacobian conditioning, alternative starting points and predictive targets not used in fitting. Reporting only the best-fitting parameter vector can conceal non-identification.

Maximum likelihood uses optimisation for statistical estimation

Maximum-likelihood estimation chooses parameters that maximise the probability of the observed data under a statistical model. Software commonly minimises the negative log-likelihood because products of small probabilities are numerically unstable. The estimated Hessian can support standard-error calculations when regularity conditions hold.

$$\hat\theta= \arg\min_{\theta} \left[-\ell(\theta;y)\right]$$

The observed information matrix is:

$$\mathcal I(\hat\theta)= -\nabla^2\ell(\hat\theta;y)$$

An approximate covariance matrix is:

$$\widehat{\operatorname{Var}}(\hat\theta) =\mathcal I(\hat\theta)^{-1}$$

Boundary solutions, weak identification, non-regular likelihoods or poor conditioning can invalidate ordinary Hessian-based intervals.

Simulation objectives contain Monte Carlo noise

When the objective is estimated by simulation, two evaluations at the same point can differ. An optimiser may chase random noise rather than genuine improvement. Increasing replications, using common random numbers or modelling the noise can make comparisons more reliable.

If estimated objective is:

$$\hat f_R(x)= \frac{1}{R}\sum_{r=1}^{R}F(x,U_r)$$

then finite (R) creates Monte Carlo error. Using the same random streams (U_r) for nearby candidate points can reduce variance in objective differences.

The final candidate should be re-evaluated with a larger independent simulation budget. Otherwise, selection among noisy estimates can make the reported optimum appear better than its true expected performance.

Uncertainty in the optimum differs from numerical convergence

Numerical convergence concerns whether the algorithm solved the specified objective accurately. Statistical or decision uncertainty concerns how the optimum changes when evidence, parameters or model structure change. Tight solver tolerances do not reduce uncertainty in the underlying evidence.

Useful analyses include:

• Re-optimisation across probabilistic parameter draws.
• Profile objectives for parameters or decisions of interest.
• Bootstrap refitting of the complete estimation procedure.
• Scenario analysis for structural assumptions.
• Robust or stochastic optimisation for uncertain constraints and outcomes.
• Value-of-information analysis when uncertainty affects the preferred decision.

Re-optimising at every parameter draw estimates adaptive decisions under each draw and may not match a decision that must be chosen before uncertainty resolves. The timing of information and decision commitment should be explicit.

Common implementation and interpretation errors

Optimisation software can return a precise parameter vector even when the objective, derivatives or stopping logic are defective. The most important errors involve local minima, scale, constraints and noisy objectives. The following failures can materially change a health-economic result.

• Accepting the first converged run as the global optimum ignores other starting points and basins.
• Interpreting a solver success code without checking gradients and constraints confuses termination with validation.
• Using inaccurate derivatives can direct the search away from the optimum or create false convergence.
• Calculating an explicit Hessian inverse unnecessarily can worsen numerical stability.
• Ignoring variable scaling can make tolerances and finite-difference steps inappropriate.
• Optimising a noisy simulation with one replication per point encourages selection of random low values.
• Rounding continuous solutions to meet integer requirements can violate feasibility or optimality.
• Treating a good calibration fit as parameter identification ignores ridges and compensating parameters.
• Reporting many decimal places because the optimiser returns them overstates input and model precision.
• Using a very small numerical tolerance as evidence of decision certainty confuses computational and substantive uncertainty.

A practical validation sequence

Validation should establish that the objective is correct before assessing the algorithm. Independent derivative, feasibility and multi-start checks should then demonstrate that the reported solution is numerically credible. The following sequence provides an auditable minimum.

1. Define every variable, bound, constraint, unit and objective component.
2. Evaluate the objective at transparent test points and reconcile each component independently.
3. Check analytic or automatic derivatives against finite differences at several interior points.
4. Scale or transform variables and constraints so numerical tolerances are meaningful.
5. Run multiple starting points and compare parameter values, objective values and termination reasons.
6. Recalculate objective and constraints outside the optimiser at the returned solution.
7. Inspect gradient, step, objective-change, feasibility and KKT residuals as applicable.
8. Use alternative credible algorithms when the problem is non-convex, noisy or poorly conditioned.
9. Profile or perturb the solution to detect flat regions, ridges and weak identification.
10. Re-evaluate stochastic objectives with a larger independent simulation sample.
11. Propagate parameter and structural uncertainty through the complete optimisation workflow.
12. Confirm that remaining numerical differences cannot change the decision at the required reporting precision.

Numerical optimisation is reliable when the objective, derivatives, constraints and stopping evidence tell a consistent story. The optimiser's output is a candidate solution that must be verified mathematically and interpreted within the health-economic decision.