\section{Introduction}
Consider the problem of maximizing an expensive-to-compute black-box objective function $f$ that depends on an \emph{optimization variable} $\mbf{x}$ and an \emph{environmental random variable} $\mbf{Z}$.
Due to the randomness in $\mbf{Z}$, the function evaluation $f(\mbf{x},\mbf{Z})$ of $f$ at $\mbf{x}$ is a random variable.
Though for such an objective function $f$, \emph{Bayesian optimization} (BO) can be naturally applied to maximize its expectation $\mbb{E}_{\mbf{Z}}[f(\mbf{x},\mbf{Z})]$ over $\mbf{Z}$ \cite{toscano2018bayesian}, this maximization objective overlooks the \emph{risks} of potentially undesirable function evaluations.
These risks can arise from either (a) the realization of an unknown distribution of $\mbf{Z}$ or (b) the realization of the random $\mbf{Z}$ given that the distribution of $f(\mbf{x},\mbf{Z})$ can be estimated well or that of $\mbf{Z}$ is known.
The issue (a) has been tackled by distributionally robust BO \cite{kirschner2020distributionally,nguyen2020distributionally} which maximizes $\mbb{E}_{\mbf{Z}}[f(\mbf{x},\mbf{Z})]$ under the worst-case realization of the distribution of $\mbf{Z}$.
To resolve the issue (b), 
the risk from the uncertainty in $\mbf{Z}$ can be controlled via the mean-variance optimization framework \cite{iwazaki2020mean}, \emph{value-at-risk} ({\textsc{VaR}}), or \emph{conditional value-at-risk} ({\textsc{CVaR}}) \cite{borisk20,torossian2020bayesian}.
The work of \citet{bogunovic2018adversarially} has considered \emph{adversarially robust BO},
where $\mbf{z}$ is controlled by an adversary deterministically.\footnote{We use upper-case letter $\mbf{Z}$ to denote the environmental random variable and lower-case letter $\mbf{z}$ to denote its realization or a (non-random) variable.}  
In this case, the objective is to find $\mbf{x}$ that maximizes the function under the worst-case realization of $\mbf{z}$, i.e., $\text{arg}\!\max_{\mbf{x}} \min_{\mbf{z}} f(\mbf{x},\mbf{z})$.


In this paper, we focus on case (b) where the distribution of 
$\mbf{Z}$ is known (or well-estimated). 
For example, in agriculture, although farmers cannot control the temperature of an outdoor farm, 
its distribution
can be estimated from historical data and controlled in an indoor environment for optimizing the plant yield.
Given the distribution of $\mbf{Z}$, the objective is to control the risk
that the function evaluation $f(\mbf{x},\mbf{z})$, 
for a $\mbf{z}$ sampled from $\mbf{Z}$,
is small.
One popular framework is to control the trade-off between the mean (viewed as reward) and the variance (viewed as risk) of the function evaluation with respect to $\mbf{Z}$ \cite{iwazaki2020mean}.
However, quantifying the risk 
using variance implies indifference between positive and negative deviations from the mean, while people often have asymmetric risk attitudes \cite{goh2012portfolio}. In our problem of maximizing the objective function, it is reasonable to assume that people  are risk-averse towards only the negative deviations 
from the mean, i.e., 
the risk of getting lower function evaluations. 
Thus, it is more appropriate to adopt
risk measures 
with this asymmetric property, such as \emph{value-at-risk} ({\textsc{VaR}})
which is a widely adopted risk measure in real-world applications (e.g., banking \cite{basel06}).
Intuitively, the risk 
that the random $f(\mbf{x},\mbf{Z})$ is less than {\textsc{VaR}} at level $\alpha \in (0,1)$ does not exceed $\alpha$, 
e.g., by specifying a small value of $\alpha$ as $0.1$, this risk is controlled to be at most $10\%$.
Therefore, to maximize the function $f$ while controlling the risk of undesirable (i.e., small) function evaluations,
we aim to maximize {\textsc{VaR}} of the random function $f(\mbf{x},\mbf{Z})$ over $\mbf{x}$.

The recent work of \citet{borisk20} has used BO to maximize {\textsc{VaR}} and has achieved state-of-the-art empirical performances.
They have assumed that we are able to select both $\mbf{x}$ and $\mbf{z}$ to query during BO,
which is motivated by fact that physical experiments can usually be studied by 
simulation \cite{williams2000sequential}. 
In the example on agriculture given above, we can control the temperature, light and water ($\mbf{z}$) in a small indoor environment to optimize the amount of fertilizer ($\mbf{x}$), which can then be used in an outdoor environment with random weather factors.
\citet{borisk20} have exploited the ability to select $\mbf{z}$ to model the function $f(\mbf{x},\mbf{z})$ as a GP, which allows them to retain the appealing closed-form posterior belief of the objective function.
To select the queries $\mbf{x}$ and $\mbf{z}$, they have designed a one-step lookahead approach based on the well-known \emph{knowledge gradient} (KG) acquisition function \cite{scott2011correlated}.
However, the one-step lookahead incurs an expensive nested optimization procedure, which is computationally expensive and hence requires approximations.
Besides, the acquisition function can only be approximated using samples of the objective function $f$ from the GP posterior and the environmental random variable $\mbf{Z}$.
While they have analysed the asymptotically unbiased and consistent estimator of the gradients, it is challenging to obtain a guarantee for the convergence of their algorithm.
Another recent work \cite{torossian2020bayesian} has also applied BO to maximize {\textsc{VaR}} using an asymmetric Laplace likelihood function and variational approximation of the posterior belief. However, in contrast to \citet{borisk20} and our work, they have focused on a different setting where the realizations of $\mbf{Z}$ are not observed.





In this paper, we adopt the setting of \citet{borisk20} which allows us to choose both $\mbf{x}$ and $\mbf{z}$ to query, and assume that the distribution of $\mbf{Z}$ is known or well-estimated. 
Our contributions include:


\textbf{Firstly}, we propose a novel BO algorithm named \emph{Value-at-risk Upper Confidence Bound} ({V-UCB}) in Section~\ref{sec:vy}.
Unlike the work of \citet{borisk20}, {V-UCB} is equipped with a no-regret convergence guarantee and is more computationally efficient.
To guide its query selection and facilitate its proof of the no-regret guarantee, the classical GP-UCB algorithm~\cite{srinivas10ucb} constructs a \emph{confidence bound} of the objective function. Similarly, to maximize the {\textsc{VaR}} of a random function, we, for the first time to the best of our knowledge, construct a confidence bound of {\textsc{VaR}} (Lemma~\ref{lemma:confbound}). 
The resulting confidence bound of {\textsc{VaR}} naturally gives rise to a strategy to select $\mbf{x}$.
However, it remains a major challenge to select $\mbf{z}$ to preserve the no-regret convergence of GP-UCB.
To this end, we firstly prove that our algorithm is no-regret as long as we ensure that at the selected $\mbf{z}$, the confidence bound of {\textsc{VaR}} \emph{lies within} the confidence bound of the objective function.
Next, we also prove that this query selection strategy is \emph{feasible}, i.e., such values of $\mbf{z}$, referred to as \emph{lacing values} (LV), exist

\textbf{Secondly}, although our theoretical no-regret property allows the selection of \emph{any} LV, we design a heuristic to select an LV such that it improves our empirical performance over random selection of LV (Section~\ref{subsec:improvedselectz}).
We also discuss the implications when $\mbf{z}$ cannot be selected by BO and is instead randomly sampled by the environment during BO (Remark~\ref{rmk:unknownZ}).
\textbf{Thirdly}, we show that adversarially robust BO \cite{bogunovic2018adversarially} can be cast as a special case of our V-UCB
when the risk level $\alpha$ of {\textsc{VaR}} approaches $0$ from the right and the domain of $\mbf{z}$ is the support of $\mbf{Z}$.
In this case, 
adversarially robust BO \cite{bogunovic2018adversarially} selects the same input queries as those selected by V-UCB since the set of LV collapse into the set of minimizers of the lower bound of the objective function (Section~\ref{subsec:vy-stableopt}).
\textbf{Lastly}, we provide practical techniques for
implementing {V-UCB} with continuous random variable $\mbf{Z}$ (Section~\ref{sec:continuousz}): we (a) introduce \emph{local neural surrogate optimization} with the \emph{pinball loss} to optimize {\textsc{VaR}}, and (b) construct an objective function to search for an LV in the continuous support of $\mbf{Z}$.

The performance of our proposed algorithm is empirically demonstrated in optimizing several synthetic benchmark functions, a portfolio optimization problem, and a simulated robot task in Section~\ref{sec:experiments}.



