Nested-Set Hessian

• Nested-Set Hessian is a derivative-free method that approximates the Hessian of a scalar function using nested difference schemes.
• It efficiently reuses function evaluations to attain minimal cost and provides explicit first-order (O(Δ)) error bounds.
• The method extends to overdetermined configurations and calculus-based variants, enhancing stability and applicability in optimization.

The nested-set Hessian is a derivative-free, first-order accurate method for approximating the Hessian matrix $\nabla^2 f(x^0)$ of a twice-differentiable scalar function $f: \mathbb{R}^n \rightarrow \mathbb{R}$, specifically designed for use in settings where explicit derivatives are unavailable. Rooted in the theory of generalized simplex gradients, the approach constructs the Hessian by combining two nested difference schemes that efficiently reutilize function evaluations. For carefully chosen sample structures, the nested-set Hessian achieves a minimal evaluation cost of $(n+1)(n+2)/2$ while providing explicit error bounds of order $O(\Delta)$ in the sampling radii. Its framework further permits calculus-based composition rules and extension to overdetermined sets, supporting flexible application within model-based derivative-free optimization routines (Hare et al., 2020).

1. Mathematical Foundations and Formal Construction

Given a base point $x^0 \in \mathbb{R}^n$, the method relies on two sets of directions:

• Inner set $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$ (full row rank)
• Outer set $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$ (full row rank)

The generalized simplex gradient at $x^0$ with respect to $T$ is defined as

$$\nabla_s f(x^0; T) = (T^\top)^\dagger \delta_f(x^0; T),$$

where

$f: \mathbb{R}^n \rightarrow \mathbb{R}$0

and $f: \mathbb{R}^n \rightarrow \mathbb{R}$1 denotes the Moore–Penrose pseudoinverse. For $f: \mathbb{R}^n \rightarrow \mathbb{R}$2 and nonsingular $f: \mathbb{R}^n \rightarrow \mathbb{R}$3, this recovers the classical simplex gradient.

The nested-set Hessian introduces a further set $f: \mathbb{R}^n \rightarrow \mathbb{R}$4 for a second layer of differencing. For each $f: \mathbb{R}^n \rightarrow \mathbb{R}$5, compute both $f: \mathbb{R}^n \rightarrow \mathbb{R}$6 and $f: \mathbb{R}^n \rightarrow \mathbb{R}$7. Construct the matrix

$f: \mathbb{R}^n \rightarrow \mathbb{R}$8

and set

$f: \mathbb{R}^n \rightarrow \mathbb{R}$9

This estimator approximates directional second derivatives via differences of approximate gradients. The error bound is given by

$(n+1)(n+2)/2$0

where $(n+1)(n+2)/2$1, and $(n+1)(n+2)/2$2 is explicit in the pseudoinverse norms and the Lipschitz constant of the true Hessian (Hare et al., 2020).

2. Evaluation Complexity and Minimal Point Sets

The generic requirement for function evaluations is at most $(n+1)(n+2)/2$3, comprising $(n+1)(n+2)/2$4, $(n+1)(n+2)/2$5, $(n+1)(n+2)/2$6, and $(n+1)(n+2)/2$7 (for all $(n+1)(n+2)/2$8). However, if $(n+1)(n+2)/2$9 and both $O(\Delta)$0 and $O(\Delta)$1 have favourable structure, numerous points coincide and only

$O(\Delta)$2

distinct function evaluations are necessary. A minimal poised set for the nested-set Hessian consists of $O(\Delta)$3 with $O(\Delta)$4, such that the collection

$O(\Delta)$5

has cardinality exactly $O(\Delta)$6. One canonical construction uses $O(\Delta)$7 and $O(\Delta)$8 (specific difference matrices), and all full-row-rank $O(\Delta)$9 can be generated via transformations from this canonical case (Hare et al., 2020).

The method readily generalizes to $x^0 \in \mathbb{R}^n$0 or $x^0 \in \mathbb{R}^n$1 (overdetermined geometry), enhancing stability via reduced pseudoinverse norms, though at the cost of additional function values.

