∎

11institutetext: Shigeng Sun 22institutetext: Department of Engineering
Sciences and Applied Mathematics, Northwestern University, Evanston, IL, USA
22email<EMAIL_ADDRESS>33institutetext: Corresponding
author: Jorge Nocedal 44institutetext: Department of Industrial Engineering
and Management Sciences, Northwestern University, Evanston, IL, USA
44email<EMAIL_ADDRESS>
# A Trust Region Method for the Optimization of Noisy Functions††thanks: Sun
was supported by NSF grant DMS-1620022. Nocedal was supported by AFOSR grant
FA95502110084 and ONR grant N00014-21-1-2675.

Shigeng Sun  Jorge Nocedal

(Received: 02 January 2022 / Accepted: )

###### Abstract

Classical trust region methods were designed to solve problems in which
function and gradient information are exact. This paper considers the case
when there are bounded errors (or noise) in the above computations and
proposes a simple modification of the trust region method to cope with these
errors. The new algorithm only requires information about the size of the
errors in the function evaluations and incurs no additional computational
expense. It is shown that, when applied to a smooth (but not necessarily
convex) objective function, the iterates of the algorithm visit a neighborhood
of stationarity infinitely often, and that the rest of the sequence cannot
stray too far away, as measured by function values. Numerical results
illustrate how the classical trust region algorithm may fail in the presence
of noise, and how the proposed algorithm ensures steady progress towards
stationarity in these cases.

###### Keywords:

Trust Region Method Nonlinear Optimization Noisy Optimization

###### MSC:

65K05 68Q25 65G99 90C30

##  1 Introduction

Trust region methods are powerful techniques for nonlinear optimization that
have the ability to incorporate second-order information, without requiring it
to be positive definite. They are endowed with strong global convergence
properties and have proven to be effective in practice. Although the design
and analysis of trust region methods are well established in the absence of
noise (or errors), this is not the case when noise is present.

In this paper, we show how to redesign the classical trust region method for
unconstrained optimization to handle problems where the objective function,
gradient, and (possibly) Hessian, are subject to bounded, non-diminishing
noise. This involves only one modification in the algorithm: the ratio of
actual/predicted reduction used for step acceptance is now relaxed by a term
proportional to the noise level. All other aspects of the classical trust
region method remain unchanged. We show that, under mild conditions, the
proposed algorithm converges to a neighborhood of stationary points, where the
size of the neighborhood is determined by the level of noise. This analysis is
more complex than for line search methods due to the effects of memory
encapsulated in the trust region update. Our convergence results do not assume
convexity of the objective function but only that it is sufficiently smooth.

Examples of practical optimization applications with bounded noise include
those that employ mixed-precision arithmetic; problems where derivatives are
approximated by finite differences; and problems in which the evaluation of
the objective function (and gradient) contain computational noise.

This investigation was motivated by numerical experiments performed by the
authors that indicated that, although the classical trust region approach
often tolerates significant levels of noise, it can fail in certain
situations. This raises the question of how to best modify the method to avoid
failures. The algorithm proposed here is inspired by work on line search
methods for unconstrained optimization berahas2019derivative ; xie2020analysis
and equality constrained optimization figen2021ConDFO . In those papers,
convergence-to-neighborhood results were derived but the analysis presented
here follows different lines, as trust region methods require different proof
techniques.

The paper is organized into 5 sections. In the rest of this section, we
provide a review of the relevant literature. In section 2, we describe the
problem setting and the proposed trust region algorithm. The main convergence
results are presented in section 3. Numerical experiments, summarized in
section 4, indicate that the new algorithm is more robust than the classical
method. Section 5 presents the final remarks on the contributions of this
work.

###  1.1 Literature Review

The study of nonlinear optimization problems with errors or noise in the
function and gradient has attracted attention in recent years, motivated by
the use of finite difference approximations to derivatives nesterov2017random
; shi2021numerical ; shi2021automatic  and by applications in machine
learning; see curtis2020adastoc  for a review of some recent work.

One of the earliest investigations of trust region methods with errors is
carter , which proved global convergence assuming that the errors in the
gradient diminish at a rate that is proportional to the norm of the true
gradient; this condition is referred to as the _norm test_ in
bollapragada2018adaptive ; exact2018ima . The importance of the norm test was
promoted in byrd2012sample , which established linear convergence and
complexity bounds for an adaptive sampling method for empirical risk
minimization, as well as in cartis2018global ; paquette , which establishes
convergence in probability for a stochastic line search method.

Prior studies of optimization methods for minimization of functions with non-
diminishing, bounded errors include berahas2019derivative , which employed a
relaxed Armijo back-tracking line search and established linear convergence to
a neighborhood of the solution for strongly convex functions. Stopping time
guarantees for the same relaxed line search is proven in berahas2019global . A
similar relaxed Armijo back-tracking line search technique is used in
katya2021stolinesearch , which considered different oracles from paquette  to
allow biased estimates, and provided complexity bounds for different noise
structures under probabilistic frameworks. Quasi-Newton methods were analyzed
in xie2020analysis , which described a noise tolerant modification of the BFGS
method; shi2020noise  showed ways to make this method robust and efficient in
practice.

For constrained optimization, berahas2021stochastic ; berahas2021sequential ;
curtis2021inexact  studied a sequential quadratic programming (SQP) method for
equality constrained optimization in the case when the objective function is
stochastic and the constraints are deterministic. Those three papers give
conditions under which convergence can be expected, giving careful attention
to the behavior of the penalty parameter. Using a relaxed Armijo line search
procedure, figen2021ConDFO  shows global convergence to a neighborhood of the
solution for an SQP method for equality constrained problems.

Analysis for trust region methods with more general (unbounded) noise is
presented in katya2018storm , which establishes almost sure global convergence
under the assumption that function and gradient information is sufficiently
accurate with high enough probability. cartis2019TRsupermartingale  views the
optimization as a generic stochastic process, and improves upon the results of
katya2018storm . The analysis presented in cartis2019TRsupermartingale
establishes convergence results for a trust region method and, under the
assumption of sufficiently accurate stochastic gradient information, derives a
stopping time result and a second order global complexity bound. A method
inspired by trust region techniques is curtis2019stochTR , which uses step
normalization techniques in the stochastic optimization setting, and
establishes conditions for linear and sublinear convergence. A series of
papers, including cartis2021TRqDA ; toint2021TRqEDAN ; toint2021ARqpEDA2 ,
analyze regularization and trust region methods with adaptive accuracy in the
function and gradient evaluations, and establish worst case complexity bounds.

The style of analysis presented in katya2018storm ;
cartis2019TRsupermartingale ; curtis2019stochTR , which is used to prove
convergence in probability, stands in contrast with the deterministic
technique employed in this paper, which assumes bounded errors. It remains to
be seen which approach is more useful for the design of noise tolerant
optimization methods—or whether the two approaches complement each other.

##  2  Problem Statement and Algorithm

Our goal is to design a trust region method to solve the unconstrained
minimization problem

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

in the case when the function $f(x)$ and gradient $g(x)=\nabla f(x)$ cannot be
evaluated exactly. Instead, we have access to noisy observations of the above
quantities, which we denote as $\tilde{f}(x)$, and $\tilde{g}(x)$. We write

$\tilde{f}(x)=f(x)+\delta_{f}(x),\quad\text{and}\quad\tilde{g}(x)=g(x)+\delta_{g}(x),$
(2)

where the error functions (or noise) $\delta_{f}(x)$, $\delta_{g}(x)$ are
assumed to be bounded, i.e.,

$|\delta_{f}(x)|\leq\epsilon_{f},\qquad\|\delta_{g}(x)\|\leq\epsilon_{g},\qquad\
\forall x\in\mathbb{R}^{n}.$ (3)

Throughout the paper $\|\cdot\|$ stands for the Euclidean norm.

Let us apply a classical trust region method to problem 1. At each iterate,
the method constructs a quadratic model

$m_{k}(p)=\tilde{f}(x_{k})+\tilde{g}(x_{k})^{T}p+\frac{1}{2}p^{T}\tilde{B}_{k}p,$
(4)

