Sequential online subsampling for thinning experimental designs

Suppose that $X$ is distributed with the probability measure $\mu$ on ${\mathscr{X}}\subseteq\mathds{R}^{d}$, a subset of $\mathds{R}^{d}$ with nonempty interior, with $d\geq 1$. For any $\xi\in{\mathscr{P}}^{+}({\mathscr{X}})$, the set of positive measure $\xi$ on ${\mathscr{X}}$ (not necessarily of mass one), we denote $\mathbf{M}(\xi)=\int_{\mathscr{X}}{\mathscr{M}}(x)\,\xi(\mbox{\rm d}x)$ where, for all $x$ in ${\mathscr{X}}$, ${\mathscr{M}}(x)\in\mathbb{M}^{\geq}$, the set (cone) of symmetric non-negative definite $p\times p$ matrices. We assume that $p>1$ in the rest of the paper (the optimal selection of information in the case $p=1$ forms a variant of the secretary problem for which an asymptotically optimal solution can be derived, see Albright and Derman, (1972); Pronzato, (2001)).

We denote by $\Phi:\mathbb{M}^{\geq}\rightarrow\mathds{R}\cup\{-\infty\}$ the design criterion we wish to maximize, and by $\lambda_{\min}(\mathbf{M})$ and $\lambda_{\max}(\mathbf{M})$ the minimum and maximum eigenvalues of $\mathbf{M}$, respectively; we shall use the $\ell_{2}$ norm for vectors and Frobenius norm for matrices, $\|\mathbf{M}\|=\operatorname{trace}^{1/2}[\mathbf{M}\mathbf{M}^{\top}]$; all vectors are column vectors. For any $t\in\mathds{R}$, we denote $[t]^{+}=\max\{t,0\}$ and, for any $t\in\mathds{R}^{+}$, $\lfloor t\rfloor$ denotes the largest integer smaller than $t$. For $0\leq\ell\leq L$ we denote by $\mathbb{M}^{\geq}_{\ell,L}$ the (convex) set defined by

$\displaystyle\mathbb{M}^{\geq}_{\ell,L}=\{\mathbf{M}\in\mathbb{M}^{\geq}:\ell<\lambda_{\min}(\mathbf{M})\mbox{ and }\lambda_{\max}(\mathbf{M})<L\}\,,$

and by $\mathbb{M}^{>}$ the open cone of symmetric positive definite $p\times p$ matrices. We make the following assumptions on $\Phi$.

H_Φ

$\Phi$ is strictly concave on $\mathbb{M}^{>}$, linearly differentiable and increasing for Loewner ordering; its gradient $\nabla_{\Phi}(\mathbf{M})$ is well defined in $\mathbb{M}^{\geq}$ for any $\mathbf{M}\in\mathbb{M}^{>}$ and satisfies $\|\nabla_{\Phi}(\mathbf{M})\|<A(\ell)$ and $\lambda_{\min}[\nabla_{\Phi}(\mathbf{M})]>a(L)$ for any $\mathbf{M}\in\mathbb{M}^{\geq}_{\ell,L}$, for some $a(L)>0$ and $A(\ell)<\infty$; moreover, $\nabla_{\Phi}$ satisfies the following Lipschitz condition: for all $\mathbf{M}_{1}$ and $\mathbf{M}_{2}$ in $\mathbb{M}^{\geq}$ such that $\lambda_{\min}(\mathbf{M}_{i})>\ell>0$, $i=1,2$, there exists $K_{\ell}<\infty$ such that $\|\nabla_{\Phi}(\mathbf{M}_{2})-\nabla_{\Phi}(\mathbf{M}_{1})\|<K_{\ell}\,\|\mathbf{M}_{2}-\mathbf{M}_{1}\|$.

