# A Stochastic Trust Region Method for Unconstrained Optimization Problems

Authors: Ningning Li; Dan Xue; Wenyu Sun; Jing Wang
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

In this paper, a stochastic trust region method is proposed to solve unconstrained minimization problems with stochastic objectives. Particularly, this method can be used to deal with nonconvex problems. At each iteration, we construct a quadratic model of the objective function. In the model, stochastic gradient is used to take the place of deterministic gradient for both the determination of descent directions and the approximation of the Hessians of the objective function. The behavior and the convergence properties of the proposed method are discussed under some reasonable conditions. Some preliminary numerical results show that our method is potentially efficient.


---

## Body

## 1. Introduction

The computational complexity of learning algorithm becomes the critical limiting factor when one envisions very large datasets. This contribution advocates many stochastic optimization algorithms for large-scale problems. To a certain extent, stochastic optimization is effective to solve the problems in stochastic systems. A number of industrial, biological, engineering, and economic problems can be considered as stochastic systems, for example, area of communication, gene, signal processing, geography, civil engineering, aerospace, and banking. Particularly, stochastic optimization algorithms are used to solve the problem of optimizing an objective function over a set of feasible values in situations where the objective function is defined as an expectation over a set of random functions. To be more precise, consider an optimization variablex∈Rn and a random variable θ∈Θ⊆Rp that determines the choice of a function f(x,θ):Rn×p→R. Stochastic optimization problems considered in this paper entail determination of the argument x∗ that minimizes the expected value F(x):=Eθ[f(x,θ)]:(1)x∗≔arg⁡minx⁡Eθfx,θ≔arg⁡min⁡Fx.We refer to f(x,θ) as the random or instantaneous functions and to F(x):=Eθ[f(x,θ)] as the average function. We assume that the instantaneous function f(x,θ) is continuously differentiable for all θ from which it follows that the average function F(x) is also continuously differentiable. Problems having the form in (1) are often used in machine learning [1–3] as well as in optimal resource allocation in wireless systems [4, 5].Many descent algorithms can be used for solving such problems as (1). However, descent methods require exact gradient of the objective function, i.e., ∇xF(x)=Eθ[∇xf(x,θ)], which is impracticable generally as a result of the very large dataset. Stochastic gradient descent (SGD) methods overcome this obstacle by using gradient approximation which is based on small data samples and are regarded as the workhorse methodology to solve large-scale stochastic optimization problems [3, 6–9]. But practical appeal of SGD methods remains limited, and the number of iterations required to approximate optimal arguments can be prohibitive in high dimensional problems. In fact, SGD inherits slow convergence from its use of gradients which is exasperated by their substitution with stochastic approximation; consecutive stochastic gradients may vary a lot or even point in opposite directions. Several accelerations of SGD have been presented, and lots of recent works have focused on reducing randomness in SGD by combining gradient and stochastic gradient or updating descent direction so that they become more close to gradient gradually. For example, the Stochastic Average Gradient (SAG) method achieves a faster convergence rate than SGD by incorporating a memory of previous gradient values, which often outperform existing SGD methods and were shown to exhibit linear convergence [10]. The Semistochastic Gradient (S2GD) method in [11] runs for one or several epochs in each of which a single full gradient and a random number of stochastic gradients are computed, following a geometric law. These algorithms end up showing a faster asymptotic convergence rate than SGD [12]; however, practical numerical experiments show that reducing randomness is of no use for problems with a challenging curvature profile.To overcome problems with the objective function’s curvature, one may consider stochastic version of Newton’s method. But computing approximation of Newton steps is difficult except in problems with specific structures because of the challenging curvature [13]. Recourse to quasi-Newton methods, regularized stochastic BFGS (RES) method, is proposed. RES utilizes stochastic gradients in lieu of deterministic gradients for determining descent directions and approximating the objective function’s curvature. As we can see in [14], RES has been shown to outperform SGD in large-dimensional problems or ill-conditioned functions. However, RES inherits the disadvantage that the Hessian matrix of objective functions must be positive definite which implies that the strongly convex objective functions are indispensable to guarantee provable convergence. Thus, RES can only be used to solve strongly convex stochastic optimization problems. But there are lots of nonconvex optimization problems in practical applications. So, in this paper, we would like to exploit an effective method for solving more general stochastic optimization problems, i.e., the objective functions can be nonconvex.The trust region methods are deemed to be invaluable tools for solving nonlinear and nonconvex optimization problems [15–17]. The main idea of typical trust region method is to construct a quadratic model and choose a step to be the approximate minimizer of the quadratic model in the trust region around the current point in which the model is trusted to be adequate to the objective function. In the iterative procedures, the trust region radii are adjusted depending on the agreement between the model functions and the objective functions.In this paper, we develop a stochastic version of trust region method, where the Hessian matrices are approximated by regularized stochastic BFGS method [14]. Similar to classical trust region methods, we use a suitable quadratic model to replace the complex objective function. In the model, deterministic gradients are replaced by stochastic gradients and exact Hessian matrices of the objective function are approximated by its stochastic approximations. Due to the fact that stochastic gradients are computable at manageable cost, stochastic trust region method is feasible in practice. Convergence theory and numerical results verify the reliability of our method to deal with both convex and nonconvex optimization problems.The rest of the paper is organized as follows. The new algorithm is illustrated in Section2. In Section 3, convergence of the algorithm is established under suitable conditions. Some preliminary numerical results are reported in Section 4. Finally, some conclusions are given in Section 5.