and solves the following trust region subproblem for the step $p_{k}$:

$\min_{p\in\mathbb{R}^{n}}m_{k}(p)\quad\text{ s.t. }\|p\|\leq\Delta_{k}.$ (5)

In 4, $\tilde{B}_{k}$ could be defined as a noisy evaluation of the Hessian, a
quasi-Newton matrix, or some other approximation. To decide if the step
$p_{k}$ should be accepted—and if the trust region radius $\Delta_{k}$ should
be modified— classical trust region methods employ the ratio of actual to
predicted reduction in the objective function, defined as

$\frac{\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)}{m_{k}(0)-m_{k}\left(p_{k}\right)}.$
(6)

This ratio is, however, not adequate in the presence of noise because if
$\Delta_{k}$ becomes very small, the numerator can be of order $\epsilon_{f}$,
while the denominator will be proportional to $\Delta_{k}$. Thus, if
$\Delta_{k}<<\epsilon_{f}$, the ratio 6 may exhibit wild oscillations that can
cause the algorithm to perform erratically; see the examples in Section 4.

To address this issue, we propose the following noise tolerant variant of 6:

$\rho_{k}=\frac{\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)+r\epsilon_{f}}{m_{k}(0)-m_{k}\left(p_{k}\right)+r\epsilon_{f}},$
(7)

where $r>2$ is a constant specified below. The reason for relaxing both the
numerator and denominator in 7 is to be consistent with the classical
narrative of trust region methods where a ratio close to 1 is an indication
that the model is adequate. An alternative approach would be to relax only the
numerator and interpret the condition $\rho_{k}>c$ (where $c>0$ is a constant)
as a relaxed Armijo condition of the type studied in berahas2019derivative ;
figen2021ConDFO . We find the first interpretation to be easier to motivate
and to yield tighter bounds in the convergence analysis. We state the
algorithm as follows.

1 Initialize $\Delta_{0}$, and chose constants $0<c_{0}\leq c_{1}<c_{2}<1$ and
$\nu>1$

2while  _a termination condition is not met_ do

3 Compute $p_{k}$ by solving 5 (exactly or approximately);

4 Evaluate $\rho_{k}$ as in 7;

5

6 if  _$\rho_{k} <c_{1}$_ then

7 $\Delta_{k+1}=\frac{1}{\nu}\Delta_{k}$;

8 else if  _$\rho_{k} >c_{2}$ _ then

9 $\Delta_{k+1}=\nu\Delta_{k}$;

10 else

11 $\Delta_{k+1}=\Delta_{k}$;

12 end if

13 if  _$\rho_{k} >c_{0}$_ then

14 $x_{k+1}=x_{k}+p_{k}$;

15 else

16 $x_{k+1}=x_{k}$;

17 end if

18 Set $k\leftarrow k+1$;

19 end while

Algorithm 1 Noisy Trust-Region Algorithm

Typical values of the parameters are $c_{0}=0.1$, $c_{1}=\frac{1}{4}$,
$c_{2}=\frac{1}{2}$, $\nu=2$, but other values can be used in practice. The
global convergence result presented in the next section holds if the constant
$r$ in 7 is chosen as

$r=2/{(1-c_{2})}.$ (8)

We assume that the step $p_{k}$ computed in step 3 yields a decrease in the
model $m_{k}$ that is at least as large as that given by the Cauchy step
(defined below). This provides much freedom in the design of the algorithm,
and includes the dogleg and Newton-CG methods, as well as the exact solution
of the trust region problem; see, e.g., mybook .

In practice it can be useful to increase the trust region radius in Step 7
only if $\rho_{k}>c_{2}$ and $\|p_{k}\|=\Delta_{k}$, as this can prevent
unnecessary oscillations in the trust region radius. The convergence result
presented in the next section can easily be extended to that case, assuming
certain technical conditions on the step computation—which are satisfied by
the dogleg and Newton-CG methods.

##  3 Global Convergence Analysis

In this section, we establish a global convergence result for Algorithm 1 that
applies to general objective functions. The proof is based on the observation
that, when the gradient is large enough, the trust region radius will
eventually become large too, ensuring sufficient descent in the objective
function despite the presence of noise. This drives the iteration toward
regions where the stationarity measure is small (i.e., comparable to the noise
level).

We begin by establishing a standard requirement on the step computation based
on the _Cauchy step_ $p_{k}^{c}$ for problem 1, which is defined as

$p_{k}^{c}=-\tau_{k}\frac{\Delta_{k}}{\left\|\tilde{g}_{k}\right\|}\tilde{g}_{k},$
(9)

where