The criterion $\Phi_{0}(\mathbf{M})=\log\det(\mathbf{M})$ and criteria $\Phi_{q}(\mathbf{M})=-\operatorname{trace}(\mathbf{M}^{-q})$, $q\in(-1,\infty)$, $q\neq 0$, with $\Phi_{q}(\mathbf{M})=-\infty$ if $\mathbf{M}$ is singular, which are often used in optimal design (in particular with $q$ a positive integer) satisfy H_Φ; see, e.g., Pukelsheim, (1993, Chap. 6). Their gradients are $\nabla_{\Phi_{0}}(\mathbf{M})=\mathbf{M}^{-1}$ and $\nabla_{\Phi_{q}}(\mathbf{M})=q\,\mathbf{M}^{-(q+1)}$, $q\neq 0$; the constants $a(L)$ and $A(\ell)$ are respectively given by $a(L)=1/L$, $A(\ell)=\sqrt{p}/\ell$ for $\Phi_{0}$ and $a(L)=q/L^{q+1}$, $A(\ell)=q\sqrt{p}/\ell^{q+1}$ for $\Phi_{q}$. The Lispchitz condition follows from the fact that the criteria are twice differentiable on $\mathbb{M}^{>}$. The positively homogeneous versions $\Phi_{0}^{+}(\mathbf{M})=\det^{1/p}(\mathbf{M})$ and $\Phi_{q}^{+}(\mathbf{M})=[(1/p)\operatorname{trace}(\mathbf{M}^{-q})]^{-1/q}$, which satisfy $\Phi^{+}(a\mathbf{M})=a\,\Phi^{+}(\mathbf{M})$ for any $a>0$ and any $\mathbf{M}\in\mathbb{M}^{\geq}$, and $\Phi^{+}(\mathbf{I}_{p})=1$, with $\mathbf{I}_{p}$ the $p\times p$ identity matrix, could be considered too; see Pukelsheim, (1993, Chaps. 5, 6). The strict concavity of $\Phi$ implies that, for any convex subset $\widehat{\mathbb{M}}$ of $\mathbb{M}^{>}$, there exists a unique matrix $\mathbf{M}^{*}$ maximizing $\Phi(\mathbf{M})$ with respect to $\mathbf{M}\in\widehat{\mathbb{M}}$.

We denote by $F_{\Phi}(\mathbf{M},\mathbf{M}^{\prime})$ the directional derivative of $\Phi$ at $\mathbf{M}$ in the direction $\mathbf{M}^{\prime}$,

$\displaystyle F_{\Phi}(\mathbf{M},\mathbf{M}^{\prime})=\lim_{\gamma\rightarrow 0^{+}}\frac{\Phi[(1-\gamma)\mathbf{M}+\gamma\mathbf{M}^{\prime}]-\Phi(\mathbf{M})}{\gamma}=\operatorname{trace}[\nabla_{\Phi}(\mathbf{M})(\mathbf{M}^{\prime}-\mathbf{M})]\,,$

and we make the following assumptions on $\mu$ and ${\mathscr{M}}$.

H_μ

$\mu$ has a bounded positive density $\varphi$ with respect to the Lebesgue measure on every open subset of ${\mathscr{X}}$.

H_M

(i) ${\mathscr{M}}$ is continuous on ${\mathscr{X}}$ and satisfies $\int_{\mathscr{X}}\|{\mathscr{M}}(x)\|^{2}\,\mu(\mbox{\rm d}x)<B<\infty$;

(ii) for any ${\mathscr{X}}_{\epsilon}\subset{\mathscr{X}}$ of measure $\mu({\mathscr{X}}_{\epsilon})=\epsilon>0$, $\lambda_{\min}\left\{\int_{{\mathscr{X}}_{\epsilon}}{\mathscr{M}}(x)\,\mu(\mbox{\rm d}x)\right\}>\ell_{\epsilon}$ for some $\ell_{\epsilon}>0$.

Since all the designs considered will be formed by points sampled from $\mu$, we shall confound ${\mathscr{X}}$ with the support of $\mu$: ${\mathscr{X}}=\{x\in\mathds{R}^{d}:\mu({\mathscr{B}}_{d}(x,\epsilon))>0\quad\forall\epsilon>0\}$, with ${\mathscr{B}}_{d}(x,\epsilon)$ the open ball with center $x$ and radius $\epsilon$. Notice that H_M-(i) implies that $\lambda_{\max}[\mathbf{M}(\mu)]<\sqrt{B}$ and $\|\mathbf{M}(\mu)\|<\sqrt{p\,B}$.

Our sequential selection procedure will rely on the estimation of the $(1-\alpha)$-quantile $C_{1-\alpha}(\mathbf{M})$ of the distribution $F_{\mathbf{M}}(z)$ of the directional derivative $Z_{\mathbf{M}}(X)=F_{\Phi}[\mathbf{M},{\mathscr{M}}(X)]$ when $X\sim\mu$, and we shall assume that H_μ,M below is satisfied. It implies in particular that $C_{1-\alpha}(\mathbf{M})$ is uniquely defined by $F_{\mathbf{M}}(C_{\mathbf{M},1-\alpha})=1-\alpha$.

H_μ,M

For all $\mathbf{M}\in\mathbb{M}^{\geq}_{\ell,L}$, $F_{\mathbf{M}}$ has a uniformly bounded density $\varphi_{\mathbf{M}}$; moreover, for any $\alpha\in(0,1)$, there exists $\epsilon_{\ell,L}>0$ such that $\varphi_{\mathbf{M}}[C_{1-\alpha}(\mathbf{M})]>\epsilon_{\ell,L}$ and $\varphi_{\mathbf{M}}$ is continuous at $C_{1-\alpha}(\mathbf{M})$.