\section{Problem Statement and Background}

Let the objective function be defined as $f: \mcl{D}_{\mbf{x}} \times \mcl{D}_{\mbf{z}} \rightarrow \mbb{R}$ where $\mcl{D}_{\mbf{x}} \subset \mbb{R}^{d_x}$ and $\mcl{D}_{\mbf{z}}\subset \mbb{R}^{d_z}$ are the bounded domain of the optimization variable $\mbf{x}$ and the support of the environmental random variable $\mbf{Z}$, respectively; $d_x$ and $d_z$ are the dimensions of $\mbf{x}$ and $\mbf{z}$, respectively. The support of $\mbf{Z}$ is defined as the smallest closed subset $\mcl{D}_{\mbf{z}}$ of $\mbb{R}^{d_z}$ such that $P(\mbf{Z} \in \mcl{D}_{\mbf{z}}) = 1$.
Let $\mbf{z} \in \mcl{D}_{\mbf{z}}$ denote a realization of the random variable $\mbf{Z}$. 
Let $f(\mbf{x},\mbf{Z})$ denote a random variable whose randomness comes from $\mbf{Z}$.
The {\textsc{VaR}} of $f(\mbf{x},\mbf{Z})$ at \emph{risk level} $\alpha \in (0,1)$ is defined as:
\begin{equation}
V_{\alpha}(f(\mbf{x},\mbf{Z})) \triangleq \inf \{\omega: P(f(\mbf{x},\mbf{Z}) \le \omega) \ge \alpha\}
\label{eq:var}
\end{equation}
which implies the risk that $f(\mbf{x},\mbf{Z})$ is less than its {\textsc{VaR}} at level $\alpha$ does not exceed $\alpha$.


Our objective is to search for $\mbf{x} \in \mcl{D}_{\mbf{x}}$ that maximizes $V_\alpha(f(\mbf{x},\mbf{Z}))$ at a user-specified risk level $\alpha \in (0,1)$. 
Intuitively, the goal is find $\mbf{x}$ where the evaluations of the objective function are as large as possible under most realizations of the environmental random variable $\mbf{Z}$ which is characterized by the probability of $1-\alpha$.

The unknown objective function $f(\mbf{x},\mbf{z})$ is modeled with a GP. That is, every finite subset of $\{f(\mbf{x},\mbf{z})\}_{(\mbf{x}, \mbf{z}) \in \mcl{D}_{\mbf{x}} \times \mcl{D}_{\mbf{z}}}$ follows a multivariate Gaussian distribution \cite{rasmussen06}.
The GP is fully specified by its \emph{prior} mean and covariance function $k_{(\mbf{x},\mbf{z}), (\mbf{x}',\mbf{z}')} \triangleq \text{cov}[f(\mbf{x},\mbf{z}), f(\mbf{x}',\mbf{z}')]$ for all $\mbf{x}, \mbf{x}'$ in $\mcl{D}_{\mbf{x}}$ and $\mbf{z}, \mbf{z}'$ in $\mcl{D}_{\mbf{z}}$. For notational simplicity (and w.l.o.g.), the former is assumed to be zero, while we use the \emph{squared exponential} (SE) kernel as its bounded maximum information gain can be used for later analysis \cite{srinivas10ucb}.

To identify the optimal $\mbf{x}_* \triangleq \text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_{\alpha}(f(\mbf{x},\mbf{Z}))$, BO algorithm selects an input query $(\mbf{x}_t,\mbf{z}_t)$ in the $t$-th iteration to obtain a noisy function evaluation $y_{(\mbf{x}_t,\mbf{z}_t)} \triangleq f(\mbf{x}_t,\mbf{z}_t) + \epsilon_t$ where $\epsilon_t \sim \mcl{N}(0,\sigma_n^2)$ are i.i.d. Gaussian noise with variance $\sigma_n^2$. Given noisy observations $\mbf{y}_{\mcl{D}_t} \triangleq (y_{(\mbf{x},\mbf{z})})_{(\mbf{x},\mbf{z}) \in \mcl{D}_t}^\top$ at observed inputs $\mcl{D}_t \triangleq \mcl{D}_{t-1} \cup \{(\mbf{x}_t, \mbf{z}_t)\}$ (and $\mcl{D}_0$ is the initial observed inputs), the GP posterior belief of function evaluation at any input $(\mbf{x},\mbf{z})$ is a Gaussian $p(f(\mbf{x},\mbf{z})|\mbf{y}_{\mcl{D}_t}) \triangleq \mcl{N}( f(\mbf{x},\mbf{z})|\mu_t(\mbf{x},\mbf{z}), \sigma_t^2(\mbf{x},\mbf{z}))$:
\begin{equation}
\begin{array}{r@{}l}
\mu_t(\mbf{x},\mbf{z}) &\triangleq \mbf{K}_{(\mbf{x},\mbf{z}),\mcl{D}_t} \bm{\Lambda}_{\mcl{D}_t\mcl{D}_t} \mbf{y}_{\mcl{D}_t}\ ,\\
\sigma_t^2(\mbf{x},\mbf{z}) &\triangleq k_{(\mbf{x},\mbf{z}),(\mbf{x},\mbf{z})}  - \mbf{K}_{(\mbf{x},\mbf{z}),\mcl{D}_t} \bm{\Lambda}_{\mcl{D}_t\mcl{D}_t} \mbf{K}_{\mcl{D}_t,(\mbf{x},\mbf{z})}
\label{eq:gppost}
\end{array}
\end{equation}
where $\bm{\Lambda}_{\mcl{D}_t\mcl{D}_t} \triangleq \left(
	\mbf{K}_{\mcl{D}_t\mcl{D}_t} + \sigma_n^2 \mbf{I}
\right)^{-1}$, 
$\mbf{K}_{(\mbf{x},\mbf{z}),\mcl{D}_t} \triangleq (k_{(\mbf{x},\mbf{z}),(\mbf{x}',\mbf{z}')})_{(\mbf{x}',\mbf{z}') \in \mcl{D}_t}$, 
$\mbf{K}_{\mcl{D}_t,(\mbf{x},\mbf{z})} \triangleq \mbf{K}_{(\mbf{x},\mbf{z}),\mcl{D}_t}^\top$,
$\mbf{K}_{\mcl{D}_t\mcl{D}_t} \triangleq (k_{(\mbf{x}',\mbf{z}'), (\mbf{x}'',\mbf{z}'')})_{(\mbf{x}',\mbf{z}'), (\mbf{x}'',\mbf{z}'') \in \mcl{D}_t}$,  $\mbf{I}$ is the identity matrix.


\section{BO of {\textsc{VaR}}}
\label{sec:vy}

Following the seminal work \cite{srinivas10ucb}, we use the \emph{cumulative regret} as the performance metric to quantify the performance of our BO algorithm. It is defined as $R_T \triangleq \sum_{t=1}^T r(\mbf{x}_t)$ where $r(\mbf{x}_t) \triangleq V_{\alpha}(f(\mbf{x}_*,\mbf{Z})) - V_{\alpha}(f(\mbf{x}_t,\mbf{Z}))$ is the \emph{instantaneous regret} and $\mbf{x}_* \triangleq \text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_{\alpha}(f(\mbf{x},\mbf{Z}))$. We would like to design a query selection strategy that incurs \emph{no regret}, i.e., $\lim_{T\rightarrow \infty} R_T / T = 0$. 
Furthermore, we have that 
$\min_{t \le T} r(\mbf{x}_t) \le R_T / T$, equivalently, $\max_{t \le T} V_{\alpha}(f(\mbf{x}_t,\mbf{Z})) \ge V_{\alpha}(f(\mbf{x}_*,\mbf{Z})) - R_T / T$.
Thus, $\text{lim}_{T \rightarrow \infty} \max_{t \le T} V_{\alpha}(f(\mbf{x}_t,\mbf{Z})) = V_{\alpha}(f(\mbf{x}_*,\mbf{Z}))$ for a no-regret algorithm.

The proof of the upper bound on the cumulative regret of GP-UCB is based on confidence bounds of the objective function \cite{srinivas10ucb}.
Similarly, in the next section, we start by constructing a confidence bound of $V_{\alpha}(f(\mbf{x},\mbf{Z}))$, which naturally leads to a query selection strategy for $\mbf{x}_t$.
\subsection{A Confidence Bound of $V_\alpha(f(\mbf{x},\mbf{Z}))$ and the Query Selection Strategy for $\mbf{x}_t$}
\label{subsec:confboundf}
Firstly, we adopt a confidence bound of the function $f(\mbf{x},\mbf{z})$ from \citet{chowdhury2017kernelized}, which assumes that $f$ belongs to a \emph{reproducing kernel Hilbert space} $\mcl{F}_k(B)$ such that its RKHS norm is bounded $\Vert f\Vert_k \le B$.