$\tau_{k}=\left\\{\begin{array}[]{ll}1&\text{ if
}\tilde{g}_{k}^{T}\tilde{B}_{k}\tilde{g}_{k}\leq 0\\\
\min\left(\left\|\tilde{g}_{k}\right\|^{3}\big{/}\left(\Delta_{k}\tilde{g}_{k}^{T}\tilde{B}_{k}\tilde{g}_{k}\right),1\right)&\text{
otherwise.}\end{array}\right.$ (10)

As is well known (see e.g. (mybook, , Lemma 4.3)), the reduction in the model
provided by the Cauchy step satisfies

$m_{k}(0)-m_{k}(p_{k}^{c})\geq\frac{1}{2}\left\|\tilde{g}_{k}\right\|\min\left(\Delta_{k},\frac{\left\|\tilde{g}_{k}\right\|}{\left\|\tilde{B}_{k}\right\|}\right).$
(11)

We assume that the step $p_{k}$ computed by Algorithm 1 yields a reduction in
the model that is not less than that produced by the Cauchy step, i.e.,

$m_{k}(0)-m_{k}(p_{k})\geq
m_{k}(0)-m_{k}(p_{k}^{c})\geq\frac{1}{2}\left\|\tilde{g}_{k}\right\|\min\left(\Delta_{k},\frac{\left\|\tilde{g}_{k}\right\|}{\left\|\tilde{B}_{k}\right\|}\right).$
(12)

We can now state the assumptions on the problem and the algorithm under which
the global convergence results are established.

######  Assumption 3.1.

The objective function $f$ is Lipschitz continuously differentiable with
constant $L$, i.e.,

$\|g(x)-g(y)\|<L\|x-y\|.$ (13)

######  Assumption 3.2.

The error in the function and gradient evaluations is bounded, i.e., 3 holds
for some constants $\epsilon_{f},\epsilon_{g}.$

We impose no other conditions on the errors, other than boundedness. Next, we
impose a minimal requirement on the Hessian approximations.

######  Assumption 3.3.

There is a constant $L_{B}>0$ such that the matrices $\tilde{B}_{k}$ satisfy

$\|\tilde{B}_{k}\|<L_{B},\ \forall k.$ (14)

There is freedom in the computation of the step $p_{k}$, but it must yield
Cauchy decrease.

######  Assumption 3.4.

The step $p_{k}$ computed by Algorithm 1 satisfies (12).

This assumption can be relaxed so as to require only a fraction of Cauchy
decrease, but we do not do so here to avoid the introduction of more
constants. The final requirement is standard.

######  Assumption 3.5.

The sequence $\\{\tilde{f}_{k}\\}$ generated by Algorithm 1 is bounded below.

We now proceed with the analysis.

###  3.1 Properties of the ratio $\rho_{k}$

We begin by establishing a bound between $\rho_{k}$ and 1. From 7, we have

$\displaystyle\left|\rho_{k}-1\right|$
$\displaystyle=\left|\frac{m_{k}\left(p_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)}{m_{k}(0)-m_{k}\left(p_{k}\right)+r\epsilon_{f}}\right|.$
(15)

From Taylor’s Theorem we have

$\displaystyle\tilde{f}(x_{k}+p_{k})$
$\displaystyle=f(x_{k}+p_{k})+\delta_{f}(x_{k}+p_{k})$
$\displaystyle=f\left(x_{k}\right)+g_{k}^{T}p_{k}+\int_{0}^{1}\left[g\left(x_{k}+tp_{k}\right)-g_{k}\right]^{T}p_{k}dt+\delta_{f}(x_{k}+p_{k}).$

With this, by 13, 14, and 3, we obtain

$\displaystyle\left|m_{k}(p_{k})-\tilde{f}(x_{k}+p_{k})\right|\leq$
$\displaystyle\
\tfrac{1}{2}(L_{B}+L)\|{p_{k}}\|^{2}+\epsilon_{g}\|{p_{k}}\|+2\epsilon_{f}$
(16) $\displaystyle\equiv$ $\displaystyle\
M\|{p_{k}}\|^{2}+\epsilon_{g}\|{p_{k}}\|+2\epsilon_{f},$

where

$M=\tfrac{1}{2}(L_{B}+L).$ (17)

By substituting 16 and 12 into 15, we establish the following result.

######  Lemma 1.

If $\rho_{k}$ is defined by 7, then for all $k$,

$\left|\rho_{k}-1\right|\leq\frac{M\Delta_{k}^{2}+\epsilon_{g}\Delta_{k}+2\epsilon_{f}}{\frac{1}{2}\|\tilde{g}_{k}\|\min(\Delta_{k},\|\tilde{g}_{k}\|/\|\tilde{B}_{k}\|)+r\epsilon_{f}}.$
(18)

This lemma suggests that $\rho_{k}$ can be made close to 1 by decreasing
$\Delta_{k}$, up until the noise term $\epsilon_{f}$ dominates. This assertion
will be made more precise below.

###  3.2 Lower Bound on Trust Region Radius

We now show that if $\Delta_{k}$ is very small and the gradient is large
compared to the noise $\epsilon_{g}$, Algorithm 1 will increase the trust
region radius. We recall that $r$ is defined in 8 and that $\nu>1$.

######  Lemma 2 (Increase of Trust Region Radius).

Suppose that, at iteration $k$,

$\|\tilde{g}_{k}\|>r\epsilon_{g}+{\gamma},$ (19)

for some constant ${\gamma}>0$. Then, if

$\Delta_{k}\leq\bar{\Delta}=:\frac{{\gamma}}{rM},$ (20)

we have that

$\Delta_{k+1}=\nu\Delta_{k}.$ (21)

###### Proof.

Since $r>2$, we have from 14, 17 and 19 that

$rM>2M>\|\tilde{B}_{k}\|\quad\text{and }\quad{\gamma}<\|\tilde{g}_{k}\|,$ (22)

and thus

$\bar{\Delta}<\|\tilde{g}_{k}\|/\|\tilde{B}_{k}\|.$ (23)

Thus, if $\Delta_{k}\leq\bar{\Delta}$, we have

$\min(\Delta_{k},\|\tilde{g}_{k}\|/\|\tilde{B}_{k}\|)=\Delta_{k}.$ (24)

In addition, if $\Delta_{k}\leq\bar{\Delta}$, we also have

$M\Delta_{k}+\epsilon_{g}\leq
M\bar{\Delta}+\epsilon_{g}=\frac{{\gamma}}{r}+\epsilon_{g}=\frac{1}{r}(r\epsilon_{g}+{\gamma}).$
(25)

Substituting 24, 19, 25 and 8 into 18, we have that for all
$\Delta_{k}\leq\bar{\Delta}$

$\displaystyle\left|\rho_{k}-1\right|$
$\displaystyle\leq\frac{M\Delta_{k}^{2}+\epsilon_{g}\Delta_{k}+2\epsilon_{f}}{\frac{1}{2}\|\tilde{g}_{k}\|\Delta_{k}+r\epsilon_{f}}$
$\displaystyle<\frac{M\Delta^{2}_{k}+\epsilon_{g}\Delta_{k}+2\epsilon_{f}}{\frac{1}{2}(r\epsilon_{g}+{\gamma})\Delta_{k}+r\epsilon_{f}}$
$\displaystyle<\frac{\frac{1}{r}(r\epsilon_{g}+{\gamma})\Delta_{k}+2\epsilon_{f}}{\frac{1}{2}(r\epsilon_{g}+{\gamma})\Delta_{k}+r\epsilon_{f}}$
$\displaystyle=\frac{2}{r}$ $\displaystyle=1-c_{2}.$ (26)

This implies that $\rho_{k}>c_{2}$, and by step 8 of Algorithm 1 we have that
$\Delta_{k+1}=\nu\Delta_{k}.~{}~{}~{}$ ∎

A consequence of this lemma is that there is a lower bound for the trust
region radius if the norm of the noisy gradient remains greater than
$r\epsilon_{g}$.

######  Corollary 1 (Lower Bound on Trust Region Radius).

Given ${\gamma}>0$, if there exist $K>0$ such that for all $k\geq K$

$\|\tilde{g}_{k}\|>r\epsilon_{g}+{\gamma},$ (27)

then there exist $K_{0}\geq K$ such that for all $k\geq K_{0}$,

$\Delta_{k}>\tfrac{1}{\nu}\bar{\Delta}=\frac{{\gamma}}{\nu rM}.$ (28)

###### Proof.

We apply Lemma 2 for each iterate after $K$ to deduce that, whenever
$\Delta_{k}\leq\,\bar{\Delta}$, the trust region radius will be increased.
Thus, there is an index $K_{0}$ for which $\Delta_{k}$ becomes greater than
$\bar{\Delta}$. On subsequent iterates, the trust region radius can never be
reduced below $\bar{\Delta}/\nu$ (by Step 6 of Algorithm 1) establishing the
bound 28. ∎

Remark. In traditional trust region analysis for deterministic (noiseless)
optimization, one shows that the trust region radius will not shrink below a
certain value that depends on the Lipschitz constant and the norm of the
current gradient. However, that analysis does not imply that the trust region
will increase beyond a certain threshold, which is required in the presence of
noise. We need to show that the trust region eventually becomes large enough
with respect to the noise level so that progress can be made. This
differentiates our analysis from classical trust region convergence theory.

###  3.3 Reduction of Noisy Function

The classical trust region algorithm is monotonic, as it requires a reduction
in the objective function when accepting a step. Due to the relaxation in 7,
Algorithm 1 can accept steps that increase the noisy function. However, when
the iterates are far from the solution, this is not the case. We now show that
when the noisy gradient and trust region radius are both large enough, the
reduction in the objective is large enough to overcome any increase allowed by
7.

######  Lemma 3 (Noisy Function Reduction).

Suppose that for some $k>0$

$\|\tilde{g}_{k}\|>r\epsilon_{g}+\gamma\quad\mbox{and}\quad\Delta_{k}\,{\geq}\,\frac{\bar{\Delta}}{\nu}=\frac{\gamma}{\nu
rM},$ (29)

where

$\gamma=\eta+\mu,$ (30)

with $\mu>0$ an arbitrarily small constant, and

$\eta=\frac{1}{2}\left(-r\epsilon_{g}+\beta\right),\quad\beta=\sqrt{(r\epsilon_{g})^{2}+8\nu
r^{2}\left(\frac{1}{c_{0}}-1\right)M\epsilon_{f}}.$ (31)

Then, if the step is accepted at iteration $k$ by Algorithm 1, we have

$\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)>\frac{c_{0}}{2\nu
rM}\left(\mu\beta+\mu^{2}\right).$ (32)

###### Proof.

As argued in 23,
$\bar{\Delta}=\frac{{\gamma}}{rM}<\frac{\left\|\tilde{g}_{k}\right\|}{\left\|\tilde{B}_{k}\right\|}$,
and therefore

$\min\left(\Delta_{k},\frac{\left\|\tilde{g}_{k}\right\|}{\left\|\tilde{B}_{k}\right\|}\right)\geq\frac{{\gamma}}{\nu
rM}.$ (33)

If the step $p_{k}$ is accepted, we have from Step 12 of Algorithm 1 that
$\rho_{k}>c_{0}$, which by 7 is equivalent to

$\frac{\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)+r\epsilon_{f}}{m_{k}(0)-m_{k}\left(p_{k}\right)+r\epsilon_{f}}>c_{0}.$
(34)

Thus by 12, 29, 33 and 30