H_μ,M is overrestricting (we only need the existence and boundedness of $\varphi_{\mathbf{M}}$, and its positiveness and continuity at $C_{1-\alpha}(\mathbf{M})$), but is satisfied is many common situations; see Section 4 for examples. Let us emphasize that H_μ and H_M are not enough to guarantee the existence of a density $\varphi_{\mathbf{M}}$, since $\operatorname{trace}[\nabla_{\Phi}(\mathbf{M}){\mathscr{M}}(x)]$ may remain constant over subsets of ${\mathscr{X}}$ having positive measure. Assuming the existence of $\varphi_{\mathbf{M}}$ and the continuity of $\varphi$ on ${\mathscr{X}}$ is also insufficient, since $\varphi_{\mathbf{M}}$ is generally not continuous when $Z_{\mathbf{M}}(x)$ is not differentiable in $x$, and $\varphi_{\mathbf{M}}$ is not necessarily bounded.

As mentioned in introduction, when the cardinality of ${\mathscr{X}}_{N}$ is very large, one may wish to select only $n$ candidates $X_{i}$ among the $N$ available, a fraction $n=\lfloor\alpha N\rfloor$ say, with $\alpha\in(0,1)$. For any $n\leq N$, we denote by $\mathbf{M}_{n,N}^{*}$ a design matrix (non necessarily unique) obtained by selecting $n$ points optimally within ${\mathscr{X}}_{N}$; that is, $\mathbf{M}_{n,N}^{*}$ gives the maximum of $\Phi(\mathbf{M}_{n})$ with respect to $\mathbf{M}_{n}=(1/n)\,\sum_{j=1}^{n}{\mathscr{M}}(X_{i_{j}})$, where the $X_{i_{j}}$ are $n$ distinct points in ${\mathscr{X}}_{N}$. Note that this forms a difficult combinatorial problem, unfeasible for large $n$ and $N$. If one assumes that the $X_{i}$ are i.i.d., with $\mu$ their probability measure on ${\mathscr{X}}$, for large $N$ the optimal selection of $n=\lfloor\alpha N\rfloor$ points amounts at constructing an optimal bounded design measure $\xi_{\alpha}^{*}$, such that $\Phi[\mathbf{M}(\xi_{\alpha}^{*})]$ is maximum and $\xi_{\alpha}\leq\mu/\alpha$ (in the sense $\xi_{\alpha}(\mathcal{A})\leq\mu(\mathcal{A})/\alpha$ for any $\mu$-measurable set $\mathcal{A}$, which makes $\xi_{\alpha}$ absolutely continuous with respect to $\mu$). Indeed, Lemma A.1 in Appendix A indicates that $\limsup_{N\rightarrow\infty}\Phi(\mathbf{M}_{\lfloor\alpha N\rfloor,N}^{*})=\Phi[\mathbf{M}(\xi_{\alpha}^{*})]$. Also, under H_Φ, $\mathsf{E}\{\Phi(\mathbf{M}_{n,N}^{*})\}\leq\Phi[\mathbf{M}(\xi_{n/N}^{*})]$ for all $N\geq n>0$; see Pronzato, (2006, Lemma 3).

A key result is that, when all subsets of ${\mathscr{X}}$ with constant $Z_{\mathbf{M}}(x)$ have zero measure, $Z_{\mathbf{M}(\xi_{\alpha}^{*})}(x)=F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]$ separates two sets ${\mathscr{X}}_{\alpha}^{*}$ and ${\mathscr{X}}\setminus{\mathscr{X}}_{\alpha}^{*}$, with $F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]\geq C_{1-\alpha}^{*}$ and $\xi_{\alpha}^{*}=\mu/\alpha$ on ${\mathscr{X}}_{\alpha}^{*}$, and $F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]\leq C_{1-\alpha}^{*}$ and $\xi_{\alpha}^{*}=0$ on ${\mathscr{X}}\setminus{\mathscr{X}}_{\alpha}^{*}$, for some constant $C_{1-\alpha}^{*}$; moreover, $\int_{\mathscr{X}}F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]\,\xi_{\alpha}^{*}(\mbox{\rm d}x)=\int_{{\mathscr{X}}_{\alpha}^{*}}F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]\,\mu(\mbox{\rm d}x)=0$; see Wynn, (1982); Fedorov, (1989) and Fedorov and Hackl, (1997, Chap. 4). (The condition mentioned in those references is that $\mu$ has no atoms, but the example in Section 4.4.2 will show that this is not sufficient; extension to arbitrary measures is considered in (Sahm and Schwabe,, 2001).)

For $\alpha\in(0,1)$, denote