\begin{lemma}[\citet{chowdhury2017kernelized}]
Pick $\delta \in (0,1)$ and set 
$\beta_t = (B + \sigma_n \sqrt{2(\gamma_{t-1} + 1 + \log1/ \delta)})^2$. 
Then,
$f(\mbf{x},\mbf{z}) \in I_{t-1}[f(\mbf{x},\mbf{z})] \triangleq [l_{t-1}(\mbf{x},\mbf{z}), u_{t-1}(\mbf{x},\mbf{z})]$
$\forall \mbf{x} \in \mcl{D}_{\mbf{x}}, \mbf{z} \in \mcl{D}_{\mbf{z}}, t \ge 1 $ holds with probability $\ge 1 - \delta$ where
\begin{equation}
\begin{array}{r@{}l}
l_{t-1}(\mbf{x},\mbf{z})\ &\triangleq \mu_{t-1}(\mbf{x},\mbf{z}) - \beta_t^{1/2} \sigma_{t-1}(\mbf{x},\mbf{z})\\
u_{t-1}(\mbf{x},\mbf{z})\ &\triangleq \mu_{t-1}(\mbf{x},\mbf{z}) + \beta_t^{1/2} \sigma_{t-1}(\mbf{x},\mbf{z})\ .
\end{array}\vspace{-1mm}
\label{eq:fbound}
\end{equation}
\label{lemma:ucb51}
\end{lemma}
As the above lemma holds for both finite and continuous $\mcl{D}_x$ and $\mcl{D}_z$, it is 
used
to analyse the regret in both cases. On the other hand, the confidence bound can be adopted to the Bayesian setting by changing only $\beta_t$ following the work of \citet{srinivas10ucb} as noted by \cite{bogunovic2018adversarially}.

Then, we exploit this confidence bound on the function evaluations (Lemma \ref{lemma:ucb51}) to formulate a confidence bound of $V_{\alpha}(f(\mbf{x},\mbf{Z}))$ as follows.
\begin{lemma}
\label{lemma:confbound}
Similar to the definition of $f(\mbf{x},\mbf{Z})$, let $l_{t-1}(\mbf{x}, \mbf{Z})$ and $u_{t-1}(\mbf{x},\mbf{Z})$ denote the random function over $\mbf{x}$ where the randomness comes from the random variable $\mbf{Z}$; $l_{t-1}$ and $u_{t-1}$ are defined in \eqref{eq:fbound}.
Then, $\forall \mbf{x} \in \mcl{D}_{\mbf{x}}$, $t \ge 1$,
\[
\begin{array}{r@{}l}
\displaystyle V_\alpha(f(\mbf{x},\mbf{Z})) \displaystyle &\displaystyle \in I_{t-1}[V_{\alpha}(f(\mbf{x},\mbf{Z}))]\\
	 &\displaystyle \triangleq [V_\alpha(l_{t-1}(\mbf{x},\mbf{Z})), V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))]
\end{array}
\]
holds with probability $\ge 1 - \delta$ for $\delta$ in Lemma~\ref{lemma:ucb51}, where $V_\alpha(l_{t-1}(\mbf{x},\mbf{Z}))$ and $V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$ are defined as \eqref{eq:var}.
\end{lemma}
The proof is in Appendix~\ref{app:proofvconfbound}.
Given the confidence bound $I_{t-1}[V_{\alpha}(f(\mbf{x},\mbf{Z}))] \triangleq [V_\alpha(l_{t-1}(\mbf{x},\mbf{Z})), V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$ in Lemma~\ref{lemma:confbound}, we follow the the well-known 
``optimism in the face of uncertainty'' 
principle to select $\mbf{x}_t = \text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$. This query selection strategy for $\mbf{x}_t$ leads to an upper bound of $r(\mbf{x}_t)$:
\begin{align}
r(\mbf{x}_t) \le V_\alpha(u_{t-1}(\mbf{x}_t,\mbf{Z})) - V_\alpha(l_{t-1}(\mbf{x}_t,\mbf{Z}))\ \forall t \ge 1
	\label{eq:iregretbound1}
\end{align}
which holds with probability $\ge 1 - \delta$ for $\delta$ in Lemma~\ref{lemma:ucb51}, and is proved in Appendix~\ref{app:iregretbound1}.


As our goal is $\lim_{T\rightarrow \infty} R_T/T = 0$, given the selected query $\mbf{x}_t$, a reasonable query selection strategy of $\mbf{z}_t$ should gather informative observations at $(\mbf{x}_t, \mbf{z}_t)$ that improves the confidence bound $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$ (i.e., $I_{t}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$ is a proper subset of $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$ if $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))] \neq \emptyset$) which can be viewed as the uncertainty of $V_{\alpha}(f(\mbf{x}_t,\mbf{Z}))$.


Assume that there exists
$\mbf{z}_l \in \mcl{D}_{\mbf{z}}$ such that  $l_{t-1}(\mbf{x}_t,\mbf{z}_l) = V_\alpha(l_{t-1}(\mbf{x}_t,\mbf{Z}))$
and
$\mbf{z}_u \in \mcl{D}_{\mbf{z}}$ such that  $u_{t-1}(\mbf{x}_t,\mbf{z}_u) = V_\alpha(u_{t-1}(\mbf{x}_t,\mbf{Z}))$.
Lemma~\ref{lemma:confbound} implies that $V_\alpha(f(\mbf{x}_t,\mbf{Z})) \in I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))] = [l_{t-1}(\mbf{x}_t,\mbf{z}_l), u_{t-1}(\mbf{x}_t,\mbf{z}_u)]$ with high probability. 
Hence, we may 
na\"ively
want to query for observations at $(\mbf{x}_t, \mbf{z}_l)$ and $(\mbf{x}_t,\mbf{z}_u)$ to reduce $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$. 
However, these observations may not always reduce $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$. 
The insight is that $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$ changes (i.e., shrinks) when either of its boundary values (i.e., $l_{t-1}(\mbf{x}_t,\mbf{z}_l)$ or $u_{t-1}(\mbf{x}_t,\mbf{z}_u)$) changes. 
Consider $u_{t-1}(\mbf{x}_t,\mbf{z}_u)$ and finite $\mcl{D}_{\mbf{z}}$ as an example, since $u_{t-1}(\mbf{x}_t,\mbf{z}_u) = V_\alpha(u_{t-1}(\mbf{x}_t,\mbf{Z}))$, 
a natural cause for the change in $u_{t-1}(\mbf{x}_t,\mbf{z}_u)$ is when $\mbf{z}_u$ changes. 
This happens if there exists $\mbf{z}' \neq \mbf{z}_u$ such that the \emph{ordering} of $u_{t-1}(\mbf{x}_t,\mbf{z}')$ relative to $u_{t-1}(\mbf{x}_t,\mbf{z}_u)$ changes given more observations. Thus, observations that are capable of reducing $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$ should be able to \emph{change the relative ordering} in this case.
We construct the following counterexample where observations at $\mbf{z}_u$ (and $\mbf{z}_l$) are not able to change the relative ordering, so they do not reduce $I_{t-1}[V_\alpha (f(\mbf{x}_t,\mbf{Z}))]$.