$\displaystyle\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)>$
$\displaystyle c_{0}\left[m_{k}(0)-m_{k}(p_{k})\right]+r(c_{0}-1)\epsilon_{f}$
$\displaystyle\geq$
$\displaystyle\frac{{c_{0}}}{2}\left\|\tilde{g}_{k}\right\|\min\left(\Delta_{k},\frac{\left\|\tilde{g}_{k}\right\|}{\left\|\tilde{B}_{k}\right\|}\right)+r(c_{0}-1)\epsilon_{f}$
$\displaystyle>$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left(r\epsilon_{g}+{\gamma}\right){\gamma}+r(c_{0}-1)\epsilon_{f}$
$\displaystyle>$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left(r\epsilon_{g}+\eta\right)\eta+r(c_{0}-1)\epsilon_{f}.$ (35)

We now chose $\eta$ so that the right hand side is positive. We obtain

$\eta\geq\frac{1}{2}\left(-r\epsilon_{g}+\beta\right)\quad\text{or}\quad\eta\leq\frac{1}{2}\left(-r\epsilon_{g}-\beta\right)$

We wish for $\eta$ to be the smallest positive value satisfying these
inequalities, yielding

$\eta=\frac{1}{2}\left(-r\epsilon_{g}+\beta\right).$ (36)

Substituting this quantity in 35, we have

$\displaystyle\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)>$
$\displaystyle\frac{{c_{0}}}{2\nu
rM}\left(r\epsilon_{g}+{\gamma}\right){\gamma}+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left(r\epsilon_{g}+\eta+\mu\right)(\eta+\mu)+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left(r\epsilon_{g}+\frac{1}{2}\left(-r\epsilon_{g}+\beta\right)+\mu\right)\left(\frac{1}{2}\left(-r\epsilon_{g}+\beta\right)+\mu\right)+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left(r\epsilon_{g}/2+\beta/2+\mu\right)\left(-r\epsilon_{g}/2+\beta/2+\mu\right)+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left[\left(\beta/2+\mu\right)^{2}-\left(r\epsilon_{g}/2\right)^{2}\right]+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left[(\beta/2)^{2}+\mu\beta+\mu^{2}-\left(r\epsilon_{g}/2\right)^{2}\right]+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left[\frac{\beta^{2}-(r\epsilon_{g})^{2}}{4}+\mu\beta+\mu^{2}\right]+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu
rM}\left[\frac{(r\epsilon_{g})^{2}+8\nu
r^{2}\left(\frac{1}{c_{0}}-1\right)M\epsilon_{f}-(r\epsilon_{g})^{2}}{4}+\mu\beta+\mu^{2}\right]+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle\frac{{c_{0}}}{2\nu rM}\left[{2\nu
r^{2}\left(\frac{1}{c_{0}}-1\right)M\epsilon_{f}}+\mu\beta+\mu^{2}\right]+r(c_{0}-1)\epsilon_{f}$
$\displaystyle=$ $\displaystyle r(1-c_{0})\epsilon_{f}+\frac{{c_{0}}}{2\nu
rM}\left(\mu\beta+\mu^{2}\right)+r(c_{0}-1)\epsilon_{f}$ $\displaystyle=$
$\displaystyle\frac{c_{0}}{2\nu rM}\left(\mu\beta+\mu^{2}\right).$

∎

The first inequality 29, together with 30, 31, identify the region where noise
does not dominate and progress in the objective function can be guaranteed.
The constant $\mu$ was introduced to ensure that our analysis is meaningful in
the case when noise is not present ($\epsilon_{f}=\epsilon_{g}=0$), as it
shows that a decrease in the objective is achieved. Nonetheless, the global
convergence results presented below are of interest only when noise is
present, so there we essentially absorb $\mu$ into $\eta$ by setting
$\mu=\epsilon_{g}/2$.

To summarize the results obtained so far, Lemma 2 states that when
$\|\tilde{g}_{k}\|$ is large enough, the trust region is either large enough
or will eventually be increased to be so. Lemma 3 states that when the
gradient and trust region are both large enough, every accepted iterate
reduces the noisy objective function by a non-vanishing amount. We show that
this drives iterations towards stationary points of the problem.

###  3.4 Global Convergence Theorems

Our global convergence results are presented in two parts. The first result
states that the iterates visit, infinitely often, a critical region
characterized by a small gradient norm. The second result states that after
visiting the above critical region for the first time, the iterates cannot
stray too far from it, as measured by the objective value.

######  Theorem 3.6 (Global Convergence to Critical Region).

Suppose that 3.1 through 3.5 are satisfied. Then, the sequence of iterates
$\\{x_{k}\\}$ generated by Algorithm 1 visits infinitely often the critical
region ${C}_{1}$ defined as

${C}_{1}=\left\\{x:\|g(x)\|\leq\left(r+1\right)\epsilon_{g}+\frac{\beta}{2}\right\\},$
(37)

where $r$ and $\beta$ are defined in 8, 30, 31, with $\mu=\epsilon_{g}/2$,
$\nu>1$ and $M$ given by 17.

###### Proof.

Assume by way of contradiction that there exist $K^{\prime}$ such that for all
$k>K^{\prime}$

$\|g(x_{k})\|>\left(r+1\right)\epsilon_{g}+\frac{\beta}{2}.$ (38)

Thus, by 3, definition 31 of $\eta$, and setting $\mu=\epsilon_{g}/2$, we have
that for all $k>K^{\prime}$

$\displaystyle\|\tilde{g}(x_{k})\|$
$\displaystyle>r\epsilon_{g}+\tfrac{1}{2}\beta$
$\displaystyle=-\tfrac{1}{2}r\epsilon_{g}+\tfrac{1}{2}\beta+\tfrac{3}{2}{r\epsilon_{g}}$
$\displaystyle=\eta+r\epsilon_{g}+\tfrac{1}{2}{r\epsilon_{g}}$
$\displaystyle>\,r\epsilon_{g}+\eta+\mu\qquad\qquad\mbox{(since $r>1$)}$
$\displaystyle=r\epsilon_{g}+{\gamma}.\qquad\quad\qquad\mbox{~{}~{}(by
\lx@cref{creftype~refnum}{eq:gammap})}$ (39)

We now apply Corollary 1 and deduce that there exist $K_{0}\geq K^{\prime}$,
such that for all $k\geq K_{0}$,

$\Delta_{k}>\frac{{\gamma}}{\nu rM}.$ (40)

When a step is not accepted, $\rho_{k}<c_{0}<c_{1}$, and Algorithm 1 will
reduce the trust region radius. If no step is accepted for all $k>K_{0}$, the
trust region radius would shrink to zero, contradicting 40. Therefore, there
must exist infinitely many accepted steps. Now, by 39, 40 the conditions of
Lemma 3 hold, and we deduce that each accepted step $k^{\prime}>K_{0}$
achieves the reduction

$\tilde{f}\left(x_{k^{\prime}}\right)-\tilde{f}\left(x_{k^{\prime}}+p_{k^{\prime}}\right)>\frac{c_{0}}{2\nu
rM}\left(\mu\beta+\mu^{2}\right)=\frac{c_{0}}{2\nu
rM}\left(\frac{\epsilon_{g}}{2}\beta+\frac{\epsilon_{g}^{2}}{4}\right).$ (41)

Since, as mentioned above, there is an infinite number of accepted steps, we
deduce that $\\{\tilde{f}(x_{k})\\}\rightarrow-\infty$, contradicting
Assumption 3.5. Therefore, the index $K^{\prime}$ defined above cannot exist
and we have that 38 is violated an infinite number of times. ∎

The achievable accuracy in the gradient guaranteed in 37 depends on
$\epsilon_{g}$ and $\sqrt{\epsilon_{f}}$, by the definition of $\beta$. The
dependence on $\epsilon_{g}$ is evident, while the dependence on
$\sqrt{\epsilon_{f}}$ is due to the combined (multiplicative) effect of the
gradient and the trust region radius bound.

Before stating our next theorem, we prove two simple technical results.

######  Proposition 1.

If Algorithm 1 takes a (nonzero) step at iteration $k$, then

$\tilde{f}_{k+1}-\tilde{f}_{k}<r(1-c_{0})\epsilon_{f}.$ (42)

###### Proof.

