Anatoli Juditsky
email: anatoli.juditsky@imag.fr

Arkadi Nemirovski
email: nemirovs@isye.gatech.edu

Claire Tauvel
email: claire.tauvel@imag.fr

Solving variational inequalities with Stochastic Mirror-Prox algorithm

Keywords: variational inequalities with monotone operators, stochastic convex-concave saddle-point problem, large scale stochastic approximation, reduced complexity algorithms for convex optimization

In this paper we consider iterative methods for stochastic variational inequalities (s.v.i.) with monotone operators. Our basic assumption is that the operator possesses both smooth and nonsmooth components. Further, only noisy observations of the problem data are available. We develop a novel Stochastic Mirror-Prox (SMP) algorithm for solving s.v.i. and show that with the convenient stepsize strategy it attains the optimal rates of convergence with respect to the problem parameters. We apply the SMP algorithm to Stochastic composite minimization and describe particular applications to Stochastic Semidefinite Feasibility problem and deterministic Eigenvalue minimization.

Introduction

Variational inequalities with monotone operators form a convenient framework for unified treatment (including algorithmic design) of problems with "convex structure," like convex minimization, convexconcave saddle point problems and convex Nash equilibrium problems. In this paper we utilize this framework to develop first order algorithms for stochastic versions of the outlined problems, where the precise first order information is replaced with its unbiased stochastic estimates. This situation arises naturally in convex Stochastic Programming, where the precise first order information is unavailable (see examples in section 4). In some situations, e.g. those considered in [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF]Section 3.3] and in Section 4. [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF], where passing from available, but relatively computationally expensive precise first order information to its cheap stochastic estimates allows to accelerate the solution process, with the gain from randomization growing progressively with problem's sizes.

Our "unifying framework" is as follows. Let Z be a convex compact set in Euclidean space E with inner product •, • , • be a norm on E (not necessarily the one associated with the inner product), and F : Z → E be a monotone mapping:

∀(z, z ′ ∈ Z) : F (z) -F (z ′ ), z -z ′ ≥ 0 (1)
We are interested to approximate a solution to the variational inequality (v.i.)

find z * ∈ Z : F (z), z * -z ≤ 0 ∀z ∈ Z (2)
associated with Z, F . Note that since F is monotone on Z, the condition in ( 2) is implied by F (z * ), zz * ≥ 0 for all z ∈ Z, which is the standard definition of a (strong) solution to the v.i. associated with Z, F . The inverse -a solution to v.i. as defined by (2) (a "weak" solution) is a strong solution as well -also is true, provided, e.g., that F is continuous. An advantage of the concept of weak solution is that such a solution always exists under our assumptions (F is well defined and monotone on a convex compact set Z).

We quantify the inaccuracy of a candidate solution z ∈ Z by the error Err vi (z) := max u∈Z F (u), z -u ; [START_REF] Grigoriadis | A sublinear-time randomized approximation algorithm for matrix games[END_REF] note that this error is always ≥ 0 and equals zero iff z is a solution to [START_REF] Ben-Tal | Non-Euclidean restricted memory level method for large-scale convex optimization[END_REF].

In what follows we impose on F , aside of the monotonicity, the requirement

∀(z, z ′ ∈ Z) : F (z) -F (z ′ ) * ≤ L z -z ′ + M (4)
with some known constants L ≥ 0, M ≥ 0. From now on,

ξ * = max z: z ≤1 ξ, z (5)
is the norm conjugate to • . We are interested in the case where ( 2) is solved by an iterative algorithm based on a stochastic oracle representation of the operator F (•). Specifically, when solving the problem, the algorithm acquires information on F via subsequent calls to a black box ("stochastic oracle", SO). At the ith call, i = 0, 1, ..., the oracle gets as input a search point z i ∈ Z (this point is generated by the algorithm on the basis of the information accumulated so far) and returns the vector Ξ(z i , ζ i ), where {ζ i ∈ R N } ∞ i=1 is a sequence of i.i.d. (and independent of the queries of the algorithm) random variables. We suppose that the Borel function Ξ(z, ζ) is such that

∀z ∈ Z : E {Ξ(z, ζ 1 )} = F (z), E Ξ(z, ζ i ) -F (z) 2 * ≤ σ 2 . (6)
We call a monotone v.i. [START_REF] Azuma | Weighted sums of certain dependent random variables[END_REF], augmented by a stochastic oracle (SO), a stochastic monotone v.i. (s.v.i.). To motivate our goal, let us start with known results [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF] on the limits of performance of iterative algorithms for solving large-scale stochastic monotone v.i.'s. To "normalize" the situation, assume that Z is the unit Euclidean ball in E = R n and that n is large. In this case, the accuracy after t steps of any algorithm for solving v.i.'s cannot be better than O(1) L t + M+σ √ t . In other words, for a properly chosen positive absolute constant C, for every number of steps t, all large enough values of n and any algorithm B for solving s.v.i.'s on the unit ball of R n , one can point out a monotone s.v.i. satisfying (4), [START_REF] Korpelevich | Extrapolation gradient methods and relation to Modified Lagrangeans[END_REF] and such that the expected error of the approximate solution zt generated by B after t steps , applied to such s.v.i., is at least c L t + M+σ √ t for some c > 0. To the best of our knowledge, no existing algorithm allows to achieve, uniformly in the dimension, this convergence rate. In fact, the "best approximations" available are given by Robust Stochastic Approximation (see [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF] and references therein) with the guaranteed rate of convergence O(1) L+M+σ √ t and extra-gradient-type algorithms for solving deterministic monotone v.i.'s with Lipschitz continuous operators (see [START_REF] Nemirovski | Prox-method with rate of convergence O(1/t) for variational inequalities with Lipschitz continuous monotone operators and smooth convex-concave saddle point problems[END_REF][START_REF] Nesterov | Smooth minimization of non-smooth functions[END_REF][START_REF] Nesterov | Excessive gap technique in nonsmooth convex minimization[END_REF][START_REF] Nesterov | Dual extrapolation and its applications to solving variational inequalities and related problems[END_REF]), which attain the accuracy O(1) L t in the case of M = σ = 0 or O(1) M √ t when L = σ = 0. The goal of this paper is to demonstrate that a specific Mirror-Prox algorithm [START_REF] Nemirovski | Prox-method with rate of convergence O(1/t) for variational inequalities with Lipschitz continuous monotone operators and smooth convex-concave saddle point problems[END_REF] for solving monotone v.i.'s with Lipschitz continuous operators can be extended onto monotone s.v.i.'s to yield, uniformly in the dimension, the optimal rate of convergence O(1) L t + M+σ √ t . We present the corresponding extension and investigate it in detail: we show how the algorithm can be "tuned" to the geometry of the s.v.i. in question and derive bounds for the probability of large deviations of the resulting error. We also present a number of applications where the specific structure of the rate of convergence indeed "makes a difference."

The main body of the paper is organized as follows: in Section 2, we describe several special cases of monotone v.i.'s we are especially interested in (convex Nash equilibria, convex-concave saddle point problems, convex minimization). We single out these special cases since here one can define a useful "functional" counterpart Err N (•) of the just defined error Err vi (•); both Err N and Err vi will participate in our subsequent efficiency estimates. Our main development -the Stochastic Mirror Prox (SMP) algorithm -is presented in Section 3. where we also provide some general results about its performance. Then in Section 4 we present SMP for Stochastic composite minimization and discuss its applications to Stochastic Semidefinite Feasibility problem and Eigenvalue minimization. All technical proofs are collected in the appendix.

Notations. In the sequel, lowercase Latin letters denote vectors (and sometimes matrices). Script capital letters, like E, Y, denote Euclidean spaces; the inner product in such a space, say, E, is denoted by •, • E (or merely •, • , when the corresponding space is clear from the context). Linear mappings from one Euclidean space to another, say, from E to F , are denoted by boldface capitals like A (there are also some reserved boldface capitals, like E for expectation, R k for the k-dimensional coordinate space, and S k for the space of k × k symmetric matrices). A * stands for the conjugate to mapping A: if A : E → F , then A * : F → E is given by the identity f, Ae F = A * f, e E for f ∈ F , e ∈ E. When both the origin and the destination space of a linear map, like A, are the standard coordinate spaces, the map is identified with its matrix A, and A * is identified with A T . For a norm • on E, • * stands for the conjugate norm, see [START_REF] Juditsky | Large Deviations of Vector-valued Martingales in 2-Smooth Normed Spaces E-print[END_REF].

For Euclidean spaces E 1 , ..., E m , E = E 1 × ... × E m denotes their Euclidean direct product, so that a vector from E is a collection u = [u 1 ; ...; u m ] ("MATLAB notation") of vectors u ℓ ∈ E ℓ , and u, v E = ℓ u ℓ , v ℓ E ℓ . Sometimes we allow ourselves to write (u 1 , ..., u m ) instead of [u 1 ; ...; u m ].

Preliminaries and Problem of interest 2.1 Nash v.i.'s and functional error

In the sequel, we shall be especially interested in a special case of v.i. (2) -in a Nash v.i. coming from a convex Nash Equilibrium problem, and in the associated functional error measure. The Nash Equilibrium problem can be described as follows: there are m players, the ith of them choosing a point z i from a given set Z i . The loss of the ith player is a given function φ i (z) of the collection z = (z 1 , ..., z m ) ∈ Z = Z 1 × ... × Z m of players' choices. With slight abuse of notation, we use for φ i (z) also the notation φ i (z i , z i ), where z i is the collection of choices of all but the ith players. Players are interested to minimize their losses, and Nash equilibrium z is a point from Z such that for every i the function φ i (z i , z i ) attains its minimum in z i ∈ Z i at z i = z i (so that in the state z no player has an incentive to change his choice, provided that the other players stick to their choices).

We call a Nash equilibrium problem convex, if for every i, Z i is a compact convex set, φ i (z i , z i ) is a Lipschitz continuous function convex in z i and concave in z i , and the function Φ(z) = m i=1 φ i (z) is convex. It is well known (see, e.g., [START_REF] Nemirovski | Accuracy certificates for computational problems with convex structure[END_REF]) that setting

F (z) = F 1 (z); . . . ; F m (z) , F i (z) ∈ ∂ zi φ i (z i , z i ), i = 1, ..., m
where ∂ zi φ i (z i , z i ) is the subdifferential of the convex function φ i (•, z i ) at a point z i , we get a monotone operator such that the solutions to the corresponding v.i. ( 2) are exactly the Nash equilibria. Note that since φ i are Lipschitz continuous, the associated operator F can be chosen to be bounded. For this v.i. one can consider, along with the v.i.-accuracy measure Err vi (z), the functional error measure

Err N (z) = m i=1 φ i (z) -min wi∈Zi φ i (w i , z i )
This accuracy measure admits a transparent justification: this is the sum, over the players, of the incentives for a player to change his choice given that other players stick to their choices.

Special case: saddle points