## 2. Stochastic Trust Region Algorithm

Considering that the objective functionF(x) is continuously twice differentiable and further assuming that the instantaneous functions have finite gradients, it follows that the gradients of F(x) are given by(2)gx≔∇Fx=Eθ∇fx,θ.When the number of functions f(x,θ) is large, as is the case in most problems of practical interest, exact evaluation of the gradient g(x) is impractical. This motivates the use of stochastic gradients in lieu of exact gradients. In addition, the Hessian matrices of F(x) are given by(3)Hx≔∇2Fx=Eθ∇2fx,θ.More precisely, consider a given set of L realizations θ~=[θ1;⋯;θL] and define the stochastic gradient of F(x) at given samples θ~ as(4)g^x,θ~≔1L∑l=1L∇fx,θl.Consider the trust region subproblem(5)minsq^ks=Fxk+g^xk,θ~Ts+12sTG^kss.t.s2≤hk.where s=x-xk, hk∈(0,h¯) is the trust region radius, and G^k is a stochastic approximation of Hessian matrix of F(x) (see Remark 2 for details). The model q^(k)(s) is minimized approximately to produce a step sk which is subject to having the norm less than or equal to hk. In other words, the solution of this subproblem sk represents a step toward minimizing the model q^(k)(s) of the objective function F(x) at xk. The radio r^k of the actual reduction of the objective function aredk=F(xk)-F(xk+sk) to the predicted reduction predk=q^(k)(0)-q^(k)(sk) is calculated. If r^k≥η1, where η1∈(0,1), then sk is accepted, which is called a successful step. Let xk+1=xk+sk; in this case, the trust region radius is increased or remains unchanged. Otherwise, sk is rejected, which is called an unsuccessful step in which xk remains unchanged and the trust region radius is reduced.Algorithm 1 (stochastic trust region).
Initialization. The variable is x0. The symmetric matrix is G^0. The initial radius is h0∈(0,h-), where h->0, and choose constant 0<γ1<1<γ2, 0<η1≤η2≤1.Step 1.
AcquireL independent samples θ~k=[θk1;⋯;θkL], and calculate (6)g^xk,θ~k≔1L∑l=1L∇fxk,θkl.Step 2.
Solve trust region subproblem (5), giving a trial step sk.Step 3.
Compute(7)r^k=aredkpredk,
wherearedk=F(xk)-F(xk+sk) and predk=q^(k)(0)-q^(k)(sk). If r^k≥η1, set xk+1=xk+sk; otherwise, set xk+1=xk.Step 4.
Update the trust region radius.
Ifr^k<η1, set hk+1=γ1hk;
Ifr^k>η2 and sk=hk, set hk+1=min⁡{γ2hk,h-};
Otherwise,hk+1=hk.Step 5.
Update Hessian approximation matrix to acquireG^k+1 when r^k≥η1; otherwise, keep the matrix unchanged. Set k≔k+1; go to Step 2.Remark 2.
The matrixG^k+1 can be updated by regularized stochastic BFGS formula [14] as follows:(8)G^k+1=G^k+r~kr~kTvkTr~k-G^kvkvkTG^kvkTG^kvk+δI,where δ>0 is a constant, vk=xk+1-xk, r~k=g^(xk+1,θ~k)-g^(xk,θ~k)-δvk denote the variable and corrected stochastic gradient variation at time k. The addition of the regularization term δI and the corrected stochastic gradient variation r~k avoid the near-singularity problems of more straightforward extensions.