If the step is taken, we have from Step 12 of Algorithm 1 that
$\rho_{k}>c_{0}$, which by 7 is equivalent to

$\frac{\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)+r\epsilon_{f}}{m_{k}(0)-m_{k}\left(p_{k}\right)+r\epsilon_{f}}>c_{0},$
(43)

and since $p_{k}$ cannot increase the model $m_{k}$, we have

$\tilde{f}\left(x_{k}\right)-\tilde{f}\left(x_{k}+p_{k}\right)>c_{0}\left[m_{k}(0)-m_{k}(p_{k})\right]+r(c_{0}-1)\epsilon_{f}>r(c_{0}-1)\epsilon_{f}.$
(44)

∎

Next, we employ Lemma 2 and obtain the following result.

######  Corollary 2 (Maintaining Lower Bound on Trust Region Radius).

Let ${\gamma}>0$ be defined by 30–31, and suppose there exist $K>0$ and
$\hat{K}>K$ such that for $k=K+1,...,\hat{K}-1$

$\|\tilde{g}_{k}\|>r\epsilon_{g}+{\gamma},$ (45)

and that

$\Delta_{K+1}\geq\frac{\gamma}{\nu rM}=\frac{\bar{\Delta}}{\nu}.$ (46)

Then for $k=K+1,...,\hat{K}-1$

$\Delta_{k}\geq\frac{\gamma}{\nu rM}=\frac{\bar{\Delta}}{\nu}.$ (47)

###### Proof.

The proof is by induction. Condition 47 holds for $k=K+1$. We show that if 47
it holds for some $k\in\\{K+1,\ldots,\hat{K}-2\\}$, then it holds for $k+1$.

Specifically, suppose that for such $k$ we have that

$\Delta_{k}\geq\frac{\gamma}{\nu rM}.$ (48)

By Lemma 2, if $\Delta_{k}\leq\frac{\gamma}{rM}$, the trust region radius is
increased, i.e.,

$\Delta_{k+1}=\nu\Delta_{k}\geq\frac{\gamma}{rM}>\frac{\gamma}{\nu rM}.$ (49)

If on the other hand $\Delta_{k}>\frac{\gamma}{rM}$, the trust region radius
could be decreased, but in that case

$\Delta_{k+1}\geq\frac{\Delta_{k}}{\nu}>\frac{\gamma}{\nu rM}.$ (50)

∎

The next theorem shows that after an iterate has entered the neighborhood
$C_{1}$ defined in Theorem 3.6, all subsequent iterates cannot stray too far
away in the sense that their function values remain within a band of the
largest function value in $C_{1}$.

######  Theorem 3.7 (Iterates Remain in the Level Set $C_{2}$).

Suppose that 3.1 through 3.5 are satisfied. Then, after the iterates $x_{k}$
generated by Algorithm 1 visit $C_{1}$ for the first time, they never leave
the set $C_{2}$ defined as

${C}_{2}=\left\\{x:f(x)\leq\sup_{y\in
C_{1}}f(y)+2\epsilon_{f}+\max[G,r(1-c_{0})\epsilon_{f}]\right\\},$ (51)

where

$G=\left[(r+1)\epsilon_{g}+\gamma+\frac{{\nu^{2}}L\gamma}{(\nu-1)rM}\right]\frac{\nu^{2}\gamma}{(\nu-1)rM},$
(52)

and $\gamma$ is defined in 30–31 with $\mu={\epsilon_{g}}/2$.

###### Proof.

The proof is based on the observation that, when the iterates leave $C_{1}$,
if the trust region is large enough, then by Lemma 3 the noisy objective
function starts decreasing immediately (Case 1); otherwise the smallness of
the trust region limits the increase in the objective function before the
trust region becomes large enough to ensure descent (Case 2). We now state
this precisely.

Suppose that the $K^{th}$ step is an exiting step, i.e., $x_{K}\in C_{1}$ and
$x_{K+1}\notin C_{1}$. We let $\hat{K}>K+1$ be the index of the first iterate
that returns to $C_{1}$. Such a $\hat{K}$ exists due to Theorem 3.6. We will
prove that all iterates $x_{k}$ with $k\in\\{K+1,\ldots,\hat{K}-1\\}$ are
contained in $C_{2}$.

Since $x_{k}\notin C_{1}$ for $k\in\\{K+1,\ldots,\hat{K}-1\\}$, we have by 37
that

$\|g_{k}\|>\left(r+1\right)\epsilon_{g}+\frac{\beta}{2},$ (53)

and we have seen in 38-39 that this implies that

$\|\tilde{g}_{k}\|>r\epsilon_{g}+{\gamma},\qquad\
k\in\\{K+1,\ldots,\hat{K}-1\\}.$ (54)

Also, we know that a step was taken at iterate $K$ since $x_{K}\in C_{1}$ and
$x_{K+1}\notin C_{1}$, and thus applying Proposition 1 yields

$\tilde{f}_{K+1}-\tilde{f}_{K}<r(1-c_{0})\epsilon_{f}.$ (55)

We divide the rest of the proof according to the size of $\Delta_{K+1}$
relative to $\bar{\Delta}$, which is defined in 20, i.e.,

$\bar{\Delta}=\frac{\gamma}{rM}.$ (56)

Case 1: Suppose $\Delta_{K+1}\geq\bar{\Delta}$. By 54 and the fact that
$\nu>1$, the conditions of Corollary 2 are satisfied and thus
$\Delta_{k}>\frac{\gamma}{\nu rM}$, for $k=K+1,\ldots,\hat{K}-1$. We can
therefore apply Lemma 3, with $\mu=\epsilon_{g}/2>0$, for each iterate
$k=K+1,\dots,\hat{K}-1$ to yield

$\tilde{f}(x_{K+1})\geq\tilde{f}(x_{K+2})\geq\cdots\geq\tilde{f}(x_{\hat{K}}).$
(57)

Combining this result with 55 we obtain

$\tilde{f}_{k}\leq\tilde{f}_{K+1}<\tilde{f}_{K}+r(1-c_{0})\epsilon_{f},\quad
k=K+1,..,\hat{K}.$ (58)

Since $x_{K}\in C_{1}$ and by 3, we conclude that for $k=K+1,\ldots,\hat{K}$,

$f_{k}<f_{K}+[2+r(1-c_{0})]\epsilon_{f}\leq\sup_{y\in
C_{1}}f(y)+[2+r(1-c_{0})]\epsilon_{f}.$ (59)

Therefore, the inequality in 51 is satisfied in this case.

Case 2: Suppose $\Delta_{K+1}<\bar{\Delta}$. We begin by considering the
increase in the function value while the trust region remains less than
$\bar{\Delta}.$ To this end, we define

$l=\left\lceil{\log_{\nu}{\frac{\bar{\Delta}}{\Delta_{K+1}}}}\right\rceil,\quad$
(60)

where $\lceil{\cdot}\rceil$ denotes the ceiling operation. Since the trust
region radius is increased by a factor of at most $\nu$, we have that $l$ is
the minimum number of steps required for the trust region radius to increase
from $\Delta_{K+1}$ to (at least) $\bar{\Delta}$. Now, if $K+l>\hat{K}$, then
the iterates return to $C_{1}$ before the trust region becomes at least
$\hat{\Delta}$. Therefore, the number of out-of-$C_{1}$ iterations taken by
the algorithm while $\Delta_{k}<\hat{\Delta}$ is

${\hat{l}=\min\\{l-1,\hat{K}-K-1\\}.}$ (61)

The increase in function values for iterations indexed by
$k=K+1,\ldots,K+\hat{l}+1$ is bounded as follows:

$\displaystyle|f(x_{k})-f(x_{K})|$
$\displaystyle\leq\sum_{i=0}^{k-K-1}|f(x_{K+1+i})-f(x_{K+i})|$
$\displaystyle\leq\sum_{i=0}^{{\hat{l}}}|f(x_{K+1+i})-f(x_{K+i})|$
$\displaystyle\leq\sum_{i=0}^{{\hat{l}}}\Delta_{K+i}\max_{x\in[x_{K+i},x_{K+1+i}]}\|g(x)\|$
$\displaystyle=\sum_{i=0}^{{\hat{l}}}\Delta_{K+i}\max_{x\in[x_{K+i},x_{K+1+i}]}\|g(x)-g(x_{K+i})+g(x_{K+i})\|$
$\displaystyle\leq\sum_{i=0}^{{\hat{l}}}\Delta_{K+i}\big{[}\|g(x_{K+i})\|+L\Delta_{K+i}\big{]}\qquad\mbox{(by
\lx@cref{creftype~refnum}{eq:lip})}.$ (62)