An important by its own right particular case of Nash Equilibrium problem is a zero sum game, where m = 2 and Φ(z) ≡ 0 (i.e., φ 2 (z) ≡ -φ 1 (z)). The convex case of this problem corresponds to the situation when φ(z 1 , z 2 ) ≡ φ 1 (z 1 , z 2 ) is a Lipschitz continuous function which is convex in z 1 ∈ Z 1 and concave in z 2 ∈ Z 2 , the Nash equilibria are exactly the saddle points (min in z 1 , max in z 2 ) of φ on Z 1 × Z 2 , and the functional error becomes

Err N (z 1 , z 2 ) = max (u1,u2)∈Z [φ(z 1 , u 1 ) -φ(u 2 , z 2 )] .
Recall that the convex-concave saddle point problem min z1∈Z1 max z2∈Z2 φ(z 1 , z 2 ) gives rise to the "primaldual" pair of convex optimization problems (P ) : min

z1∈Z1 φ(z 1 ), (D) : max z2∈Z2 φ(z 2 ), where φ(z 1 ) = max z2∈Z2 φ(z 1 , z 2 ), φ(z 2 ) = min z1∈Z1 φ(z 1 , z 2 ).
The optimal values Opt(P ) and Opt(D) in these problems are equal, the set of saddle points of φ (i.e., the set of Nash equilibria of the underlying convex Nash problem) is exactly the direct product of the optimal sets of (P ) and (D), and Err N (z 1 , z 2 ) is nothing but the sum of non-optimalities of z 1 , z 2 considered as approximate solutions to respective optimization problems:

Err N (z 1 , z 2 ) = φ(z 1 ) -Opt(P ) + Opt(D) -φ(z 2 ) .
In the sequel, we refer to the v.i. ( 2) coming from a convex Nash Equilibrium problem as Nash v.i., and to the just outlined particular case as the Saddle Point v.i. It is easy to verify that in the Saddle Point case the functional error Err N (z) is ≤ Err vi (z); this is not necessary so for a general Nash v.i.

Composite Optimization problem and its saddle point reformulation

While the algorithm we intend to develop is applicable to a general-type stochastic v.i. with monotone operator, the applications to be considered in this paper deal with (saddle point reformulation of) convex composite optimization problem (cf. [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF]). As the simplest motivating example, one can keep in mind the minimax problem

min x∈X max 1≤i≤m φ ℓ (x), (7)
where X ⊂ R n is a convex compact set and φ ℓ (x) are Lipschitz continuous convex functions on X. This problem can be rewritten as the saddle point problem

min x∈X max y∈Y φ(x, y) := m ell=1 y ℓ φ ℓ (x), (8)
where Y = {y ∈ R m + : m ℓ=1 y ℓ = 1} is the standard simplex. The advantages of the saddle point reformulation are twofold. First, when all φ ℓ are smooth, so is φ, in contrast to the objective in (7) which typically is nonsmooth; this makes the saddle point reformulation better suited for processing by first order algorithms. Starting with the breakthrough paper of Nesterov [START_REF] Nesterov | Smooth minimization of non-smooth functions[END_REF], this phenomenon, in its general form, is utilized in the fastest known so far first order algorithms for "well-structured" nonsmooth convex programs. Second, in the stochastic case, stochastic oracles providing unbiased estimates of the first order information on φ i oracles, while not induce a similar oracle for the objective of [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF], do induce such an oracle for the v.i. associated with [START_REF] Nemirovski | Prox-method with rate of convergence O(1/t) for variational inequalities with Lipschitz continuous monotone operators and smooth convex-concave saddle point problems[END_REF] and thus make the problem amenable to first order algorithms.

Composite minimization problem.

In this paper, we focus on a substantial extension of the minimax problem [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF], namely, on a Composite minimization problem min x∈X φ(x) := Φ(φ 1 (x), ..., φ m (x)), [START_REF] Lu | Large-Scale Semidefinite Programming via Saddle Point Mirror-Prox Algorithm[END_REF] where the inner functions φ ℓ (•) are vector-valued, and the outer function Φ is real-valued. We are about to impose structural restrictions which allow to reformulate the problem as a "good" convex-concave saddle point problem, specifically, as follows:

A. X ⊂ X is a convex compact;

B. φ ℓ (x) : X → E ℓ , 1 ≤ ℓ ≤ m, are Lipschitz continuous mappings taking values in Euclidean spaces E ℓ equipped with closed convex cones K ℓ . We assume φ ℓ to be K ℓ -convex, meaning that for any

x, x ′ ∈ X, λ ∈ [0, 1], φ ℓ (λx + (1 -λ)x ′ ) ≤ K ℓ λφ ℓ (x) + (1 -λ)φ ℓ (x ′ ),
where the notation a

≤ K b ⇔ b ≥ K a means that b -a ∈ K. C. Φ(•) is a convex function on E = E 1 × ... × E m
given by the Fenchel-type representation

Φ(u 1 , ..., u m ) = max y∈Y m ℓ=1 u ℓ , A ℓ y + b ℓ E ℓ -Φ * (y) , (10)
for

u ℓ ∈ E ℓ , 1 ≤ ℓ ≤ m. Here -Y ⊂ Y is a convex compact set, -the affine mappings y → A ℓ y + b ℓ : Y → E ℓ are such that A ℓ y + b ℓ ∈ K * ℓ for
φ ℓ (x), A ℓ y + b ℓ E ℓ -Φ * (y) (11)
and the cost function in the latter problem is Lipschitz continuous and convex-concave due to the convexity of Φ * , K ℓ -convexity of φ ℓ (•) and the condition A ℓ y + b ℓ ∈ K * ℓ whenever y ∈ Y . The associated Nash v.i. is given by the domain Z = X × Y and the monotone mapping

F (z) ≡ F (x, y) = m ℓ=1 [φ ′ ℓ (x)] * [A ℓ y + b ℓ ]; - m ℓ=1 A * ℓ φ ℓ (x) + Φ ′ * (y) . (12)
Same as in the case of minimax problem [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF], the advantage of the saddle point reformulation (11) of ( 9) is that, independently of whether Φ is smooth, φ is smooth whenever all φ ℓ are so. Another advantage, instrumental in the stochastic case, is that F is linear in φ ℓ (•), so that stochastic oracles providing unbiased estimates of the first order information on φ ℓ induce straightforwardly an unbiased SO for F .

Example: Matrix Minimax problem

For 1 ≤ ℓ ≤ m, let E ℓ = S p ℓ be the space of symmetric p ℓ × p ℓ matrices equipped with the Frobenius inner product A, B F = Tr(AB), and let K ℓ be the cone S p ℓ + of symmetric positive semidefinite p ℓ × p ℓ matrices. Now let X ⊂ X be a convex compact set, and φ ℓ : X → E ℓ be S where P jℓ are given p ℓ × q j matrices, and λ max (A) is the maximal eigenvalue of a symmetric matrix A.

Observing that for a symmetric q × q matrix A one has

λ max (A) = max S∈Sq Tr(AS)
where S q = {S ∈ S q + : Tr(S) = 1}. When denoting by Y the set of all symmetric positive semidefinite block-diagonal matrices y = Diag{y 1 , ..., y k } with unit trace and diagonal blocks y j of sizes q j × q j , we can represent (P ) in the form of ( 9), [START_REF] Nemirovski | Accuracy certificates for computational problems with convex structure[END_REF] Observe that in the simplest case of k = m, p j = q j , 1 ≤ j ≤ m and P jℓ equal to I p for j = ℓ and to 0 otherwise, the problem becomes

min x∈X max 1≤ℓ≤m λ max (φ ℓ (x)) . (13)
If, in addition, p j = q j = 1 for all j, we arrive at the convex minimax problem [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF].

Illustration: Semidefinite Feasibility problem. With X and φ ℓ as above, consider the Semidefinite Feasibility problem find x ∈ X :

ψ ℓ (x) 0, 1 ≤ ℓ ≤ m. (S)
Choosing somehow scaling factors β ℓ > 0 and setting φ ℓ (x) = β ℓ ψ ℓ (x), we can pose (S) as the Matrix Minimax problem min

x∈X max 1≤ℓ≤m λ max (φ ℓ (x)); (S)

is solvable if and only if the optimal value in the Matrix

Minimax problem is ≤ 0.

Stochastic Mirror-Prox algorithm

We are about to present the stochastic version of the deterministic Mirror-Prox algorithm proposed in [START_REF] Nemirovski | Prox-method with rate of convergence O(1/t) for variational inequalities with Lipschitz continuous monotone operators and smooth convex-concave saddle point problems[END_REF]. The method is aimed at solving v.i. ( 2) associated with a convex compact set Z ⊂ E and a bounded monotone operator F : Z → E. In contrast to the original version of the method, below we allow for errors when computing the values of F -we assume that given a point z ∈ Z, we can compute an approximation (perhaps random) F (z) ∈ E of F (z).

Algorithm's setup

The setup for SMP (Stochastic Mirror Prox algorithm) is given by 1. a norm • on E; • * stands for the conjugate norm, see (5);

2. a distance-generating function (d.-g.f.) for Z, that is, a continuous convex function ω(•) : Z → R such that
(a) with Z o being the set of all points z ∈ Z such that the subdifferential ∂ω(z) of ω(•) at z is nonempty, ∂ω(•) admits a continuous selection on Z o : there exists a continuous on

Z o vector-valued function ω ′ (z) such that ω ′ (z) ∈ ∂ω(z) for all z ∈ Z o ; (b) ω(•) is strongly convex, modulus 1, w.r.t. the norm • : ∀(z, z ′ ∈ Z o ) : ω ′ (z) -ω ′ (z ′ ), z -z ′ ≥ z -z ′ 2 . ( 14
)
In order for the SMP associated with the outlined setup to be practical, ω(•) and Z should "fit" each other, meaning that one can easily solve problems of the form

min z∈Z [ω(z) + e, z ] , e ∈ E. (15)
The prox-function associated with a setup for SMP is defined as

V (z, u) = ω(u) -ω(z) -ω ′ (z), u -z : Z o × Z → R + . We set (a) Θ(z) = max u∈Z V (z, u) [z ∈ Z o ]; (b) z c = argmin Z ω(z); (c) Ω = 2Θ(z c ). ( 16
)
Note that z c is well defined (since Z is a convex compact set and ω(•) is continuous and strongly convex on Z) and belongs to Z o (since 0 ∈ ∂ω(z c )). Note also that due to the strong convexity of ω and the origin of z c we have

∀(u ∈ Z) : 1 2 u -z c 2 ≤ Θ(z c ) ≤ max z∈Z ω(z) -ω(z c ); (17)
in particular we see that

Z ⊂ {z : z -z c ≤ Ω}. ( 18
)
Prox-mapping. Given a setup for SMP and a point z ∈ Z o , we define the associated prox-mapping as

P z(ξ) = argmin u∈Z {ω(u) + ξ -ω ′ (z), u } ≡ argmin u∈Z {V (z, u) + ξ, u } : E → Z o .
Since S is compact and ω(•) is continuous and strongly convex on ZThis mapping is clearly well defined.

Basic SMP setups

We illustrate the just-defined notions with three basic examples.