## 3. Convergence Analysis

In this section, we prove that the iterative sequence generated by Algorithm1 is convergent. For the subsequent analysis, we define the instantaneous objective function associated with samples θ~=[θ1;⋯;θL] as(9)f^x,θ~≔1L∑l=1Lfx,θl.The definition of the instantaneous objective function f^(x,θ~) in association with the fact that Fx≔Eθ[f(x,θ)] implies(10)Fx=Eθf^x,θ~.In order to prove the global convergence, we make the following assumptions.A1 The instantaneous objective function f^(x,θ~) is continuously twice differentiable.A2 The level set L={xF(x)≤F(x0)} is bounded. Moreover, the function F(x) is bounded below in L.A3 Since the stochastic gradient g^(xk,θ~k) is an unbiased estimator of ∇F(xk) in the sense of Eθ[g^(xk,θ~k)xk]=∇F(xk), there exists a positive constant N such that, for all k, it holds(11)∇Fxk-g^xk,θ~k≤N.A4 There exists a ρ>0, such that G^k≤ρ for all k.As a consequence of A1, the functionF(x) is also continuously twice differentiable owing to the linearity of the expected operation and the expression in (10).Lemma 3.
Assume that A1, A3, and A4 hold; ifsk is a solution or approximate solution of subproblem (5), then we have(12)Eθq^k0-q^ksk∣xk≥C1Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk,where C1>0 is a constant and Eθ[q^(k)(0)-q^(k)(sk)xk] denotes the conditional expectation of predk=q^(k)(0)-q^(k)(sk) with giving xk.Proof.
Considering subproblem (5), the Cauchy point can be expressed by(13)skc=-Γkhkg^xk,θ~kg^xk,θ~k,where(14)Γk=1,else.min⁡g^xk,θ~k3hkg^xk,θ~kTG^kg^xk,θ~k,1,ifg^xk,θ~kTG^kg^xk,θ~k>0.
Case  1. When g^(xk,θ~k)TG^kg^(xk,θ~k)≤0, it follows from (14) that Γk=1; thus, we have(15)Eθq^k0-q^kskc∣xk=q^k0-Eθq^k-hkg^xk,θ~kg^xk,θ~k∣xk=Fxk-Fxk+Eθhkg^xk,θ~k-hk22g^xk,θ~k2g^xk,θ~kTG^kg^xk,θ~k∣xk=Eθhkg^xk,θ~k∣xk-Eθhk22g^xk,θ~k2g^xk,θ~kTG^kg^xk,θ~k∣xk≥Eθhkg^xk,θ~k∣xk≥Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.
Case  2. When g^(xk,θ~k)TG^kg^(xk,θ~k)>0 and g^(xk,θ~k)3/hkg^(xk,θ~k)TG^kg^(xk,θ~k)≤1, it follows from (14) that Γk=g^(xk,θ~k)3/hkg^(xk,θ~k)TG^kg^(xk,θ~k); then, we have(16)Eθq^k0-q^kskc∣xk=q^k0-Eθq^k-g^xk,θ~k3hkg^xk,θ~kTG^kg^xk,θ~khkg^xk,θ~k∣xk=Eθg^xk,θ~k4g^xk,θ~kTG^kg^xk,θ~k-g^xk,θ~kTG^kg^xk,θ~k2g^xk,θ~k4g^xk,θ~kTG^kg^xk,θ~k2∣xk=Eθg^xk,θ~k42g^xk,θ~kTG^kg^xk,θ~k∣xk≥Eθg^xk,θ~k42g^xk,θ~kG^kg^xk,θ~k∣xk=Eθg^xk,θ~k22G^k∣xk≥12Eθg^xk,θ~kmin⁡hk,g^xk,θ~kG^k∣xk≥12Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.
Case  3. When g^(xk,θ~k)TG^kg^(xk,θ~k)>0 and g^(xk,θ~k)3/hkg^(xk,θ~k)TG^kg^(xk,θ~k)>1, it follows from (14) that Γk=1, and we have(17)Eθq^k0-q^kskc∣xk=q^k0-Eθq^k-hkg^xk,θ~kg^xk,θ~k∣xk=Eθhkg^xk,θ~k-hk22g^xk,θ~k2g^xk,θ~kTG^kg^xk,θ~k∣xk≥Eθhkg^xk,θ~k-hk22g^xk,θ~k2g^xk,θ~k3hk∣xk=12Eθhkg^xk,θ~k∣xk≥12Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.
With the observation of these cases above, we can come to the conclusion that, for all, the Cauchy pointskc of (5) satisfies(18)Eθq^k0-q^kskc∣xk=q^k0-Eθq^kskc∣xk≥12Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.Assume that ske is the exact solution of subproblem (5); we can draw the following inequalities from the property of Cauchy point that(19)Eθq^kske∣xk≤Eθq^kskc∣xk.In combination with inequality (18) above, we can write(20)Eθq^k0-q^kske∣xk=q^k0-Eθq^kske∣xk≥q^k0-Eθq^kskc∣xk≥12Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.Further assume that sk is the approximate solution; then, there exists a constant C2∈(0,1] satisfying(21)Eθq^k0-q^ksk∣xk≥C2Eθq^k0-q^kske∣xk.If we set C1=1/2C2, then we have(22)Eθq^k0-q^ksk∣xk≥C2Eθq^k0-q^kskc∣xk≥C1Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.Lemma 3 implies that a sufficient decrement of the model q^(k)(s) on average is guaranteed under assumptions A1 and A4.Lemma 4.
If A1, A3, and A4 hold true, we have(23)EθFxk+sk-q^ksk∣xk≤NEθsk∣xk+12ρEθsk2∣xk+Eθcsksk∣xk,where c(sk) is a function of sk and decreases with the decrease of sk.Proof.
From Taylor theorem, we have(24)Fxk+sk=Fxk+∇FxkTsk+∫01∇Fxk+tsk-∇FxkTskdt.Then, observing (24) and the model q^(k)(s), we have(25)Fxk+sk-q^ksk=∇FxkTsk-g^xk,θ~kTsk-12skTG^ksk+∫01∇Fxk+tsk-∇FxkTskdt≤∇Fxk-g^xk,θ~ksk+12ρsk2+∫01∇Fxk+tsk-∇FxkTskdt.Taking expectation with xk given in both sides of (25) and applying (11), we can write(26)EθFxk+sk-q^ksk∣xk≤Eθ∇Fxk-g^xk,θ~ksk+12ρsk2+∫01∇Fxk+tsk-∇FxkTskdt∣xk≤NEθsk∣xk+12ρEθsk2∣xk+Eθcsksk∣xk.Lemma 5.
If the same assumptions as in Lemma3 hold, g^(xk,θ~k)≠0 and Δ~ is a small tolerance. Then, there exist infinite number of k which satisfy hk+1≥hk.Proof.
From the definition ofr^k=F(xk)-F(xk+sk)/q^(k)(0)-q^(k)(sk) and properties of expectation, we can define(27)Eθr^k∣xk=EθFxk-Fxk+sk∣xkEθq^k0-q^ksk∣xk.Further assume g^(xk,θ~k)≥ε>0 and there exists a positive constant ω∈[0,1] such that sk≤ωhk, and it follows from assumptions and Lemma 3 and Lemma 4 that(28)Eθr^k-1∣xk=Eθr^k∣xk-1=EθFxk-Fxk+sk∣xk-Eθq^k0-q^ksk∣xkEθq^k0-q^ksk∣xk≤EθFxk+sk-q^ksk∣xkEθq^k0-q^ksk∣xk≤NEθsk∣xk+1/2ρEθsk2∣xk+Eθcsksk∣xkC1Eθg^xk,θ~kmin⁡hk,g^xk,θ~k/ρ∣xk≤ωhkN+1/2ρωhk+cskC1εmin⁡hk,ε/ρ.We can select Δ~ small enough to satisfy hk≤Δ~≤ε/ρ, and(29)N+12ρωhk+csk≤1-η1C1εω.Thus, we have Eθ[r^k-1∣xk]≤1-η1; this result is shown that r^k is closer to 1 than η1 on average which implies by Algorithm 1 that there are infinite numbers of k which satisfy hk+1≥hk.Theorem 6.
Considering the new algorithm defined above, suppose that A1-A4 hold, and the sequence of iterates generated by Algorithm1 satisfies(30)limk→∞inf⁡g^xk,θ~k=0,over realizations of the random samples θ~kk=0∞.Proof.
Here we use a contradiction with Lemma5 to prove (30). For that purpose, assume that there exist ε>0 and a positive index set K such that for every k∈K we have g^(xk,θ~k)≥ε>0. At the same time, it is assumed that there are infinite number of successful iterations, and we can obtain the following inequality from Algorithm 1:(31)Fxk-Fxk+sk≥η1q^k0-q^ksk,taking expectation with xk given in both sides of (31) we have(32)EθFxk-Fxk+sk∣xk≥η1Eθq^k0-q^ksk∣xk≥η1C1Eθg^xk,θ~kmin⁡hk,g^xk,θ~kρ∣xk.≥η1C1εmin⁡hk,ερ.
Considering assumption A3, ask approaches infinity, the function F(x) bounded below implies that the right side of (32) is tending to zero; thus, we have(33)limk→∞⁡hk=0.But the result in (33) is in contradiction with Lemma 5. Further assume there are finite number of successful iterations; that is, the iterations are unsuccessful for large k; then, the trust region radius is reduced, i.e. limk→∞hk=0. In this case, we still get a contradiction with Lemma 5. Hence, we obtain (30).Remark 7.
We note that the stochastic gradientg^(xk,θ~k) is an unbiased estimator of ∇F(xk), i.e., Eθ[g^(xk,θ~k)∣xk]=∇F(xk). Thus, we can also obtain limk→∞inf⁡∇F(xk)=0 from the result in Theorem 6.