\begin{example}
\label{example:counter}
This example is described by Fig.~\ref{fig:orderinguncertainty}. 
We reduce notational clutter by removing $\mbf{x}_t$ and $t$ since they are fixed in this example, i.e., we use $f(\mbf{z})$, $f(\mbf{Z})$, and $l(\mbf{z})$ to denote $f(\mbf{x}_t, \mbf{z})$, $f(\mbf{x}_t, \mbf{Z})$, and $l_{t-1}(\mbf{x}_t, \mbf{z})$ respectively. 
We condition on the event $f(\mbf{z}) \in I[f(\mbf{z})] \triangleq [l(\mbf{z}), u(\mbf{z})]$ for all $\mbf{z} \in \mcl{D}_{\mbf{z}}$ which occurs with probability $\ge 1 - \delta$ in Lemma~\ref{lemma:ucb51}. In this example, $\mbf{z}_l = \mbf{z}_1$ and $l(\mbf{z}_1) = u(\mbf{z}_1)$, so there is no uncertainty in $f(\mbf{z}_l) = f(\mbf{z}_1)$. Similarly, there is no uncertainty in $f(\mbf{z}_u) = f(\mbf{z}_2)$. Thus, new observations at $\mbf{z}_l$ and $\mbf{z}_u$ change neither $l(\mbf{z}_l)$ nor $u(\mbf{z}_u)$, so these observations do not reduce the confidence bound $I[V_{\alpha=0.4}(f(\mbf{Z}))] = [l(\mbf{z}_l), u(\mbf{z}_u)]$ (plotted as the double-headed arrow in Fig.~\ref{fig:orderinguncertainty}b).
In fact, to reduce $I[V_{\alpha=0.4}(f(\mbf{Z}))]$, we should gather new observations at $\mbf{z}_0$ which potentially change the ordering of $u(\mbf{z}_0)$ relative to $u(\mbf{z}_2)$ (which is $u(\mbf{z}_u)$ without new observations). For example, after 
getting
new observations at $\mbf{z}_0$, if $u(\mbf{z}_0)$ is improved to be in the white region between A and B ($u(\mbf{z}_0) > u(\mbf{z}_2)$ in Fig.~\ref{fig:orderinguncertainty}b changes to $u(\mbf{z}_0) < u(\mbf{z}_2)$ in Fig.~\ref{fig:orderinguncertainty}c), then $I[V_{\alpha=0.4}(f(\mbf{Z}))]$ is reduced to $[l(\mbf{z}_1), u(\mbf{z}_0)]$ 
because now $\mbf{z}_u=\mbf{z}_0$.
Thus, as the confidence bound $I[f(\mbf{z}_0)]$ is shortened with more and more observations at $\mbf{z}_0$, the confidence bound $I[V_{\alpha=0.4}(f(\mbf{Z}))]$ reduces  (the white region in Fig.~\ref{fig:orderinguncertainty} is `laced up').
\end{example}
\begin{figure*}[h!]
\centering
\begin{tabular}{@{}ccc@{}}
\includegraphics[width=0.25\textwidth]{img/uncertainty/ordering-uncertainty-top}
&
\includegraphics[width=0.25\textwidth]{img/uncertainty/ordering-uncertainty-bottom}
&
\includegraphics[width=0.25\textwidth]{img/uncertainty/ordering-uncertainty-extra}
\\
(a) & (b) & (c)
\end{tabular}
\caption{
A counterexample against selecting $\mbf{z}_u$ and $\mbf{z}_l$ as input queries.
Here $\mbf{z}$ follows a discrete uniform distribution over $\mcl{D}_{\mbf{z}} \triangleq \{\mbf{z}_0, \mbf{z}_1, \mbf{z}_2\}$.
(a) shows the mappings of $\mbf{z}$ to the upper bound $u(\mbf{z})$ and lower bound $l(\mbf{z})$.
The {\textsc{VaR}} at $\alpha=0.4$ of $u(\mbf{Z})$ and $l(\mbf{Z})$ are $u(\mbf{z}_2)$ and $l(\mbf{z}_1)$, respectively, i.e., $\mbf{z}_u = \mbf{z}_2$ and $\mbf{z}_l = \mbf{z}_1$.
(b) shows the values of $l(\mbf{z})$ and $u(\mbf{z})$ for all $\mbf{z}$ on the same axis, as well as the confidence bounds of $f(\mbf{z})$ and $V_{\alpha}(f(\mbf{Z}))$. The gray areas A and B indicate the intervals of values $\omega \in \mbb{R}$ where $\omega \le l(\mbf{z}_l) = l(\mbf{z}_1)$ and $\omega \ge u(\mbf{z}_u) = u(\mbf{z}_2)$, respectively. 
(c) shows a hypothetical scenario when $I[f(\mbf{z}_0)]$ is shortened with more observations at $\mbf{z}_0$.
}
\label{fig:orderinguncertainty}
\end{figure*}

In the next section, we define a property of $\mbf{z}_0$ in Example~\ref{example:counter} 
and prove the existence of $\mbf{z}$'s with this property.
Then, we prove that along with the optimistic selection of $\mbf{x}_t$, the selection of $\mbf{z}_t$ such that it satisfies this property leads to a no-regret algorithm.

\subsection{Lacing Value (LV) and the Query Selection Strategy for $\mbf{z}_t$}
\label{subsec:lv}

We note that in Example~\ref{example:counter}, as long as the confidence bound of the function evaluation at $\mbf{z}_0$ contains 
the confidence bound
of {\textsc{VaR}}, observations at $\mbf{z}_0$ can 
reduce the confidence bound of {\textsc{VaR}}. We name the values of $\mbf{z}$ satisfying this property as \emph{lacing values} (LV):






\begin{definition}[Lacing values]
\label{definition:lv}
\emph{Lacing values} (LV) with respect to $\mbf{x} \in \mcl{D}_{\mbf{x}}$ and $t \ge 1$ are $\mbf{z}_{\text{LV}} \in \mcl{D}_{\mbf{z}}$ that satisfies $l_{t-1}(\mbf{x},\mbf{z}_{\text{LV}}) \le V_{\alpha}(l_{t-1}(\mbf{x}, \mbf{Z})) \le V_{\alpha}(u_{t-1}(\mbf{x}, \mbf{Z})) \le u_{t-1}(\mbf{x},\mbf{z}_{\text{LV}})$, equivalently,
$I_{t-1}[V_{\alpha}(f(\mbf{x},\mbf{Z}))] \subset [l_{t-1}(\mbf{x},\mbf{z}_{\text{LV}}), u_{t-1}(\mbf{x},\mbf{z}_{\text{LV}})]\ .$
\end{definition}

Recall that the support $\mcl{D}_{\mbf{z}}$ of $\mbf{Z}$ is defined as the smallest closed subset $\mcl{D}_{\mbf{z}}$ of $\mbb{R}^{d_z}$ such that $P(\mbf{Z} \in \mcl{D}_{\mbf{z}}) = 1$ (e.g., $\mcl{D}_{\mbf{z}}$ is a finite subset of $\mbb{R}^{d_z}$ and $\mcl{D}_{\mbf{z}} = \mbb{R}^{d_z}$). 
The following theorem guarantees the existence of lacing values and is proved in Appendix~\ref{app:prooflv}.

\begin{theorem}[Existence of lacing values]
\label{theorem:lv}
$\forall \alpha \in (0,1)$, $\forall \mbf{x} \in \mcl{D}_{\mbf{x}}$, $\forall t \ge 1$, there exists a lacing value in $\mcl{D}_{\mbf{z}}$ with respect to $\mbf{x}$ and $t$.
\end{theorem}

\begin{corollary}
\label{corollary:loclv}
Lacing values with respect to $\mbf{x} \in \mcl{D}_{\mbf{x}}$ and $t \ge 1$ are in $\mcl{Z}_l^{\le} \cap \mcl{Z}_u^{\ge}$ where $\mcl{Z}_l^{\le} \triangleq \{\mbf{z} \in \mcl{D}_{\mbf{z}}: l_{t-1}(\mbf{x},\mbf{z}) \le V_{\alpha}(l_{t-1}(\mbf{x}, \mbf{Z}))\}$ and $\mcl{Z}_u^{\ge} \triangleq \{\mbf{z} \in \mcl{D}_{\mbf{z}}: u_{t-1}(\mbf{x},\mbf{z}) \ge V_{\alpha}(u_{t-1}(\mbf{x}, \mbf{Z}))\}$.
\end{corollary}

As a special case, when $\mbf{z}_l = \mbf{z}_u$,  $I_{t-1}[V_\alpha(f(\mbf{x},\mbf{Z}))] = I_{t-1}[f(\mbf{x},\mbf{z}_l)]$ which means $\mbf{z}_l = \mbf{z}_u$ is an LV.
Based on Theorem~\ref{theorem:lv}, we can always select $\mbf{z}_t$ as an LV w.r.t $\mbf{x}_t$ defined in Definition~\ref{definition:lv}.
This strategy for the selection of $\mbf{z}_t$, together with the selection of $\mbf{x}_t = \text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$ (Section \ref{subsec:confboundf}), completes our algorithm: \emph{{\textsc{VaR}} Upper Confidence Bound} ({V-UCB}) (Algorithm~\ref{alg:v-ucb}).

\textbf{Upper Bound on Regret.}
As a result of the selection strategies for $\mbf{x}_t$ and $\mbf{z}_t$, our V-UCB algorithm enjoys the following upper bound on its instantaneous regret (proven in Appendix~\ref{app:iregretbound2}):
\begin{lemma}
\label{lemma:iregretbound2}
By selecting $\mbf{x}_t$ as a maximizer of $V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$ and selecting $\mbf{z}_t$ as an LV w.r.t $\mbf{x}_t$, the instantaneous regret is upper-bounded by:
\begin{equation*}
r(\mbf{x}_t) \le 2 \beta_t^{1/2} \sigma_{t-1}(\mbf{x}_t, \mbf{z}_t)\ \forall t \ge 1
\end{equation*}
with probability $\ge 1 - \delta$ for $\delta$ in Lemma~\ref{lemma:ucb51}.
\end{lemma}
Lemma~\ref{lemma:iregretbound2}, together with Lemma 5.4 from \citet{srinivas10ucb}, implies that the cumulative regret of our algorithm is bounded (Appendix~\ref{app:rtbound}): $R_T \le \sqrt{C_1 T \beta_T \gamma_T}$ 
where $C_1 \triangleq 8/\log(1 + \sigma_n^{-2})$, and $\gamma_T$ is the maximum information gain about $f$ that can be obtained from any set of $T$ observations.
\citet{srinivas10ucb} have analyzed $\gamma_T$ for several commonly used kernels such as SE and Mat\'ern kernels, and have shown that for these kernels, the upper bound on $R_T$ grows sub-linearly.
This implies that our algorithm is \emph{no-regret} because $\lim_{T\rightarrow \infty} R_T / T = 0$.



\begin{algorithm}[tb]
   \caption{The {V-UCB} Algorithm}
\begin{algorithmic}[1]
   \STATE {\bfseries Input:} $\mcl{D}_{\mbf{x}}$, $\mcl{D}_{\mbf{z}}$, prior $\mu_0=0, \sigma_0, k$
   \FOR{$i=1,2,\dots$}
   \STATE Select $\mbf{x}_t = \text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$
   \STATE Select $\mbf{z}_t$ as a \emph{lacing value} w.r.t. $\mbf{x}_t$ (Definition~\ref{definition:lv})
   \STATE Obtain observation $y_t \triangleq f(\mbf{x}_t,\mbf{z}_t) + \epsilon_t$
   \STATE Update the GP posterior belief to obtain $\mu_t$ and $\sigma_t$
   \ENDFOR
\end{algorithmic}
\label{alg:v-ucb}
\end{algorithm}

Inspired by \citet{bogunovic2018adversarially}, at the {$T$-th} iteration of {V-UCB}, we can recommend $\mbf{x}_{t_*(T)}$ as an estimate of the maximizer $\mbf{x}_*$ of $V_\alpha(f(\mbf{x},\mbf{Z}))$, where $t_*(T) \triangleq  \text{arg}\!\max_{t \in \{1,\dots,T\}} V_{\alpha}(l_{t-1}(\mbf{x}_t,\mbf{Z}))$. Then, the instantaneous regret $r(\mbf{x}_{t_*(T)})$ is upper-bounded by $\sqrt{C_1 \beta_T \gamma_T / T}$ with high probability as we show in Appendix~\ref{app:recommendxbound}.
In our experiments in Section~\ref{sec:experiments}, 
we recommend $\text{arg}\!\max_{\mbf{x} \in \mcl{D}_T} V_{\alpha}(\mu_{t-1}(\mbf{x}, \mbf{Z}))$ (where $\mu_{t-1}(\mbf{x}, \mbf{Z})$ is a random function defined in the same manner as $f(\mbf{x},\mbf{Z})$) as an estimate of $\mbf{x}_*$ due to its empirical convergence.

\textbf{Computational Complexity.}
To compare our computational complexity
with that of the {$\rho \text{KG}$} algorithm from \citet{borisk20}, we exclude the common part of training the GP model (line 6) and assume that $\mcl{D}_{\mbf{z}}$ is finite. 
Then, the time complexity of {V-UCB} is dominated by that of the selection of $\mbf{x}_t$ (line 3) which includes the time complexity $\mcl{O}(|\mcl{D}_{\mbf{z}}| |\mcl{D}_{t-1}|^2)$ for the GP prediction at $\{\mbf{x}\} \times \mcl{D}_{\mbf{z}}$, and $\mcl{O}(|\mcl{D}_{\mbf{z}}|\log|\mcl{D}_{\mbf{z}}|)$ for the sorting of $u_{t-1}(\mbf{x},\mcl{D}_{\mbf{z}})$ and searching of {\textsc{VaR}}.
Hence, our overall complexity is $\mcl{O}(n |\mcl{D}_{\mbf{z}}|\ (|\mcl{D}_{t-1}|^2 + \log |\mcl{D}_{\mbf{z}}|))$, where $n$ is the number of iterations to maximize $V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$ (line 3). 
Therefore, our V-UCB is more computationally efficient than {$\rho \text{KG}$} and its variant with approximation {$\rho \text{KG}^{apx}$}, whose complexities are $\mcl{O}(n_{\text{out}} n_{\text{in}} K |\mcl{D}_{\mbf{z}}|\ (|\mcl{D}_{t-1}|^2 + |\mcl{D}_{\mbf{z}}| |\mcl{D}_{t-1}| + |\mcl{D}_{\mbf{z}}|^2 + M |\mcl{D}_{\mbf{z}}|))$ of {$\rho \text{KG}$} and $\mcl{O}(n_{\text{out}} |\mcl{D}_{t-1}| K |\mcl{D}_{\mbf{z}}|\ (|\mcl{D}_{t-1}|^2 + |\mcl{D}_{\mbf{z}}| |\mcl{D}_{t-1}| + |\mcl{D}_{\mbf{z}}|^2 + M |\mcl{D}_{\mbf{z}}|))$, respectively.\footnote{$n_{\text{out}}$ and $n_{\text{in}}$ are the numbers of iterations for the outer and inner optimization respectively, $K$ is the number of fantasy GP models required for their one-step lookahead, and $M$ is the number of functions sampled from the GP posterior \cite{borisk20}.}



\subsection{On the Selection of $\mbf{z}_t$}
\label{subsec:improvedselectz}

Although Algorithm~\ref{alg:v-ucb} is guaranteed to be no-regret with any choice of LV as $\mbf{z}_t$, we would like to select the LV that can reduce a large amount of the uncertainty of $V_{\alpha}(f(\mbf{x}_t,\mbf{Z}))$. 
However, relying on the information gain measure or the knowledge gradient method often incurs the expensive one-step lookahead.
Therefore, we use a simple heuristic by choosing the LV $\mbf{z}_{\text{LV}}$ with the maximum probability mass (or probability density if $\mbf{Z}$ is continuous) of $\mbf{z}_{\text{LV}}$.
We motivate this heuristic using an example with $\alpha = 0.2$, i.e., $V_{\alpha=0.2}(f(\mbf{x}_t,\mbf{Z})) = \inf \{\omega: P(f(\mbf{x}_t,\mbf{Z}) \le \omega) \ge 0.2\}$.
Suppose $\mcl{D}_{\mbf{z}}$ is finite and there are $2$ LV's $\mbf{z}_0$ and $\mbf{z}_1$ with $P(\mbf{z}_0) \ge 0.2$ and $P(\mbf{z}_1) = 0.01$. 
Then, because $P(f(\mbf{x}_t,\mbf{Z}) \le f(\mbf{x}_t,\mbf{z}_0)) \ge P(\mbf{z}_0) \ge 0.2$, it follows that $V_{\alpha=0.2}(f(\mbf{x}_t,\mbf{Z})) \le f(\mbf{x}_t,\mbf{z}_0)$, i.e., querying $\mbf{z}_0$ at $\mbf{x}_t$ gives us information about an explicit upper bound on $V_{\alpha=0.2}(f(\mbf{x}_t,\mbf{Z}))$ to reduce its uncertainty.
In contrast, this cannot be achieved by querying $\mbf{z}_1$ due to its low probability mass.
This simple heuristic can also be implemented when $\mbf{Z}$ is a continuous random variable which we will introduce in Section~\ref{sec:continuousz}.
\begin{remark}
\label{rmk:unknownZ}
Although we assume that we can select both $\mbf{x}_t$ and $\mbf{z}_t$ during our
algorithm, Corollary~\ref{corollary:loclv} also gives us some insights about the scenario where we cannot select $\mbf{z}_t$.
In this case, in each iteration $t$, we select $\mbf{x}_t$ while $\mbf{z}_t$ is randomly sampled by the environment following the distribution of $\mbf{Z}$. 
Next, we observe both $\mbf{z}_t$ and $y_{(\mbf{x}_t,\mbf{z}_t)}$ and then update the GP posterior belief of $f$.
Of note, Corollary~\ref{corollary:loclv} has shown that all LV lie in the set $\mcl{Z}_l^{\le} \cap \mcl{Z}_u^{\ge}$. 
However, the probability of this set is usually small, because $P(\mbf{Z} \in \mcl{Z}_l^{\le} \cap \mcl{Z}_u^{\ge}) \le P(\mbf{Z} \in \mcl{Z}_l^{\le}) \le \alpha$ and small values of $\alpha$ are often used by real-world applications to manage risks.
Thus, the probability that the sampled $\mbf{z}_t$ is an LV is small. 
As a result, we suggest sampling a large number of $\mbf{z}_t$'s from the environment to increase the chance that an LV is sampled.
On the other hand, the small probability of sampling an LV 
motivates the need for us to select $\mbf{z}_t$.
\end{remark}

\subsection{{V-UCB} Approaches {\textsc{StableOpt}} as $\alpha \rightarrow 0^+$}
\label{subsec:vy-stableopt}

Recall that the objective of adversarially robust BO is to find $\mbf{x} \in \mcl{D}_{\mbf{x}}$ that maximizes
$\min_{\mbf{z} \in \mcl{D}_{\mbf{z}}} f(\mbf{x}, \mbf{z})$ \cite{bogunovic2018adversarially} by iteratively specifying input query $(\mbf{x}_t, \mbf{z}_t)$ to collect noisy observations $y_{\mbf{x}_t,\mbf{z}_t}$. 
It is different from BO of {\textsc{VaR}} since its $\mbf{z}$ is not random but selected by an adversary who aims to minimize the function evaluation. 
The work of \citet{bogunovic2018adversarially} has proposed a no-regret algorithm for this setting named {\textsc{StableOpt}}, which selects
\begin{equation}
\begin{array}{r@{}l}
\mbf{x}_t &= \text{arg}\!\max_{\mbf{x} \in \mcl{D}_x} \min_{\mbf{z} \in \mcl{D}_z} u_{t-1}(\mbf{x},\mbf{z})\ ,\\
\mbf{z}_t &= \text{arg}\!\min_{\mbf{z} \in \mcl{D}_z} l_{t-1}(\mbf{x}_t,\mbf{z})
\end{array}
\label{eq:stableopt}
\end{equation}
where $u_{t-1}$ and $l_{t-1}$ are defined in \eqref{eq:fbound}.

At first glance, BO of {\textsc{VaR}} and adversarially robust BO are seemingly different problems because $\mbf{Z}$ is a random variable in the former but not in the latter. 
However, based on our key observation on the connection between the minimum value of a continuous function $h(\mbf{w})$ and the {\textsc{VaR}} of the random variable $h(\mbf{W})$ in the following theorem, these two problems and their solutions are connected as shown in Corollary~\ref{corollary:alpha0lv}, and~\ref{corollary:stableopt}.

\begin{theorem}
\label{theorem:0plus}
Let $\mbf{W}$ be a random variable with the support $\mcl{D}_w \subset \mbb{R}^{d_w}$ and dimension $d_w$.
Let $h$ be a continuous function mapping from $\mbf{w} \in \mcl{D}_w$ to $\mbb{R}$. Then, $h(\mbf{W})$ denotes the random variable whose realization is the function $h$ evaluation at a realization $\mbf{w}$ of $\mbf{W}$. Suppose $h(\mbf{w})$ has a minimizer $\mbf{w}_{\min} \in \mcl{D}_w$, then $\lim_{\alpha \rightarrow 0^+} V_{\alpha}(h(\mbf{W})) = h(\mbf{w}_{\min})\ .$
\end{theorem}

\begin{corollary}
\label{corollary:alpha0lv}
Adversarially robust BO which finds $\text{arg}\!\max_{\mbf{x}} \min_{\mbf{z}} f(\mbf{x},\mbf{z})$ can be cast as BO of {\textsc{VaR}} by letting (a) $\alpha$ approach $0$ from the right and (b) $\mcl{D}_{\mbf{z}}$ be the support of the environmental random variable $\mbf{Z}$, i.e., $\text{arg}\!\max_{\mbf{x}} \lim_{\alpha \rightarrow 0^+} V_\alpha(f(\mbf{x},\mbf{Z}))$.
\end{corollary}

Interestingly, from Theorem~\ref{theorem:0plus}, we observe that $\mcl{Z}_l^\le$ in Corollary~\ref{corollary:loclv} 
approaches 
the set of minimizers  $\text{arg}\!\min_{\mbf{z} \in \mcl{D}_{\mbf{z}}} l_{t-1}(\mbf{x}_t,\mbf{z})$ as $\alpha \rightarrow 0^+$. 
Corollary \ref{corollary:stableopt} below shows that LV w.r.t $\mbf{x}_t$ becomes a minimizer of $l_{t-1}(\mbf{x}_t,\mbf{z})$ which is same as the selected $\mbf{z}_t$ of {\textsc{StableOpt}} in \eqref{eq:stableopt}.

\begin{corollary}
\label{corollary:stableopt}
The {\textsc{StableOpt}} solution to adversarially robust BO selects the same input query as that selected by the {V-UCB} solution to the corresponding BO of {\textsc{VaR}} in Corollary~\ref{corollary:alpha0lv}.
\end{corollary}

The proof of Theorem~\ref{theorem:0plus} and its two corollaries are shown in Appendix~\ref{app:0plus}.
We note that {V-UCB} is also applicable to the optimization of $V_{\alpha}(f(\mbf{x},\mbf{Z}))$ where the distribution of $\mbf{Z}$ is a conditional distribution given $\mbf{x}$. For example,
in robotics, if there exists noise/perturbation in the control, an optimization problem of interest is
$V_{\alpha}(f(\mbf{x} + \bm{\xi}(\mbf{x})))$ where $\bm{\xi}(\mbf{x})$ is the random perturbation that depends on $\mbf{x}$.



\subsection{Implementation of V-UCB with Continuous Random Variable $\mbf{Z}$}
\label{sec:continuousz}

The V-UCB algorithm involves two steps: selecting $\mbf{x}_t = \text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ and selecting $\mbf{z}_t$ as the LV $\mbf{z}_{\text{LV}}$ with the largest probability mass (or probability density). When $|\mcl{D}_{\mbf{z}}|$ is finite, given $\mbf{x}$,  $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ can be computed exactly. 
The gradient of $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ with respect to $\mbf{x}$ can be obtained 
easily (e.g., using
automatic differentiation provided in the Tensorflow library \cite{tensorflow2015-whitepaper}) 
to aid the selection of $\mbf{x}_t$.
In this case, the latter step can also be performed by constructing the set of all LV (checking the condition in the Definition~\ref{definition:lv} for all $\mbf{z} \in \mcl{D}_{\mbf{z}}$) and choosing the LV $\mbf{z}_{\text{LV}}$ with the largest probability mass.

{\bf{Estimation of {\textsc{VaR}}.}} When $\mbf{Z}$ is a continuous random variable, estimating {\textsc{VaR}} by an ordered set of samples (e.g., in \citet{borisk20}) may require a prohibitively large number of samples, especially for small values of $\alpha$. Thus, we employ the following popular pinball (or tilted absolute value) function in quantile regression \cite{koenker1978regression} to estimate {\textsc{VaR}} as a lower $\alpha$-quantile:
\begin{equation*}
\rho_{\alpha}(w) \triangleq \begin{cases}
    \alpha w &\text{if } w \ge 0\ ,\\
    (\alpha - 1) w &\text{if } w < 0
\end{cases}
\end{equation*}
where $w \in \mbb{R}$. In particular, to estimate $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ as $\nu \in \mbb{R}$, we find $\nu$ that minimizes:
\begin{equation}
    \mbb{E}_{\mbf{z} \sim p(\mbf{Z})}[\rho_{\alpha}(u_{t-1}(\mbf{x},\mbf{z}) - \nu)]\ .
    \label{eq:lossfindvar}
\end{equation}
A well-known example is when $\alpha = 0.5$ and $\rho_{\alpha}$ is the absolute value function, then the optimal $\nu$ is the median. 
The loss in \eqref{eq:lossfindvar} can be optimized using stochastic gradient descent with a random batch of samples of $\mbf{Z}$ at each optimization iteration. 

{\bf{Maximization of $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$.}} Unfortunately, to maximize $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ over $\mbf{x} \in \mcl{D}_{\mbf{x}}$, there is no gradient of $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ with respect to $\mbf{x}$ under the above approach. This situation resembles BO where there is no gradient information, but only noisy observations at input queries. Unlike BO, the observation (samples of $u_{t-1}(\mbf{x},\mbf{Z})$ at $\mbf{x}$) is not costly. Therefore, we propose the \emph{local neural surrogate optimization} (LNSO) algorithm to find $\text{arg}\!\max_{\mbf{x} \in \mcl{D}_{\mbf{x}}} V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ which is visualized in Fig.~\ref{fig:lnso}. Suppose the optimization is initialized at $\mbf{x} = \mbf{x}^{(0)}$, instead of maximizing $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ (whose gradient w.r.t. $\mbf{x}$ is unknown), we maximize a surrogate function $g(\mbf{x},\bm{\theta}^{(0)})$ (modeled by a neural network) that approximates $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ well in a local region of $\mbf{x}^{(0)}$, e.g., a ball $\mcl{B}(\mbf{x}^{(0)}, r)$ of radius $r$ in Fig.~\ref{fig:lnso}. 
We obtain such parameters $\bm{\theta}^{(0)}$ by minimizing the following loss function:
\begin{equation}
\begin{array}{l}
\mcl{L}_g(\bm{\theta},\mbf{x}^{(0)})\\ 
\ \displaystyle\triangleq \mbb{E}_{\mbf{x} \in \mcl{B}(\mbf{x}^{(0)}, r)} \mbb{E}_{\mbf{Z} \sim p(\mbf{Z})} \left[ \rho_\alpha(u_{t-1}(\mbf{x},\mbf{z}) - g(\mbf{x};\bm{\theta}))\right]
\end{array}
\label{eq:losslnso}
\end{equation}
where the expectation $\mbb{E}_{\mbf{x} \in \mcl{B}(\mbf{x}^{(0)}, r)}$ is taken over a uniformly distributed $\mbf{X}$ in $\mcl{B}(\mbf{x}^{(0)}, r)$. Minimizing \eqref{eq:losslnso} can be performed with stochastic gradient descent.
If maximizing $g(\mbf{x},\bm{\theta}^{(0)})$ leads to a value $\mbf{x}^{(i)} \notin \mcl{B}(\mbf{x}^{(0)}, r)$ (Fig.~\ref{fig:lnso}), we update the local region to be centered at $\mbf{x}^{(i)}$ ($\mcl{B}(\mbf{x}^{(i)}, r)$) and find  $\bm{\theta}^{(i)} = \text{arg}\!\min_{\bm{\theta}} \mcl{L}_g(\bm{\theta}, \mbf{x}^{(i)})$ such that $g(\mbf{x},\bm{\theta}^{(i)})$ approximates $V_{\alpha}(u_{t-1}(\mbf{x},\mbf{Z}))$ well $\forall \mbf{x} \in \mcl{B}(\mbf{x}^{(i)},r)$. Then, $\mbf{x}^{(i)}$ is updated by maximizing $g(\mbf{x},\bm{\theta}^{(i)})$ for $\mbf{x} \in \mcl{B}(\mbf{x}^{(i)},r)$. The complete algorithm is described in Appendix~\ref{app:lnso}.

We prefer a small value of $r$ so that the ball $\mcl{B}$ is small. In such case, $V_\alpha(u_{t-1}(\mbf{x},\mbf{Z}))$ for $\mbf{x} \in \mcl{B}$ can be estimated well with a small neural network $g(\mbf{x},\bm{\theta})$ whose training requires a small number of iterations. 

\begin{figure}
    \centering
    \includegraphics[width=0.27\textwidth]{img/lnso.pdf}
    \caption{Plot of a hypothetical optimization path (as arrows) of LNSO initialized at $\mbf{x}^{(0)}$. Input $\mbf{x}$ is $2$-dimensional. The boundary of a ball $\mcl{B}$ of radius $r$ is plotted as a dotted circle. When the updated $\mbf{x}$ moves out of $\mcl{B}$, the center of $\mcl{B}$ and $\bm{\theta}$ are updated.}
    \label{fig:lnso}
\end{figure}

{\bf{Search of Lacing Values.}} 
Given a continuous random variable $\mbf{Z}$, to find an LV w.r.t $\mbf{x}_t$ in line 4 of Algorithm~\ref{alg:v-ucb}, i.e., to find a $\mbf{z}$ satisfying $d_u(\mbf{z}) \triangleq u_{t-1}(\mbf{x}_t, \mbf{z}) - V_\alpha(u_{t-1}(\mbf{x}_t,\mbf{Z})) \ge 0$ and $d_l(\mbf{z}) \triangleq V_\alpha(l_t(\mbf{x}_t,\mbf{Z})) - l_{t-1}(\mbf{x}_t, \mbf{z}) \ge 0$, we choose a $\mbf{z}$ that minimizes
\begin{equation}
\mcl{L}_{\text{LV}}(\mbf{z}) \triangleq \text{ReLU}(-d_u(\mbf{z})) + \text{ReLU}(-d_l(\mbf{z}))
\label{eq:llv}
\end{equation}
where $\text{ReLU}(\omega) = \max(\omega,0)$ is the rectified linear unit function ($\omega \in \mbb{R}$). To include the heuristic in Section~\ref{subsec:improvedselectz} which selects the LV with the highest probability density, we find $\mbf{z}$ that minimizes
\begin{equation*}
\mcl{L}_{\text{LV-P}}(\mbf{z}) \triangleq \mcl{L}_{\text{LV}}(\mbf{z}) - \mbb{I}_{d_u(\mbf{z}) \ge 0} \mbb{I}_{ d_l(\mbf{z}) \ge 0}\ p(\mbf{z})
\end{equation*}
where $\mcl{L}_{\text{LV}}(\mbf{z})$ is defined in \eqref{eq:llv}; $p(\mbf{z})$ is the probability density; $\mbb{I}_{d_u(\mbf{z}) \ge 0}$ and $\mbb{I}_{ d_l(\mbf{z}) \ge 0}$ are indicator functions.

\section{Experiments}
\label{sec:experiments}

In this section, we empirically evaluate the performance of {V-UCB}. The work of \citet{borisk20} has motivated the use of the approximated variant of their algorithm {$\rho \text{KG}^{apx}$} over its original version {$\rho \text{KG}$} by showing that {$\rho \text{KG}^{apx}$} achieves comparable empirical performances to {$\rho \text{KG}$} while incurring much less computational cost. Furthermore, {$\rho \text{KG}^{apx}$} has been shown to significantly outperform 
other competing algorithms \cite{borisk20}.
Therefore, we choose {$\rho \text{KG}^{apx}$} as the main baseline to empirically compare with {V-UCB}. The experiments using {$\rho \text{KG}^{apx}$} is performed by 
adding
new objective functions to the existing implementation of \citet{borisk20} at \texttt{https://github.com/saitcakmak/BoRisk}.

Regarding {V-UCB}, when $\mcl{D}_{\mbf{z}}$ is finite and the distribution of $\mbf{Z}$ is not uniform, we perform {V-UCB} by selecting $\mbf{z}_t$ as an LV at random, labeled as \emph{{V-UCB} Unif}, and by selecting $\mbf{z}_t$ as the LV with the maximum probability mass, labeled as \emph{{V-UCB} Prob}.

The performance metric 
is defined as $V_\alpha(f(\mbf{x}_*,\mbf{Z})) - V_\alpha(f(\tilde{\mbf{x}},\mbf{Z}))$ where $\tilde{\mbf{x}}$ is the recommended input. 
The evaluation of {\textsc{VaR}} is described in Section~\ref{sec:continuousz}.
The recommended input is $\text{arg}\!\max_{\mbf{x} \in \mcl{D}_T} V_{\alpha}(\mu_{t-1}(\mbf{x}, \mbf{Z}))$ for {V-UCB}, and $\text{arg}\!\min_{\mbf{x} \in \mcl{D}_{\mbf{x}}} \mbb{E}_{t-1}[V_{\alpha}(f(\mbf{x}, \mbf{Z}))]$ for {$\rho \text{KG}^{apx}$} \cite{borisk20}, where $\mbb{E}_{t-1}$ is the conditional expectation over the unknown $f$ given the observations $\mcl{D}_{t-1}$ (approximated by a finite set of functions sampled from the GP posterior belief).\footnote{While the work of \citet{borisk20} considers a minimization problem of {\textsc{VaR}}, our work considers a maximization problem of {\textsc{VaR}}. Therefore, the objective functions for {$\rho \text{KG}^{apx}$} are the negation of those for {V-UCB}. For {V-UCB} at risk level $\alpha$, 
the risk level for {$\rho \text{KG}^{apx}$} is $1-\alpha$.}
We repeat each experiment $10$ times with different random initial observations $\mbf{y}_{\mcl{D}_0}$ and plot both the mean (as lines) and the $70\%$ confidence interval (as shaded areas) of the $\log10$ of the performance metric.
The detailed descriptions of experiments are deferred to Appendix~\ref{app:experiment}.

\subsection{Synthetic Benchmark Functions}
\label{subsec:syn}

We use $3$ synthetic benchmark functions: Branin-Hoo, Goldstein-Price, and Hartmann-3D functions to construct $4$ optimization problems: Branin-Hoo-$(1,1)$, Goldstein-Price-$(1,1)$, Hartmann-$(1,2)$, and Hartmann-$(2,1)$.
The tuples represent $(d_x, d_z)$ corresponding to the 
dimensions of $\mbf{x}$ and $\mbf{z}$. 
The noise variance $\sigma_n^2$ is set to $0.01$. The risk level $\alpha$ is $0.1$. 
There are $2$ different settings: finite $\mcl{D}_{\mbf{z}}$ ($|\mcl{D}_{\mbf{z}}| = 64$ for Hartmann-$(1,2)$ and $|\mcl{D}_{\mbf{z}}|=100$ for the others) and continuous $\mcl{D}_{\mbf{z}}$. In the latter setting, $r=0.1$ and the surrogate function is a neural network of $2$ hidden layers with $30$ hidden neurons each, and sigmoid activation functions.

The results are shown in Fig.~\ref{fig:synfinite} and Fig.~\ref{fig:syncont} for the settings of discrete $\mcl{D}_{\mbf{z}}$ and continuous $\mcl{D}_{\mbf{z}}$, respectively. When $\mcl{D}_{\mbf{z}}$ is discrete (Fig.~\ref{fig:synfinite}), {V-UCB} Unif is on par with {$\rho \text{KG}^{apx}$} in optimizing Branin-Hoo-$(1,1)$ and Goldstein-Price-$(1,1)$, and outperforms {$\rho \text{KG}^{apx}$} in optimizing Hartmann-$(1,2)$ and Hartmann-$(2,1)$. {V-UCB} Prob is also able to exploit the probability distribution of $\mbf{Z}$ to outperform {V-UCB} Unif. When $\mcl{D}_{\mbf{z}}$ is continuous (Fig.~\ref{fig:syncont}), {V-UCB} Prob outperforms {$\rho \text{KG}^{apx}$}.
The unsatisfactory performance of {$\rho \text{KG}^{apx}$} in some experiments may be attributed to its approximation of the inner optimization problem in the acquisition function \cite{borisk20}, and the approximation of {\textsc{VaR}} using samples of $\mbf{Z}$ and the GP posterior belief.
\newcommand\figheight{0.201}
\begin{figure}[ht]
\centering
\begin{tabular}{@{}cc@{}}
\includegraphics[height=\figheight\textwidth]{img/negative_branin_discrete}
&
\includegraphics[height=\figheight\textwidth]{img/negative_goldstein_discrete}
\\
(a) Branin-Hoo-$(1,1)$
&
(b) Goldstein-Price-$(1,1)$
\\
\\
\includegraphics[height=\figheight\textwidth]{img/negative_hartmann3d_12_discrete}
&
\includegraphics[height=\figheight\textwidth]{img/negative_hartmann3d_21_discrete}
\\
(c) Hartmann-$(1,2)$
&
(d) Hartmann-$(2,1)$
\end{tabular}
\caption{Synthetic benchmark functions with finite $\mcl{D}_{\mbf{z}}$.}
\label{fig:synfinite}
\end{figure}
\begin{figure}[ht]
\centering
\begin{tabular}{@{}cc@{}}
\includegraphics[height=\figheight\textwidth]{img/negative_branin_continuous.pdf}
&
\includegraphics[height=\figheight\textwidth]{img/negative_goldstein_continuous.pdf}
\\
(a) Branin-Hoo-$(1,1)$
&
(b) Goldstein-Price-$(1,1)$
\\
\\
\includegraphics[height=\figheight\textwidth]{img/negative_hartmann3d_12_continuous.pdf}
&
\includegraphics[height=\figheight\textwidth]{img/negative_hartmann3d_21_continuous.pdf}
\\
(c) Hartmann-$(1,2)$
&
(d) Hartmann-$(2,1)$
\end{tabular}
\caption{Synthetic benchmark functions with continuous $\mcl{D}_{\mbf{z}}$.}
\label{fig:syncont}
\end{figure}
\subsection{Simulated Optimization Problems}
\label{subsec:simulated}


The first problem is portfolio optimization adopted by \cite{borisk20}. 
There are $d_x = 3$ optimization variables (risk and trade aversion parameters, and the holding cost multiplier) and $d_z = 2$ environmental random variables (bid-ask spread and borrow cost). The variable $\mbf{Z}$ follows a discrete uniform distribution with $|\mcl{D}_{\mbf{z}}| = 100$. Hence, there is no difference between {V-UCB} Unif and {V-UCB} Prob. Thus, we only report the results of the latter.
The objective function is the posterior mean of a trained GP on the dataset in \citet{borisk20} of size $3000$ generated from CVXPortfolio. The noise variance $\sigma_n^2$ is set to $0.01$. The risk level $\alpha$ is set to $0.2$.

The second problem is a simulated robot pushing task for which we use the implementation from the work of \citet{wang17mes}.
The simulation is viewed as a $3$-dimensional function $\mbf{h}(r_x, r_y, t_p)$ returning the 2-D location of the pushed object, where $r_x,r_y \in [-5,5]$ are the robot location and $t_p \in [1,30]$ is the pushing duration. 
The objective is to minimize the distance to a fixed goal location $\mbf{g}=(g_x,g_y)^\top$, i.e., the objective function of the maximization problem is $f_0(r_x, r_y, t_p) = -\Vert \mbf{h}(r_x, r_y, t_p) - \mbf{g}\Vert$. We assume that there are perturbations in specifying the robot location $W_x, W_y$ whose support $\mcl{D}_{\mbf{z}}$ includes $64$ equi-distant points in $[-1,1]^2$ and whose probability mass is proportional to $\exp(-(W_x^2 + W_y^2) / 0.4^2)$. It leads to a random objective function $f(r_x, r_y, t_p, W_x, W_y) \triangleq f_0(r_x + W_x, r_y + W_y, t_p)$. 
We aim to maximize the {\textsc{VaR}} of $f$ which is more difficult than maximizing that of $f_0$. 
Moreover, a random noise following $\mcl{N}(0,0.01)$ is added to the evaluation of $f$. 
The risk level $\alpha$ is set to $0.1$.

The results are shown in Fig.~\ref{fig:portrobot}. 
We observe that {V-UCB} outperforms {$\rho \text{KG}^{apx}$} in both problems. Furthermore, in comparison to our synthetic experiments, the difference between {V-UCB} Unif and {V-UCB} Prob is not significant in the robot pushing experiment. This is because the chance that a uniform sample of LV has a large probability mass is higher in the robot pushing experiment due to a larger region of $\mcl{D}_{\mbf{z}}$ having high probabilities.

\begin{figure}[ht]
\centering
\begin{tabular}{@{}cc@{}}
\includegraphics[height=\figheight\textwidth]{img/negative_portfolio_optimization_discrete}
&
\includegraphics[height=\figheight\textwidth]{img/robot_pushing_optimization_discrete}
\\
(a) Portfolio optimization $(3,2)$
&
(b) Robot pushing $(3,2)$
\end{tabular}
\caption{Simulated experiments.}
\label{fig:portrobot}
\end{figure}




\section{Conclusion}
\label{sec:conclusion}

To tackle the BO of {\textsc{VaR}} problem, we construct a no-regret algorithm, namely \emph{{\textsc{VaR}} upper confidence bound} (V-UCB), through the design of a confidence bound of {\textsc{VaR}} and a set of \emph{lacing values} (LV) that is guaranteed to exist. Besides, we introduce a heuristic to select an LV that improves the emprical performance of {V-UCB} over random selection of LV. We also draw 
an elegant
connection between BO of {\textsc{VaR}} and adversarially robust BO in terms of both problem formulation and solutions. Lastly, we provide practical techniques for implementing {\textsc{VaR}} with continuous $\mbf{Z}$. 
While {V-UCB} is more computationally efficient than the the state-of-the-art {$\rho \text{KG}^{apx}$} algorithm for BO of {\textsc{VaR}}, it also demonstrates competitive empirical performances in our experiments.