$\displaystyle\mathbb{M}(\alpha)=\left\{\mathbf{M}(\xi_{\alpha})=\int_{\mathscr{X}}{\mathscr{M}}(x)\,\xi_{\alpha}(\mbox{\rm d}x):\,\xi_{\alpha}\in{\mathscr{P}}^{+}({\mathscr{X}}),\,\xi_{\alpha}\leq\frac{\mu}{\alpha},\,\int_{\mathscr{X}}\xi_{\alpha}(\mbox{\rm d}x)=1\right\}\,.$

In (Pronzato,, 2006), it is shown that, for any $\mathbf{M}\in\mathbb{M}^{>}$,

$\displaystyle\mathbf{M}^{+}(\mathbf{M},\alpha)=\arg\max_{\mathbf{M}^{\prime}\in\mathbb{M}(\alpha)}F_{\Phi}(\mathbf{M},\mathbf{M}^{\prime})=\frac{1}{\alpha}\,\int_{\mathscr{X}}\mathbb{I}_{\{F_{\Phi}[\mathbf{M},{\mathscr{M}}(x)]\geq C_{1-\alpha}\}}\,{\mathscr{M}}(x)\,\mu(\mbox{\rm d}x)\,,$

(2.1)

where, for any proposition $\mathcal{A}$, $\mathbb{I}_{\{\mathcal{A}\}}=1$ if $\mathcal{A}$ is true and is zero otherwise, and $C_{1-\alpha}=C_{1-\alpha}(\mathbf{M})$ is an $(1-\alpha)$-quantile of $F_{\Phi}[\mathbf{M},{\mathscr{M}}(X)]$ when $X\sim\mu$ and satisfies

$\displaystyle\int_{\mathscr{X}}\mathbb{I}_{\{F_{\Phi}[\mathbf{M},{\mathscr{M}}(x)]\geq C_{1-\alpha}(\mathbf{M})\}}\,\mu(\mbox{\rm d}x)=\alpha\,.$

(2.2)

Therefore, $\mathbf{M}_{\alpha}^{*}=\mathbf{M}(\xi_{\alpha}^{*})$ is the optimum information matrix in $\mathbb{M}(\alpha)$ (unique since $\Phi$ is strictly concave) if and only if it satisfies $\max_{\mathbf{M}^{\prime}\in\mathbb{M}(\alpha)}F_{\Phi}(\mathbf{M}_{\alpha}^{*},\mathbf{M}^{\prime})=0$, or equivalently $\mathbf{M}_{\alpha}^{*}=\mathbf{M}^{+}(\mathbf{M}_{\alpha}^{*},\alpha)$, and the constant $C_{1-\alpha}^{*}$ equals $C_{1-\alpha}(\mathbf{M}_{\alpha}^{*})$; see (Pronzato,, 2006, Th. 5); see also Pronzato, (2004). Note that $C_{1-\alpha}^{*}\leq 0$ since $\int_{\mathscr{X}}F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]\,\xi_{\alpha}^{*}(\mbox{\rm d}x)=0$ and $F_{\Phi}[\mathbf{M}(\xi_{\alpha}^{*}),{\mathscr{M}}(x)]\geq C_{1-\alpha}^{*}$ on the support of $\xi_{\alpha}^{*}$.

Suppose that the $X_{i}$ are i.i.d. with $\mu$. The solution of $\mathbf{M}=\mathbf{M}^{+}(\mathbf{M},\alpha)$, $\alpha\in(0,1)$, with respect to $\mathbf{M}$ by stochastic approximation yields the iterations

$\displaystyle\begin{array}[]{rcl}n_{k+1}&=&n_{k}+\mathbb{I}_{\{F_{\Phi}[\mathbf{M}_{n_{k}},{\mathscr{M}}(X_{k+1})]\geq C_{1-\alpha}(\mathbf{M}_{n_{k}})\}}\,,\\ \\ \mathbf{M}_{n_{k+1}}&=&\mathbf{M}_{n_{k}}+\frac{1}{n_{k}+1}\,\mathbb{I}_{\{F_{\Phi}[\mathbf{M}_{n_{k}},{\mathscr{M}}(X_{k+1})]\geq C_{1-\alpha}(\mathbf{M}_{n_{k}})\}}\,\left[{\mathscr{M}}(X_{k+1})-\mathbf{M}_{n_{k}}\right]\,.\end{array}$

(3.4)

Note that $\mathsf{E}\left\{\mathbb{I}_{\{F_{\Phi}[\mathbf{M},{\mathscr{M}}(X)]\geq C_{1-\alpha}(\mathbf{M})\}}\,\left[{\mathscr{M}}(X)-\mathbf{M}\right]\right\}=\alpha\,[\mathbf{M}^{+}(\mathbf{M},\alpha)-\mathbf{M}]$. The almost sure (a.s.) convergence of $\mathbf{M}_{n_{k}}$ in (3.4) to $\mathbf{M}(\xi_{\alpha}^{*})$ that maximizes $\Phi(\mathbf{M})$ with respect $\mathbf{M}\in\mathbb{M}(\alpha)$ is proved in (Pronzato,, 2006) under rather weak assumptions on $\Phi$, ${\mathscr{M}}$ and $\mu$.