Configuration	Function Evaluations	Accuracy
General ($x^0 \in \mathbb{R}^n$2)	$x^0 \in \mathbb{R}^n$3	$x^0 \in \mathbb{R}^n$4
Minimal $x^0 \in \mathbb{R}^n$5	$x^0 \in \mathbb{R}^n$6	$x^0 \in \mathbb{R}^n$7
Quadratic Interpolation	$x^0 \in \mathbb{R}^n$8	$x^0 \in \mathbb{R}^n$9 gradient / exact Hessian for quadratics

3. Error Bounds and Stability Considerations

Assuming $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$0 near $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$1 with Lipschitz constant $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$2, the error in the nested-set Hessian is bounded by

$$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$3

where $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$4. For $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$5,

$$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$6

Stability is governed by the conditioning of $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$7 and $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$8 (through the pseudoinverse norms). For highly nonlinear $$T = [t^1\, \ldots\, t^k] \in \mathbb{R}^{n \times k}$$9, $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$0 must be small to guarantee control of $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$1 error (Hare et al., 2020).

4. Algorithmic Procedure

The computation proceeds through the following steps:

Input: x0 ∈ ℝⁿ, S=[s¹…sᵐ] ∈ ℝⁿ×ᵐ, T=[t¹…tᵏ] ∈ ℝⁿ×ᵏ

1. Evaluate f0 = f(x0)

2. For j = 1…k, evaluate fT[j] = f(x0 + t^j)

3. Compute δ0 = fT – f0 * 1ₖ, then ∇s0 = (Tᵀ)† δ0

4. For i = 1…m:
     Evaluate fS[i] = f(x0 + s^i)
     For j = 1…k, evaluate fST[i,j] = f(x0 + s^i + t^j)
     Compute δi[j] = fST[i,j] – fS[i]
     Compute ∇s_i = (Tᵀ)† δi
     Set row i of Δ_∇ = (∇s_i – ∇s0)ᵀ

5. Compute H = (Sᵀ)† Δ_∇
Output: Approximate Hessian H ≈ ∇²f(x⁰)

This nested structure enables efficient reuse of function values, particularly in minimal sets (Hare et al., 2020).

5. Calculus-Based Approximation Variants

Two calculus-based Hessian approximation schemes are developed to address composite expressions:

• Simplex calculus Hessian replaces gradient and Hessian terms in standard product, quotient, and power rules with the corresponding nested-set estimators, e.g.,

$$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$2

This variant exhibits $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$3 error and exactness for affine components.

• Quadratic calculus Hessian constructs quadratic interpolation models $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$4 for $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$5 and uses these for gradient terms, while retaining nested-set Hessians in second-derivative terms. For composite $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$6:

$$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$7

This formula yields $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$8 accuracy and is exact whenever all factors are quadratic (Hare et al., 2020).

6. Comparative Assessment and Applicability

Compared to alternative approaches, the nested-set Hessian matches the minimal evaluation cost of quadratic interpolation while simplifying implementation by relying on nested simplex-gradient calculations rather than solving a global system. In contrast, finite-difference Hessians using central difference of gradients require $$S = [s^1\, \ldots\, s^m] \in \mathbb{R}^{n \times m}$$9 gradient evaluations—translating to $x^0$0 function values due to the cost of gradient estimation.

The nested-set Hessian offers:

• Explicit $x^0$1 error for the Hessian estimator,
• Minimal $x^0$2 cost for well-structured sets,
• Flexible extension to overdetermined (stabilized) configurations,
• Compatibility with calculus-based rules supporting combinations of estimated Hessians.

Its limitations include only first-order ($x^0$3) accuracy, sensitivity to ill-conditioning in $x^0$4 and $x^0$5, and the need for sufficiently small $x^0$6 to control error for highly nonlinear functions. Quadratic interpolation, by contrast, achieves higher gradient accuracy ($x^0$7) and exactness for quadratics but is more complex to implement (Hare et al., 2020).

Overall, the nested-set Hessian supplies an efficient, robust, and implementable framework for derivative-free second-order information that is readily integrated into model-based optimization routines and advanced calculus operations on estimated models.