Example 1: Euclidean setup. Here E is R N with the standard inner product, • 2 is the standard Euclidean norm on R N (so that • * = • ) and ω(z) = 1 2 z T z (i.e., Z o = Z). Assume for the sake of simplicity that 0 ∈ Z. Then z c = 0 and Ω = max z∈Z z 2 2 . The prox-function and the prox-mapping are given by

V (z, u) = 1 2 z -u 2 2 , P z (ξ) = argmin u∈Z (z -ξ) -u 2 . Example 2: Simplex setup. Here E is R N , N > 1, with the standard inner product, z = z 1 := N j=1 |z j | (so that ξ * = max j |ξ j |), Z is a closed convex subset of the standard simplex D N = {z ∈ R N : z ≥ 0, N j=1 z j = 1}
containing its barycenter, and ω(z) = N j=1 z j ln z j is the entropy. Then

Z o = {z ∈ Z : z > 0} and ω ′ (z) = [1 + ln z 1 ; ...; 1 + ln z N ], z ∈ Z o .
It is easily seen (see, e.g., [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF]) that here

z c = [1/N ; ...; 1/N ], Ω ≤ 2 ln(N )
(the latter inequality becomes equality when Z contains a vertex of D N ). The prox-function is

V (z, u) = N j=1 u j ln(u j /z j ),
and the prox-mapping is easy to compute when Z = D N :

(P z (ξ)) j = N i=1 z i exp{-ξ i } -1 z j exp{-ξ j }.
Example 3: Spectahedron setup. This is the "matrix analogy" of the Simplex setup. Specifically, now E is the space of N ×N block-diagonal symmetric matrices, N > 1, of a given block-diagonal structure equipped with the Frobenius inner product a, b F = Tr(ab) and the trace norm

|a| 1 = N i=1 |λ i (a)|
, where λ 1 (a) ≥ ... ≥ λ N (a) are the eigenvalues of a symmetric N × N matrix a; the conjugate norm |a| ∞ is the usual spectral norm (the largest singular value) of a. Z is assumed to be a closed convex subset of the spectahedron S = {z ∈ E : z 0, Tr(z) = 1} containing the matrix N -1 I N . The d.-g.f. is twice the matrix entropy

ω(z) = 2 N j=1 λ j (z) ln λ j (z), so that Z o = {z ∈ Z : z ≻ 0} and ω ′ (z) = 2 ln(z) + 2I N . This setup, similarly to the Simplex one, results in z c = N -1 I N and Ω ≤ 2 √ ln N [2]
. When Z = S, it is relatively easy to compute the prox-mapping (see [START_REF] Ben-Tal | Non-Euclidean restricted memory level method for large-scale convex optimization[END_REF][START_REF] Nemirovski | Prox-method with rate of convergence O(1/t) for variational inequalities with Lipschitz continuous monotone operators and smooth convex-concave saddle point problems[END_REF]); this task reduces to the singular value decomposition of a matrix from E. It should be added that the matrices from S are exactly the matrices of the form a = H(b) ≡ (Tr(exp{b})) -1 exp{b} with b ∈ E. Note also that when Z = S, the prox-mapping becomes "linear in matrix logarithm": if z = H(a), then P z (ξ) = H(a -ξ/2).

Algorithm: the construction

The t-step SMP algorithm is applied to the v.i. ( 2), works as follows:

Algorithm 1.

1. Initialization: Choose r 0 ∈ Z o and stepsizes γ τ > 0, 1 ≤ τ ≤ t.

2.

Step τ , τ = 1, 2, ..., t: Given r τ -1 ∈ Z o , set

w τ = P rτ-1 (γ τ F (r τ -1 )), r τ = P rτ-1 (γ τ F (w τ )) (19)
When τ < t, loop to step t + 1.

At step t, output

z t = t τ =1 γ τ -1 t τ =1 γ τ w τ . ( 20
)
Here F is the approximation of F (•) available to the algorithm, so that F (z) ∈ E is the output of the "black box" -the oracle -representing F , the input to the oracle being z ∈ Z. In what follows we assume that F is a bounded monotone operator represented by a Stochastic Oracle.

Stochastic Oracle (SO). At the ith call to the SO, the input being z ∈ Z, the oracle returns the vector

F = Ξ(z, ζ i ), where {ζ i ∈ R N } ∞ i=1 is a sequence of i.i.d. random variables, and Ξ(z, ζ) : Z × R N → E is a Borel function satisfying the following Assumption I: With some µ ∈ [0, ∞), for all z ∈ Z we have (a) E {Ξ(z, ζ i ) -F (z)} * ≤ µ (b) E Ξ(z, ζ i ) -F (z) 2 * ≤ σ 2 . ( 21
)
which is slightly milder than [START_REF] Korpelevich | Extrapolation gradient methods and relation to Modified Lagrangeans[END_REF]. The associated version of Algorithm 1 will be referred to as Stochastic Mirror Prox (SMP) algorithm.

In some cases, we augment Assumption I by the following Assumption II: For all z ∈ Z and all i we have

E exp{ Ξ(z, ζ i ) -F (z) 2 * /σ 2 } ≤ exp{1}. ( 22
)
Note that Assumption II implies (21.b), since

exp{E Ξ(z, ζ i ) -F (z) 2 * /σ 2 } ≤ E exp{ Ξ(z, ζ i ) -F (z) 2 * /σ 2 }
by the Jensen inequality.

Algorithm: Main result

From now on, assume that the starting point r 0 in Algorithm 1 is the minimizer z c of ω(•) on Z. Further, to avoid unnecessarily complicated formulas (and with no harm to the efficiency estimates) we stick to the constant stepsize policy γ τ ≡ γ, 1 ≤ τ ≤ t, where t is a fixed in advance number of iterations of the algorithm. Our main result is as follows:

Theorem 1. Let v.i.
(2) with monotone operator F satisfying (4) be solved by t-step Algorithm 1 using a SO, and let the stepsizes

γ τ ≡ γ, 1 ≤ τ ≤ t, satisfy 0 < γ ≤ 1 √ 3L . Then (i) Under Assumption I, one has E Err vi ( z t ) ≤ K 0 (t) ≡ Ω 2 tγ + 7γ 2 [M 2 + 2σ 2 ] + 2µΩ, ( 23
)
where M is the constant from (4) and Ω is given by (16).

(ii) Under Assumptions I, II, one has, in addition to (23), for any Λ > 0,

Prob Err vi ( z t ) > K 0 (t) + ΛK 1 (t) ≤ exp{-Λ 2 /3} + exp{-Λt}, (24)
where

K 1 (t) = 7σ 2 γ 2 + 2σΩ √ t .
In the case of a Nash v.i., Err vi (•) in (23), (24) can be replaced with Err N (•).

When optimizing the bound (23) in γ, we get the following Corollary 1. In the situation of Theorem 1, let the stepsizes γ τ ≡ γ be chosen according to

γ = min 1 √ 3L , Ω 2 7t(M 2 + 2σ 2 ) . ( 25
)
Then under Assumption I one has

E Err vi ( z t ) ≤ K * 0 (t) ≡ max 7 4 Ω 2 L t , 7Ω M 2 +2σ 2 3t + 2µΩ, (26)
(see (16)). Under Assumptions I, II, one has, in addition to (26), for any Λ > 0,

Prob Err vi ( z t ) > K * 0 (t) + ΛK * 1 (t) ≤ exp{-Λ 2 /3} + exp{-Λt} (27) with K * 1 (t) = 7 2 Ωσ √ t . ( 28
)
In the case of a Nash v.i., Err vi (•) in (26), ( 27) can be replaced with Err N (•).

Remark 1. Observe that the upper bound (26) for the error of Algorithm 1 with stepsize strategy (25), in agreement with the lower bound of [START_REF] Nemirovski | Problem complexity and method efficiency in Optimization[END_REF], depends in the same way on the "size" σ of the perturbation Ξ(z, ζ i ) -F (z) and on the bound M for the non-Lipschitz component of F . From now on to simplify the presentation, with slight abuse of notation, we denote M the maximum of these quantities. Clearly, the latter implies that the bounds (21.b) and (22), and thus the bounds (26) -(28) of Corollary 1 hold with M substituted for σ.

Comparison with Robust Mirror SA Algorithm

Consider the case of a Nash s.v.i. with operator F satisfying (4) with L = 0, and let the SO be unbiased (i.e., µ = 0). In this case, the bound (26) reads

E {Err N ( z t )} ≤ 7ΩM √ t , (29)
where

M 2 = max sup z,z ′ ∈Z F (z) -F (z ′ ) 2 * , sup z∈Z E Ξ(z, ζ i ) -F (z) 2 *
The bound (29) looks very much like the efficiency estimate

E {Err N (z t )} ≤ O(1) ΩM √ t (30)
(from now on, all O(1)'s are appropriate absolute positive constants) for the approximate solution zt of the t-step Robust Mirror SA (RMSA) algorithm [4] 1) . In the latter estimate, Ω is exactly the same as in (29), and M is given by

M 2 = max sup z F (z) 2 * ; sup z∈Z E Ξ(z, ζ i ) -F (z) 2 *
.

Note that we always have M ≤ 2M , and typically M and M are of the same order of magnitude; it may happen, however (think of the case when F is "almost constant"), that M ≪ M . Thus, the bound (29) never is worse, and sometimes can be much better than the SA bound (30). It should be added that as far as implementation is concerned, the SMP algorithm is not more complicated than the RMSA (cf. the description of Algorithm 1 with the description

r t = P rt-1 ( F (r t-1 )), z t = t τ =1 γ τ -1 t τ =1 γ τ r τ ,
of the RMSA).

The just outlined advantage of SMP as compared to the usual Stochastic Approximation is not that important, since "typically" M and M are of the same order. We believe that the most interesting feature of the SMP algorithm is its ability to take advantage of a specific structure of a stochastic optimization problem, namely, insensitivity to the presence in the objective of large, but smooth and well-observable components.

We are about to consider several less straightforward applications of the outlined insensitivity of the SMP algorithm to smooth well-observed components in the objective.

Application to Stochastic Composite minimization

Our present goal is to apply the SMP algorithm to the Composite minimization problem (9) in the case when the associated monotone operator [START_REF] Nesterov | Excessive gap technique in nonsmooth convex minimization[END_REF] is given by a Stochastic Oracle.

Throughout this section, the structural assumptions A -C from Section 2.2 are in force.

Assumptions

We start with augmenting the description of the problem of interest [START_REF] Lu | Large-Scale Semidefinite Programming via Saddle Point Mirror-Prox Algorithm[END_REF], see Section 2.2, with additional assumptions specifying the SO and the SMP setup for the v.i. reformulation

find z * ∈ Z := X × Y : F (z), z -z * ≥ 0 ∀z ∈ Z (31)
of ( 9); here F is the monotone operator [START_REF] Nesterov | Excessive gap technique in nonsmooth convex minimization[END_REF]. Specifically, we assume that D. The embedding space X of X is equipped with a norm • x , and X itself -with a d.-g.f. ω x (x), the associated parameter (16.c) being some Ω x ;

E. The spaces E ℓ , 1 ≤ ℓ ≤ m, where the functions φ ℓ take their values, are equipped with norms (not necessarily the Euclidean ones)

• (ℓ) with conjugates • (ℓ, * ) such that ∀v, v ′ ∈ X : (a) [φ ′ ℓ (v) -φ ′ ℓ (v ′ )]h (ℓ) ≤ [L x v -v ′ x + M x ] h x (b) [φ ′ ℓ (v)]h (ℓ) ≤ Ω x L x h x (32) for certain selections φ ′ ℓ (v) ∈ ∂ K ℓ φ ℓ (v), v ∈ X 2)
and certain nonnegative constants L x , M x .