To estimate the right hand side, we need to bound the total displacement made
by the algorithm during those iterations. It follows from 60 that

$\bar{\Delta}/\nu\leq\nu^{l-1}\Delta_{K+1}<\bar{\Delta}\leq\nu^{l}\Delta_{K+1},$
(63)

and thus for $i=0,...,{\hat{l}}$,

$\Delta_{K+1+i}\leq\nu^{i}\Delta_{K+1}\leq\nu^{\hat{l}}\Delta_{K+1}\leq\nu^{l-1}\Delta_{K+1}<\bar{\Delta}.$
(64)

By 54, 64, we can apply Lemma 2 to each iterate $i=0,...,\hat{l}$, and obtain

$\Delta_{i+1}=\nu\Delta_{i}.$ (65)

Thus for $i=0,...,\hat{l}$,

$\Delta_{K+1+i}=\nu^{i}\Delta_{K+1}\leq\nu^{\hat{l}}\Delta_{K+1}\leq\nu^{l-1}\Delta_{K+1}<\bar{\Delta}.$
(66)

Summing from $i=0$ to ${\hat{l}}$, we have

$\sum_{i=0}^{{\hat{l}}}\Delta_{K+1+i}=\sum_{i=0}^{{\hat{l}}}\nu^{i}\Delta_{K+1}<\frac{\bar{\Delta}}{\nu^{\hat{l}}}\sum_{i=0}^{{\hat{l}}}\nu^{i}=\frac{\bar{\Delta}}{\nu^{{\hat{l}}}}\frac{\nu^{{\hat{l}}+1}-1}{\nu-1}<\frac{\bar{\Delta}}{\nu^{{\hat{l}}}}\frac{\nu^{{\hat{l}}+1}}{\nu-1}=\frac{\nu}{\nu-1}\bar{\Delta}.$
(67)

By assumption, $\Delta_{K+1}<\bar{\Delta}$, which implies
$\Delta_{K}<\nu\bar{\Delta}$; adding this to 67 we obtain

$\sum_{i=0}^{{\hat{l}}+1}\Delta_{K+i}<\frac{\nu^{2}}{\nu-1}\bar{\Delta}.$ (68)

Therefore, for $i=0,\ldots,{\hat{l}}$,

$\displaystyle\|g(x_{K+i})\|+L\Delta_{K+i}$
$\displaystyle=\|g(x_{K})+\sum_{j=0}^{i-1}\left[g(x_{K+j+1})-g(x_{K+j})\right]\|+L\Delta_{K+i}$
$\displaystyle\leq\left\|g(x_{K})\right\|+\sum_{j=0}^{i-1}\left\|g(x_{K+j+1})-g(x_{K+j})\right\|+L\Delta_{K+i}$
$\displaystyle\leq\left\|g(x_{K})\right\|+\left(\sum_{j=0}^{i-1}L\Delta_{K+j}\right)+L\Delta_{K+i}$
$\displaystyle<\left\|g(x_{K})\right\|+L\sum_{j=0}^{\hat{l}+1}\Delta_{K+j}\quad(\mbox{since}\
i<{\hat{l}}+1)$
$\displaystyle<\left\|g(x_{K})\right\|+\frac{\nu^{2}}{\nu-1}L\bar{\Delta}\qquad(\mbox{by
\lx@cref{creftype~refnum}{eq:sumdel}}).$ (69)

Substituting this inequality into 62, we obtain for any
$k=K+1,...,K+\hat{l}+1$,

$\displaystyle|f(x_{k})-f(x_{K})|$
$\displaystyle\leq\sum_{i=0}^{{\hat{l}}}\Delta_{K+i}\left[\left\|g(x_{K})\right\|+\frac{\nu^{2}}{\nu-1}L\bar{\Delta}\right]$
$\displaystyle<\left[\|g(x_{K})\|+\frac{\nu^{2}}{\nu-1}L\bar{\Delta}\right]\frac{\nu^{2}}{\nu-1}\bar{\Delta}$
$\displaystyle\leq\left[(r+1)\epsilon_{g}+\gamma+\frac{\nu^{2}}{\nu-1}L\bar{\Delta}\right]\frac{\nu^{2}}{\nu-1}\bar{\Delta}\qquad(\mbox{since
$x_{K}\in C_{1}$})$
$\displaystyle=\left[(r+1)\epsilon_{g}+\gamma+\frac{{\nu^{2}}L\gamma}{(\nu-1)rM}\right]\frac{\nu^{2}\gamma}{(\nu-1)rM}\qquad\mbox{(by
\lx@cref{creftype~refnum}{eq:dhat+})}$ $\displaystyle=G.$ (70)

Therefore, for $k=K+1,\ldots,K+\hat{l}+1$,

$f(x_{k})<f(x_{K})+G\leq\sup_{y\in C_{1}}f(y)+G.$ (71)

We now consider two possibilities.

Case 2a): Suppose $K+1+l>\hat{K}$. Then, $\hat{K}-K-1\leq l-1$ and by 61 we
have that $\hat{l}=\hat{K}-K-1$. Condition 71, thus reads

$f(x_{k})<f(x_{K})+G\leq\sup_{y\in C_{1}}f(y)+G,\qquad k=K+1,\ldots,\hat{K},$
(72)

and thus the inequality in 51 is satisfied for $k=K+1,\ldots,\hat{K}-1$.

Case 2b): suppose $K+1+l\leq\hat{K}$. Then, by 60 we have that $\hat{l}=l-1$,
and 71 reads

$f(x_{k})<f(x_{K})+G\leq\sup_{y\in C_{1}}f(y)+G\qquad k=K+1,\ldots,K+l.$ (73)

Let us now consider the iterates following $K+l$ that are outside $C_{1}$,
i.e., those indexed by $k=K+l+1,...,\hat{K}-1$. Letting $i={\hat{l}}=l-1$ in
66 and recalling the first inequality in 63,

$\Delta_{K+l}=\nu^{l-1}\Delta_{K+1}\geq\frac{\bar{\Delta}}{\nu}.$ (74)

We can therefore apply Corollary 2 to iterates indexed by
$k=K+l+1,...,\hat{K}-1$ and deduce that

$\Delta_{k}\geq\frac{\bar{\Delta}}{\nu},\quad k=K+l+1,...,\hat{K}-1.$

This fact, together with 54, allow us to invoke Lemma 3, for
$k=K+l,...,\hat{K}-1$, to yield

$\tilde{f}(x_{K+l})\geq\tilde{f}(x_{K+1+l})\geq\tilde{f}(x_{K+2+l})\geq...\geq\tilde{f}(x_{\hat{K}}).$
(75)

Recalling 73 with $k=K+l$ and using 3 we obtain

$\tilde{f}(x_{K+l})<\sup_{y\in C_{1}}f(y)+G+\epsilon_{f}.$ (76)

This condition together with 75 yields

$f(x_{k})\leq\tilde{f}(x_{K+l})+\epsilon_{f}<\sup_{y\in
C_{1}}f(y)+G+2\epsilon_{f}\qquad k=K+l,...,\hat{K}-1.$ (77)

Combining this bound with 73 we conclude

$f(x_{k})<\sup_{y\in C_{1}}f(y)+G+2\epsilon_{f},\qquad k=K+1,...,\hat{K}-1,$
(78)

and thus the inequality in 51 is satisfied.

∎

The constant $G$ defined in 52 is proportional to
$\epsilon_{g}^{2},\epsilon_{g}\sqrt{\epsilon_{f}},\epsilon_{f}$. Since that
$G$ characterizes the function value bounds, the dependence on $\epsilon_{f}$
is expected; the dependence on $\epsilon_{g}$ and
$\epsilon_{g}\sqrt{\epsilon_{f}}$ arises from the combined effect of the trust
region radius and gradient norm.