## 4. Numerical Experiments

In this section, we apply the proposed stochastic trust region method to solve convex and nonconvex problems with stochastic objectives and also compare it with SGD in terms of convergence time and central processing unit (CPU) runtime. Both methods will be tested on the following problems: problem 1 and 2 are convex stochastic optimization problems while problem 3 is nonconvex.


### 4.1. Example 1: Standard Quadratic Function

We use a stochastic quadratic objective function as a test case. In particular, consider a positive definite diagonal matrixA∈Rn×n, a vector b∈Rn, a random vector θ∈Rn, and diagonal matrix diag(θ) defined by θ. The function F(x) is defined as(34)Fx≔Eθfx,θ≔Eθ12xTA+Adiagθx+bTx.In (34), the vector θ is chosen uniformly at random from the n-dimensional box Θ=[-θ0,θ0]n for some given constant θ0<1. The linear term bTx is added so that the instantaneous functions f(x,θ) have different minima which are (almost surely) different from the minimum of the average function F(x). The quadratic term is chosen so that the condition number of F(x) is the condition number of A. Indeed, since Eθ[θ]=0, the average function in (34) can be written as F(x)=1/2xTAx+bTx. The parameter θ0 controls the variability of the instantaneous functions f(x,θ). For small θ0≈0, the instantaneous functions are close to each other and to the average function. For large θ0≈1, the instantaneous functions vary over a large range. Observe that we can write the optimum argument as x∗=A-1b for comparison against iterates xk. Further consider a given σ and study the convergence metric(35)Γ≔mink⁡k:xk-x∗x∗≤σ,which represents the time needed to achieve a given relative distance to optimality. To study the effect of the problem’s condition number, we generate instances of (34) by choosing b uniformly at random from the box [0,1]n and the matrix A as diagonal with elements aii uniformly drawn from the discrete set 1,10-1,…,10-ξ. This choice of A yields problems with condition number 10ξ. Representative runs of stochastic trust region method and SGD for n=50, θ0=0.5 and ξ=3 are shown in Figure 1. For the stochastic trust region method and SGD run, the stochastic gradients g^(x,θ~) in (4) are computed as an average of L=10 realizations. Moreover, we set γ1=0.95, γ2=2, and h¯=50 for our method due to the randomness of instantaneous functions f(x,θ~) and the large condition number.Figure 1
Convergence of the stochastic trust region method and SGD for function (35). Relative distance to opotimality xk-x∗/x∗ is shown with respect to the number of the stochastic function processed. For our method, the number of iterations required to achieve a certain accuracy is smaller than the corresponding number for SGD.
(a)
1 (b)
2As expected for a problem with large condition number, the condition number ofF(x) is 103 since we are using ξ=3, and the stochastic trust region method is much faster than SGD. After k=1892, the distance to optimality for the SGD iterate is xk-x∗/x∗=4.2×10-3. Comparable accuracy for stochastic trust region method is achieved after k=105 iterations. Conversely, upon processing k=2×103 random functions, our method achieves accuracy xk-x∗/x∗=1.7×10-3.