F. Functions φ ℓ (•) are represented by an unbiased SO. At the ith call to the oracle, x ∈ X being the input, the oracle returns vectors

f ℓ (x, ζ i ) ∈ E ℓ and linear mappings G ℓ (x, ζ i ) from X to E ℓ , 1 ≤ ℓ ≤ m ({ζ i } are i.i.d.
random "oracle noises") such that for any x ∈ X and i = 1, 2, ...,

(a) E {f ℓ (x, ζ i )} = φ ℓ (x), 1 ≤ ℓ ≤ m (b) E max 1≤ℓ≤m f ℓ (x, ζ i ) -φ ℓ (x) 2 (ℓ) ≤ M 2 x Ω 2 x ; (c) E {G ℓ (x, ζ i )} = φ ′ ℓ (x), 1 ≤ ℓ ≤ m, (d) E max h∈X h x≤1 [G ℓ (x, ζ i ) -φ ′ ℓ (x)]h 2 (ℓ) ≤ M 2 x , 1 ≤ ℓ ≤ m. (33)
G. The data participating in the Fenchel-type representation [START_REF] Nemirovski | Accuracy certificates for computational problems with convex structure[END_REF] of Φ are such that (a) the embedding space Y of Y is equipped with a norm • y , and Y itself --with a d.-g.f. ω y (y), the associated parameter (16.c) being some Ω y ;

(b) we have y y ≤ 2Ω y for all y ∈ Y ;3 (c) The convex function Φ * (y) is given by the precise deterministic first order oracle, and

Φ ′ * (y) -Φ ′ * (y ′ ) y, * ≤ L y y -y ′ y + M y (34)
for certain selection Φ ′ * (y) ∈ ∂Φ * (y), y ∈ Y , and some nonnegative L y , M y .

Stochastic Oracle for (31).

The assumptions E -G induce an unbiased SO for the operator F in (31), specifically, the oracle

Ξ(x, y, ζ i ) = m ℓ=1 G * ℓ (x, ζ i )[A ℓ y + b ℓ ]; - m ℓ=1 A * ℓ f ℓ (x, ζ i ) + Φ ′ * (y) , ( 35
)
and this is the oracle we will use when solving (9) by the SMP algorithm.

2) For a K-convex function φ :

X → E (X ⊂ X is convex, K ⊂ E is a closed convex cone) and x ∈ X, the K-subdifferential ∂ K φ(x) is comprised of all linear mappings h → Ph : X → E such that φ(u) ≥ K φ(x) + P(u -x) for all u ∈ X. When φ is Lipschitz continuous on X, ∂ K φ(x) = ∅ for all x ∈ X; if φ is differentiable at x ∈ int X (
as it is the case almost everywhere on int X), one has ∂φ(x) ∂x ∈ ∂ K φ(x).

4.2 Setup for the SMP as applied to (31), [START_REF] Nesterov | Excessive gap technique in nonsmooth convex minimization[END_REF] In retrospect, the setup for SMP we are about to present is kind of the best -resulting in the best possible efficiency estimate (26) -we can build from the entities participating in the description of the problem (9) as given by the assumptions A -G. Specifically, we equip the space E = X × Y with the norm

(x, y) ≡ x 2 x /Ω 2 x + y 2 y /Ω 2 y ⇒ (ξ, η) * = Ω 2 x ξ 2 x, * + Ω 2 y η 2 y, *
and equip Z = X × Y with the d.-g.f.

ω(x, y) = 1 Ω 2 x ω x (x) + 1 Ω 2 y ω y (y)
(it is immediately seen that ω(•) indeed is a d.-g.f. w.r.t. X, • ). The SMP-related properties of our setup are summarized in the following

Lemma 1. Let A = max y∈Y: y y ≤1 m ℓ=1 A ℓ y (ℓ, * ) , B = m ℓ=1 b ℓ (ℓ, * ) . ( 36
) (i) The parameter Ω associated with ω(•), • , Z by (16), is ≤ √ 2. (ii) One has ∀(z, z ′ ∈ Z) : F (z) -F (z ′ ) * ≤ L z -z ′ + M, (37)
where

L = 4AΩ 2 x Ω y L x + Ω x B + Ω 2 y L y , M = [3AΩ y + B] Ω x M x + Ω y M y . Besides this, ∀(z ∈ Z, i) : E {Ξ(z, ζ i )} = F (z); E Ξ(z, ζ i ) -F (z) 2 * ≤ M 2 . ( 38
)
Finally, when relations (33.b,d) are strengthened to

E exp max 1≤ℓ≤m f ℓ (x, ζ i ) -φ ℓ (x) 2 (ℓ) /(Ω x M ) 2 ≤ exp{1}, E exp max h∈X , h x≤1 [G ℓ (x) -φ ′ ℓ (x)]h 2 (ℓ) /M 2 ≤ exp{1}, 1 ≤ ℓ ≤ m, (39)
then E exp{ Ξ(z, ζ i ) -F (z) 2 * /M 2 } ≤ exp{1}. ( 40
)
Combining Lemma 1 with Corollary 1, we get explicit efficiency estimates for the SMP algorithm as applied to the Stochastic composite minimization problem (9); these are nothing than estimates (26) with σ, Ω replaced with M , √ 2, respectively.

Application to Matrix Minimax problem

m ℓ=1 φ ℓ (x), A ℓ y F , A ℓ y = k j=1 P jℓ y j P T jℓ , ( 42
)
where Y is the spectahedron in the space Y of block-diagonal symmetric matrices y = Diag{y 1 , ..., y k } with k diagonal blocks of sizes q 1 , ..., q k . Note that we are in the situation described by assumptions A -C from Section 2.2, with Φ * (•) ≡ 0, and K ℓ := S p ℓ + ⊂ E ℓ := S p ℓ . Using the spectahedron setup for Y and setting M y = L y = 0, Ω y = 2 ln(q 1 + ... + q k ), we meet all assumptions in G. Let us specify the norms • (ℓ) on the spaces E ℓ = S p ℓ as the standard matrix norms | • | ∞ (maximal singular value), and let assumptions D -F from Section 4.1 take place. Note that in the case in question (36) reads

A = max 1≤j≤k max ξ∈R q j : ξ 2=1 m ℓ=1 P jℓ ξ 2 2 = max 1≤j≤k | m ℓ=1 P T jℓ P jℓ | ∞ , B = 0.
(look what are the extreme points of Y ), so that the quantities L, M as given by Lemma 1 become

L = O(1)A ln(q 1 + ... + q k )Ω 2 x L x , M = O(1)A ln(q 1 + ... + q k )Ω x M x . ( 43
)
Note that in the case of problem ( 13) one has A = 1.

Application to Stochastic Semidefinite Feasibility problem

Now consider the stochastic version of the Semidefinite Feasibility problem (S):

find x ∈ X : ψ ℓ (x) 0, 1 ≤ ℓ ≤ m. (44)
Assuming form now on that the latter problem is feasible, we can rewrite it as the Matrix Minimax problem Opt = min

x∈X max 1≤ℓ≤m λ max (φ ℓ (x)), φ ℓ (x) = β ℓ ψ ℓ (x), (45)
where β ℓ > 0 are "scale factors" we are free to choose. We are about to show how to use this freedom in order to improve the SMP efficiency estimates. Assuming, same as in the case of a general-type Matrix Minimax problem, that D takes place, let us modify assumptions E, F as follows:

E ′ : the -convex Lipschitz functions ψ ℓ : X → S p ℓ are such that max

h∈X , h x ≤1 |[ψ ′ ℓ (x) -ψ ′ ℓ (x ′ )]h| ∞ ≤ L ℓ x -x ′ x + M ℓ , max h∈X , h x ≤1 |ψ ′ ℓ (x)h| ∞ ≤ Ω x L ℓ (46)
for certain selections ψ ′ ℓ (x) ∈ ∂ K ℓ ψ ℓ (x), x ∈ X, with some known nonnegative constants L ℓ , M ℓ . F ′ : ψ ℓ (•) are represented by an SO which at the ith call, the input being x ∈ X, returns the matrices f ℓ (x, ζ i ) ∈ S p ℓ and the linear maps G ℓ (x, ζ i ) from X to S p ℓ ({ζ i } are i.i.d. random "oracle noises") such that for any x ∈ X it holds

(a) E f ℓ (x, ζ i ) = ψ ℓ (x), E G ℓ (x, ζ i ) = ψ ′ ℓ (x), 1 ≤ ℓ ≤ m (b) E max 1≤ℓ≤m | f ℓ (x, ζ i ) -ψ ℓ (x)| 2 ∞ /(Ω x M ℓ ) 2 ≤ 1 (c) E max h∈X , h x ≤1 |[ G ℓ (x, ζ i ) -ψ ′ ℓ (x)]h| 2 ∞ /M 2 ℓ ≤ 1, 1 ≤ ℓ ≤ m. ( 47
)
Given a number t of steps of the SMP algorithm, let us act as follows.

(I): We compute the m quantities µ ℓ = ΩxL ℓ √ t + M ℓ , ℓ = 1, ..., m, and set

µ = max 1≤ℓ≤m µ ℓ , β ℓ = µ µ ℓ , φ ℓ (•) = β ℓ ψ ℓ (•), L x = Ω -1 x µ √ t, M x = µ. (48)
Note that by construction β ℓ ≥ 1 and L x /L ℓ ≥ β ℓ , M x /M ℓ ≥ β ℓ for all ℓ, so that the functions φ ℓ satisfy (32) with the just defined L x , M x . Further, the SO for ψ ℓ (•)'s can be converted into an SO for φ ℓ (•)'s by setting

f ℓ (x, ζ) = β ℓ f ℓ (x, ζ), G ℓ (x, ζ) = β ℓ G ℓ (x, ζ).
By (47) and due to M x , and the fact that A = 1, B = 0, we arrive at the following result:

L x /L ℓ ≥ β ℓ , M x /M ℓ ≥ β ℓ ,
Proposition 1. With the outlined construction, the t-step SMP algorithm with the setup presented in Section 4.2 (where one uses A = 1, B = L y = M y = 0 and the just defined L x , M x ) and constant stepsizes γ τ ≡ γ defined by (25), yields an approximate solution z t = ( x t , y t ) such that