##  4 Numerical Experiments

To illustrate the performance of the proposed Algorithm 1, we coded it in
matlab and applied it to a small selection of unconstrained optimization
problems. We injected uniformly distributed noise in the evaluations of the
function and gradient. Specifically, we let (c.f. 2)

$\delta_{f}=X_{f}\in\mathbb{R},\quad X_{f}\sim
U(-\epsilon_{f},\epsilon_{f}),\quad\mbox{and}\quad\delta_{g}=X_{g}\in\mathbb{R}^{n},\quad
X_{g}\sim\mathbb{B}_{n}(0,\epsilon_{g}),$ (79)

where $U(-a,a)$ denotes the uniform distribution from $-a$ to $a$, and
$\mathbb{B}_{n}(0,a)$ denotes the $n$ dimensional ball centered at $0$ with
radius $a$. By generating noise in this way we satisfy 3.2.

We set the parameters in Algorithm 1 as follows:
$c_{0}=0.1,c_{1}=1/4,c_{2}=1/2$ and $\nu=2$. The solution of the trust region
subproblem (Step 3 of Algorithm 1) was computed using the standard Newton-CG
method described e.g. in mybook , with termination accuracy $10^{-8}$. In
order to better illustrate the performance of the algorithm in the presence of
noise, we did not include a stop test and simply ran it for 200 iterations,
which was sufficient to observe its asymptotic behavior.

###  4.1 Failure of the Classical Trust Region Algorithm

We present two examples showing failure of the classical trust region
algorithm, in contrast with Algorithm 1. First, we consider the simple
quadratic function

$f=x^{T}Dx,$ (80)

where $x\in\mathbb{R}^{8}$ and $D$ is the diagonal matrix

$D={\rm diag}(1e-5,1e-4.75,1e-4.5,....,1e-3.25).$ (81)

The condition number of $D$ is roughly $56$. We set $\epsilon_{f}=10^{-1}$ and
$\epsilon_{g}=10^{-5}$ in 79. The Hessian of the quadratic model 4 was defined
as $B_{k}=\nabla^{2}f(x_{k})$; i.e., we did not inject noise in this
experiment. We started both algorithms from $x_{0}=(1000,0,0,....,0)$, with an
initial trust region radius $\Delta_{0}=1$. The results are displayed Figure
1.

Figure 1: New and classical trust region algorithms applied to a simple
quadratic problem.

The four panels in Figure 1 compare the performance of the classical algorithm
(red dashed line) and Algorithm 1 (blue solid line). The horizontal axis in
each panel records the iteration number. In the upper left panel (a) we report
the norm of the (noiseless) gradient $\|\nabla f(x_{k})\|$, along with the
injected noise level $\epsilon_{g}$ (solid black line); the light blue dashed
line plots the lowest value generated by Algorithm 1 in the past 25
iterations. In the upper right panel (b) we report the trust region radius; in
the lower-left panel (c) the distance to solution; and in the lower right
panel (d), the computed actual-to-predicted reduction ratio $\rho_{k}$; for
graphical clarity, ratios greater than $5$ or less than $-5$ were plotted as
$+/-5$ in panel (d).

We observe that the classical algorithm exhibits large oscillations in
$\rho_{k}$, which causes the trust region radius to shrink so much that
significant progress cannot be made. In contrast, $\rho_{k}$ is controlled
well in Algorithm 1. In this test, initial the trust region radius
$\Delta_{0}$ is not small.

In the next experiment, we illustrate the damaging effect that a very small
$\Delta_{0}$ can have on the classical algorithm, but not on the proposed
algorithm. We applied the two algorithms to the following tri-diagonal
function

$f(x)=\frac{1}{2}\left(x^{(1)}-1\right)^{2}+\frac{1}{2}\sum_{i=1}^{N-1}\left(x^{(i)}-2x^{(i+1)}\right)^{4},\quad
N=200.$ (82)

The results are reported in Figure 2. In the upper left panel, we additionally
plot in purple the size of the critical region $C_{1}$, i.e. the value of the
right-hand side in 37. (The latter requires knowledge of the constant $M$,
which we approximate by the norm of the Hessian at the solution.) This panel
shows that the theoretical prediction given in Theorem 3.6 is pessimistic when
compared to the final achieved accuracy in the gradient, as is to be expected
of convergence results that assume that the largest possible error occurs at
every iteration. The upper right hand panel illustrates that Algorithm 1 is
able to quickly increase the trust region radius an allow progress, unlike the
classical algorithm.

Figure 2: New and classical trust region algorithms initialized with small
trust region radius.

###  4.2 General Performance of the Proposed Algorithm

We also tested the two algorithms on a subset of problems from
schittkowski1987more ; the results are presented in the supplementary
material. As a representative of these runs, we report the results for the
tri-diagonal objective function 82. This time, the Hessian $B_{k}$ of the
quadratic model 4 is obtained by injecting noise in the true Hessian matrix.
We define

$B_{k}=\nabla^{2}f(x_{k})+\delta_{B},$ (83) $\delta_{B}=\frac{A^{T}\Lambda
A}{\|A\|^{2}},\quad A_{ij}\sim U(0,1),\quad(\Lambda)_{ii}\sim
U(-\epsilon_{B},\epsilon_{B}),$ (84)

where $\Lambda$ is a diagonal matrix. Thus, the matrices $B_{k}$ are symmetric
but not necessarily positive definite. We employed larger noise levels than in
the previous experiments: $\epsilon_{f}=10$, $\epsilon_{g}=100$, and
$\epsilon_{B}=1000$. This simulates the situation that may occur when
employing finite difference approximations, where the error increases with the
order of differentiation. Both algorithms were initialized from the same
starting point $x_{0}$, which was generated such that each entry in $x_{0}$ is
sampled uniformly from $-50$ to $50$. To ensure a fair comparison, at each
iterate we inject exactly the same noise into both algorithms.

We report the results in Figure 3, which displays the same information as in
Figure 2. We observe that both algorithms perform similarly before entering
the noisy regime. Algorithm 1 exhibits larger oscillations in the gradient
norm due to the larger trust region radius, but achieves a lower objective
function value. Whereas the large reduction in the trust region radius led to
failures of the classical algorithm in the examples reported above, in many
test runs such as that given in Figure 3, it can be beneficial by producing
increasingly smaller steps that yield milder oscillations in the gradient norm
than Algorithm 1. We cannot, however, recommend this type of trust region
reduction as a general procedure for handling noise since failures can happen
unexpectedly.

Figure 3: Comparison of the new and classical trust region algorithms when
solving problem 82 with uniform noise given by 79 83.

###  4.3 Evaluating the Theoretical Results

We have seen that the critical region $C_{1}$ gives a pessimistic estimate of
the achievable accuracy in the gradient because the analysis assumes worst-
case behavior at each iteration, rather than providing estimates in high
probability. Nevertheless, Theorem 3.6 identifies the functional relationship
between the achievable accuracy and the noise level: the right hand side in 37
scales as a function of $\epsilon_{g}$ and $\sqrt{\epsilon_{f}}$. We performed
numerical tests to measure if the accuracy achieved in practice scales in that
manner.

We employed the tridiagonal function 82, for which we can estimate the
constant $M$, as mentioned above. For given $\epsilon_{f}$ and $\epsilon_{g}$,
we compute the right hand side in 37, which we denote as
$C(\epsilon_{f},\epsilon_{g})$, and ran Algorithm 1 as in the previous test.
We repeated the run 10 times using different seeds, $s=1,\ldots,10$, to
generate noise. For each run, we track the smallest value of
$\|\tilde{g}_{k}\|$ during the most recent 25 iterations and record the
smallest such value observed during the run, which we denote as
$\|\tilde{g}^{*}({\epsilon_{g},\epsilon_{f},s})\|$, where $s$ denotes the
seed. In Figure 4, we report the quantity