### 4.2. Example 2: Extended Powell Singular Function

Powell singular function (PSF) is also known as Powell Quartic Function. The Hessian matrix at minimizer is doubly singular; thus, PSF is a severe test problem. To analyse the numerical performance of our algorithm, for a random vectorθ∈Rn/4, a variable x∈Rn, we develop stochastic Extended PSF as the test function as follows [18]:(36)Fx=Eθfx,θ=Eθ∑i=1n/41+θix4i-3+10x4i-22+5x4i-1-x4i2+x4i-2-2x4i-14+10x4i-3-x4i4,where θ is selected uniformly at random from the box Θ=[-θ0,θ0]n/4 with θ0=0.5. Indeed, observe that since Eθ[θ]=0, the average function F(x) in (36) can be written as(37)Fx=∑i=1n/4x4i-3+10x4i-22+5x4i-1-x4i2+x4i-2-2x4i-14+10x4i-3-x4i4.Standard starting point is x0=(3,-1,0,1,…,3,-1,0,1)T and the Hessian matrix at the standard starting point is nonsingular. The average function F(x) in (36) is convex and the unique unconstrained minimizer is x∗=(0,0,…,0)T with F(x∗)=0. In Figure 2, for n=40, L=5, specially γ1=0.5, γ2=2, and h¯=5, since the minimizer is x∗=(0,0,…,0)T, we plot the variation of distance to optimality instead of relative distance with the increase of the number of iterations.Figure 2
Convergence of the stochastic trust region method and SGD for function (37). Distance to opotimality xk-x∗ is shown with respect to the number k of stochastic function processed.
(a)
1 (b)
2As the number of iterations increases, we can observe that SGD achieves accuracyxk=8.9×10-3 after k=2000 iterations and has a slow convergence to optimality. However, our method needs only k=82 iterations to achieve the same accuracy. Furthermore, after k=285 iterations, stochastic trust region method has xk=5.2×10-4 with F(xk)=4.8532×10-7. The result indicates that our method has adequate accuracy and a faster convergence rate than SGD for severe convex problems.