The construction (3.4) requires the calculation of the $(1-\alpha)$-quantile $C_{1-\alpha}(\mathbf{M}_{n_{k}})$ for all $n_{k}$, see (2.2), which is not feasible when $\mu$ is unknown and has a prohibitive computational cost when we know $\mu$. For that reason, it is proposed in (Pronzato,, 2006) to replace $C_{1-\alpha}(\mathbf{M}_{n_{k}})$ by the empirical quantile $\widetilde{C}_{\alpha,k}(\mathbf{M}_{n_{k}})$ that uses the empirical measure $\mu_{k}=(1/k)\,\sum_{i=1}^{k}\delta_{X_{i}}$ of the $X_{i}$ that have been observed up to stage $k$. This construction preserves the a.s. convergence of $\mathbf{M}_{n_{k}}$ to $\mathbf{M}(\xi_{\alpha}^{*})$ in (3.4), but its computational cost and storage requirement increase with $k$, which makes it unadapted to situations with very large $N$. The next section considers the recursive estimation of $C_{1-\alpha}(\mathbf{M}_{n_{k}})$ and contains the main result of the paper.

The idea is to plug a recursive estimator of the $(1-\alpha)$-quantile $C_{1-\alpha}(\mathbf{M}_{n_{k}})$ in (3.4). Under mild assumptions, for random variables $Z_{i}$ that are i.i.d. with distribution function $F$ such that the solution of the equation $F(z)=1-\alpha$ is unique, the recursion

$\displaystyle\widehat{C}_{k+1}=\widehat{C}_{k}+\frac{\beta}{k+1}\,\left(\mathbb{I}_{\{Z_{k+1}\geq\widehat{C}_{k}\}}-\alpha\right)$

(3.5)

with $\beta>0$ converges a.s. to the quantile $C_{1-\alpha}$ such that $F(C_{1-\alpha})=1-\alpha$. Here, we shall use a construction based on (Tierney,, 1983). In that paper, a clever dynamical choice of $\beta=\beta_{k}$ is shown to provide the optimal asymptotic rate of convergence of $\widehat{C}_{k}$ towards $C_{1-\alpha}$, with $\sqrt{k}(\widehat{C}_{k}-C_{1-\alpha})\stackrel{{\scriptstyle\rm d}}{{\rightarrow}}{\mathscr{N}}(0,\alpha(1-\alpha)/f^{2}(C_{1-\alpha}))$ as $k\rightarrow\infty$, where $f(z)=\mbox{\rm d}F(z)/\mbox{\rm d}z$ is the p.d.f. of the $Z_{i}$ — note that it coincides with the asymptotic behavior of the sample (empirical) quantile. The only conditions on $F$ are that $f(z)$ exists for all $z$ and is uniformly bounded, and that $f$ is continuous and positive at the unique root $C_{1-\alpha}$ of $F(z)=1-\alpha$.

There is a noticeable difference, however, with the estimation of $C_{1-\alpha}(\mathbf{M}_{n_{k}})$: in our case we need to estimate a quantile of $Z_{k}(X)=F_{\Phi}[\mathbf{M}_{n_{k}},{\mathscr{M}}(X)]$ for $X\sim\mu$, with the distribution of $Z_{k}(X)$ evolving with $k$. For that reason, we shall impose a faster dynamic to the evolution of $\widehat{C}_{k}$, and replace (3.5) by

$\displaystyle\widehat{C}_{k+1}=\widehat{C}_{k}+\frac{\beta_{k}}{(k+1)^{q}}\,\left(\mathbb{I}_{\{Z_{k}(X_{k+1})\geq\widehat{C}_{k}\}}-\alpha\right)$

(3.6)

for some $q\in(0,1)$. The combination of (3.6) with (3.4) yields a particular nonlinear two-time-scale stochastic approximation scheme. There exist advanced results on the convergence of linear two-time-scale stochastic approximation, see Konda and Tsitsiklis, (2004); Dalal et al., (2018). To the best of our knowledge, however, there are few results on convergence for nonlinear schemes. Convergence is shown in (Borkar,, 1997) under the assumption of boundedness of the iterates using the ODE method of Ljung, (1977); sufficient conditions for stability are provided in (Lakshminarayanan and Bhatnagar,, 2017), also using the ODE approach. In the proof of Theorem 3.1 we provide justifications for our construction, based on the analyses and results in the references mentioned above.