$R(\epsilon_{f},\epsilon_{g})=\log_{10}\frac{C(\epsilon_{f},\epsilon_{g})}{\sum_{s=1}^{10}\|\tilde{g}^{*}({\epsilon_{g},\epsilon_{f},s})\|}$
(85)

as we vary $\epsilon_{f}$ and $\epsilon_{g}$ from $10^{-2}$ to $10^{2}$.

Figure 4: $R(\epsilon_{f},\epsilon_{g})$ given in 85: Log10 of the ratio
between predicted and actual accuracy in the gradient, as a function of these
noise level $\epsilon_{f},\epsilon_{g}$. The small variation in these numbers
suggests that Theorem 3.6 gives the correct dependence on the noise levels.

The fact that the ratio between the theoretical bound and the smallest
gradient norm measured in practice remained roughly constant gives numerical
support to the claim that the achievable gradient norm is proportional to
$\epsilon_{g}$ and $\sqrt{\epsilon_{f}}$. We should note that these
observations are valid only when averaging multiple runs with different seeds,
as one can observe significant variations among individual runs of Algorithm
1.

##  5 Final Remarks

In this paper, we proposed a noise-tolerant trust region algorithm that avoids
the pitfall of the classical algorithm, which can shrink the trust region
prematurely, preventing progress toward a stationary point. Robustness is
achieved by relaxing the ratio test used in the step acceptance, so as to
account for errors in the function.

We showed that when the noise in the function and gradient evaluations is
bounded by the constants $\epsilon_{f},\epsilon_{g}$, an infinite subsequence
of iterates satisfies

$\|g_{k}\|=O({\sqrt{\epsilon_{f}}},\epsilon_{g}).$ (86)

When noise is not present, our results yield the limit
$\\{\|g_{k}\|\\}\rightarrow 0$ (the sets $C_{1}$ and $C_{2}$ in Theorem 3.6
and Theorem 3.7 coincide in this case).

The technique and analysis presented here are relevant to the case when noise
can be diminished as needed, as assumed e.g. in katya2018storm ;
cartis2019TRsupermartingale ; bollapragada2018adaptive . Algorithm 1 can be
run until it ceases to make significant progress, at which point the accuracy
in the function and gradient is increased (i.e., $\epsilon_{f},\epsilon_{g}$
are reduced) and the algorithm is restarted with the new value of
$\epsilon_{f}$ in 7; this process can then be repeated. This provides a
disciplined approach for achieving high accuracy in the solution using a
noise-tolerant trust region algorithm.

## Acknowledgments

We thank Richard Byrd for many valuable suggestions and insights, as well as
Figen Öztoprack, Andreas Wächter and Melody Xuan for proofreading the paper
and providing valuable feedback.

## References

  * (1) Bellavia, S., Gurioli, G., Morini, B., Toint, P.: The impact of noise on evaluation complexity: The deterministic trust-region case.  arXiv preprint arXiv:2104.02519 (2021)
  * (2) Berahas, A.S., Byrd, R.H., Nocedal, J.: Derivative-free optimization of noisy functions via quasi-Newton methods.  SIAM Journal on Optimization 29(2), 965–993 (2019)
  * (3) Berahas, A.S., Cao, L., Scheinberg, K.: Global convergence rate analysis of a generic line search algorithm with noise.  SIAM Journal on Optimization 31((2)), 1489–1518 (2021)
  * (4) Berahas, A.S., Curtis, F.E., O’Neill, M.J., Robinson, D.P.: A stochastic sequential quadratic optimization algorithm for nonlinear equality constrained optimization with rank-deficient Jacobians.  arXiv preprint arXiv:2106.13015 (2021)
  * (5) Berahas, A.S., Curtis, F.E., Robinson, D., Zhou, B.: Sequential quadratic optimization for nonlinear equality constrained stochastic optimization.  SIAM Journal on Optimization 31(2), 1352–1379 (2021)
  * (6) Blanchet, J., Cartis, C., Menickelly, M., Scheinberg, K.: Convergence rate analysis of a stochastic trust region method via submartingales.  INFORMS Journal on Optimization 1((2)), 92–119 (2019)
  * (7) Bollapragada, R., Byrd, R., Nocedal, J.: Adaptive sampling strategies for stochastic optimization.  SIAM Journal on Optimization 28(4), 3312–3343 (2018)
  * (8) Bollapragada, R., Byrd, R.H., Nocedal, J.: Exact and inexact subsampled newton methods for optimization.  IMA Journal of Numerical Analysis (2018).  DOI 10.1093/imanum/dry009.  URL http://dx.doi.org/10.1093/imanum/dry009
  * (9) Byrd, R.H., Chin, G.M., Nocedal, J., Wu, Y.: Sample size selection in optimization methods for machine learning.  Mathematical Programming 134(1), 127–155 (2012)
  * (10) Carter, R.G.: On the global convergence of trust region algorithms using inexact gradient information.  SIAM Journal on Numerical Analysis 28(1), 251–265 (1991)
  * (11) Cartis, C., Gould, N.I.M., Toint, P.: Strong evaluation complexity of an inexact trust-region algorithm with for arbitrary-order unconstrained nonconvex optimization.  arXiv preprint arXiv:2001.10802 (2021)
  * (12) Cartis, C., Scheinberg, K.: Global convergence rate analysis of unconstrained optimization methods based on probabilistic models.  Mathematical Programming 169(2), 337–375 (2018)
  * (13) Chen, R., Menickelly, M., Scheinberg, K.: Stochastic optimization using a trust-region method and random models.  Mathematical Programming 169(2), 447–487 (2018)
  * (14) Curtis, F.E., Robinson, D.P., Zhou, B.: Inexact sequential quadratic optimization for minimizing a stochastic objective function subject to deterministic nonlinear equality constraints.  arXiv preprint arXiv:2107.03512 (2021)
  * (15) Curtis, F.E., Scheinberg, K.: Adaptive stochastic optimization: A framework for analyzing stochastic optimization algorithms.  IEEE Signal Processing Magazine 37(5), 32–42 (2020)
  * (16) Curtis, F.E., Scheinberg, K., Shi, R.: A stochastic trust region algorithm based on careful step normalization.  INFORMS Journal on Optimization 1(3), 200–220 (2019)
  * (17) Gould, N.I.M., Toint, P.: An adaptive regularization algorithm for unconstrained optimization with inexact function and derivatives values.  arXiv preprint arXiv:2111.14098 (2021)
  * (18) Jin, B., Scheinberg, K., Xie, M.: High probability complexity bounds for line search based on stochastic oracles.  arXiv preprint arXiv:2106.06454 (2021)
  * (19) Nesterov, Y., Spokoiny, V.: Random gradient-free minimization of convex functions.  Foundations of Computational Mathematics 17(2), 527–566 (2017)
  * (20) Nocedal, J., Wright, S.: Numerical Optimization, 2 edn.  Springer New York (1999)
  * (21) Öztoprak, F., Byrd, R., Nocedal, J.: Constrained optimization in the presence of noise.  arXiv preprint arXiv:2110.04355 (2021)
  * (22) Paquette, C., Scheinberg, K.: A stochastic line search method with convergence rate analysis.  arXiv preprint arXiv:1807.07994 (2018)
  * (23) Schittkowski, K.: More test examples for nonlinear programming codes.  Lecture Notes in Econom. and Math. Systems 282 (1987)
  * (24) Shi, H.J.M., Xie, Y., Byrd, R., Nocedal, J.: A noise-tolerant quasi-newton algorithm for unconstrained optimization.  arXiv preprint arXiv:2010.04352 (2020)
  * (25) Shi, H.J.M., Xie, Y., Xuan, M.Q., Nocedal, J.: Adaptive finite-difference interval estimation for noisy derivative-free optimization.  arXiv preprint arXiv:2110.06380 (2021)
  * (26) Shi, H.J.M., Xuan, M.Q., Oztoprak, F., Nocedal, J.: On the numerical performance of derivative-free optimization methods based on finite-difference approximations.  arXiv preprint arXiv:2102.09762 (2021)
  * (27) Xie, Y., Byrd, R.H., Nocedal, J.: Analysis of the BFGS method with errors.  SIAM Journal on Optimization 30(1), 182–209 (2020)