### 4.3. Example 3: Extended Rosenbrock Function

In mathematical optimization, the Rosenbrock function is frequently used in nonconvex optimization problems as a performance test problem for optimization algorithms, introduced by Howard H. Rosenbrock in 1960. The global minimum lies in a narrow, long, parabolic shaped flat valley. It is trivial to find the valley; however, convergence to the global minimum is difficult. Here we use a stochastic version of the Extended Rosenbrock function as a test case [18]. Specifically, consider a random vector θ∈Rn/2, a variable x∈Rn. The function F(x) is defined as(38)Fx≔Eθfx,θ≔Eθ∑i=1n/21+θi100x2i-12-x2i2+1-x2i-12.Similarly, we choose θ uniformly at random from the box Θ=[-θ0,θ0]n/2 and set θ0=0.5. In fact, from the observation Eθ[θ]=0, we can write the average function as(39)Fx=∑i=1n/2100x2i-12-x2i2+1-x2i-12.The function in (39) has a minimizer x∗=[1,1,…,1]T with F(x∗)=0. For the purpose of this test problem, we choose n=50, L=5. In particular, set γ1=0.5, γ2=2, and h¯=5 for our method. The relative distances in (38) after running of stochastic trust region method and SGD are provided in Figure 3. After k=99 iterations, stochastic trust region method achieves accuracy xk-x∗/x∗=1.94×10-4, and correspondingly SGD achieves such accuracy after k=1523 iterations.Figure 3
Convergence of the stochastic trust region method and SGD for function (39). Relative distance to opotimality xk-x∗/x∗ is shown with respect to the number k of the stochastic function processed. For our method, the number of iterations required to achieve a certain accuracy is smaller than the corresponding number for SGD.
(a)
1 (b)
2Since each iteration of stochastic trust region method is more complex than SGD, we also compare the performances in terms of central processing unit (CPU) runtime required to achieve accuracy10-2. We report the average runtime of stochastic trust region method and SGD for problems above in Table 1; the value of parameters is the same as above. In particular, we call stochastic trust region method STR for short in the table. As we can see that our method enjoys a significant improvement in CPU runtime than SGD for both convex and nonconvex problems.Table 1
CPU runtime of STR and SGD.
STR SGD Standard Quadratic function 0.22s 0.7s  Extended Powell Singular function 0.18s 0.64s  Extended Rosenbrock function 0.25s 0.73s