The construction is summarized in Algorithm 1 below. The presence of the small number $\epsilon_{1}$ is only due to technical reasons: setting $z_{k+1}=+\infty$ when $n_{k}/k<\epsilon_{1}$ in (3.9) has the effect of always selecting $X_{k+1}$ when less than $\epsilon_{1}\,k$ points have been selected previously; it ensures that $n_{k+1}/k>\epsilon_{1}$ for all $k$ and thus that $\mathbf{M}_{n_{k}}$ always belongs to $\mathbb{M}^{\geq}_{\ell,L}$ for some $\ell>0$ and $L<\infty$; see Lemma B.1 in Appendix.

• Algorithm 1: sequential selection ($\alpha$ given).

• 0)

  Choose $k_{0}\geq p$, $q\in(1/2,1)$, $\gamma\in(0,q-1/2)$, and $0<\epsilon_{1}\ll\alpha$.

• 1)

  Initialization: select $X_{1},\ldots,X_{k_{0}}$, compute $\mathbf{M}_{n_{k_{0}}}=(1/k_{0})\,\sum_{i=1}^{k_{0}}{\mathscr{M}}(X_{i})$. If $\mathbf{M}_{n_{k_{0}}}$ is singular, increase $k_{0}$ and select the next points until $\mathbf{M}_{n_{k_{0}}}$ has full rank. Set $k=n_{k}=k_{0}$, the number of points selected.

  Compute $\zeta_{i}=Z_{k_{0}}(X_{i})$, for $i=1,\ldots,k_{0}$ and order the $\zeta_{i}$ as $\zeta_{1:k_{0}}\leq\zeta_{2:k_{0}}\leq\cdots\leq\zeta_{k_{0}:k_{0}}$; denote $k_{0}^{+}=\lceil(1-\alpha/2)\,k_{0}\rceil$ and $k_{0}^{-}=\max\{\lfloor(1-3\,\alpha/2)\,k_{0}\rfloor,1\}$.

  Initialize $\widehat{C}_{k_{0}}$ at $\zeta_{\lceil(1-\alpha)\,k_{0}\rceil:k_{0}}$; set $\beta_{0}=k_{0}/(k_{0}^{+}-k_{0}^{-})$, $h=(\zeta_{k_{0}^{+}:k_{0}}-\zeta_{k_{0}^{-}:k_{0}})$, $h_{k_{0}}=h/k_{0}^{\gamma}$ and $\widehat{f}_{k_{0}}=\left[\sum_{i=1}^{k_{0}}\mathbb{I}_{\{|\zeta_{i}-\widehat{C}_{k_{0}}|\leq h_{k_{0}}\}}\right]/(2\,k_{0}\,h_{k_{0}})$.

• 2)

  Iteration $k+1$: collect $X_{k+1}$ and compute $Z_{k}(X_{k+1})=F_{\Phi}[\mathbf{M}_{n_{k}},{\mathscr{M}}(X_{k+1})]$.

  $\displaystyle\begin{array}[]{ll}\mbox{If }n_{k}/k<\epsilon_{1}&\mbox{ set }z_{k+1}=+\infty\,;\\ \mbox{otherwise }&\mbox{ set }z_{k+1}=Z_{k}(X_{k+1})\,.\end{array}$

  (3.9)

  If $z_{k+1}\geq\widehat{C}_{k}$, update $n_{k}$ into $n_{k+1}=n_{k}+1$ and $\mathbf{M}_{n_{k}}$ into

  $\displaystyle\mathbf{M}_{n_{k+1}}=\mathbf{M}_{n_{k}}+\frac{1}{n_{k}+1}\,\left[{\mathscr{M}}(X_{k+1})-\mathbf{M}_{n_{k}}\right]\,;$

  (3.10)

  otherwise, set $n_{k+1}=n_{k}$.

• 3)

  Compute $\beta_{k}=\min\{1/\widehat{f}_{k},\beta_{0}\,k^{\gamma}\}$; update $\widehat{C}_{k}$ using (3.6).

  Set $h_{k+1}=h/(k+1)^{\gamma}$ and update $\widehat{f}_{k}$ into

  $\displaystyle\widehat{f}_{k+1}=\widehat{f}_{k}+\frac{1}{(k+1)^{q}}\,\left[\frac{1}{2\,h_{k+1}}\,\mathbb{I}_{\{|Z_{k}(X_{k+1})-\widehat{C}_{k}|\leq h_{k+1}\}}-\widehat{f}_{k}\right]\,.$

• 4)

  $k\leftarrow k+1$, return to Step 2.