E max 1≤ℓ≤m max[β ℓ λ max (ψ ℓ ( x t ), 0] ≤ E max 1≤ℓ≤m β ℓ λ max (ψ ℓ ( x t )) -Opt ≤ K 0 (t) ≡ 80 Ωxµ √ ln( m ℓ=1 p ℓ ) √ t , (50)
(cf. (26) and take into account that we are in the case of Ω = √ 2, while the optimal value in (49) is nonpositive, since (44) is feasible).

Furthermore, if assumptions (47.b,c) are strengthened to

E max 1≤ℓ≤m exp{| f ℓ (x, ζ i ) -ψ ℓ (x)| 2 ∞ /(Ω x M ℓ ) 2 } ≤ exp{1}, E exp{ max h∈X , h x ≤1 |[ G ℓ (x, ζ i ) -ψ ′ ℓ (x)]h| 2 ∞ /M 2 ℓ } ≤ exp{1}, 1 ≤ ℓ ≤ m,
then, in addition to (50), we have for any Λ > 0:

Prob max 1≤ℓ≤m max[β ℓ λ max (ψ ℓ ( x t )), 0] > K 0 (t) + ΛK 1 (t) ≤ exp{-Λ 2 /3} + exp{-Λt},
where

K 1 (t) = 15Ω x µ ln( m ℓ=1 p ℓ ) √ t .
Discussion Imagine that instead of solving the system of matrix inequalities (44), we were interested to solve just a single matrix inequality ψ ℓ (x) 0, x ∈ X. When solving this inequality by the SMP algorithm as explained above, the efficiency estimate would be

E max[λ max (ψ ℓ ( x ℓ t )), 0] ≤ O(1) ln(p ℓ + 1)Ω x Ω x L ℓ t + M ℓ √ t = O(1) ln(p ℓ + 1)β -1 ℓ Ω x µ √ t ,
(recall that the matrix inequality in question is feasible), where x ℓ t is the resulting approximate solution. Looking at (50), we see that the expected accuracy of the SMP as applied, in the aforementioned manner, to (44) is only by a logarithmic in ℓ p ℓ factor worse:

E {max[λ max (ψ ℓ ( x t ), 0]} ≤ O(1) ln( m ℓ=1 p ℓ )β -1 ℓ Ωxµ √ t = O(1) ln( m ℓ=1 p ℓ ) Ωxµ ℓ √ t . (51)
Thus, as far as the quality of the SPM-generated solution is concerned, passing from solving a single matrix inequality to solving a system of m inequalities is "nearly costless". As an illustration, consider the case where some of ψ ℓ are "easy" -smooth and easy-to-observe (M ℓ = 0), while the remaining ψ ℓ are "difficult", i.e., might be non-smooth and/or difficult-to-observe (Ω x L ℓ / √ t ≤ M ℓ ). In this case, (51) reads

E {ψ ℓ ( x t )} ≤ O(1) ln( m ℓ=1 p ℓ ) • Ω 2 x L ℓ t , ψ ℓ is easy, ΩxM ℓ √
t , ψ ℓ is difficult. In other words, the violations of the easy and the difficult constraints in (44) converge to 0 as t → ∞ with the rates O(1/t) and O(1/ √ t), respectively. It should be added that when X is the unit Euclidean ball in X = R n and X, X are equipped with the Euclidean setup, the rates of convergence O(1/t) and O(1/ √ t) are the best rates one can achieve without imposing bounds on n and/or imposing additional restrictions on ψ ℓ 's.

Eigenvalue optimization via SMP

The problem we are interested in now is

Opt = min x∈X f (x) := λ max (A 0 + x 1 A 1 + ... + x n A n ), X = {x ∈ R n : x ≥ 0, n i=1 x i = 1}, ( 52
)
where A 0 , A 1 , ..., A n , n > 1, belong to the space S of symmetric matrices with block-diagonal structure (p 1 , ..., p m ) (i.e., a matrix A ∈ S is block-diagonal with p ℓ × p ℓ diagonal blocks A ℓ , 1 ≤ ℓ ≤ m). We set

p (κ) = m ℓ=1 p κ ℓ , κ = 1, 2, ...; p max = max ℓ p ℓ ; A ∞ = max 1≤j≤n |A j | ∞ . Setting φ ℓ : X → E ℓ = S p ℓ , φ ℓ (x) = A ℓ 0 + n j=1 x j A ℓ j , 1 ≤ ℓ ≤ m,
we represent (52) as a particular case of the Matrix Minimax problem [START_REF] Nesterov | Dual extrapolation and its applications to solving variational inequalities and related problems[END_REF], with all functions φ ℓ (x) being affine and X being the standard simplex in X = R n . Now, since A j are known in advance, there is nothing stochastic in our problem, and it can be solved either by interior point methods, or by "computationally cheap" gradient-type methods which are preferable when the problem is large-scale and medium accuracy solutions are sought. For instance, one can apply the t-step (deterministic) Mirror Prox algorithm DMP from [START_REF] Nemirovski | Prox-method with rate of convergence O(1/t) for variational inequalities with Lipschitz continuous monotone operators and smooth convex-concave saddle point problems[END_REF] to the saddle point reformulation [START_REF] Nesterov | Smooth minimization of non-smooth functions[END_REF]

f (x t ) -Opt ≤ O(1) ln(n) ln(p (1) )A ∞ t . ( 54
)
This efficiency estimate is the best known so far among those attainable with "computationally cheap" deterministic methods. On the other hand, the complexity of one step of the algorithm is dominated, up to an absolute constant factor, by the necessity, given x ∈ X and y ∈ Y , 1. to compute A 0 + n j=1 x j A j and [Tr(Y A 1 ); ...; Tr(Y A n )]; 2. to compute the eigenvalue decomposition of an y ∈ S.

When using the standard Linear Algebra, the computational effort per step is (2) + p (3) ] (55) arithmetic operations. We are about to demonstrate that one can equip the deterministic problem in question by an "artificial" SO in such a way that the associated SMP algorithm, under certain circumstances, exhibits better performance than deterministic algorithms. Let us consider the following construction of the SO for F (different from the SO (35)!). Observe that the monotone operator associated with the saddle point problem (53) is F (x, y) = [Tr(yA 1 ); ...; Tr(yA n )]

C det = O(1)[np
F x (x,y) ; -A 0 - n j=1 x j A j F y (x,y) . ( 56
)
Given x ∈ X, y = Diag{y 1 , ..., y m } ∈ Y , we build a random estimate Ξ = [Ξ x ; Ξ y ] of F (x, y) = [F x (x, y); F y (x, y)] as follows:

1. we generate a realization  of a random variable taking values 1, ..., n with probabilities x 1 , ..., x n (recall that x ∈ X, the standard simplex, so that x indeed can be seen as a probability distribution), and set

Ξ y = A 0 + A  ; (57
The just defined random estimate Ξ of F (x, y) can be expressed as a deterministic function Ξ(x, y, η) of (x, y) and random variable η uniformly distributed on [0, 1]. Assuming all matrices A j directly available (so that it takes O(1) arithmetic operations to extract a particular entry of A j given j and indexes of the entry) and given x, y and η, the value Ξ(x, y, ξ) can be computed with the arithmetic cost O(1)(n(p max ) 2 + p (2) ) (indeed, O(1)(n + p (1) ) operations are needed to convert η into ı and , O(1)p (2) operations are used to write down the y-component -A 0 -A  of Ξ, and O(1)n(p max ) 2 operations are needed to compute Ξ x ). Now consider the SO's Ξ k (k is a positive integer) obtained by averaging the outputs of k calls to our basic oracle Ξ. Specifically, at the ith call to the oracle Ξ k , z = (x, y) ∈ Z = X × Y being the input, the oracle returns the vector

Ξ k (z, ζ i ) = 1 k k s=1 Ξ(z, η is ),
where ζ i = [η i1 ; ...; η ik ] and {η is } 1≤i, 1≤s≤k are independent random variables uniformly distributed on [0, 1]. Note that the arithmetic cost of a single call to Ξ k is

C k = O(1)k(n(p max ) 2 + p (2)
).

The Nash v.i. associated with (53) and the stochastic oracle Ξ k (k is the first parameter of our construction) specify a Nash s.v.i. on the domain Z = X × Y . We equip the standard simplex X and its embedding space X = R n with the Simplex setup, and the spectahedron Y and its embedding space S with the Spectahedron setup (see Section 3.2). Let us next combine the x-and the y-setups, exactly as explained in the beginning of Section 4.2, into an SMP setup for the domain Z = X × Y -a d.-g.f. ω(•) and a norm • on the embedding space R n × (S p1 × ... × S p ℓ ) of Z. The SMP-related properties of the resulting setup are summarized in the following statement.

Lemma 2. Let n ≥ 3, p (1) ≥ 3. Then (i) The parameter of the just defined d.-g.f. ω w.r.t. the just defined norm

• is Ω = √ 2.
(ii) For any z, z ′ ∈ Z one has

F (z) -F (z ′ ) * ≤ L z -z ′ , L = √ 2 ln(n) + ln(p (1) ) A ∞ . (59)
Besides this, for any (z ∈ Z, i = 1, 2, ...,

(a) E {Ξ k (z, ζ i )} = F (z); (b) E exp{ Ξ(z, ζ i ) -F (z) 2 * /M 2 } ≤ exp{1}, M = 27[ln(n) + ln(p (1) )]A ∞ / √ k. (60)
Combining Lemma 2 and Corollary 1, we arrive at the following Proposition 2. With properly chosen positive absolute constants O(1), the t-step SMP algorithm with constant stepsizes

γ τ = O(1) min[1, k/t] ln(np (1) )A ∞ , 1 ≤ τ ≤ t [A ∞ = max 1≤j≤n |A j | ∞ ]
as applied to the saddle point reformulation of problem (52), the stochastic oracle being Ξ k , produces a random feasible approximate solution x t to the problem with the error

ǫ( x t ) = λ max A 0 + n j=1 [ x t ] j A j -Opt satisfying E {ǫ( x t )} ≤ O(1) ln(np (1) )A ∞ 1 t + 1 √ kt , (61)
and for any Λ > 0:

Prob ǫ( x t ) > O(1) ln(np (1) )A ∞ 1 t + 1 + Λ √ kt ≤ exp{-Λ 2 /3} + exp{-Λt}.
Further, assuming that all matrices A j are directly available, the overall computational effort to compute

x t is C = O(1)t k(n(p max ) 2 + p (2) ) + p (3) (62)
arithmetic operations.

To justify the bound (62) it suffices to note that O(1)k(n(p max ) 2 +p (2) ) operations per step is the price of two calls to the stochastic oracle Ξ k and O(1)(n + p (3) ) operations per step is the price of computing two prox mappings.

Discussion Let us find out whether randomization can help when solving a large-scale problem (52), that is, whether, given quality of the resulting approximate solution, the computational effort to build such a solution with the Stochastic Mirror Prox algorithm SMP can be essentially less than the one for the deterministic Mirror Prox algorithm DMP. To simplify our considerations, assume from now on that p ℓ = p, 1 ≤ ℓ ≤ m, and that ln(n) = O(1) ln(mp). Assume also that we are interested in a (perhaps, random) solution x t which with probability ≥ 1 -δ satisfies ǫ( x t ) ≤ ǫ. We fix a tolerance δ ≪ 1 and the relative accuracy ν = ǫ ln(mnp)A∞ ≤ 1 and look what happens when (some of) the sizes m, n, p of the problem become large.

Observe first of all that the overall computational effort to solve (52) within relative accuracy ν with the DMP algorithm is

C DMP (ν) = O(1)m(n + p)p 2 ν -1
operations (see (54), (55)). As for the SMP algorithm, let us choose k which balances the per step computational effort O(1)k(n(p max ) 2 + p (2) ) = O(1)k(n + m)p 2 to produce the answers of the stochastic oracle and the per step cost of prox mappings O(1)(n + mp 3 ), that is, let us set k = Ceil mp m+n 4 . With this choice of k, Proposition 2 says that to get a solution of the required quality, it suffices to carry out

t = O(1) ν -1 + ln(1/δ)k -1 ν -2 (63)
steps of the method, provided that this number of steps is ≥ ln(2/δ). The latter assumption is automatically satisfied when the absolute constant factor in (63) is ≥ 1 and ν ln(2/δ) ≤ 1, which we assume from now on. Combining (63) and the upper bounds on the arithmetic cost of an SMP step stated in Proposition 2, we conclude that the overall computational effort to produce a solution of the required quality with the SMP algorithm is

C SMP (ν, δ) = O(1)k(n + m)p 2 ν -1 + ln(1/δ)k -1 ν -2
operations, so that

R := C DMP (ν) C SMP (ν, δ) = O(1) m(n + p) (m + n)(k + ln(1/δ)ν -1 ) [k = Ceil mp m+n ]
We see that when ν, δ are fixed, m ≥ n/p and n/p is large, then R is large as well, that is, the randomized algorithm significantly outperforms its deterministic counterpart. Another interesting observation is as follows. In order to produce, with probability ≥ 1 -δ, an approximate solution to (52) with relative accuracy ν, the just defined SMP algorithm requires t steps, with t given by ( 63), and at every one of these steps it "visits" O(1)k(m+ n)p 2 randomly chosen entries in the data matrices A 0 , A 1 , ..., A n . The overall number of data entries visited by the algorithm is therefore

N SMP = O(1)tk(m + n)p 2 = O(1) kν -1 + ln(1/δ)ν -2 (m + n)p 2 .
At the same time, the total number of data entries is N tot = m(n + 1)p 2 . Therefore

ϑ := N SMP N tot = O(1) kν -1 + ln(1/δ)ν -2 1 m + 1 n [k = Ceil mp m+n ]
We see that when δ, ν are fixed, m ≥ n/p and n/p is large, ϑ is small, i.e., the approximate solution of the required quality is built when inspecting a tiny fraction of the data. This sublinear time behavior [START_REF] Rubinfeld | Sublinear time algorithms[END_REF] was already observed in [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF] for the Robust Mirror Descent Stochastic Approximation as applied to a matrix game (the latter problem is the particular case of (52) with p 1 = ... = p m = 1 and A 0 = 0). Note also that an "ad hoc" sublinear time algorithm for a matrix game, in retrospect close to the one from [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF], was discovered in [START_REF] Grigoriadis | A sublinear-time randomized approximation algorithm for matrix games[END_REF] as early as in 1995.

Then

Err vi ( z t ) ≤ t τ =1 γ τ -1 Γ(t), (66)
where Err vi ( z t ) is defined in (3),

Γ(t) = 2Θ(r 0 ) + t τ =1 3γ 2 τ 2 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 + ǫ 2 wτ 3 (67) + t τ =1 γ τ ∆ τ , w τ -y τ -1
and Θ(•) is defined by (16).

Finally, when (2) is a Nash v.i., one can replace Err vi ( z t ) in (66) with Err N ( z t ).

Proof of Theorem 2 1 0 . We start with the following simple observation: if r e is a solution to (15), then ∂ Z ω(r e ) contains -e and thus is nonempty, so that r e ∈ Z o . Moreover, one has ω ′ (r e ) + e, u -r e ≥ 0 ∀u ∈ Z.

Indeed, by continuity argument, it suffices to verify the inequality in the case when u ∈ rint(Z) ⊂ Z o .

For such an u, the convex function

f (t) = ω(r e + t(u -r e )) + r e + t(u -r e ), e , t ∈ [0, 1]
is continuous on [0, 1] and has a continuous on [0, 1] field of subgradients

g(t) = ω ′ (r e + t(u -r e )) + e, u -r e .
It follows that f is continuously differentiable on [0, 1] with the derivative g(t). Since the function attains its minimum on [0, 1] at t = 0, we have g(0) ≥ 0, which is exactly (68). 2 0 . At least the first statement of the following Lemma is well-known:

Lemma 3. For every z ∈ Z o , the mapping ξ → P z (ξ) is a single-valued mapping of E onto Z o , and this mapping is Lipschitz continuous, specifically,

P z (ζ) -P z (η) ≤ ζ -η * ∀ζ, η ∈ E. (69)
Besides this, for all u ∈ Z,

(a) V (P z (ζ), u) ≤ V (z, u) + ζ, u -P z (ζ) -V z (z, P z (ζ)) (b) ≤ V (z, u) + ζ, u -z + ζ 2 * 2 .
(70)

Proof : Let v ∈ P z (ζ), w ∈ P z (η). As V ′ u (z, u) = ω ′ (u) -ω ′ (z), invoking (68), we have v, w ∈ Z o and ω ′ (v) -ω ′ (z) + ζ, v -u ≤ 0 ∀u ∈ Z. ( 71
) ω ′ (w) -ω ′ (z) + η, w -u ≤ 0 ∀u ∈ Z. ( 72
)
Setting u = w in (71) and u = v in (72), we get

ω ′ (v) -ω ′ (z) + ζ, v -w ≤ 0, ω ′ (w) -ω ′ (z) + η, v -w ≥ 0, whence ω ′ (w) -ω ′ (v) + [η -ζ], v -w ≥ 0, or η -ζ * v -w ≥ η -ζ, v -w ≥ ω ′ (v) -ω ′ (w), v -w ≥ v -w 2 ,
Proof : (a): this is nothing but (69). (b): Using (70.a) in Lemma 3 we can write for u = r + :

V (w, r + ) ≤ V (z, r + ) + ζ, r + -w -V (z, w).
This results in

V (z, r + ) ≥ V (w, r + ) + V (z, w) + ζ, w -r + . (77)
Now using (70.a) with η substituted for ζ we get

V (r + , u) ≤ V (z, u) + η, u -r + -V (z, r + ) = V (z, u) + η, u -w + η, w -r + -V (z, r + ) ≤ V (z, u) + η, u -w + η -ζ, w -r + -V (z, w) -V (w, r + ) [by (77)] ≤ V (z, u) + η, u -w + η -ζ, w -r + -1 2 [ w -z 2 + w -r + 2 ],
where the concluding inequalities are due to the strong convexity of ω(•). To conclude the bound (b) of (76) it suffices to note that by the Young inequality,

η -ζ, w -r + ≤ η -ζ 2 * 2 + 1 2 w -r + 2 .
5 0 . We are able now to prove Theorem 2. By (4) we have

F (w τ ) -F (r τ -1 ) 2 * ≤ (L r τ -1 -w τ + M + ǫ rτ-1 + ǫ wτ ) 2 ≤ 3L 2 w τ -r τ -1 2 + 3M 2 + 3(ǫ rτ-1 + ǫ wτ ) 2 . (78)
Applying Lemma 4 with z = r τ -1 , ζ = γ τ F (r τ -1 ), η = γ τ F (w τ ) (so that w = w τ and r + = r τ ), we have for any u ∈ Z

γ τ F (w τ ), w τ -u + V (r τ , u) -V (r τ -1 , u) ≤ γ 2 τ 2 F (w τ ) -F (r τ -1 ) 2 -1 2 w τ -r τ -1 2 ≤ 3γ 2 τ 2 L 2 w τ -r τ -1 2 + M 2 + (ǫ rτ-1 + ǫ wτ ) 2 -1 2 w τ -r τ -1 2 [by (78)] ≤ 3γ 2 τ 2 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 [by (65)]
When summing up from τ = 1 to τ = t we obtain

t τ =1 γ τ F (w τ ), w τ -u ≤ V (r 0 , u) -V (r t , u) + t τ =1 3γ 2 τ 2 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 ≤ Θ(r 0 ) + t τ =1 3γ 2 τ 2 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 .
Hence, for all u ∈ Z,

t τ =1 γ τ F (w τ ), w τ -u ≤ Θ(r 0 ) + t τ =1 3γ 2 τ 2 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 + t τ =1 γ τ ∆ τ , w τ -u = Θ(r 0 ) + t τ =1 3γ 2 τ 2 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 + t τ =1 γ τ ∆ τ , w τ -y τ -1 + t τ =1 γ τ ∆ τ , y τ -1 -u (79)
where y τ are given by (64). Since the sequences {y τ }, {ξ τ = γ τ ∆ τ } satisfy the premise of Corollary 2, we have (∀u ∈ Z) :

t τ =1 γ τ ∆ τ , y τ -1 -u ≤ V (r 0 , u) + t τ =1 γ 2 τ 2 ∆ τ 2 * ≤ Θ(r 0 ) + t τ =1 γ 2 τ 2 ǫ 2 wτ ,
and thus (79) implies that for any u ∈ Z

t τ =1 γ τ F (w τ ), w τ -u ≤ Γ(t) (80)
with Γ(t) defined in (67). To complete the proof of (66) in the general case, note that since F is monotone, (80) implies that for all u ∈ Z,

t τ =1 γ τ F (u), w τ -u ≤ Γ(t), whence ∀(u ∈ Z) : F (u), z t -u ≤ t τ =1 γ τ -1 Γ(t).
When taking the supremum over u ∈ Z, we arrive at (66).

In the case of a Nash v.i., setting w τ = (w τ,1 , ..., w τ,m ) and u = (u 1 , ..., u m ) and recalling the origin of F , due to the convexity of φ i (z i , z i ) in z i , for all u ∈ Z we get from (80):

t τ =1 γ τ m i=1 [φ i (w τ ) -φ i (u i , (w τ ) i )] ≤ t τ =1 γ τ m i=1 F i (w τ ), (w τ ) i -u i ≤ Γ(t). Setting φ(z) = m i=1 φ i (z), we get t τ =1 γ τ φ(w τ ) - m i=1 φ i (u i , (w τ ) i ) ≤ Γ(t).
Recalling that φ(•) is convex and φ i (u i , •) are concave, i = 1, ..., m, the latter inequality implies that

t τ =1 γ τ φ( z t ) - m i=1 φ i (u i , ( z t ) i ) ≤ Γ(t),
or, which is the same,

m i=1 φ i ( z t ) - m i=1 φ i (u i , ( z t ) i ) ≤ t τ =1 γ τ -1 Γ(t).
This relation holds true for all u = (u 1 , ..., u m ) ∈ Z; taking maximum of both sides in u, we get

Err N ( z t ) ≤ t τ =1 γ τ -1 Γ(t).

Proof of Theorem 1

In what follows, we use the notation from Theorem 2. By this theorem, in the case of constant stepsizes

γ τ ≡ γ we have Err vi ( z t ) ≤ [tγ] -1 Γ(t), (81)
where

Γ(t) = Ω 2 + 3γ 2 2 t τ =1 M 2 + (ǫ rτ-1 + ǫ wτ ) 2 + ǫ 2 wτ 3 + γ t τ =1 ∆ τ , w τ -y τ -1 ≤ Ω 2 + 7γ 2 2 t τ =1 M 2 + ǫ 2 rτ-1 + ǫ 2 wτ + γ t τ =1 ∆ τ , w τ -y τ -1 . (82)
For a Nash v.i., Err vi in this relation can be replaced with Err N . Let us suppose that the random vectors ζ i are defined on the probability space (Ω, F , P). We define two nested families of σ-fields

F i = σ(r 0 , ζ 1 , ζ 2 , . . . , ζ 2i-1 ) and G i = σ(r 0 , ζ 1 , ζ 2 , . . . , ζ 2i ), i = 1, 2, ..., so that F 1 ⊂ ... G i-1 ⊂ F i ⊂ G i ⊂ ....
Then by description of the algorithm r τ -1 is G τ -1 -measurable and w τ is F τ -measurable. Therefore ǫ rτ-1 is F τ -measurable, and ǫ wτ and ∆ τ are G τ -measurable. We conclude that under Assumption I we have

E ǫ 2 rτ-1 |G τ -1 ≤ σ 2 , E ǫ 2 wτ |F τ ≤ σ 2 , E {∆ τ |F τ } * ≤ µ, (83)
and under Assumption II, in addition,

E exp{ǫ 2 rτ-1 σ -2 } G τ -1 ≤ exp{1}, E exp{ǫ 2 wτ σ -2 } F τ ≤ exp{1}. (84) Now, let Γ 0 (t) = 7γ 2 2 t τ =1 M 2 + ǫ 2 rτ-1 + ǫ 2 wτ .
We conclude by ( 83) that

E {Γ 0 (t)} ≤ 7γ 2 t 2 [M 2 + 2σ 2 ]. ( 85
) Further, y τ -1 clearly is F τ -1 -measurable, whence w τ -y τ -1 is F τ -measurable. Therefore E ∆ τ , w τ -y τ -1 F τ = E {∆ τ |F τ } , w τ -y τ -1 ≤ µ w τ -y τ -1 ≤ 2µΩ, (86)
where the concluding inequality follows from the fact that Z is contained in the • -ball of radius Ω centered at z c , see (18). From (86) it follows that

E γ t τ =1 ∆ τ , w τ -y τ -1 ≤ 2µγtΩ.
Combining the latter relation, (81), ( 82) and (85), we arrive at (23). (i) is proved.

To prove (ii), observe, first, that setting

J t = t τ =1 σ -2 ǫ 2 rτ-1 + σ -2 ǫ 2 wτ , we get Γ 0 (t) = 7γ 2 M 2 t 2 + 7γ 2 σ 2 2 J t . (87)
At the same time, setting H j = σ(r 0 , ξ 1 , ..., ξ j }, we can write

J t = 2t j=1 ξ j ,
where ξ j ≥ 0 is H j -measurable, and

E {exp{ξ j }|H j-1 } ≤ exp{1}, see (84). It follows that E exp{ k+1 j=1 ξ j } = E E exp{ k j=1 ξ j } exp{ξ k+1 } H k = E exp{ k j=1 ξ j }E exp{ξ k+1 } H k ≤ exp{1}E exp{ k j=1 ξ j } . (88)
Whence E[exp{J}] ≤ exp{2t}, and applying the Tchebychev inequality, we get

∀Λ > 0 : Prob {J > 2t + Λt} ≤ exp{-Λt}.
Along with (87) it implies that

∀Λ ≥ 0 : Prob Γ 0 (t) > 7γ 2 t 2 [M 2 + 2σ 2 ] + Λ 7γ 2 σ 2 t 2 ≤ exp{-Λt}. (89)
Let now ξ τ = ∆ τ , w τ -y τ -1 . Recall that w τ -y τ -1 is F τ -measurable. Besides this, we have seen that w τ -y τ -1 ≤ D ≡ 2Ω. Taking into account (84) and (86), we get

(a) E {ξ τ |F τ } ≤ ρ ≡ µD, (b) E exp{ξ 2 τ R -2 }|F τ ≤ exp{1}, with R = σD. (90)
Observe that exp{x} ≤ x + exp{9x 2 /16} for all x. Thus (90.b) implies for 0

≤ s ≤ 4 3R E {exp{sξ τ }|F τ } ≤ E{sξ τ |F τ } + E exp 9s 2 ξ 2 τ 16 |F τ ≤ sρ + exp 9s 2 R 2 16 ≤ exp sρ + 9s 2 R 2 16 . (91)
Further, we have

sξ τ ≤ 3 8 s 2 R 2 + 2 3 ξ 2 τ R -2
, hence for all s ≥ 0,

E {exp{sξ τ }|F τ } ≤ exp{3s 2 R 2 /8}E exp 2ξ 2 τ 3R 2 |F τ ≤ exp 3s 2 R 2 8 + 2 3 .
When s ≥ 4 3R , the latter quantity is ≤ exp{3s 2 R 2 /4}, which combines with (91) to imply that for s ≥ 0,

E {exp{sξ τ }|F τ } ≤ exp{sρ + 3s 2 R 2 /4}. (92)
Acting as in (88), we derive from (92) that

s ≥ 0 ⇒ E exp{s t τ =1 ξ τ } ≤ exp{stρ + 3s 2 tR 2 /4},
and by the Tchebychev inequality, for all Λ > 0,

Prob t τ =1 ξ τ > tρ + ΛR √ t ≤ inf s≥0 exp{3s 2 tR 2 /4 -sΛR √ t} = exp{-Λ 2 /3}.
Finally, we arrive at

Prob γ t τ =1 ∆ τ , w τ -y τ -1 > 2γ µt + Λσ √ t Ω ≤ exp{-Λ 2 /3}. (93)
for all Λ > 0. Combining (81), (82), ( 89) and (93), we get (24).

Proof of Lemma 1

Proof of (i) We clearly have

Z o = X o × Y o , and ω(•) is indeed continuously differentiable on this set. Let z = (x, y) and z ′ = (x ′ , y ′ ), z, z ′ ∈ Z. Then ω ′ (z) -ω ′ (z ′ ), z -z ′ = 1 Ω 2 x ω ′ x (x) -ω ′ x (x ′ ), x -x ′ + 1 Ω 2 y ω ′ y (y) -ω ′ y (y ′ ), y -y ′ ≥ 1 Ω 2 x x -x ′ 2 x + 1 Ω 2 y y -y ′ 2 y ≥ [x ′ -x; y ′ -y] 2 .
Thus, ω(•) is strongly convex on Z, modulus 1, w.r.t. the norm • . Further, the minimizer of ω(•) on Z clearly is z c = (x c , y c ), and it is immediately seen that max z∈Z

V (z c , z) = 1, whence Ω = √ 2.
Proof of (ii) 1 0 . Let z = (x, y) and z ′ = (x ′ , y ′ ) with z, z ′ ∈ Z. Note that by assumption G in Section 4.1 we have y ′ y ≤ 2Ω y .

Further, we have from ( 12)

F (z ′ ) -F (z) = [∆ x ; ∆ y ],
where

∆ x = m ℓ=1 [φ ′ ℓ (x ′ ) -φ ′ ℓ (x)] * [A ℓ y ′ + b ℓ ] + m ℓ=1 [φ ′ ℓ (x)] * A ℓ [y ′ -y], ∆ y = - m ℓ=1 A * ℓ [φ ℓ (x) -φ ℓ (x ′ )] + Φ ′ * (y ′ ) -Φ ′ * (y).
We have

∆x x, * = max h∈X h x ≤1 h, m ℓ=1 [[φ ′ ℓ (x ′ ) -φ ′ ℓ (x)] * [A ℓ y ′ + b ℓ ] + [φ ′ ℓ (x)] * A ℓ [y ′ -y]] X ≤ m ℓ=1   max h∈X h x≤1 h, [φ ′ ℓ (x ′ ) -φ ′ ℓ (x)] * [A ℓ y ′ + b ℓ ] X + max h∈X , h x≤1 h, [φ ′ ℓ (x)] * A ℓ [y ′ -y] X   = m ℓ=1 max h∈X h x≤1 [φ ′ ℓ (x ′ ) -φ ′ ℓ (x)]h, A ℓ y ′ + b ℓ X + max h∈X , h x≤1 [φ ′ ℓ (x)]h, A ℓ [y ′ -y] X ≤ m ℓ=1 max h∈X h x ≤1 [φ ′ ℓ (x ′ ) -φ ′ ℓ (x)]h (ℓ) A ℓ y ′ + b ℓ (ℓ, * ) + max h∈X h x≤1 φ ′ ℓ (x)h (ℓ) A ℓ [y ′ -y] (ℓ, * ) .
Then by (32),

where all ρ ℓ ≥ 0, ℓ ρ ℓ ≤ 2AΩ y + B and

ξ ℓ = ξ ℓ (ζ i ) = max h∈X , h x ≤1 [φ ′ ℓ (x) -G ℓ (x, ζ i )]h (ℓ)
Denoting by p 2 (η) the second moment of a scalar random variable η, observe that p(•) is a norm on the space of square summable random variables representable as deterministic functions of ζ i , and that and the latter quantity is ≤ M , see (37). We have established the second relation in (38).

3 0 . It remains to prove that in the case of (39), relation (40) takes place. To this end, one can repeat word by word the reasoning from item 2 0 with the function p e (η) = inf t > 0 : E exp{η 2 /t 2 } ≤ exp{1} in the role of p(η). Note that similarly to p(•), p e (•) is a norm on the space of random variables η which are deterministic functions of ζ i and are such that p e (η) < ∞.

Proof of Lemma 2

Item (i) can be verified exactly as in the case of Lemma 1; the facts expressed in (i) depend solely on the construction from Section 4.2 preceding the latter Lemma, and are independent of what are the setups for X, X and Y, Y.

Let us verify item (ii). Note that we are in the situation (x, y) = Combining the latter bounds with (102) we get (60.b).

  p ℓ + -convex Lipschitz continuous mappings. These data induce the Matrix Minimax problem min x∈X max 1≤j≤k λ max m ℓ=1 P T jℓ φ ℓ (x)P jℓ , (P )

  jℓ y j P jℓ = max y∈Y m ℓ=1 u ℓ , A ℓ y F , A ℓ y = k j=1 P jℓ y j P T jℓ

  of our specific Matrix Minimax problem, i.e., to the saddle point problem min x∈X max y∈Y y, A 0 + n j=1 x j A j F , Y = y = Diag{y 1 , ..., y m } : y ℓ ∈ S p ℓ + , 1 ≤ ℓ ≤ m, Tr(y) = 1 .(53)The accuracy of the approximate solution xt of the DMP algorithm is [8,Example 2]

) 2 .

 2 we compute the quantities ν ℓ = Tr(y ℓ ), 1 ≤ ℓ ≤ m. Since y ∈ Y , we have ν ℓ ≥ 0 and m ℓ=1 ν ℓ = 1. We further generate a realization ı of random variable taking values 1, ..., m with probabilities ν 1 , ..., ν m , and set Ξ x = [Tr(A ı 1 ȳı ); ...; Tr(A ı n ȳı )], ȳı = (Tr(y ı )) -1 y ı .

[[

  Lx xx ′ x + Mx][ A ℓ y ′ (ℓ, * ) + b ℓ (ℓ, * ) ] + ΩxLx A ℓ [yy ′ ] (ℓ, * ) = [Lx xx ′ x + Mx] m ℓ=1 A ℓ y ′ (ℓ, * ) + b ℓ (ℓ, * ) ] + ΩxLx m ℓ=1 A ℓ [yy ′ ] (ℓ, * ) ≤ [Lx xx ′ x + Mx][A y ′ y + B] + ΩxLxA yy ′ y ,by definition of A and B. Next, due to (94) we get by definition of• ∆x x, * ≤ [Lx xx ′ x + Mx][2AΩy + B] + ΩxLxA yy ′ y ≤ [ΩxLx zz ′ + Mx][2AΩy + B] + ΩxLxAΩy zz ′ , which implies ∆ x x, * ≤ [3AΩ x L x Ω y + B] z -z ′ + [2AΩ y + B]M x . ℓ (x)φ ℓ (x ′ )] + Φ ′ * (y ′ ) -Φ ′ * (y) Y ≤ max η∈Y, η y ≤1 m ℓ=1 η, A * ℓ [φ ℓ (x)φ ℓ (x ′ )] Y + Φ ′ * (y ′ ) -Φ ′ * (y) y, * = max η∈Y, η y ≤1 m ℓ=1 A ℓ η, φ ℓ (x)φ ℓ (x ′ ) E ℓ + Φ ′ * (y ′ ) -Φ ′ * (y) y, * ≤ max η∈Y, η y ≤1 m ℓ=1 A ℓ η (ℓ, * ) φ ℓ (x)φ ℓ (x ′ ) (ℓ) + Φ ′ * (y ′ ) -Φ ′ * (y) y, * ≤ max η∈Y, η y ≤1 m ℓ=1 A ℓ η (ℓ, * ) ΩxLx xx ′ x + [Ly yy ′ y + My],by (32.b) and (34). Therefore∆ y y, * ≤ AΩ x L x x -x ′ x + L y y -y ′ y + M y , whence ∆ y y, * ≤ [AΩ 2 x L x + Ω y L y ] z -z ′ + M y .(96)From (95), (96) it follows thatF (z) -F (z ′ ) * ≤ Ω x ∆ x x, * + Ω y ∆ y y, * ≤ 4AL x Ω 2 x Ω y + Ω x B + Ω 2 y L y z -z ′ + [2AΩ x Ω y + Ω x B] M x + Ω y M y .We have justified (37) 2 0 . Let us verify (38). The first relation in (38) is readily given by (33.a,c). Let us fix z = (x, y) ∈ Z and i, and let∆ = F (z) -Ξ(z, ζ i ) ℓ (x) -G ℓ (x, ζ i )] * ψ ℓ [A ℓ y + b ℓ ] ℓ (x) -f ℓ (x, ζ i )] ℓ, * ) ≤ 2AΩ y + B. (98)Besides this, for u ℓ ∈ E ℓ we havem ℓ=1 A * ℓ u ℓ y, * = max η∈Y, η y ≤1 m ℓ=1 A * ℓ u ℓ , η Y = max η∈Y, η y ≤1 m ℓ=1 u ℓ , A ℓ η Y ≤ max η∈Y, η y ≤1 1≤ℓ≤m u ℓ (ℓ) A ℓ η (ℓ, * ) ≤ max η∈Y, η y ≤1 max 1≤ℓ≤m u ℓ (ℓ) 1≤ℓ≤m A ℓ η (ℓ, * ) = A max 1≤ℓ≤m u ℓ (ℓ) .(99)Hence, settingu ℓ = φ ℓ (x) -f ℓ (x, ζ i ) we obtain ∆ y y, * = m ℓ=1 A * ℓ [φ ℓ (x) -f ℓ (x, ζ i )] y, * ≤ A max 1≤ℓ≤m φ ℓ (x) -f ℓ (x, ζ i ) (ℓ) ℓ (x) -G ℓ (x, ζ i )] * ψ ℓ X = max h∈X , h x ≤1 m ℓ=1 [φ ′ ℓ (x) -G ℓ (x, ζ i )]h, ψ ℓ X ≤ max h∈X , h x ≤1 m ℓ=1 [φ ′ ℓ (x) -G ℓ (x, ζ i )]h(ℓ) ψ ℓ (ℓ, * ) ≤ m ℓ=1 max h∈X , h x ≤1 [φ ′ ℓ (x) -G ℓ (x, ζ i )]h (ℓ) ξ ℓ =ξ ℓ (ζi) ψ ℓ (ℓ, * ) ρ ℓ Invoking (98), we conclude that ∆ x x, * ≤ m ℓ=1 ρ ℓ ξ ℓ ,