## 4.1. Example 1: Standard Quadratic Function

We use a stochastic quadratic objective function as a test case. In particular, consider a positive definite diagonal matrixA∈Rn×n, a vector b∈Rn, a random vector θ∈Rn, and diagonal matrix diag(θ) defined by θ. The function F(x) is defined as(34)Fx≔Eθfx,θ≔Eθ12xTA+Adiagθx+bTx.In (34), the vector θ is chosen uniformly at random from the n-dimensional box Θ=[-θ0,θ0]n for some given constant θ0<1. The linear term bTx is added so that the instantaneous functions f(x,θ) have different minima which are (almost surely) different from the minimum of the average function F(x). The quadratic term is chosen so that the condition number of F(x) is the condition number of A. Indeed, since Eθ[θ]=0, the average function in (34) can be written as F(x)=1/2xTAx+bTx. The parameter θ0 controls the variability of the instantaneous functions f(x,θ). For small θ0≈0, the instantaneous functions are close to each other and to the average function. For large θ0≈1, the instantaneous functions vary over a large range. Observe that we can write the optimum argument as x∗=A-1b for comparison against iterates xk. Further consider a given σ and study the convergence metric(35)Γ≔mink⁡k:xk-x∗x∗≤σ,which represents the time needed to achieve a given relative distance to optimality. To study the effect of the problem’s condition number, we generate instances of (34) by choosing b uniformly at random from the box [0,1]n and the matrix A as diagonal with elements aii uniformly drawn from the discrete set 1,10-1,…,10-ξ. This choice of A yields problems with condition number 10ξ. Representative runs of stochastic trust region method and SGD for n=50, θ0=0.5 and ξ=3 are shown in Figure 1. For the stochastic trust region method and SGD run, the stochastic gradients g^(x,θ~) in (4) are computed as an average of L=10 realizations. Moreover, we set γ1=0.95, γ2=2, and h¯=50 for our method due to the randomness of instantaneous functions f(x,θ~) and the large condition number.Figure 1
Convergence of the stochastic trust region method and SGD for function (35). Relative distance to opotimality xk-x∗/x∗ is shown with respect to the number of the stochastic function processed. For our method, the number of iterations required to achieve a certain accuracy is smaller than the corresponding number for SGD.
(a)
1 (b)
2As expected for a problem with large condition number, the condition number ofF(x) is 103 since we are using ξ=3, and the stochastic trust region method is much faster than SGD. After k=1892, the distance to optimality for the SGD iterate is xk-x∗/x∗=4.2×10-3. Comparable accuracy for stochastic trust region method is achieved after k=105 iterations. Conversely, upon processing k=2×103 random functions, our method achieves accuracy xk-x∗/x∗=1.7×10-3.


## 4.2. Example 2: Extended Powell Singular Function