Note that $\widehat{C}_{k}$ is updated whatever the value of $Z_{k}(X_{k+1})$. Recursive quantile estimation by (3.6) follows (Tierney,, 1983). To ensure a faster dynamic for the evolution of $\widehat{C}_{k}$ than for $\mathbf{M}_{n_{k}}$, we take $q<1$ instead of $q=1$ in (Tierney,, 1983), and the construction of $\widehat{f}_{k}$ and the choices of $\beta_{k}$ and $h_{k}$ are modified accordingly. Following the same arguments as in the proof of Proposition 1 of (Tierney,, 1983), the a.s. convergence of $\widehat{C}_{k}$ to $C_{1-\alpha}$ in the modified version of (3.5) is proved in Theorem C.1 (Appendix C).

The next theorem establishes the convergence of the combined stochastic approximation schemes with two time-scales.

The o.d.e. associated with (3.6), for a fixed matrix $\mathbf{M}$ and thus a fixed $Z(\cdot)$, such that $Z(X)=F_{\Phi}[\mathbf{M},{\mathscr{M}}(X)]$ has the distribution function $F$ and density $f$, is

$\displaystyle\frac{\mbox{\rm d}C(t)}{\mbox{\rm d}t}=1-F[C(t)]-\alpha=F(C_{1-\alpha})-F[C(t)]\,,$

where $C_{1-\alpha}=C_{1-\alpha}(\mathbf{M})$ satisfies $F(C_{1-\alpha})=1-\alpha$. Consider the Lyapunov function $L(C)=[F(C)-F(C_{1-\alpha})]^{2}$. It satisfies $\mbox{\rm d}L[C(t)]/\mbox{\rm d}t=-2\,f[C(t)]\,L[C(t)]\leq 0$, with $\mbox{\rm d}L[C(t)]/\mbox{\rm d}t=0$ if and only if $C=C_{1-\alpha}$. Moreover, $C_{1-\alpha}$ is Lipschitz continuous in $\mathbf{M}$; see Lemma D.1 in Appendix D. The conditions for Theorem 1.1 in (Borkar,, 1997) are thus satisfied concerning the iterations for $\widehat{C}_{k}$.

Denote $\widehat{\mathbf{M}}_{k}=\mathbf{M}_{n_{k}}$ and $\rho_{k}=n_{k}/k$, so that (3.9) implies $k\,\rho_{k}\geq\epsilon_{1}\,(k-1)$ for all $k$; see Lemma B.1 in Appendix. They satisfy

$\displaystyle\rho_{k+1}=\rho_{k}+\frac{R_{k+1}}{k+1}\mbox{ and }\widehat{\mathbf{M}}_{k+1}=\widehat{\mathbf{M}}_{k}+\frac{\boldsymbol{\Omega}_{k+1}}{k+1}\,,$

(3.11)

where $R_{k+1}=\mathbb{I}_{\{Z_{k}(X_{k+1})\geq\widehat{C}_{k}\}}-\rho_{k}$, and $\boldsymbol{\Omega}_{k+1}=(1/\rho_{k+1})\,\mathbb{I}_{\{Z_{k}(X_{k+1})\geq\widehat{C}_{k}\}}\,\left[{\mathscr{M}}(X_{k+1})-\widehat{\mathbf{M}}_{k}\right]$. We have $\mathsf{E}\{R_{k+1}|{\mathscr{F}}_{k}\}=\int_{\mathscr{X}}\mathbb{I}_{\{Z_{k}(x)\geq\widehat{C}_{k}\}}\,\mu(\mbox{\rm d}x)-\rho_{k}$ and $\mathsf{var}\{R_{k+1}|{\mathscr{F}}_{k}\}=F_{k}(\widehat{C}_{k})[1-F_{k}(\widehat{C}_{k})]$, with $F_{k}$ the distribution function of $Z_{k}(X)$, which, from H_μ,M, has a well defined density $f_{k}$ for all $k$. Also,

$\displaystyle\mathsf{E}\{\boldsymbol{\Omega}_{k+1}|{\mathscr{F}}_{k}\}$

$\displaystyle=$

$\displaystyle\frac{{\mathcal{I}}_{k}}{\rho_{k}+\frac{1-\rho_{k}}{k+1}}=\frac{1+1/k}{\rho_{k}+1/k}\,{\mathcal{I}}_{k}\,,$

where ${\mathcal{I}}_{k}=\int_{\mathscr{X}}\mathbb{I}_{\{Z_{k}(x)\geq\widehat{C}_{k}\}}\,\left[{\mathscr{M}}(x)-\widehat{\mathbf{M}}_{k}\right]\,\mu(\mbox{\rm d}x)$. Denote $\Delta_{k+1}=\boldsymbol{\Omega}_{k+1}-{\mathcal{I}}_{k}/\rho_{k}$, so that