E ∆ 2 * 1 2 = E Ω 2 x ∆ x 2 x 2 ≤

 2222 p(ξ) ≤ Ω x M x , p(ξ ℓ ) ≤ M x by (33.b,d). Now by (100), (101), p (Ω x ∆ x x, * + Ω y ∆ y y, * ) ≤ p Ω x m ℓ=1 ρ ℓ ξ ℓ + Ω y Aξ ≤ Ω x ℓ ρ ℓ max ℓ p(ξ ℓ ) + Ω y Ap(ξ) ≤ Ω x [2AΩ y + B]M x + Ω y AΩ x M x ,

x 2 1 /( 2

 12 ln(n)) + |y| 2 1 /(4 ln(p(1) )), (ξ, η) * = 2 ln(n) ξ 2 ∞ + 4 ln(p(1) )|η| 2 ∞ .(102)For z = (x, y), z ′ = (x ′ , y ′ ) ∈ Z we have F (z) -F (z ′ ) = [Tr((y -y ′ )A 1 ); ...; Tr((y -y ′ )A n )] ≤ |y -y ′ | 1 max 1≤j≤n |A j | ∞ ≤ 2 ln(p (1) )A ∞ z -z ′ , |∆ y | ∞ ≤ x -x ′ ∞ max 1≤j≤n |A j | ∞ ≤ 2 ln(n)A ∞ z -z ′ , and (∆ x , ∆ y ) * ≤ 2 2 ln(n) ln(p (1) )A ∞ z -z ′ ,as required in (59). Further, relation (60.a) is clear from the construction of Ξ k . To prove (60.b), observe that when (x, y) ∈ Z, we have Ξx (x, y, η) ∞ ≤≤ A ∞ , |Ξ y (x, y, η)| ∞ ≤ A ∞ (see (57), (58)), whence Ξ x (x, y, η) -F x (x, y) ∞ ≤ 2A ∞ , |Ξ y (x, y, η) -F y (x, y)| ∞ ≤ 2A ∞ (103)due to F (x, y) = E η {Ξ(x, y, η)}, Applying [5, Theorem 2.1(iii), Example 3.2, Lemma 1], we derive from (103) that for every (x, y) ∈ Z and every i = 1, 2, ... it holdsE exp{ Ξ x k (x, y, ζ i ) -F x (x, y) 2 ∞ /N 2 k,x } ≤ exp{1}, N k,x = 2A ∞ 2 exp{1/2} ln(n) + 3 k -1/2 and E exp{ Ξ y k (x, y, ζ i ) -F y (x, y) 2 ∞ /N 2 k,y } ≤ exp{1},N k,y = 2A ∞ 2 exp{1/2} ln(p(1) ) + 3 k -1/2 .

  X → S p ℓ are -convex Lipschitz continuous mappings and P jℓ ∈ R p ℓ ×qj . As it was shown in Section 2.2.2, (41) admits saddle point representation min

	Consider Matrix Minimax problem from Section 2.2.2, that is, the problem
	min x∈X	max 1≤j≤k	λ max	m ℓ=1	P T jℓ φ ℓ (x)P jℓ ,	(41)
	where φ ℓ (•) :