Powell singular function (PSF) is also known as Powell Quartic Function. The Hessian matrix at minimizer is doubly singular; thus, PSF is a severe test problem. To analyse the numerical performance of our algorithm, for a random vectorθ∈Rn/4, a variable x∈Rn, we develop stochastic Extended PSF as the test function as follows [18]:(36)Fx=Eθfx,θ=Eθ∑i=1n/41+θix4i-3+10x4i-22+5x4i-1-x4i2+x4i-2-2x4i-14+10x4i-3-x4i4,where θ is selected uniformly at random from the box Θ=[-θ0,θ0]n/4 with θ0=0.5. Indeed, observe that since Eθ[θ]=0, the average function F(x) in (36) can be written as(37)Fx=∑i=1n/4x4i-3+10x4i-22+5x4i-1-x4i2+x4i-2-2x4i-14+10x4i-3-x4i4.Standard starting point is x0=(3,-1,0,1,…,3,-1,0,1)T and the Hessian matrix at the standard starting point is nonsingular. The average function F(x) in (36) is convex and the unique unconstrained minimizer is x∗=(0,0,…,0)T with F(x∗)=0. In Figure 2, for n=40, L=5, specially γ1=0.5, γ2=2, and h¯=5, since the minimizer is x∗=(0,0,…,0)T, we plot the variation of distance to optimality instead of relative distance with the increase of the number of iterations.Figure 2
Convergence of the stochastic trust region method and SGD for function (37). Distance to opotimality xk-x∗ is shown with respect to the number k of stochastic function processed.
(a)
1 (b)
2As the number of iterations increases, we can observe that SGD achieves accuracyxk=8.9×10-3 after k=2000 iterations and has a slow convergence to optimality. However, our method needs only k=82 iterations to achieve the same accuracy. Furthermore, after k=285 iterations, stochastic trust region method has xk=5.2×10-4 with F(xk)=4.8532×10-7. The result indicates that our method has adequate accuracy and a faster convergence rate than SGD for severe convex problems.


## 4.3. Example 3: Extended Rosenbrock Function

In mathematical optimization, the Rosenbrock function is frequently used in nonconvex optimization problems as a performance test problem for optimization algorithms, introduced by Howard H. Rosenbrock in 1960. The global minimum lies in a narrow, long, parabolic shaped flat valley. It is trivial to find the valley; however, convergence to the global minimum is difficult. Here we use a stochastic version of the Extended Rosenbrock function as a test case [18]. Specifically, consider a random vector θ∈Rn/2, a variable x∈Rn. The function F(x) is defined as(38)Fx≔Eθfx,θ≔Eθ∑i=1n/21+θi100x2i-12-x2i2+1-x2i-12.Similarly, we choose θ uniformly at random from the box Θ=[-θ0,θ0]n/2 and set θ0=0.5. In fact, from the observation Eθ[θ]=0, we can write the average function as(39)Fx=∑i=1n/2100x2i-12-x2i2+1-x2i-12.The function in (39) has a minimizer x∗=[1,1,…,1]T with F(x∗)=0. For the purpose of this test problem, we choose n=50, L=5. In particular, set γ1=0.5, γ2=2, and h¯=5 for our method. The relative distances in (38) after running of stochastic trust region method and SGD are provided in Figure 3. After k=99 iterations, stochastic trust region method achieves accuracy xk-x∗/x∗=1.94×10-4, and correspondingly SGD achieves such accuracy after k=1523 iterations.Figure 3
Convergence of the stochastic trust region method and SGD for function (39). Relative distance to opotimality xk-x∗/x∗ is shown with respect to the number k of the stochastic function processed. For our method, the number of iterations required to achieve a certain accuracy is smaller than the corresponding number for SGD.
(a)
1 (b)
2Since each iteration of stochastic trust region method is more complex than SGD, we also compare the performances in terms of central processing unit (CPU) runtime required to achieve accuracy10-2. We report the average runtime of stochastic trust region method and SGD for problems above in Table 1; the value of parameters is the same as above. In particular, we call stochastic trust region method STR for short in the table. As we can see that our method enjoys a significant improvement in CPU runtime than SGD for both convex and nonconvex problems.Table 1
CPU runtime of STR and SGD.
STR SGD Standard Quadratic function 0.22s 0.7s  Extended Powell Singular function 0.18s 0.64s  Extended Rosenbrock function 0.25s 0.73s


## 5. Conclusions

In this paper, we propose a stochastic trust region method and show that the new algorithm is convergent for solving unconstrained minimization problems with stochastic objectives. Based on the trust region framework and the BFGS update, our method can deal with convex optimization problems with ill-conditioned objective functions as well as nonconvex optimization problems. With careful analysis, we are able to show that our method is convergent. Numerical results illustrate that the method can efficiently solve the given test problems. Therefore, the new method is potentially efficient and thus paves the way towards developing concrete algorithms for specific applications.



---

*Source: 8095054-2019-02-03.xml*