$\displaystyle\widehat{\mathbf{M}}_{k+1}=\widehat{\mathbf{M}}_{k}+\frac{1}{k+1}\,\frac{{\mathcal{I}}_{k}}{\rho_{k}}+\frac{\Delta_{k+1}}{k+1}\,.$

(3.12)

We get $\mathsf{E}\{\Delta_{k+1}|{\mathscr{F}}_{k}\}=(\rho_{k}-1)/[\rho_{k}\,(k\,\rho_{k}+1)]\,{\mathcal{I}}_{k}$ and

	$\displaystyle\mathsf{var}\{\{\Delta_{k+1}\}_{i,j}|{\mathscr{F}}_{k}\}$	$\displaystyle=$	$\displaystyle\mathsf{var}\{\{\boldsymbol{\Omega}_{k+1}\}_{i,j}|{\mathscr{F}}_{k}\}$
		$\displaystyle=$	$\displaystyle\frac{(k+1)^{2}}{(k\,\rho_{k}+1)^{2}}\,\left[\int_{\mathscr{X}}\mathbb{I}_{\{Z_{k}(x)\geq\widehat{C}_{k}\}}\,\{{\mathscr{M}}(x)-\widehat{\mathbf{M}}_{k}\}_{i,j}^{2}\,\mu(\mbox{\rm d}x)-\{{\mathcal{I}}_{k}\}_{i,j}^{2}\right]\,,$

where (3.9) implies that $\rho_{k}>\epsilon_{1}/2$, and therefore $(k+1)\,(1-\rho_{k})/[\rho_{k}\,(k\,\rho_{k}+1)]<4/\epsilon_{1}^{2}$, and $\mathsf{var}\{\{\Delta_{k+1}\}_{i,j}|{\mathscr{F}}_{k}\}$ is a.s. bounded from H_M-(i). This implies that $\sum_{k}\Delta_{k+1}/(k+1)<\infty$ a.s. The limiting o.d.e. associated with (3.11) and (3.12) are

	$\displaystyle\frac{\mbox{\rm d}\rho(t)}{\mbox{\rm d}t}$	$\displaystyle=$	$\displaystyle\int_{\mathscr{X}}\mathbb{I}_{\{F_{\Phi}[\widehat{\mathbf{M}}(t),{\mathscr{M}}(x)]\geq C_{1-\alpha}[\widehat{\mathbf{M}}(t)]\}}\,\mu(\mbox{\rm d}x)-\rho(t)=\alpha-\rho(t)\,,$
	$\displaystyle\frac{\mbox{\rm d}\widehat{\mathbf{M}}(t)}{\mbox{\rm d}t}$	$\displaystyle=$	$\displaystyle\frac{1}{\rho(t)}\,\int_{\mathscr{X}}\mathbb{I}_{\{F_{\Phi}[\widehat{\mathbf{M}}(t),{\mathscr{M}}(x)]\geq C_{1-\alpha}[\widehat{\mathbf{M}}(t)]\}}\,\left[{\mathscr{M}}(x)-\widehat{\mathbf{M}}(t)\right]\,\mu(\mbox{\rm d}x)$
		$\displaystyle=$	$\displaystyle\frac{\alpha}{\rho(t)}\,\left\{\mathbf{M}^{+}[\widehat{\mathbf{M}}(t),\alpha]-\widehat{\mathbf{M}}(t)\right\}\,,$

where $\mathbf{M}^{+}[\widehat{\mathbf{M}}(t),\alpha]$ is defined by (2.1). The first equation implies that $\rho(t)$ converges exponentially fast to $\alpha$, with $\rho(t)=\alpha+[\rho(0)-\alpha]\,\exp(-t)$; the second equation gives

$\displaystyle\frac{\mbox{\rm d}\Phi[\widehat{\mathbf{M}}(t)]}{\mbox{\rm d}t}=\operatorname{trace}\left[\nabla_{\Phi}[\widehat{\mathbf{M}}(t)]\,\frac{\mbox{\rm d}\widehat{\mathbf{M}}(t)}{\mbox{\rm d}t}\right]=\frac{\alpha}{\rho(t)}\,\max_{\mathbf{M}^{\prime}\in\mathbb{M}(\alpha)}F_{\Phi}[\widehat{\mathbf{M}}(t),\mathbf{M}^{\prime}]\geq 0\,,$

with a strict inequality if $\widehat{\mathbf{M}}(t)\neq\mathbf{M}_{\alpha}^{*}$, the optimal matrix in $\mathbb{M}(\alpha)$. The conditions of Theorem 1.1 in (Borkar,, 1997) are thus satisfied, and $\widehat{\mathbf{M}}_{k}$ converges to $\mathbf{M}_{\alpha}^{*}$ a.s.