x∈X max y=Diag{y1,...,y k }∈Y

  , ..., φ m and solve this Stochastic composite problem by t-step SMP algorithm. Combining Lemma 1, Corollary 1 and taking into account the origin of the quantities L x ,

	this oracle satisfies (33).
	(II) We then build the Stochastic Matrix Minimax problem
	Opt = min x∈X	max 1≤ℓ≤m	λ max (φ ℓ (x)),	(49)
	associated with the just defined φ 1

[START_REF] Azuma | Weighted sums of certain dependent random variables[END_REF] In this reference, only the Minimization and the Saddle Point problems are considered. However, the results of[START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF] can be easily extended to s.v.i.'s.

This requirement can be ensured by shifting Y to include the origin and the associated shift in ωy(•), see (18).

The rationale behind balancing is clear: with the just defined k, the arithmetic cost of an iteration still is of the same order as when k = 1, while the left hand side in the efficiency estimate (61) becomes better than for k = 1.

Research of this author was partly supported by the Office of Naval Research grant # N000140811104 and the NSF grant DMS-0914785.

Appendix

Preliminaries

We need the following technical result about the algorithm (19), (20): Theorem 2. Consider t-step algorithm 1 as applied to a v.i. [START_REF] Ben-Tal | Non-Euclidean restricted memory level method for large-scale convex optimization[END_REF] with a monotone operator F satisfying [START_REF] Juditsky | Stochastic Approximation Approach to Stochastic Programming[END_REF]. For τ = 1, 2, ..., let us set

for z belonging to the trajectory {r 0 , w 1 , r 1 , ..., w t , r t } of the algorithm, let

and let {y τ ∈ Z o } t τ =0 be the sequence given by the recurrence

Assume that

and (69) follows. This relation, as a byproduct, implies that P z (•) is single-valued.

To prove (70), let v = P z (ζ). We have

Indeed, observe that by definition, V (z, •) is strongly convex modulus 1 w.r.t. • , and

3 0 . We have the following simple corollary of Lemma 3:

.. be a sequence of elements of E. Define the sequence {y τ } ∞ τ =0 in Z o as follows:

Then y τ is a measurable function of y 0 and ξ 1 , ..., ξ τ such that

with |ζ τ | ≤ r ξ τ * (here r = max u∈Z u ). Further,

Proof : Using the bound (70.a) with ζ = ξ τ and z = y τ -1 , so that y τ = P yτ-1 (ξ τ ), we obtain for any u ∈ Z:

summing up these inequalities over τ we get (73). Further, by definition of P z (ξ) we have

so that ζ τ ≤ r ξ τ * due to V ≥ 0, and

Summing up these inequalities over τ , we get (74).

4 0 . We also need the following result. (76)