\section{Introduction} \label{sec: intro}
Combinatorial optimization problems have been studied with respect to various application contexts including network routing problems, facility location, minimum spanning trees and etc.; see \cite{Korte2011, Conforti2014} and the references therein. In a standard deterministic formulation a linear cost function is minimized over a binary set of feasible decisions. For notational convenience, we introduce an instance of a combinatorial optimization problem as a triplet $(A, \mathbf{c}, X)$. Specifically, $A$ is a ground set (following a game theoretic notation the elements $a \in A$ are referred to as \textit{actions}), $\mathbf{c} = \{c_a, a \in A\}$ is a vector of associated costs and $X \subseteq \{0, 1\}^{|A|}$ is a set of feasible decisions. 

In this study we assume that the cost coefficients $\mathbf{c} \in \mathbb{R}^{|A|}$ form a random vector governed by some unknown probability distribution $\mathbb{Q}^*$, which is only observable through a finite data set. Ideally, if the nominal distribution $\mathbb{Q}^*$ is precisely known, then we endeavor to solve the following stochastic programming problem: 

\begin{equation} \label{stochastic programming problem}
\min_{\mathbf{x} \in X}  f(\mathbf{x}, \mathbb{Q}^*) := \mathbb{E}_{\mathbb{Q}^*} \{\gamma(\mathbf{c}, \mathbf{x})\},
\end{equation}        
where we set $\gamma(\mathbf{c}, \mathbf{x}) = \mathbf{c}^\top \mathbf{x}$ and hence, $f(\mathbf{x}, \mathbb{Q}^*)$ is a \textit{nominal expected loss} incurred by the decision-maker. 

At the same time, it is argued in \cite{Esfahani2018, VanParys2020} that usually the decision-maker cannot observe the nominal distribution $\mathbb{Q}^*$, but has access only to a finite number of independent samples from this distribution. From a computational perspective, if $\mathbb{Q}^*$ is not available to the decision-maker, then we cannot even calculate the objective function value in (\ref{stochastic programming problem}) for a fixed decision $\mathbf{x} \in X$. For this reason, the related solution approaches attempt to construct a fairly good approximation of (\ref{stochastic programming problem}) based on the available data.

\subsection{Related literature}
In this section we briefly discuss how a finite set of random observations can be transformed to a feasible decision in (\ref{stochastic programming problem}). In particular, this problem can be addressed withing the framework of \textit{data-driven stochastic programming} and \textit{distributionally robust optimization}. 

In the former approach the expected value in (\ref{stochastic programming problem}) is approximated, e.g., by leveraging Monte Carlo simulation and the resulting approximation is optimized over the set $X$ of feasible decisions~\cite{Kleywegt2002}. Sample average approximation (SAA) methods provide deterministic formulations that enjoy strong \textit{asymptotic performance guarantees} due to the central limit theorem. In other words, the solutions obtained from SAA converge (in some sense) to the nominal solution of (\ref{stochastic programming problem}) as the sample size grows to infinity. Furthermore, a robust modification of SAA, which satisfies some performance guarantees both in an asymptotic and a finite sample sense, can be applied; see, e.g., \cite{Bertsimas2018saa}.   

Another approach to approximate (\ref{stochastic programming problem}) is to construct an ambiguity set (or a family) of probability distributions $\widehat{\mathcal{Q}}$ that are deemed possible under the existing information. In particular, we must guarantee that the nominal distribution $\mathbb{Q}^*$ belongs to $\widehat{\mathcal{Q}}$ with a sufficiently high probability. Then the decision-maker may hedge against uncertainty by solving a distributionally robust optimization (DRO) problem, where the expected value in (\ref{stochastic programming problem}) is optimized under the worst-case possible distribution within $\widehat{\mathcal{Q}}$; see, e.g., \cite{Goh2010, Delage2010, Wiesemann2014}. 

In particular, the problem of constructing an ambiguity set from data observations is addressed by relatively many authors. For example, Delage and Ye \cite{Delage2010} design confidence sets for the support, mean and covariance matrix of uncertain parameters based on a finite set of independent samples. Alternatively, several studies explore a distance metric in the space of probability distributions ``centered'' at the empirical distribution of the training samples; we refer to the studies in \cite{Esfahani2018, VanParys2020} and the references therein.      

Henceforth, a procedure that transforms the data to an estimate of the expected loss in (\ref{stochastic programming problem}) for a fixed decision is referred to as a \textit{prediction rule}. On the other hand, a procedure that computes an optimal decision by minimizing the predicted loss is referred to as a \textit{prescription rule}; see, e.g., the related definitions in \cite{VanParys2020}. The major drawback of the aforementioned solution techniques is that they decouple the prediction and subsequent prescription phases. In a recent study Van Parys~et~al.~\cite{VanParys2020} propose a meta-optimization problem, which aims at finding the \textit{least conservative} prediction and prescription rules that satisfy specified asymptotic guarantees.    

More specifically, the authors in \cite{VanParys2020} seek prediction and prescription rules, whose \textit{out-of-sample disappointment} decays exponentially at some constant rate, $r > 0$, with the increase of the sample size, irrespective of a data-generating distribution. The out-of-sample disappointment measures the probability that the nominal expected loss, i.e., the expected loss under $\mathbb{Q}^*$, exceeds a predicted expected loss for some fixed decision. Van Parys~et~al.~\cite{VanParys2020} demonstrate that the least conservative prediction and prescription rules can be obtained by solving a DRO problem over all distributions within a ball with respect to the relative entropy distance. In particular, the radius of the ball coincides with the decay rate $r$ and the center is at the empirical distribution of the data.  

\subsection{Our approach and contribution}
The solution techniques proposed in \cite{VanParys2020} are applicable whenever the decision-maker has access to complete observations of the cost vector. That is, the data set in \cite{VanParys2020} is comprised of $T \in \mathbb{Z}_{>0}$ vectors (of the same size), which are drawn independently from the nominal distribution~$\mathbb{Q}^*$. Furthermore, in the recent e-print Sutter~et~al.~\cite{Sutter2020} apply a similar methodology and seek the least conservative prediction and prescription rules in the case when the data set is non-i.i.d. and the original stochastic optimization problem is not necessarily convex in $\mathbb{Q}^*$. 

In this paper we propose another extension of the problem setting in \cite{VanParys2020} with respect to the construction of the data set. Specifically, we assume that the number of random observations for particular components of the cost vector $\mathbf{c}$ may vary, i.e., the associated data set is \textit{incomplete}. This assumption seems to be relevant in a combinatorial optimization setting, where the historical data can be collected individually for each action in the ground set $A$. 

For instance, one may consider the shortest path problem in a network, where $A$ is a set of directed arcs and $\mathbf{c}$ is a vector of associated arc costs/travel times. The components $c_a$, $a \in A$, can be estimated, e.g., using Bluetooth sensors \cite{Haghani2010, Asudegi2013} or GPS-equipped floating vehicles \cite{Vo2011}. In this situation travel time estimates are subject to measurement errors and can be viewed as random observations, but not the ground truth arc travel times. Furthermore, the sensors can be placed at predefined nodes in the network and thereby the number of random observations for particular arcs is different by construction.         

Another application of incomplete data sets is motivated in the context of online learning; see, e.g., the study in \cite{Audibert2014} and the references therein. More specifically, the decision-maker may collect the data by taking a decision $\mathbf{x} \in X$ and observing the costs associated with this decision, i.e., the costs $c_a$ such that $x_a = 1$. In view of the exploration-exploitation trade-off~\cite{Sutton2018}, some of the actions in the resulting data set are observed more often than the others.  

Overall, the goal of this paper is to construct prediction and prescription rules that
\begin{itemize}
	\item[(\textit{i})] are the least conservative in some predefined sense;
	\item[(\textit{ii})] obey both asymptotic and finite sample performance guarantees in terms of the out-of-sample disappointment; \item[(\textit{iii})] can be computed rather effectively using off-the-shelf mixed-integer programming (MIP) solvers. 
\end{itemize}

Following the study of Van Parys~et~al.~\cite{VanParys2020}, we formulate a problem of finding the least conservative prediction (prescription) rules under the assumption that the data set is incomplete. In contrast to \cite{VanParys2020}, it turns out that the prediction and prescription problems (which are in fact vector optimization problems over a functional space) may not admit a unique \textit{strongly optimal} solution. However, we prove that a \textit{weakly optimal} solution can be obtained by solving a DRO problem over all probability distributions with given component-wise relative entropy distances from the empirical marginal distributions. Moreover, the obtained prediction and prescription rules are shown to obey both asymptotic and finite sample performance guarantees; recall our point (\textit{ii}) above. 

In view of the construction of our data set, we construct the family of distributions as a composition of univariate sets for each component of the cost vector. This construction enables to decouple the solution procedure for the DRO problem into two consecutive stages. At the first stage we retrieve the worst-case expected costs by solving univariate convex optimization problems; at the second stage we solve a deterministic version of the combinatorial optimization problem. This observation satisfies the requirement in (\textit{iii}). 

At the same time, even if the data set is complete, then a straightforward implementation of the DRO approach in \cite{VanParys2020} to a combinatorial optimization problem setting is computationally difficult. Namely, it results in a non-linear MIP problem. Nevertheless, some modification of the approach in \cite{VanParys2020} is examined within our numerical experiments.



One natural limitation of our approach is that we exploit the decision-maker's \textit{expected loss} as the objective criterion in (\ref{stochastic programming problem}), i.e., $\gamma(\mathbf{c}, \mathbf{x}) = \mathbf{c}^\top \mathbf{x}$. First, we note that routing decisions with the least expected cost are used in intelligent transportation and in-vehicle route guidance systems \cite{Fu1998}; we also refer to \cite{Ketkov2021} for some motivation behind the expected loss criterion in a data-driven framework. 

Importantly, the expected loss is completely defined by the mean of the cost vector $\mathbf{c}$. This observation simplifies the form of stochastic optimization problem (\ref{stochastic programming problem}), i.e., the objective function $f(\mathbf{x}, \mathbb{Q}^*)$ can be partitioned into $|A|$ terms with respect to individual actions $a \in A$. Some other possible objective criteria are briefly discussed in Section \ref{sec: conclusion}.   

Finally, we conduct computational experiments where the proposed DRO approach is applied to a class of data-driven shortest path problems. It turns out that under some parameter settings our approach outperforms several other benchmark approaches. Specifically, we consider prediction and prescription rules induced by standard measure concentration inequalities for the expected loss and by some modification of the DRO formulation in \cite{VanParys2020}.  


The remaining structure of the paper can be summarized as follows:
\begin{itemize}
	\item In Section \ref{sec: problem} we introduce prediction and prescription problems (that are motivated and similar to the related meta-optimization problems in \cite{VanParys2020}) under the assumption that the data set is incomplete. We also define problem-specific asymptotic and finite sample performance guarantees.
	\item Section \ref{sec: solution techniques} provides some analysis of the prediction and prescription problems. By leveraging the inclusion–exclusion principle and large deviations theory we establish that the outlined problems admit weakly optimal solutions. One of these solutions can be obtained by solving a DRO problem over all probability distributions with given component-wise relative entropy distances from the empirical marginal distributions.   
	\item In Section \ref{sec: comp study} we provide numerical experiments with applications to a data-driven shortest path problem. The proposed DRO approach is compared with several benchmark approaches in terms of their average out-of-sample performance.    
	\item Finally, Section~\ref{sec: conclusion} concludes the paper and outlines possible directions for future research.
\end{itemize}

\noindent \textbf{Notation.} All vectors and matrices are labelled by bold letters. A vector of all ones is referred to as $\mathbf{1}$. The natural logarithm for some $q \in \mathbb{R}_{>0}$ is denoted as $\ln(q)$. Furthermore, we adapt the conventions $0 \ln(\frac{0}{q}) = 0$ and $q' \ln(\frac{q'}{0}) = \infty$ for any $q, q' \in \mathbb{R}_{>0}$.  

We use $\mathbbm{1}\{Z\}$ as an indicator of event~$Z$. The uniform distribution on an interval $[l,u]$ is referred to as $\mathcal{U}(l,u)$. Finally, we denote by
$\mathcal{Q}_0(\mathbb{R}^{k})$ the space of all probability distributions on $\mathbb{R}^{k}$ for some $k \in \mathbb{Z}_{>0}$.

\section{Problem formulation} \label{sec: problem} 
\subsection{Modeling assumptions and terminology} \label{subsec: terminology}
As briefly outlined in Section \ref{sec: intro}, we focus on combinatorial optimization problems given by a ground set of possible actions $A$, a nonnegative cost vector $\mathbf{c} \in \mathbb{R}^{|A|}$ and a set of feasible decisions $X \subseteq \{0, 1\}^{|A|}$.  We assume that the cost vector $\mathbf{c}$ is governed by some unknown probability distribution $\mathbb{Q}^* \in \mathcal{Q}_0(\mathbb{R}^{|A|})$, which is only observable through a finite data set $\widehat{\mathbf{C}}$. 

The data set $\widehat{\mathbf{C}}$ is constructed as follows. For each $a \in A$ we form a vector of $T_a \in \mathbb{Z}_{>0}$ independent identically distributed (i.i.d.) observations of the random cost $c_a$. Formally, for each $a \in A$ we have:
\begin{equation} \nonumber
\widehat{\mathbf{c}}_a := (\widehat{c}_{a,1}, \ldots, \widehat{c}_{a,T_a})^\top
\end{equation} 
and $\widehat{\mathbf{C}} := \{\widehat{\mathbf{c}}_a, a \in A\}$.

Throughout the paper we make the following modeling assumptions:\\
\textbf{A1.} Each component of the cost vector $\mathbf{c}$ is strictly positive and has a finite discrete support, i.e., $c_a \in \{z_{a,1}, \ldots, z_{a,d_a}\}$ for $d_a \in \mathbb{Z}_{>0}$ and $a \in A$. Furthermore, $c_a = z_{a,i}$ for each $i \in \{1, \ldots, d_a\}$ holds with a nonzero probability. \\ 
\textbf{A2.} We assume that $\gamma(\mathbf{c}, \mathbf{x}) = \mathbf{c}^\top \mathbf{x}$ in the definition of stochastic optimization problem (\ref{stochastic programming problem}). 

The intuition behind Assumption \textbf{A2} is discussed in the previous section.
Then Assumption \textbf{A1} coincides with the related assumption in \cite{VanParys2020} and simplifies our further theoretical observations. Without loss of generality, elements of the support with zero probabilities can be dropped. Furthermore, we leave the case of a continuous support as a future research direction; see Section \ref{sec: conclusion} for a more detailed discussion.  

We observe that under Assumptions \textbf{A1} and \textbf{A2} the objective function in (\ref{stochastic programming problem}) can be expressed as:
\begin{equation} \label{eq: objective function linearity}
f(\mathbf{x}, \mathbb{Q}^*) := \mathbb{E}_{\mathbb{Q}^*} \{\mathbf{c}^\top \mathbf{x}\} = \sum_{a \in A} \mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} x_a = \sum_{a \in A} \sum_{i = 1}^{d_a} (z_{a,i} q^*_{a,i})x_a  
\end{equation} 
Taking into account the form of (\ref{eq: objective function linearity}), the ambiguity set of candidate data-generating distributions can be described as follows:
\begin{equation} \label{eq: set of all possible distributions}
\begin{gathered}
\mathcal{Q} :=  \Big\{ \mathbb{Q} \in \mathcal{Q}_0(\mathbb{R}^{|A|}): \; 0 \leq q_{a,i} \leq 1 \quad \forall i \in \{1, \ldots, d_a\}, \; \forall a \in A; \; \sum_{i = 1}^{d_a} q_{a,i} = 1 \quad \forall a \in A\Big\}
\end{gathered}
\end{equation}
Here, for a fixed $a \in A$ we denote by $q_{a,i}$ the probability that $c_a$ equals to $z_{a,i}$, $i \in \{1, \ldots, d_a\}$.

 Next, we recast some definitions described in \cite{VanParys2020}  under the assumption that the data set is given by $\widehat{\mathbf{C}}$. Without loss of generality, the ground actions $a \in A$ are indexed by the integers $1, \ldots, m$, where $m = |A|$. Also, for any possible data-generating distribution $\mathbb{Q} \in \mathcal{Q}$, we denote by $\mathbb{Q}_a$, $a \in A$, the marginal distributions induced by $\mathbb{Q}$ and by $\mathbf{Q}$ a vector of marginal distributions, i.e, $\mathbf{Q} = (\mathbb{Q}_1, \ldots, \mathbb{Q}_m)^\top$. 




First, we note that the value of $f(\mathbf{x}, \mathbb{Q}^*)$ for a given $\mathbf{x} \in X$ is completely defined by the marginal probability distributions $\mathbb{Q}^*_a$, $a \in A$; recall~(\ref{eq: objective function linearity}). In particular, 
\begin{equation} \label{marginal ambiguity set}
\mathbb{Q}^*_a \in \mathcal{Q}_a := \Big\{ \mathbb{Q} \in \mathcal{Q}_0(\mathbb{R}): \; 0 \leq q_{a,i} \leq 1 \quad \forall i \in \{1, \ldots, d_a\}, \; \sum_{i = 1}^{d_a} q_{a,i} = 1 \Big\}
\end{equation}
By leveraging the available data we construct empirical estimators of $\mathbb{Q}^*_a$ as follows:
\begin{equation}
\widehat{\mathbb{Q}}_{a,T_a}\{c_a = i\} = \widehat{q}_{a,i} = \frac{1}{T_a}\sum_{j = 1}^{T_a}\mathbbm{1}\{\widehat{c}_{a,j} = z_{a,i}\} \quad \forall i \in \{1, \ldots, d_a\}
\end{equation} 
The distributions $\widehat{\mathbb{Q}}_{a,T_a}$, $a \in A$, are referred to as \textit{empirical marginal distributions}. 

As we outline in Section \ref{sec: intro}, stochastic optimization problem (\ref{stochastic programming problem}) cannot be solved exactly and thus, we attempt to approximate the nominal expected loss $f(\mathbf{x}, \mathbb{Q}^*)$ and the optimal decision $\mathbf{x}^*(\mathbb{Q}^*) \in \argmin_{\mathbf{x} \in X} f(\mathbf{x}, \mathbb{Q}^*)$ by some functions of the empirical marginal distributions. More precisely, we define \textit{data-driven prediction} and \textit{prescription rules}.  

\begin{definition} [\textbf{Prediction and prescription rules}] \upshape \label{def: data-driven predictor and prescriptor} 
An arbitrary function $\hat{f}: X \times \mathcal{Q}_1 \times \ldots \times~\mathcal{Q}_m \rightarrow \mathbb{R}$ is called a \textit{prediction rule}. A function $\hat{x}: \mathcal{Q}_1 \times \ldots \times \mathcal{Q}_m \rightarrow \mathbb{R}$ is called an associated \textit{prescription rule}, if  $\hat{x}(\mathbf{Q}) \in \argmin_{\mathbf{x} \in X} \hat{f}(\mathbf{x}, \mathbf{Q})$ for any fixed $\mathbb{Q} \in \mathcal{Q}$. \vspace{-8.75mm}\flushright$\square$  
\end{definition} 

Let $\widehat{\mathbf{Q}}_{\mathbf{T}}$ be a vector of empirical marginal distributions, i.e., 
$$\widehat{\mathbf{Q}}_{\mathbf{T}} = (\widehat{\mathbb{Q}}_{1,T_1}, \ldots, \widehat{\mathbb{Q}}_{m,T_m})^\top$$ 
Then by leveraging Definition \ref{def: data-driven predictor and prescriptor} one may approximate the nominal expected loss $f(\mathbf{x}, \mathbb{Q}^*)$ and the optimal decision $\mathbf{x}^*(\mathbb{Q}^*)$ with the estimates $\hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}})$ and $\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}})$, respectively. 
Furthermore, following \cite{VanParys2020}, we assess the quality of approximation by using an \textit{out-of-sample disappointment}. 

\begin{definition} [\textbf{Out-of-sample disappointment}] \upshape \label{def: out-of-sample prediction disappointment} 
For any prediction rule $\hat{f}$ and feasible decision $\mathbf{x} \in X$ we introduce an \textit{out-of-sample prediction disappointment} as:
\begin{equation} \label{eq: out-of-sample prediction}
\Pr\Big\{f(\mathbf{x}, \mathbb{Q}^*) > \hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \},
\end{equation}
where the probability is taken with respect to the nominal distribution $\mathbb{Q}^*$. Similarly, for a given pair $(\hat{f}, \hat{x})$  
\begin{equation} \label{eq: out-of-sample prescription}
\Pr\Big\{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}^*) > \hat{f}(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \}
\end{equation}
is reffered to as an \textit{out-of-sample prescription disappointment}.  \vspace{-10mm}\flushright$\square$  
\end{definition} 

The out-of-sample disappointment for a given decision $\mathbf{x} \in X$ quantifies the probability that the nominal expected loss $f(\mathbf{x}, \mathbb{Q}^*)$ is underestimated. It is pointed out in \cite{VanParys2020} that underestimated losses are more harmful for the decision-maker than overestimated losses (as it strives to minimize its loss over a set of feasible decisions.)  

In this paper we focus on a class of prediction rules that cater an \textit{exponential decay rate} for the out-of-sample prediction disappointment irrespective of a decision $\mathbf{x} \in X$ and a data-generating distribution $\mathbb{Q}^* \in \mathcal{Q}$. In other words, we exploit the following \textit{asymptotic} and \textit{finite sample performance guarantees}: 
\begin{itemize} 
\item \textbf{Asymptotic guarantee.} Let $T_{min} = \min_{a \in A} T_a$. Then the out-of-sample prediction disappointment (\ref{eq: out-of-sample prediction}) decays exponentially at some rate $r > 0$ when $T_{min}$ goes to infinity, i.e., 
\begin{align} \label{eq: asymptotic guarantee} \tag{$\mathcal{AG}$}
& \lim_{T_{min} \rightarrow +\infty} \sup  \frac{1}{T_{min}} \ln  \Big(\Pr\Big\{f(\mathbf{x}, \mathbb{Q}) > \hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \} \Big) \leq -r \quad \forall \mathbf{x} \in X, \; \forall \mathbb{Q} \in \mathcal{Q} 
\end{align}  
\item \textbf{Finite sample guarantee.}  
Let $g(T)$ be a given function that decays exponentially at some rate $r > 0$ as $T$ goes to infinity. Then the out-of-sample prediction disappointment is bouve by $g(T_{min})$, i.e.,  
\begin{align} \label{eq: finite sample guarantee} \tag{$\mathcal{FG}$}
& \Pr\Big\{f(\mathbf{x}, \mathbb{Q}) > \hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \} \leq g(T_{min}) \quad \forall \mathbf{x} \in X, \; \forall \mathbb{Q} \in \mathcal{Q}, \; \forall T_{min} \in \mathbb{Z}_{>0}
\end{align}  
\end{itemize}

In the following example we demonstrate that the finite sample guarantee (\ref{eq: finite sample guarantee}) holds for a rather simple class of prediction rules. That is, we approximate the expected value in (\ref{stochastic programming problem}) by the sum of empirical mean of the objective function $\mathbf{c}^\top \mathbf{x}$ and any positive constant. 

\begin{example} \upshape \label{example 1}
We assume that the prediction rule $\hat{f}$ for a fixed $\mathbf{x} \in X$ is defined as:
\begin{equation} \label{eq: sample-average predictor + epsilon}
\hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) = \sum_{a \in A} \Big( \frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} \Big) x_a + \varepsilon = \sum_{a \in A} \Big(\sum_{i = 1}^{d_a} z_{a,i} \widehat{q}_{a,i} \Big) x_a + \varepsilon,
\end{equation}	
where $\varepsilon > 0$ is a given positive constant. Also, let $A_\mathbf{x} = \{a \in A: x_a = 1\}$, $z_{min} = \min\{z_{a, i}: i \in \{1, \ldots, d_a\}, a \in A\}$ and $z_{max} = \max\{z_{a, i}: i \in \{1, \ldots, d_a\}, a \in A\}$. Using (\ref{eq: objective function linearity}) we observe that the out-of-sample prediction disappointment (\ref{eq: out-of-sample prediction}) can be bounded from above as:
\allowdisplaybreaks 
\begin{align} \label{ex 1: proof} 
\centering
\Pr\Big\{f(\mathbf{x}, \mathbb{Q}^*) > \hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \} = & \Pr\Big\{\sum_{a \in A} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} \Big) x_a > \varepsilon \Big\}  \leq \nonumber \\ \leq
\Pr\Big\{\bigcup_{a \in A_\mathbf{x}} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} > & \frac{\varepsilon}{|A_\mathbf{x}|} \Big) \Big\} \leq \sum_{a \in A_\mathbf{x}} \Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} > \frac{\varepsilon}{|A_\mathbf{x}|}\Big\} \leq \\ 
\leq \sum_{a \in A_\mathbf{x}} \exp\Big(\frac{-2 T_a \; \varepsilon^2}{|A_\mathbf{x}|^2(z_{max} - z_{min})^2} \Big) & \leq |A| \exp\Big(\frac{-2 T_{min} \;\varepsilon^2 }{|A|^2(z_{max} - z_{min})^2} \Big) = g(T_{min}) \nonumber
\end{align}
The first inequality can be checked by contradiction, while the second inequality is implied by the standard union bound. The third inequality in (\ref{ex 1: proof}) exploits Hoeffding's inequality for the sum of bounded i.i.d. random variables; see \cite{Boucheron2003}. The upper bounding function $g(T_{min})$ has an exponential decay rate given by:
\begin{equation} \nonumber
r(\varepsilon) = -\lim_{T_{min} \to \infty} \sup \Big(\frac{1}{T_{min}} \ln g(T_{min}) \Big) = \frac{2 \varepsilon^2}{|A|^2(z_{max} - z_{min})^2} ,
\end{equation}  
Furthermore, $r(\varepsilon)$ can be made arbitrarily large by an appropriate choice of parameter $\varepsilon$. We conclude that the finite sample guarantee (\ref{eq: finite sample guarantee}) for the constructed prediction rule holds.  
\vspace{-8.94mm}\flushright$\square$
\end{example}
In the next section we introduce prediction and prescription problems, which are motivated by the related formulations in \cite{VanParys2020}. 

\subsection{Prediction and prescription problems} \label{subsec: prediction and prescription problems}
First, we impose a partial order on the space of prediction rules. That is, a prediction rule $\hat{f}_1$ weakly outperforms a prediction rule $\hat{f}_2$, i.e., $\hat{f}_1 \preceq \hat{f}_2$, if and only if
\begin{equation} \label{eq: partial order prediction}
\hat{f}_1(\mathbf{x}, \mathbb{Q}_1, \ldots, \mathbb{Q}_m) \leq \hat{f}_2(\mathbf{x}, \mathbb{Q}_1, \ldots, \mathbb{Q}_m) \quad \forall \mathbf{x} \in X, \; \forall \mathbb{Q} \in \mathcal{Q}
\end{equation}  


\begin{table}
	\footnotesize
	\centering
	\onehalfspacing
	\begin{tabular}{c | c | c}
		\hline
		Properties & Our formulation & Formulation of Van~Parys~et~al.~\cite{VanParys2020}\\
		\hline
		\multirow{2}{*}{Sample size} & \multirow{2}{*}{$T_a \in \mathbb{Z}_+$ for each $a \in A$} & \multirow{2}{*}{the same size $T \in \mathbb{Z}_+$ for each $a \in A$}\\
		& & \\\hline
		\multirow{2}{*}{Set of feasible decisions} & \multirow{2}{*}{$X \subseteq \{0, 1\}^{|A|}$} & \multirow{2}{*}{$X \subseteq \mathbb{R}^{n}$, $X$ is compact} \\
		& & \\\hline
		\multirow{2}{*}{Objective function in (\ref{stochastic programming problem})} & \multirow{2}{*}{$\gamma(\mathbf{c}, \mathbf{x}) = \mathbf{c}^\top \mathbf{x}$} & $\gamma(\mathbf{c}, \mathbf{x})$ is any continuous function \\
		& & with respect to $\mathbf{x} \in X$ \\\hline
		\multirow{2}{*}{Prediction rule} & an arbitrary function & a continuous function \\
		& of decisions and marginal distributions & of decisions and joint distributions \\
		\hline
	\end{tabular}
	\caption{\scriptsize The table summarizes the difference between our formulations (\ref{prediction problem}) and (\ref{prescription problem}) and the related formulations in \cite{VanParys2020}.}
	\label{tab: difference}
\end{table}

We denote by $\mathcal{F}$ a set of all real-valued functions on the domain $X \times \mathcal{Q}_1 \times \ldots \times \mathcal{Q}_m$. Then a \textit{prediction problem} can be expressed as:
\allowdisplaybreaks\begin{subequations} 
	\begin{align} \label{prediction problem} \tag{$\mathcal{P}$}
	& \min_{\hat{f}(\cdot) \in \mathcal{F}} \hat{f}(\cdot) \quad (\mbox{ w.r.t. }\preceq) \\
	\mbox{\upshape s.t.}
	& \lim_{T_{min} \rightarrow +\infty} \sup  \frac{1}{T_{min}} \ln\Big(\Pr\Big\{f(\mathbf{x}, \mathbb{Q}) > \hat{f}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \}\Big) \leq -r \quad \forall \mathbf{x} \in X, \; \forall \mathbb{Q} \in \mathcal{Q} \nonumber
	\end{align}
\end{subequations}
Formally, vector optimization problem (\ref{prediction problem}) aims at finding the least conservative prediction rule (with respect to the partial order defined above) that satisfies asymptotic guarantee (\ref{eq: asymptotic guarantee}). 

Analogously, we define a partial order on the space of prescription rules and introduce a prescription problem. In particular, we assume that a prediction-prescription pair $(\hat{f}_1, \hat{x}_1)$ weakly outperforms a pair $(\hat{f}_2, \hat{x}_2)$, if the following set of inequalities hold:
\begin{equation} \label{eq: partial order prescription}
\hat{f}_1(\hat{x}_1(\mathbf{Q}), \mathbf{Q}) \leq \hat{f}_2(\hat{x}_2(\mathbf{Q}), \mathbf{Q}) \quad \forall \mathbb{Q} \in \mathcal{Q}
\end{equation}  
(recall that $\mathbf{Q} = (\mathbb{Q}_1, \ldots, \mathbb{Q}_m)^\top$ is a vector of marginal distributions induced by some $\mathbb{Q} \in \mathcal{Q}$.)

Let $\mathcal{X}$ be a set of all pairs $(\hat{f}, \hat{x})$, where $\hat{f} \in \mathcal{F}$ and $\hat{x}$ is the prescription rule induced by $\hat{f}$. Then a \textit{prescription problem} can be formulated as follows:
\allowdisplaybreaks\begin{subequations} 
	\begin{align} \label{prescription problem} \tag{$\mathcal{P}'$}
	& \min_{(\hat{f}, \hat{x}) \in \mathcal{X}} (\hat{f}, \hat{x}) \quad (\mbox{ w.r.t. }\preceq) \\
	\mbox{\upshape s.t.}
	& \lim_{T_{min} \rightarrow +\infty} \sup  \frac{1}{T_{min}}\ln\Big( \Pr\Big\{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}) > \hat{f}(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \}\Big) \leq -r \quad \forall \mathbb{Q} \in \mathcal{Q} \nonumber
	\end{align}
\end{subequations}


Table \ref{tab: difference} provides some comparison of prediction and prescription problems (\ref{prediction problem}) and (\ref{prescription problem}) with the related meta-optimization problems in \cite{VanParys2020}. In fact, we consider another class of prediction rules that are some functions of empirical marginal distributions; recall Definition \ref{def: data-driven predictor and prescriptor} and the form of objective criterion in (\ref{stochastic programming problem}). Also, in our subsequent theoretical results (and in contrast to Theorem 7 in \cite{VanParys2020}) we do not restrict the prediction rules to be continuous functions. In particular, the feasible set $X$ is assumed to be discrete and thus, continuity with respect to $\mathbf{x} \in X$ is not well-defined.


In contrast to \cite{VanParys2020}, it turns out that vector optimization problems (\ref{prediction problem}) and~(\ref{prescription problem}) may not admit a unique \textit{strongly optimal solution} (we operate with the standard notion of optimality for vector optimization problems \cite{Jahn2009}.) At the same time, a \textit{weakly optimal solution} can be derived by solving a particular distributionally robust optimization problem. The intuition behind these observations is made precise in the next section. 

\section{Solution techniques} \label{sec: solution techniques}
In this section we analyze a particular class of distributionally robust prediction and prescription rules based on the relative entropy distance from empirical marginal distributions. More precisely, for any fixed action $a \in A$ a relative entropy distance or Kullback–Leibler divergence (see, e.g., \cite{VanParys2020} and the references therein) between an empirical distribution $\widehat{\mathbb{Q}}_{a, T_a}$ and some other distribution $\mathbb{Q}_a \in \mathcal{Q}_a$ is defined as follows:  
\begin{equation} \label{eq: relative entropy}
D_{KL}(\widehat{\mathbb{Q}}_{a, T_a}\; \Vert \; \mathbb{Q}_a) = \sum_{i = 1}^{d_a} \widehat{q}_{a,i} \ln\frac{\widehat{q}_{a,i}}{q_{a,i}}
\end{equation}

Then we consider the following ambiguity set of probability distributions:
\begin{equation} \label{ambiguity set}
\begin{gathered}
\widehat{\mathcal{Q}} :=  \Big\{ \mathbb{Q} \in \mathcal{Q}_0(\mathbb{R}^{|A|}): \; 0 \leq q_{a,i} \leq 1 \quad \forall i \in \{1, \ldots, d_a\}, \; \forall a \in A;  \\ \quad \qquad \;
\sum_{i = 1}^{d_a} q_{a,i} = 1 \quad \forall a \in A; \; D_{KL}(\widehat{\mathbb{Q}}_{a, T_a}\; \Vert \; \mathbb{Q}_a) \leq r_a \quad \forall a \in A \Big\} 
\end{gathered}
\end{equation}
Formally, $\widehat{\mathcal{Q}}$ accounts for univariate probability distributions within balls centered at the empirical marginal distributions. Naturally, the radius $r_a$, $a \in A$, depends on the number of random observations, $T_a$. An analytical form of $r_a$ is described via Theorem \ref{theorem 1} in Section \ref{subsec: asymptotic and finite sample guarantees}.  

Note that the distributional constraints in (\ref{ambiguity set}) are partitioned into $|A|$ non-overlapping subsets with respect to each action in the ground set. In this regard, for each $a \in A$ we define univariate ambiguity sets given by:  
\begin{equation} \nonumber
\begin{gathered}
\widehat{\mathcal{Q}}_a :=  \Big\{ \mathbb{Q}_a \in \mathcal{Q}_0(\mathbb{R}): \; 0 \leq q_{a,i} \leq 1 \quad \forall i \in \{1, \ldots, d_a\}, \; \sum_{i = 1}^{d_a} q_{a,i} = 1, \; D_{KL}(\widehat{\mathbb{Q}}_{a, T_a}\; \Vert \; \mathbb{Q}_a) \leq r_a \Big\} 
\end{gathered}
\end{equation}
Then for any fixed decision $\mathbf{x} \in X$ and empirical distributions $\widehat{\mathbf{Q}}_{\mathbf{T}}$ distributionally robust prediction and prescription rules are introduced as:
\begin{subequations} 
\begin{align}
& \hat{f}_{\mathbf{r}}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) = \max_{\mathbb{Q} \in \widehat{\mathcal{Q}}}\mathbb{E}_{\mathbb{Q}}\{\mathbf{c}^\top \mathbf{x}\} = \sum_{a \in A} \Big(\max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a}\mathbb{E}_{\mathbb{Q}_a}\{c_a\}\Big) x_a \label{distributionally robust predictor} \\
& \hat{x}_{\mathbf{r}}(\widehat{\mathbf{Q}}_{\mathbf{T}}) \in \argmin_{\mathbf{x} \in X} \hat{f}_{\mathbf{r}}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) 
\label{distributionally robust prescriptor}
\end{align}
\end{subequations}
In the sequel, we analyze theoretical properties of (\ref{distributionally robust predictor}) and (\ref{distributionally robust prescriptor}).

\subsection{Properties of the distributionally robust prediction rule}
The section is organized as follows. At first, we consider a dual reformulation of (\ref{distributionally robust predictor}). Then we establish that (\ref{distributionally robust predictor}) and (\ref{distributionally robust prescriptor}) can be treated as weakly optimal solutions of the prediction problem (\ref{prediction problem}) and the prescription problem (\ref{prescription problem}), respectively. 
\subsubsection{Dual reformulation}
We demonstrate that optimization problem (\ref{distributionally robust predictor}) admits a dual reformulation. Furthermore, the related DRO problem 
\begin{equation} \label{DRO problem} \tag{$\mathcal{DRO}$}
\min_{\mathbf{x} \in X} \hat{f}_{\mathbf{r}}(\mathbf{x},  \widehat{\mathbf{Q}}_{\mathbf{T}}) = \min_{\mathbf{x} \in X} \max_{\mathbb{Q} \in \widehat{\mathcal{Q}}} \mathbb{E}_{\mathbb{Q}^*}\{\mathbf{c}^\top \mathbf{x}\} 
\end{equation}
can be solved by leveraging $|A|$ one-dimensional convex problems and a single deterministic mixed-integer programming (MIP) problem. 
\begin{proposition} \upshape \label{proposition 1}
Distributionally robust optimization problem (\ref{DRO problem}) is equivalent to the following minimization problem:
\begin{subequations}
\begin{align}
& \min_{\boldsymbol{\beta}, \mathbf{x}} \sum_{a \in A} \Big( \beta_a - e^{-r_a}\prod_{i = 1}^{d_a} (\beta_a - z_{a,i})^{\widehat{q}_{a,i}} \Big) x_a \\
\mbox{s.t. } & \mathbf{x} \in X \\
& \beta_a \geq d_a \quad \forall a \in A 
\end{align}
\end{subequations}	
\begin{proof}
Note that a single-level reformulation of (\ref{DRO problem}) can be obtained by dualizing the worst-case expectation problems for each particular $a \in A$, that is,
\begin{subequations} \label{worst-case expectation problem}
\begin{align}
& \max_{\mathbf{q}_a}\mathbb{E}_{\mathbb{Q}_a}\{c_a\} = \sum_{i = 1}^{d_a} z_{a,i}q_{a,i} \\
\mbox{s.t. } & 0 \leq q_{a,i} \leq 1 \quad \forall i \in \{1, \ldots, d_a\} \\
& \sum_{i = 1}^{d_a} q_{a,i} = 1 \\
& \sum_{i = 1}^{d_a} \widehat{q}_{a,i} \ln\frac{\widehat{q}_{a,i}}{q_{a,i}} \leq r_a 
\end{align}
\end{subequations}	
It is rather straightforward to verify that the feasibility set in (\ref{worst-case expectation problem}) is convex and the objective function is linear. The dual problem can be formulated as a one-dimensional convex optimization problem of the form (see Proposition 2 in \cite{VanParys2020}):
\begin{subequations} \label{individual dual problems}
\begin{align}
& \min_{\beta_a} \Big( \beta_a - e^{-r_a}\prod_{i = 1}^{d_a} (\beta_a - z_{a,i})^{\widehat{q}_{a,i}} \Big) \\
\mbox{s.t. } & \beta_a \geq d_a 
\end{align}
\end{subequations}	
If we combine the minimization over $\mathbf{x} \in X$ and the minimization over the dual variables $\beta_a$, $a \in A$, then the result follows.   		
\end{proof}
\end{proposition}
Thus, (\ref{DRO problem}) can be recast as a single deterministic MIP problem with the worst-case expected costs induced by (\ref{individual dual problems}). Furthermore, despite the feasible set in (\ref{individual dual problems}) is unbounded, the objective function is bounded from below; see the proof of Proposition 2 in \cite{VanParys2020}.
\subsubsection{Asymptotic and finite sample guarantees} \label{subsec: asymptotic and finite sample guarantees}
It turns out that the distributionally robust prediction rule (\ref{distributionally robust predictor}) with an appropriate choice of the parameters $r_a$, $a \in A$, satisfies both the finite sample guarantee (\ref{eq: finite sample guarantee}) and the asymptotic guarantee (\ref{eq: asymptotic guarantee}). In particular, this assertion follows from a strong \textit{large deviation principle}; see, e.g., \cite{VanParys2020, Groeneboom1979}, which is formulated and briefly discussed below.       

\begin{proposition} [\textbf{Strong LDP, univariate case}] If for some fixed $a \in A$ the samples $\widehat{c}_{a,1}, \ldots, \widehat{c}_{a,T_a}$ are drawn independently from some $\mathbb{Q}^*_a \in \mathcal{Q}_a$, then for every Borel set $\mathcal{D}_a \subseteq \mathcal{Q}_a$ and any $T_a \in \mathbb{N}$ we have:
\begin{equation} \label{eq: strong LDP 1}
\mbox{\upshape Pr} \Big\{\widehat{\mathbb{Q}}_{a,T_a} \in \mathcal{D}_a \Big\} \leq (T_a + 1)^{d_a}\exp \Big(-T_a \inf_{\mathbb{Q}_a \in \mathcal{D}_a} D_{KL}(\mathbb{Q}_a \; \Vert \; \mathbb{Q}^*_a)\Big) 
\end{equation}
Additionally, for any $\mathbb{Q}_a \in \mathcal{Q}_a$ 
\begin{equation} \label{eq: strong LDP 2}
 \mbox{\upshape Pr} \Big\{\widehat{\mathbb{Q}}_{a,T_a} = \mathbb{Q}_a \Big\} \geq (T_a + 1)^{-d_a}\exp \Big(-T_a D_{KL}(\mathbb{Q}_a \; \Vert \; \mathbb{Q}^*_a)\Big) 
\end{equation}
\begin{proof} The former inequality follows from Theorem 2 in \cite{VanParys2020}. The latter inequality is implied by the proof of Theorem 1 in \cite{VanParys2020}; see page 28 in the Appendix section in \cite{VanParys2020}.
\end{proof}
\end{proposition}

In other words, if $\mathbb{Q}^*_a \notin \mathcal{D}$, then the empirical distribution $\widehat{\mathbb{Q}}_{a,T_a}$ belongs to $\mathcal{D}$ with a probability that decays exponentially in $T_a$. In the next result we mimic the ideas behind Example~\ref{example 1}. That is, we bound from above the out-of-sample prediction disappointment with respect to (\ref{distributionally robust predictor}) and cater a fixed exponential decay rate.   

\begin{theorem} \label{theorem 1}
Let $r > 0$ be a required exponential decay rate and 
\begin{equation} \label{eq: radius} 
r_a = \frac{1}{T_a}\Big(\ln|A| + d_a \ln(T_a + 1) + T_{min} r \Big) \quad \forall a \in A
\end{equation}
Then for any $T_a \in \mathbb{Z}_{>0}$, $a \in A$, we have:
\begin{equation} \label{eq: finite sample prediction}
\mbox{\upshape Pr}\Big\{f(\mathbf{x}, \mathbb{Q}^*) > \hat{f}_{\mathbf{r}}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \} \leq e^{-T_{min} r} \quad \forall \mathbf{x} \in X, \; \forall \mathbb{Q}^* \in \mathcal{Q} 
\end{equation}
\begin{proof}
First, we fix $\mathbf{x} \in X$ and let $A_{\mathbf{x}} = \{a \in A: x_a = 1\}$. Note that
\begin{equation} \label{eq: feasibility proof 1} \nonumber
\begin{gathered}
\Pr\Big\{f(\mathbf{x}, \mathbb{Q}^*) > \hat{f}_{\mathbf{r}}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) \Big \} = \Pr\Big\{\mathbb{E}_{\mathbb{Q}^*}\{\mathbf{c}^\top \mathbf{x}\} -  \max_{\mathbb{Q} \in \widehat{\mathcal{Q}}}\mathbb{E}_{\mathbb{Q}_a}\{\mathbf{c}^\top \mathbf{x}\} > 0\Big\} = \\ = \Pr\Big\{\sum_{a \in A} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a}\mathbb{E}_{\mathbb{Q}_a}\{c_a\}\Big)x_a  > 0 \Big\}  \leq \\ \leq \Pr\Big\{\bigcup_{a \in A_{\mathbf{x}}} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a}\mathbb{E}_{\mathbb{Q}_a}\{c_a\} > 0 \Big) \Big\}
\end{gathered}
\end{equation}
Furthermore, each inequality 
$$\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a}\mathbb{E}_{\mathbb{Q}_a}\{c_a\} > 0$$
implies that $\mathbb{Q}^*_a \notin \widehat{\mathcal{Q}}_a$. Hence, by leveraging the union bound and the strong LDP principle we observe that: 
\begin{equation} \label{eq: feasibility proof 2} \nonumber
\begin{gathered}
\Pr\Big\{\bigcup_{a \in A_{\mathbf{x}}} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a}\mathbb{E}_{\mathbb{Q}_a}\{c_a\} > 0 \Big) \Big\} \leq \\ \leq \sum_{a \in A_{\mathbf{x}}} \Pr\Big\{\mathbb{Q}^*_a \notin \widehat{\mathcal{Q}}_a\Big\} = \sum_{a \in A_{\mathbf{x}}} \Pr\Big\{D_{KL}(\widehat{\mathbb{Q}}_{a,T_a} \; \Vert \; \mathbb{Q}^*_a) > r_a\Big\} \\ \leq \sum_{a \in A_{\mathbf{x}}}(T_a + 1)^{d_a} e^{-T_a r_a} \leq e^{-T_{min} r} \frac{|A_{\mathbf{x}}|}{|A|} \leq e^{-T_{min} r} = g(T_{min})
\end{gathered}
\end{equation}
 The last two inequalities follow from (\ref{eq: radius}) and the fact that $|A_{\mathbf{x}}| \leq |A|$ for any $\mathbf{x} \in X$.  This observation concludes the proof. 
\end{proof}
\end{theorem}

We conclude that the finite sample and asymptotic guarantees hold for the distributionlly robust predictor (\ref{distributionally robust predictor}) with $r_a$, $a \in A$, given by (\ref{eq: radius}). Furthermore, for a given confidence level $\alpha \in (0, 1)$, the corresponding value of parameter $r$ is given by $-\frac{1}{T_{min}}\ln \alpha$. 

\subsubsection{Weak optimality} \label{subsubsec: weak optimality}
In this section we demonstrate that (\ref{distributionally robust predictor}) provides a weakly optimal solution of the prediction problem (\ref{prediction problem}). To this end, we introduce a standard definition of weak optimality.
\begin{definition}[\textbf{Weak optimality, prediction}] \label{def: weak optimality} \upshape
A prediction rule $\hat{f}$ is called \textit{weakly optimal} for (\ref{prediction problem}), if it is feasible in (\ref{prediction problem}) and there exist no other feasible prediction rule that is less conservative than $\hat{f}$ uniformly over all the decisions $\mathbf{x} \in X$ and the distributions $\mathbb{Q} \in \mathcal{Q}$.  
\end{definition} 

\begin{theorem} \label{theorem 2}
Let $r > 0$ be a required exponential decay rate. The distributionlly robust predictor (\ref{distributionally robust predictor}) with $\mathbf{r} = \{r_a, a \in A\}$ given by (\ref{eq: radius}) is weakly optimal for (\ref{prediction problem}).
\begin{proof}
Assume that there exists another prediction rule $\hat{f}$, which is less conservative than $\hat{f}_{\mathbf{r}}$ and satisfies the asymptotic guarantee (\ref{eq: asymptotic guarantee}). Therefore, for any decision $\mathbf{x} \in X$ and any probability distribution $\mathbb{Q} \in \mathcal{Q}$ there exists $\varepsilon \in \mathbb{R}_{>0}$ such that:
\begin{equation} \label{eq: weak optimality proof 1}
\hat{f}(\mathbf{x}, \mathbf{Q}) \leq \hat{f}_{\mathbf{r}}(\mathbf{x}, \mathbf{Q}) - \varepsilon
\end{equation} 
In the remainder of the proof we show that $\hat{f}$ is infeasible in (\ref{prediction problem}). Let $b \in A$ be an element of $A$ with the minimal number of random observations, i.e., $T_b = T_{min}$. We pick a feasible decision $\widetilde{\mathbf{x}} \in X$ such that $\widetilde{x}_b = 1$ (if there are no feasible decisions that contain $b$, then the element $b$ can be removed.) Also, set $A_{\widetilde{\mathbf{x}}} = \{a \in A: \widetilde{x}_a = 1\}$. 

For the outlined decision $\widetilde{\mathbf{x}} \in X$ and some data-generating distribution $\mathbb{Q}^* \in \mathcal{Q}$ (the form of $\mathbb{Q}^*$ is specified further) we bound the out-of-sample disappointment (\ref{eq: out-of-sample prediction}) from below as:
\begin{equation} \label{eq: weak optimality proof 2} \tag{$LB_1$}
\begin{gathered}
\Pr\Big\{f(\widetilde{\mathbf{x}}, \mathbb{Q}^*) > \hat{f}(\widetilde{\mathbf{x}}, \widehat{\mathbf{Q}}_{\mathbf{T}})\Big\} \geq \Pr\Big\{f(\widetilde{\mathbf{x}}, \mathbb{Q}^*) - \hat{f}_{\mathbf{r}}(\widetilde{\mathbf{x}}, \widehat{\mathbf{Q}}_{\mathbf{T}}) > -\varepsilon \Big\} = \\ = \Pr\Big\{\sum_{a \in A} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a} \mathbb{E}_{\mathbb{Q}_a}\{c_a\} \Big) \widetilde{x}_a > -\varepsilon \Big\} 
\end{gathered}
\end{equation}
It turns out that for a particular choice of the data-generating distribution $\mathbb{Q}^*$ the inequality 
\begin{equation} \label{eq: weak optimality proof joint}
\sum_{a \in A} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a} \mathbb{E}_{\mathbb{Q}_a}\{c_a\} \Big) \widetilde{x}_a > -\varepsilon
\end{equation}  
is implied by the individual inequality  
\begin{equation} \label{eq: weak optimality proof individual} 
\begin{gathered}
\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\varepsilon' 
\end{gathered}
\end{equation}
for some $\varepsilon' > 0$.

\textbf{Step 1.} Our first step is to bound the probability of (\ref{eq: weak optimality proof individual}) from below by some function that decays exponentially at some rate $r' < r$ as $T_{min}$ tends to infinity. First, we prove that: 
\begin{equation}  \label{eq: weak optimality proof 3} \tag{$LB_2$}
\mbox{\upshape Pr} \Big\{\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\varepsilon' \Big\} \geq \Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} > 0 \Big\},
\end{equation}
where the relative entropy ball $\widehat{\mathcal{Q}}_b$ of the radius $r_b$ is replaced by another ball $\widehat{\mathcal{Q}}'_b$ with the same center, $\widehat{\mathbb{Q}}_{b,T_b}$, but of a smaller radius $r'_b < r_b$. 

Indeed, let
\begin{equation} \nonumber
\mathbb{Q}_b^{(w)} \in \argmax_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\}
\end{equation}
be the worst-case marginal distribution within the relative entropy ball $\widehat{\mathcal{Q}}_b$ and 
$$\mathbb{Q}_b(\lambda) := \lambda \widehat{\mathbb{Q}}_{b,T_b} + (1 - \lambda) \mathbb{Q}_b^{(w)}$$
for some $\lambda \in (0, 1)$. Then by convexity of the relative entropy we have:
\begin{equation} \nonumber
D_{KL}(\widehat{\mathbb{Q}}_{b,T_b} \; \Vert \; \mathbb{Q}_b(\lambda)) \leq \lambda D_{KL}(\widehat{\mathbb{Q}}_{b,T_b} \; \Vert \; \widehat{\mathbb{Q}}_{b,T_b}) + (1 - \lambda) D_{KL}(\widehat{\mathbb{Q}}_{b,T_b} \; \Vert \; \mathbb{Q}_b^{(w)})  \leq (1 - \lambda)r_b = r'_b < r_b
\end{equation} 
Furthermore, we observe that:
\begin{equation} \label{eq: weak optimality proof 4}
\begin{gathered}
\max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} - \varepsilon' = \mathbb{E}_{\mathbb{Q}^{(w)}_b}\{c_b\} - \varepsilon' = \sum_{i = 1}^{d_b} q^{(w)}_{b,i} z_{b,i}  - \varepsilon' < \\ < \sum_{i = 1}^{d_b} \Big((1 - \lambda)q^{(w)}_{b,i} + \lambda \widehat{q}_{b,i}) z_{b,i} \Big) = \mathbb{E}_{\mathbb{Q}_b(\lambda)}\{c_b\} \leq \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\}, 
\end{gathered}
\end{equation}
where the ambiguity set $\widehat{\mathcal{Q}}'_b$ is obtained from $\widehat{\mathcal{Q}}_b$ by substituting $r'_b$ instead of $r_b$.
It is rather straightforward to verify that the strict inequality in (\ref{eq: weak optimality proof 4}) holds for a sufficiently small value of parameter $\lambda$. Hence, the inequality in (\ref{eq: weak optimality proof 3}) follows.

Next, we consider a situation where $\widehat{\mathbb{Q}}_{b, T_b} = \widetilde{\mathbb{Q}}_b$ for some $\widetilde{\mathbb{Q}}_b \in \mathcal{Q}_b$. Our goal is to pick the nominal distribution $\mathbb{Q}^*_b$ so as to satisfy the inequality
\begin{equation} \label{eq: weak optimality proof 3a}
\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} > 0
\end{equation}
in the right-hand side of (\ref{eq: weak optimality proof 3}). To this end, we assume that 
\begin{equation} \nonumber
\widetilde{\mathbb{Q}}^{(w)}_b \in \argmax_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} 
\end{equation}
and set $\mathbb{Q}^*_b = \lambda' \widetilde{\mathbb{Q}}_b + (1 - \lambda') \widetilde{\mathbb{Q}}_b^{(w)}$ for some $\lambda' < 0$. 
We observe that: 
\begin{equation} \label{eq: weak optimality proof 5} 
\begin{gathered}
\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} = \lambda' \mathbb{E}_{\widetilde{\mathbb{Q}}_b}\{c_b\} + (1 - \lambda') \mathbb{E}_{\widetilde{\mathbb{Q}}_b^{(w)}}\{c_b\} - \mathbb{E}_{\widetilde{\mathbb{Q}}^{(w)}_b}\{c_b\} = \\ = \lambda'(\mathbb{E}_{\widetilde{\mathbb{Q}}_b}\{c_b\} - \mathbb{E}_{\widetilde{\mathbb{Q}}_b^{(w)}}\{c_b\}) = \lambda'(\mathbb{E}_{\widetilde{\mathbb{Q}}_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\}) > 0 
\end{gathered}
\end{equation}
The last strict inequality is implied by the fact that a linear function $\mathbb{E}_{\mathbb{Q}_b} \{c_b\}$ is maximized over a convex set $\widehat{\mathcal{Q}}'_b$ and $\widetilde{\mathbb{Q}}_b \in \mbox{int}(\widehat{\mathcal{Q}}'_b)$. Thus, with the aforementioned choice of the empirical distribution of $\widehat{\mathbb{Q}}_{b, T_b}$ and the nominal distribution $\mathbb{Q}^*_b$, the inequality (\ref{eq: weak optimality proof 3a}) holds; recall the lower bound in (\ref{eq: weak optimality proof 3}).  

Furthermore, by an appropriate choice of $\lambda'$ we may guarantee that $\mathbb{Q}^*_b$ belongs to the initial relative entropy ball, i.e., the ball centered at $\widehat{\mathbb{Q}}_{b, T_b}$ with the radius $r_b$.  
More specifically, by the continuity of the logarithm for any $\delta > 0$ we have: 
\begin{equation} \nonumber
\begin{gathered}
D_{KL}(\widetilde{\mathbb{Q}}_b \; \Vert \; \mathbb{Q}^*_b) - D_{KL}(\widetilde{\mathbb{Q}}_b \; \Vert \; \widetilde{\mathbb{Q}}_b^{(w)}) =
\sum_{i=1}^{d_b} \widetilde{q}_{b,i} \Big(\ln \frac{\widetilde{q}_{b,i}}{\lambda' \widetilde{q}_{b,i} + (1 - \lambda')q^{(w)}_{b,i}} - \ln \frac{\widetilde{q}_{b,i}}{q^{(w)}_{b,i}}\Big) = \\ = \sum_{i=1}^{d_b} \widetilde{q}_{b,i} \ln \frac{q^{(w)}_{b,i}}{\lambda' \widetilde{q}_{b,i} + (1 - \lambda')q^{(w)}_{b,i}} \leq \delta,
\end{gathered}
\end{equation} 
if $\lambda' < 0$ is sufficiently close to zero. Furthermore, $D_{KL}(\widetilde{\mathbb{Q}}_b \; \Vert \; \widetilde{\mathbb{Q}}_b^{(w)}) \leq r'_b < r_b$ and hence,
\begin{equation} \label{eq: weak optimality proof 6}
D_{KL}(\widetilde{\mathbb{Q}}_b \; \Vert \; \mathbb{Q}^*_b)  \leq r'_b + \delta = \widetilde{r}_b < r_b
\end{equation} 
Finally, by leveraging (\ref{eq: weak optimality proof 5}) and (\ref{eq: weak optimality proof 6}) we observe that:  
\begin{equation} \label{eq: weak optimality proof 7} \tag{$LB_3$}
\begin{gathered}
\Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}'_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} > 0 \Big\} \geq \Pr\{\widehat{\mathbb{Q}}_{b,T_b} = \widetilde{\mathbb{Q}}_b\} \geq \\ \geq (T_{min} + 1)^{-d_b} e^{-T_{min} D_{KL}(\widetilde{\mathbb{Q}}_b \; \Vert \; \mathbb{Q}^*_b)} \geq (T_{min} + 1)^{-d_b} e^{-T_{min} \widetilde{r}_b}
\end{gathered}
\end{equation}
The second inequality in (\ref{eq: weak optimality proof 7}) follows from the lower bound (\ref{eq: strong LDP 2}). Let $\delta' = r_b - \widetilde{r}_b > 0$. Then bringing the lower bounds (\ref{eq: weak optimality proof 3}) and (\ref{eq: weak optimality proof 7}) together yields
\begin{equation} \label{eq: weak optimality result of step 1!} \tag{$LB_4$}
\begin{gathered}
\lim_{T_{min} \to +\infty} \sup \frac{1}{T_{min}} \ln \Big( \mbox{\upshape Pr} \Big\{\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b}\mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\varepsilon' \Big\} \Big) \geq \\
\geq \lim_{T_{min} \to +\infty} \sup \frac{1}{T_{min}} \ln \Big((T_{min} + 1)^{-d_b} e^{-T_{min} (r_b - \delta')}\Big) = -r + \delta' > -r
\end{gathered}
\end{equation}
Here, we exploit that 
$$r_b = \frac{1}{T_{min}}\Big(\ln|A| + d_b \ln(T_{min} + 1) + T_{min} r \Big) 
$$
by definition of $b$. 
This observation concludes the first step of the proof. 

\textbf{Step 2.}  In this part of the proof we pick the nominal distribution $\mathbb{Q}^* \in \mathcal{Q}$ (which marginal distribution with respect to the component $c_{b}$ is precisely defined at the first step) so as to demonstrate that $\hat{f}$ is not feasible in (\ref{prediction problem}).  Our choice of $\mathbb{Q}^*$ suggests a \textit{linear functional dependence} among the components $c_a$, $a \in A_{\widetilde{\mathbf{x}}}$. Formally, we assume that:
\begin{equation} \label{eq: linear dependence}
c_a = w_a c_b + u_a \quad \forall a \in A_{\widetilde{\mathbf{x}}}, \; a \neq b,
\end{equation}
where $0 < w_a < W$ for some $W \in \mathbb{R}_{>0}$. In particular, for $a \in A_{\widetilde{\mathbf{x}}}$ we have $q_{a,i} = q_{b,i}$ and the support of $c_a$ is given by:
$$z_{a,i} = w_a z_{b,i} + u_a, \quad i \in \{1, \ldots d_a\}$$
with $d_a = d_b$. 
Taking into account (\ref{eq: linear dependence}) we calculate the sum in (\ref{eq: weak optimality proof 2}), i.e.,  
\begin{equation} \label{eq: weak optimality proof 8} \tag{$LB_5$}
\begin{gathered}
\sum_{a \in A} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a} \mathbb{E}_{\mathbb{Q}_a}\{c_a\} \Big) \widetilde{x}_a = \sum_{a \in A_{\widetilde{\mathbf{x}}}} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a} \mathbb{E}_{\mathbb{Q}_a}\{c_a\} \Big) = \\ 
= \sum_{a \in A_{\widetilde{\mathbf{x}}}} \Big( w_a \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} + u_a - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_a} (w_a \mathbb{E}_{\mathbb{Q}_b}\{c_b\} + u_a) \Big) = \\
= \sum_{a \in A_{\widetilde{\mathbf{x}}}} w_a \Big(\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_a} \mathbb{E}_{\mathbb{Q}_b}\{c_b\}\Big) \geq \sum_{a \in A_{\widetilde{\mathbf{x}}}} w_a \Big(\mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b} \mathbb{E}_{\mathbb{Q}_b}\{c_b\}\Big) 
\end{gathered}
\end{equation}
The last inequality in (\ref{eq: weak optimality proof 8}) holds whenever $T_{min} \geq T_0$ for some $T_0 \in \mathbb{Z}_{>0}$ and follows from the fact that 
$$\widehat{\mathcal{Q}}_a \subseteq \widehat{\mathcal{Q}}_b \quad \forall a \in A_{\widetilde{\mathbf{x}}}, \; a \neq b$$
More precisely, if $T_{min}$ is sufficiently large, then $r_a \leq r_b$; see (\ref{eq: radius}). 

In conclusion, using (\ref{eq: weak optimality proof 2}), (\ref{eq: weak optimality result of step 1!}) and (\ref{eq: weak optimality proof 8}) for $T_{min} \geq T_0$ and $\varepsilon' = \frac{\varepsilon}{\sum_{a \in A_{\widetilde{\mathbf{x}}}} w_a}$ results in: 
\begin{equation} \label{eq: weak optimality proof 9} \tag{$LB_6$}
\begin{gathered}
\Pr\Big\{f(\widetilde{\mathbf{x}}, \mathbb{Q}^*) > \hat{f}(\widetilde{\mathbf{x}}, \widehat{\mathbf{Q}}_{\mathbf{T}})\Big\} \geq \Pr\Big\{\sum_{a \in A} \Big(\mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} - \max_{\mathbb{Q}_a \in \widehat{\mathcal{Q}}_a} \mathbb{E}_{\mathbb{Q}_a}\{c_a\} \Big) \widetilde{x}_a > -\varepsilon \Big\} \geq \\ 
\geq \Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b} \mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\frac{\varepsilon}{\sum_{a \in A_{\widetilde{\mathbf{x}}}} w_a} \Big\} \geq \\
\geq \Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b} \mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\frac{\varepsilon}{|A| W} \Big\} \geq (T_{min} + 1)^{-d_b} e^{-T_{min} \widetilde{r}_b}
\end{gathered}
\end{equation}
Thus, the prediction rule $\hat{f}$ does not satisfy the asymptotic guarantee (\ref{eq: asymptotic guarantee}) due to (\ref{eq: weak optimality result of step 1!}) and the result follows. 
\end{proof}
\end{theorem}

\begin{remark} \upshape
We conclude that the distributionally robust prediction rule (\ref{distributionally robust predictor}) caters a weakly optimal solution of (\ref{prediction problem}). However, this fact does not indicate whether (\ref{distributionally robust predictor}) is a strongly optimal solution for (\ref{prediction problem}) or not. We did not manage to answer this question by using the methodology from the proof of Theorem \ref{theorem 3} and thus, it is a possible direction for future research.  \vspace{-9mm}\flushright$\square$
\end{remark} 


\subsection{Properties of the distributionally robust prescriptor}
In this section we describe theoretical properties of the distributionally robust prescription rule (\ref{distributionally robust prescriptor}), which is defined by the pair of functions $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$. Recall that $\hat{f}_{\mathbf{r}}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}})$ serves an approximation of the expected cost in (\ref{stochastic programming problem}) for a fixed decision and $\hat{x}_{\mathbf{r}}(\widehat{\mathbf{Q}}_{\mathbf{T}})$ provides an estimate of the optimal decision. 

First, we claim that the pair $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ is feasible in (\ref{prescription problem}). Indeed, by Theorem \ref{theorem 1} the asymptotic guarantee (\ref{eq: asymptotic guarantee}) holds for any feasible decision $\mathbf{x} \in X$ and, in particular, for $\mathbf{x} = \hat{x}_{\mathbf{r}}(\widehat{\mathbf{Q}}_{\mathbf{T}})$. Furthermore, we demonstrate that $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ is a weakly optimal solution of (\ref{prescription problem}). In this regard, we slightly modify the definition of weak optimality.  

\begin{definition}[\textbf{Weak optimality, prescription}] \label{def: weak optimality 2} \upshape
A prediction-prescription pair $(\hat{f}, \hat{x})$ is called \textit{weakly optimal} for (\ref{prescription problem}), if it is feasible in (\ref{prescription problem}) and there exist no other feasible prediction-prescription pair that is less conservative than $\hat{f}$ uniformly over all distributions $\mathbb{Q} \in \mathcal{Q}$. \vspace{-9.5mm}\flushright$\square$  
\end{definition} 

\begin{theorem} \label{theorem 3}
Let $r > 0$ be a required exponential decay rate. 
Then the prediction-prescription pair $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ with $\mathbf{r} = \{r_a, a \in A\}$ given by (\ref{eq: radius}) is weakly optimal for (\ref{prescription problem}).
\begin{proof}
Assume that there exists another prediction-prescription pair $(\hat{f}, \hat{x})$, which is less conservative than $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ and is feasible in (\ref{prescription problem}). Therefore, for any decision $\mathbb{Q} \in \mathcal{Q}$ there exists $\varepsilon \in \mathbb{R}_{>0}$ such that:
\begin{equation} \label{eq: weak optimality proof pescription 1}
\hat{f}(\hat{x}(\mathbf{Q}), \mathbf{Q}) \leq \hat{f}_{\mathbf{r}}(\hat{x}_{\mathbf{r}}(\mathbf{Q}), \mathbf{Q}) - \varepsilon
\end{equation} 
In the remainder of the proof we show that $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ is infeasible in (\ref{prescription problem}).

To this end, we bound from below the out-of-sample disappointment of $\hat{f}$ at $\mathbf{x} = \hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}})$ for some data-generating distribution $\mathbb{Q}^* \in \mathcal{Q}$, that is,
\begin{equation} \label{eq: weak optimality proof prescription 2}
\begin{gathered}
\Pr\Big\{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}^*) > \hat{f}(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \widehat{\mathbf{Q}}_{\mathbf{T}})\Big\} \geq \Pr\Big\{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}^*) - \hat{f}_{\mathbf{r}}(\hat{x}_{\mathbf{r}}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \widehat{\mathbf{Q}}_{\mathbf{T}}) > -\varepsilon \Big\} \geq \\ 
\geq \Pr\Big\{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}^*) - \hat{f}_{\mathbf{r}}(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \widehat{\mathbf{Q}}_{\mathbf{T}}) > -\varepsilon \Big\}
\end{gathered}
\end{equation}
Here, the first inequality is implied by (\ref{eq: weak optimality proof pescription 1}) and the second inequality follows from the definition of $\hat{x}_{\mathbf{r}}$; recall (\ref{distributionally robust prescriptor}).

In contrast to Theorem \ref{theorem 2}, $\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}})$ is a random vector, which depends on the empirical marginal distributions. Actually, it suffices to reiterate the proof of Theorem \ref{theorem 2} with some minor changes concerning the choice of the nominal distribution $\mathbb{Q}^*$. Namely, we set $b = \argmin_{a \in A} T_a$ and impose a linear functional dependence over all the actions $a \in A$, that is,
\begin{equation} \label{eq: linear dependence prescription}
c_a = w_a c_b + u_a \quad \forall a \in A, \; a \neq b,
\end{equation}
where $0 < w_a < W$ for some $W \in \mathbb{R}_{>0}$ (the marginal distribution $\mathbb{Q}_b^*$ is constructed exactly as in the proof of Theorem \ref{theorem 2} .) 
This prerequisite allows us to bound the out-of-sample performance uniformly over all the decisions $\mathbf{x} \in X$.  
In fact, by using (\ref{eq: weak optimality proof 2}), (\ref{eq: weak optimality result of step 1!}) and (\ref{eq: weak optimality proof 8}) for $T_b \geq T_0$, $\widetilde{x} = \hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}})$ and $\varepsilon' = \frac{\varepsilon}{|A|W}$ observe that:
\begin{equation} \label{eq: weak optimality proof prescription 3}
\begin{gathered}
 \Pr\Big\{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}^*) - \hat{f}_{\mathbf{r}}(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \widehat{\mathbf{Q}}_{\mathbf{T}}) > -\varepsilon \Big\} \geq \\ \geq \Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b} \mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\frac{\varepsilon}{\sum_{a \in A_{\widetilde{\mathbf{x}}}} w_a} \Big\} \geq \\ \geq \Pr\Big\{ \mathbb{E}_{\mathbb{Q}^*_b}\{c_b\} - \max_{\mathbb{Q}_b \in \widehat{\mathcal{Q}}_b} \mathbb{E}_{\mathbb{Q}_b}\{c_b\} > -\frac{\varepsilon}{|A|W} \Big\} \geq (T_{min} + 1)^{-d_b} e^{-T_{min} (r_b - \delta')},
\end{gathered}
\end{equation}
where $\delta' > 0$. Following the proof of Theorem \ref{theorem 2} we claim that the out-of-sample prescription disappointment for $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ decays at the rate of at most $r - \delta'$. 
Hence, the pair $(\hat{f}_{\mathbf{r}}, \hat{x}_{\mathbf{r}})$ is infeasible in (\ref{prescription problem}) and the result follows. 
\end{proof}
\end{theorem}



In view of the discussion in Section \ref{subsec: prediction and prescription problems}, our formulations (\ref{prediction problem}) and (\ref{prescription problem}) are slightly different from those in \cite{VanParys2020}; recall Table \ref{tab: difference}. As a result, prediction and prescription rules (\ref{distributionally robust predictor}) and (\ref{distributionally robust prescriptor}) can be computed rather effectively whenever there exist effective algorithms for the underlying combinatorial optimization problem; recall Proposition \ref{proposition 1}. Meanwhile, a naive implementation of the DRO approach from \cite{VanParys2020} with a discrete set of feasible decisions, $X$, results in a non-linear MIP problem; see Section \ref{subsec: benchmark approaches} for details. As a remark, we claim that Theorems \ref{theorem 1} - \ref{theorem 3} remain valid, if all feasible decisions are merely nonnegative and bounded.

Second, we allow the number of random observations for individual actions $a \in A$ to vary. In this setting, we propose a new class of prediction rules (recall the last row of Table \ref{tab: difference}), where (\ref{distributionally robust predictor}) and (\ref{distributionally robust prescriptor}) serve as a weakly optimal solution for (\ref{prediction problem}) and (\ref{prescription problem}), respectively. On the other hand, the DRO approach in \cite{VanParys2020} stipulates a strongly optimal solution for their prediction and prescription problems; see Theorems 4 and 7 in \cite{VanParys2020}.


\subsection{Improved finite sample guarantees}
We note that most of the previous theoretical results have an asymptotic flavor, i.e., we limit the out-of-sample disappointment (recall Definition \ref{def: out-of-sample prediction disappointment}) whenever $T_{min}$ tends to infinity. In practice we need to exploit the finite sample guarantee (\ref{eq: finite sample guarantee}), which is relevant for any fixed sample size.  

At the same time, the finite sample guarantee of Theorem \ref{theorem 1} is implied by the strong LDP bound (\ref{eq: strong LDP 1}), which contains a polynomial term $(T_a + 1)^{d_a}$. If the number of observations, $T_a$, $a \in A$, is sufficiently small, then this polynomial has a rather negative effect on the quality of distributionally robust decisions. Recently Agrawal \cite{Agrawal2020} proposed a tighter upper bound for the probability in (\ref{eq: strong LDP 1}) induced by a thorough analysis of the moment generating function of the relative entropy distance. More specifically, for any $r_a > \frac{d_a-1}{T_a}$ the following inequality holds (see Theorem 2.1 in \cite{Agrawal2020}):
\begin{equation} \label{eq: strong LDP 3} \nonumber
\mbox{\upshape Pr} \Big\{D_{KL}(\mathbb{Q}_a \; \Vert \; \mathbb{Q}^*_a) > r_a \Big\} \leq \Big( \frac{e}{d_a-1} r_a T_a  \Big)^{d_a-1} e^{-r_a T_a}
\end{equation}        

For this reason, in lieu of the formula (\ref{eq: radius}) the radius of the relative entropy ball in (\ref{ambiguity set}) can be derived by solving the following equation with respect to $r_a > \frac{d_a-1}{T_a}$:
\begin{equation} \label{eq: radius improved} 
\Big( \frac{e}{d-1} r_a T_a  \Big)^{d_a-1} e^{-r_a T_a} = \alpha_a,
\end{equation}
where $\alpha_a \in (0, 1)$ is a given confidence level.
It is rather straightforward to verify that this equation has at most one feasible solution. 

Another class of upper bounds for the probability in (\ref{eq: strong LDP 1}) is proposed by Mardia~et~al.~\cite{Mardia2020}. In order to provide consistency and completeness of the current manuscript we reformulate the key theoretical result of \cite{Mardia2020} by using our notation; see Theorem 3 in \cite{Mardia2020}. 

Assume that $u_0 = \pi$, $u_1 = 2$ and $K_{-1} = 1$. Define
\begin{equation} \nonumber
u_i = \begin{cases}\pi \times \frac{1 \times 3 \times \ldots \times (i-1)}{2 \times 4 \times \ldots \times i}, \mbox{ if } i \mbox{ is even and } i \geq 2 \\
2 \times \frac{2 \times 4 \times \ldots \times (i-1)}{1 \times 3 \times \ldots \times i}, \mbox{ if } i \mbox{ is odd and } i \geq 3
\end{cases} \mbox{ and } K_{j} = \prod_{i = 0}^{j}u_i \quad \forall j \geq 1
\end{equation}
Then for any $d_a \geq 2$, $T_a \geq 2$ and $a \in A$ Theorem 3 in \cite{Mardia2020} provides the following upper bound:
\begin{equation} \label{eq: radius improved 2}
\mbox{\upshape Pr} \Big\{D_{KL}(\mathbb{Q}_a \; \Vert \; \mathbb{Q}^*_a) > r_a \Big\} \leq \Big(\frac{3u_1}{u_2} \sum_{j = 0}^{d_a-2} K_{j-1} (\frac{e\sqrt{T_a}}{2\pi})^j\Big) e^{-T_a r_a}
\end{equation}
Eventually, the radius $r_a$ can be obtained by setting the right-hand side of (\ref{eq: radius improved 2}) equal to $\alpha_a$. 


It is rather straightforward to verify that the estimate of the radius (\ref{eq: radius improved 2}) is typically smaller than the other estimates, (\ref{eq: radius}) and (\ref{eq: radius improved}). Nevertheless, in our computational experiments we compute $r_a$, $a \in A$, by leveraging all the proposed techniques and adhere the minimal found estimate of the radius.
 

\section{Computational study} \label{sec: comp study}
This section is focused on a class of data-driven shortest path problems (SPP), where $A$ is set of directed arcs, $\mathbf{c}$ is a vector of associated costs/travel times and $X$ is a set of all simple paths between two fixed nodes in the network. We provide a numerical comparison of our prediction and prescription rules, namely, (\ref{distributionally robust predictor}) and (\ref{distributionally robust prescriptor}), with some standard confidence bounds for the mean of the cost vector. In addition, the proposed DRO approach is compared with some modification of the DRO approach described in \cite{VanParys2020}. In general, we examine different parameter settings with respect to both the form of the nominal distribution $\mathbb{Q}^*$ and the sample size.    

For some given $\mathbb{Q}^* \in \mathcal{Q}$ the quality of a prescription rule $\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}})$ is measured by using a \textit{nominal relative loss}, that is,
\begin{eqnarray} \label{eq: nominal relative loss}
\rho(\mathbb{Q}^*, \widehat{\mathbf{Q}}_{\mathbf{T}}) = \frac{f(\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}), \mathbb{Q}^*)}{\min_{\mathbf{x} \in X} f(\mathbf{x}, \mathbb{Q}^*)} 
\end{eqnarray} 
In fact, $\rho \geq 1$ evaluates the out-of-sample performance of $\hat{x}(\widehat{\mathbf{Q}}_{\mathbf{T}}) \in X$ under the nominal distribution $\mathbb{Q}^*$. Ideally, we have $\rho = 1$, while it is not the case in general as long as the decision-maker operates with a limited information about $\mathbb{Q}^*$.  

All experiments are performed on a PC with \textit{CPU i7-9700} and \textit{RAM 32 GB}. The deterministic version of the shortest path problem is solved in Python with \textit{CPLEX 12.10.0.0}. The dual formulations in (\ref{individual dual problems}) are tackled using the function \textit{scipy.optimize.minimize}() and the method of Nelder-Mead with default parameters. In particular, we verify that the strong duality holds by solving the primal optimization problems (\ref{worst-case expectation problem}) with \textit{CVX 1.0.31}. Finally, the integrals are computed numerically using the function \textit{scipy.integrate.quad()} in Python (we refer to the definition of a discretized normal distribution in Table \ref{tab: data-generating distributions}.) 

\subsection{Benchmark approaches} \label{subsec: benchmark approaches}
The idea of improving the sample-average approximation for the expected value in (\ref{stochastic programming problem}) is not novel. For example, Anderson and Philpott \cite{Anderson2019} analyze different mechanisms that hedge against uncertainty in a more sophisticated way, e.g., a CVaR-based risk measure, phi-divergence using total variation and a Wasserstein metric. At the same time, we are not aware about DRO formulations (except for a formulation in \cite{Ketkov2021}) that can address incomplete data sets. 

We consider two possible approximations of the stochastic programming problem in~(\ref{stochastic programming problem}). Henceforth, instead of the exponential decay rate, $r > 0$, we focus on the confidence level $\alpha \in (0, 1)$. Recall that for any given $\alpha$ we have: 
$$r = -\frac{1}{T_{min}}\ln \alpha$$   
Also, for simplicity of exposition we assume that the support of $c_a$ for each $a \in A$ is given by the integers $\{1, \ldots, d\}$, i.e., $d_a = d$ and $z_{a,i} = i$ for $i \in \{1, \ldots, d\}$.  

\textbf{Confidence bounds for the mean.} Following the idea of Example \ref{example 1} we construct an upper confidence bound for the mean of $c_a$, $a \in A$, via Hoeffding's inequality. Formally, for any $a \in A$ and $\varepsilon_a > 0$ we have:
\begin{equation} \label{hoeffding bound} 
\begin{gathered} 
\Pr\Big\{\frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} - \mathbb{E}_{\mathbb{Q}^*_a}\{c_a\} > \varepsilon_a\Big\} \leq \exp\Big(-2 T_a \frac{\varepsilon_a^2}{(d - 1)^2} \Big),
\end{gathered}
\end{equation}
where value of parameter $\varepsilon_a$ can be obtained by setting the right-hand side of (\ref{hoeffding bound}) equal to $\alpha_a \in (0, 1)$, where $\sum_{a \in A}{\alpha_a} = \alpha$.

Then we consider the following prediction and prescription rules, respectively:
\begin{subequations} 
	\begin{align}
	& \hat{f}_{hoef}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) = \sum_{a \in A} \min\Big \{\frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} + \varepsilon_a; d \Big \} x_a \label{hoeffding predictor} \\
	& \hat{x}_{hoef}(\widehat{\mathbf{Q}}_{\mathbf{T}}) \in \argmin_{\mathbf{x} \in X} \hat{f}_{hoef}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) 
	\label{hoeffding prescriptor}
	\end{align}
\end{subequations}
In fact, we solve a deterministic version of the shortest path problem with the worst-case expected costs induced by (\ref{hoeffding bound}).  Let $A' \subseteq A$ be a set of actions for which the minimum in the right-hand side of (\ref{hoeffding predictor}) equals $d$. By using the union bound for some fixed $\mathbf{x} \in X$ we bound the out-of-sample disappointment of (\ref{hoeffding predictor}) as follows:
\begin{equation} \nonumber
\begin{gathered}
\Pr\Big\{f(\mathbf{x}, \mathbb{Q}^*) - \hat{f}_{hoef}(\mathbf{x}, \widehat{\mathbf{Q}}_{\mathbf{T}}) > 0 \Big\} = \\ = 
\Pr\Big\{\sum_{a \in A} \Big( \mathbb{E}_{\mathbb{Q}^*_a} \{c_a\} - \min \{\frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} + \varepsilon_a; d \} \Big)x_a > 0 \Big\} \leq \\
\leq \Pr\Big\{\sum_{a \in A \setminus A'} \Big( \mathbb{E}_{\mathbb{Q}^*_a} \{c_a\} - \frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j} - \varepsilon_a \Big)x_a > 0 \Big\} \leq \\ \leq \Pr\Big\{\bigcup_{a \in A} \Big( \mathbb{E}_{\mathbb{Q}^*_a} \{c_a\} - \frac{1}{T_a}\sum_{j = 1}^{T_a} \widehat{c}_{a,j}  > \varepsilon_a \Big)\Big\} \leq \sum_{a \in A} \alpha_a = \alpha
\end{gathered}
\end{equation}
Here, we exploit that $\mathbb{E}_{\mathbb{Q}^*_a} \{c_a\} \leq d$ by construction of the support. 
Hence, the out-of-sample disappointment of (\ref{hoeffding predictor}) is bounded from above by $\alpha$ for any feasible $\mathbf{x} \in X$. 

\textbf{Truncated DRO methods (DRO$_1$ and DRO$_2$).}  Here, we consider several alternative approaches based on a truncation of the data set $\widehat{\mathbf{C}}$. Namely, the truncated DRO$_1$ approach inherits the idea of Van Parys~et~al. \cite{VanParys2020} and the truncated DRO$_2$ approach applies our prediction and prescriptions rules (\ref{distributionally robust predictor}) and (\ref{distributionally robust prescriptor}) to a truncated data set.  

First, in view of the discussion in Section \ref{sec: intro}, we examine some modification of the optimization approach described in \cite{VanParys2020}.  Since the study in \cite{VanParys2020} operates with complete data sets we store only $T_{min}$ random observations for each particular action $a \in A$.  By leveraging the truncated data set we can construct a unique ball with respect to the relative entropy distance centered at the joint empirical distribution of the data. The outlined prediction and prescription rules can be described by their dual formulations as follows  (we omit some minor technical details for brevity, see Proposition 2 in \cite{VanParys2020}): 
\begin{subequations} 
	\begin{align}
	& \hat{f}_{trunc}(\mathbf{x}, \widetilde{\mathbb{Q}}_{T_{min}}) = \min_{\beta}\Big\{\beta - e^{-r}\prod_{i = 1}^{d^{|A|}} (\beta - (\widetilde{\mathbf{c}}^{(i)})^\top \mathbf{x})^{\widetilde{q}_{i}}: \; \beta \geq d \sum_{a \in A}x_a\Big\} \label{truncated predictor} \\
	& \hat{x}_{trunc}(\widetilde{\mathbb{Q}}_{T_{min}}) \in \argmin_{\mathbf{x} \in X} \hat{f}_{trunc}(\mathbf{x}, \widetilde{\mathbb{Q}}_{T_{min}}),
	\label{truncated prescriptor}
	\end{align}
\end{subequations}
where $\widetilde{\mathbb{Q}}_{T_{min}}$ is an empirical distribution of the cost vector $\mathbf{c}$ with respect to the truncated data set, $d^{|A|}$ is a number of possible realizations of the cost vector and $\widetilde{q}_{i}$ is an empirical probability that $\mathbf{c} = \widetilde{\mathbf{c}}^{(i)}$, $i \in \{1, \ldots, d^{|A|}\}$. Eventually, $r$ is a radius of the relative entropy ball centered at the joint empirical distribution $\widetilde{\mathbb{Q}}_{T_{min}}$. Henceforth, the outlined solution approach is referred to as a truncated DRO$_1$ approach.   

Despite the prediction and prescription rules (\ref{truncated predictor}) and (\ref{truncated prescriptor}) are strongly optimal in the case of continuous $X$ (among all continuous functions of the empirical \textit{joint distribution}, see Theorems 4 and 7 in \cite{VanParys2020}), we argue that effectiveness of the outlined solutions in a combinatorial optimization setting is substantially limited.  That is, the objective function in (\ref{truncated prescriptor}) is non-linear and the feasible set is mixed-integer. Hence, computation of (\ref{truncated prescriptor}) results in a non-linear MIP problem.  

In our computational experiments we simply enumerate all the decisions $\mathbf{x} \in X$ and solve the resulting univariate convex optimization problems in the parameter $\beta$. Therefore, a further implementation of the truncated DRO$_1$ approach is restricted to comparatively small instances of the shortest path problem.       

\begin{remark} \upshape
In order to provide a fair comparison between the truncated DRO problem induced by (\ref{truncated prescriptor}) and the baseline problem (\ref{DRO problem}) we need to compute the radius $r$ in three different ways. Namely, we use (\ref{eq: radius}), (\ref{eq: radius improved}) and (\ref{eq: radius improved 2}) assuming that the number of random observations and possible realizations of the cost vector are given by $T_{min}$ and $d^{|A|}$, respectively. In particular, the partial sum in (\ref{eq: radius improved 2}) is computed by eliminating zero terms in the sense of a floating-point precision. \vspace{-9mm}\flushright$\square$
\end{remark} 

As an alternative, we also examine the baseline DRO approach induced by (\ref{distributionally robust prescriptor}) with the truncated data set (this approach is referred to as a truncated DRO$_2$ approach). Put differently, we assume that $T_a = T_{min}$ for every $a \in A$. Similar to the truncated DRO$_1$ method we do not take into account a part of the data, but additionally preserve the complexity of (\ref{distributionally robust prescriptor}).     

\subsection{Test instances}
In our experiments we consider \textit{a fully-connected layered graph} with $h$ intermediate layers and $w$ nodes at each layer $i \in \{1, \ldots, h\}$. The first and the last layer consist of unique nodes, which are the source and the destination nodes, respectively. For example, a network with $h = 3$ and $w = 3$ is depicted in Figure \ref{fig: layered graph}.

\begin{figure}
	\centering
	\begin{tikzpicture}[scale=0.75,transform shape]
	\Vertex[x=0,y=0]{1}
	\Vertex[x=3,y=3]{2}
	\Vertex[x=3,y=0]{3}
	\Vertex[x=3,y=-3]{4}
	\Vertex[x=6,y=3]{5}
	\Vertex[x=6,y=0]{6}
	\Vertex[x=6,y=-3]{7}
	\Vertex[x=9,y=3]{8}
	\Vertex[x=9,y=0]{9}
	\Vertex[x=9,y=-3]{10}
	\Vertex[x=12,y=0]{11}
	\tikzstyle{LabelStyle}=[fill=white,sloped]
	\tikzstyle{EdgeStyle}=[post]
	\Edge(1)(2)
	\Edge(1)(3)
	\Edge(1)(4)
	\Edge(8)(11)
	\Edge(9)(11)
	\Edge(10)(11)
	\Edge(2)(5)
	\Edge(2)(6)
	\Edge(2)(7)
	\Edge(3)(5)
	\Edge(3)(6)
	\Edge(3)(7)
	\Edge(4)(5)
	\Edge(4)(6)
	\Edge(4)(7)
	\Edge(5)(8)
	\Edge(5)(9)
	\Edge(5)(10)
	\Edge(6)(8)
	\Edge(6)(9)
	\Edge(6)(10)
	\Edge(7)(8)
	\Edge(7)(9)
	\Edge(7)(10)
	\end{tikzpicture}
	\caption{\footnotesize A fully-connected layered graph with $h = 3$ intermediate layers and $w = 3$ nodes at each layer.}
	\label{fig: layered graph}
\end{figure}

\begin{table}
	
	\footnotesize
	\centering
	\onehalfspacing
	\begin{tabular}{c c c c c}
		\hline
		Joint distribution & Marginal distributions & Support & Mean & Variance\\
		\hline
		\multirow{2}{*}{Product of marginals} & \multirow{2}{*}{Binomial($p_a$, $d - 1$)} & \multirow{2}{*}{$c_a \in \{1, \ldots, d\}$} & \multirow{2}{*}{$(d - 1)p_a + 1$} & \multirow{2}{*}{$p_a(1 - p_a)$} \\
		& & & & \\\hline
		\multirow{2}{*}{Multinomial($\mathbf{p}$, $\mathbf{d - 1}$)} & Binomial($p_a$, $d - 1$) & $c_a \in \{1, \ldots, d\}$, & \multirow{2}{*}{$(d - 1)p_a + 1$} & \multirow{2}{*}{$p_a(1 - p_a)$}\\
		& with $\sum_{a \in A} p_a = 1$ & $\sum_{a \in A} c_a = d - 1 + |A|$ & & \\\hline
	    \multirow{2}{*}{Product of marginals} & Discretization of normal & \multirow{2}{*}{$c_a \in \{1, \ldots, d\}$} & \multirow{2}{*}{$\sum_{i = 1}^d i q^*_{a, i}$} & $\sum_{i = 1}^d i^2 q^*_{a, i} - $\\
	    & with mean $\mu_a$, variance $\sigma_a$ & & & $\Big(\sum_{i = 1}^d i q^*_{a, i}\Big)^2$\\
		\hline
	\end{tabular}
	\caption{\scriptsize The table summarizes discrete distributions that are used to model $\mathbb{Q}^*$. The values of binomial and multinomial distributions are shifted to the right by one. The support, mean and variance are component-wise. Parameters $q^*_{a, i}$, $a \in A$, of the normal distributions are computed using (\ref{eq: prob normal discretized}).}
	\label{tab: data-generating distributions}
\end{table}

In addition, we propose several types of data-generating distributions that are restricted to take values in $\{1, \ldots, d \}$; see Table \ref{tab: data-generating distributions}. The binomial and multinomial distributions are standard, but shifted to the right by one. As long as the variance of these two distributions cannot be controlled directly, we propose a discretization of univariate normal distributions for each $a \in A$. That is, for given $\mu_a$, $\sigma_a$ and a probability density function of the normal distribution $$\theta_a(x) = \frac{1}{\sqrt{2\pi} \sigma_a} e^{\frac{-(x - \mu_a)^2}{2\sigma_a^2}}$$ 
the nominal probabilities $q^*_{a, i}$, $i \in \{1, \ldots, d\}$, are computed as follows:
\begin{equation} \label{eq: prob normal discretized}
q^*_{a, i} = \Pr\{c_a = i\} = \frac{1}{C}\int_{i - 0.5}^{i + 0.5} \theta_a(x)\mbox{d}x,
\end{equation}
where $C > 0$ is a normalization constant, which guarantees that the sum of probabilities is equal to $1$. 

Next, in order to pick the values of $T_a$, $a \in A$, we fix some $\widetilde{T}_{min} \in \mathbb{Z}_{>0}$  and set  
\begin{equation} \label{eq: T max}
\widetilde{T}_{max} = \widetilde{T}_{min} + \Delta,
\end{equation}
where $\Delta > 0$ is some positive integer constant. By leveraging the values of $\widetilde{T}_{min}$ and $\widetilde{T}_{max}$ we set the values of $T_a$, $a \in A$, in three different ways; see Table \ref{tab: choice of T}.

 The intuition behind the choice of $T_a$ is as follows.  The Binomial$_1$ distribution sorts the mean values of $T_a$, $a \in A$, in an increasing order with respect to their nominal expected costs. Oppositely, the Binomial$_2$ distribution sorts the mean values of $T_a$, $a \in A$, in the decreasing order. In other words, in the former situation the actions with higher expected costs can be observed sufficiently often; in the latter situation the same holds for the actions with smaller expected costs.  Finally, the uniform distribution is somewhat in the middle between the Binomial$_1$ and Binomial$_2$ distributions. 

In the next section we demonstrate that the distribution of $T_a$, $a \in A$, substantially affects the out-of-sample tests both for Hoeffding bounds and the DRO approach induced by (\ref{distributionally robust prescriptor}).   

\begin{table}	
	\footnotesize
	\centering
	\onehalfspacing
	\begin{tabular}{c c c c c}
		\hline
		Distribution of $T_a$, $a \in A$ & Support & Parameters \\
		\hline
			\multirow{2}{*}{Uniform} & 	\multirow{2}{*}{$[\widetilde{T}_{min}, \widetilde{T}_{max}]$} & 	\multirow{2}{*}{-} \\
			& & \\\hline
		\multirow{2}{*}{Binomial$_1$($p_a, \widetilde{T}_{max} - \widetilde{T}_{min}$)} & \multirow{2}{*}{$[\widetilde{T}_{min}, \widetilde{T}_{max}]$} & \multirow{2}{*}{$p_a = \frac{\mathbb{E}_{\mathbb{Q}_a^*}\{c_a\} - \min_{b \in A}\mathbb{E}_{\mathbb{Q}_b^*}\{c_b\}}{\max_{b \in A}\mathbb{E}_{\mathbb{Q}_b^*}\{c_b\} - \min_{b \in A}\mathbb{E}_{\mathbb{Q}_b^*}\{c_b\}}$} \\
		& & \\\hline
		\multirow{2}{*}{Binomial$_2$($p_a, \widetilde{T}_{max} - \widetilde{T}_{min}$)} & \multirow{2}{*}{$[\widetilde{T}_{min}, \widetilde{T}_{max}]$} & \multirow{2}{*}{$p_a = 1 - \frac{\mathbb{E}_{\mathbb{Q}_a^*}\{c_a\} - \min_{b \in A}\mathbb{E}_{\mathbb{Q}_b^*}\{c_b\}}{\max_{b \in A}\mathbb{E}_{\mathbb{Q}_b^*}\{c_b\} - \min_{b \in A}\mathbb{E}_{\mathbb{Q}_b^*}\{c_b\}}$} \\
		& & \\\hline
	\end{tabular}
	\caption{\scriptsize The table summarizes the ways to pick $T_a$, $a \in A$. The binomial distributions are shifted to the right by $\widetilde{T}_{min}$.}
	\label{tab: choice of T}
\end{table}
\subsection{Results and discussion} \label{subsec: results and discussion}
In view of the discussion above, we compare the aforementioned solution approaches in terms of the nominal relative loss (\ref{eq: nominal relative loss}). Specifically, we compute the average relative loss and median absolute deviations around the mean (MADs) over $N_0 = 200$ randomly generated test instances. Also, in all the experiments we set $d = 50$ and $\alpha = 0.05$. The confidence levels $\alpha_a$, $a \in A$, for individual distributional constraints are supposed to be in an inverse ratio with the numbers of observations, $T_a$, that is,
$$\alpha_a = \frac{\alpha}{T_a} \Big(\sum_{\widetilde{a} \in A} \frac{1}{T_{\widetilde{a}}}\Big)^{-1}$$
In particular, according to the union bound we have $\sum_{a \in A} \alpha_a = \alpha$. Finally, we assume that $w = 4$ and $h = 7$, if other is not specified. 

\subsubsection{Dependence on the form of the nominal distribution}
Throughout this section we assume that the sample sizes $T_a$, $a \in A$, are distributed uniformly on the interval $[\widetilde{T}_{min}, \widetilde{T}_{max}]$ and $\Delta = 5$; recall the first row of Table \ref{tab: choice of T} and (\ref{eq: T max}). The parameters of binomial (multinomial) and normal distributions for each $a \in A$ are set as $p_a \sim \mathcal{U}(0,1)$ and $\mu_a \sim \mathcal{U}(1,d)$, respectively; recall Table \ref{tab: data-generating distributions}. As a preliminary step, we focus on a comparison of the baseline DRO approach (\ref{distributionally robust prescriptor}) with standard confidence bounds for the expected loss (\ref{hoeffding prescriptor}). 

In our first experiment (see Figures \ref{fig: binomial, Tmin} and \ref{fig: multinomial, Tmin}) we consider the nominal relative loss (\ref{eq: nominal relative loss}) as a function of $\widetilde{T}_{min}$, $\widetilde{T}_{min} \in \{5, 7, \ldots, 35\}$, under the binomial and multinomial distributions, respectively. It turns out that the baseline DRO approach provides better out-of-sample tests for the binomial distribution, but demonstrates a rather poor out-of-sample performance in the case of multinomial distribution.

The intuition behind this fact can be described in terms of the in-sample performance. That is, in Figures \ref{fig: binomial, costs} and \ref{fig: multinomial, costs} we consider the nominal costs (sorted in an increasing order) and associated worst-case expected costs averaged over $N_0 = 200$ test instances with $\widetilde{T}_{min} = 25$. The key observation is that the smaller expected costs are better estimated using Hoeffding bounds, while the larger costs are better estimated using the baseline DRO approach. By construction the multinomial distribution has an additional  constraint $\sum_{a \in A} p_a = 1$ (see Table \ref{tab: data-generating distributions}), which alludes smaller expected costs than those for the binomial distribution. This observation provides a practical evidence behind the plots in Figures \ref{fig: binomial, Tmin} and \ref{fig: multinomial, Tmin}. 

Eventually, we set $\widetilde{T}_{min} = 25$ and consider the nominal relative loss as a function of variance $\sigma_a = \sigma \in \{1, 3, \ldots, d - 1\}$, $a \in A$, under the discretized normal distribution; see Figure \ref{fig: normal distribution, variance}. Observe that the baseline DRO approach outperforms Hoeffding bounds whenever the variance is sufficiently small. 

Our hypothesis is that for a sufficiently large variance the estimates of larger costs become more conservative and approach to the upper bound, $d$. Therefore, the advantages of the baseline DRO approach illustrated in Figure \ref{fig: binomial, costs} become less pronounced for marginal distributions with a relatively high variance.  

\begin{figure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}
		\includegraphics[width = \linewidth]{graphs_clear/exp2/graph_loss.jpg}
		\caption{Binomial distribution.}	
		\label{fig: binomial, Tmin}
	\end{subfigure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}	
		\includegraphics[width = \linewidth]{graphs_clear/exp3/graph_loss.jpg}
		\caption{Multinomial distribution.}
		\label{fig: multinomial, Tmin}
	\end{subfigure}
	\caption{Average relative loss (\ref{eq: nominal relative loss}) and MADs as a function of $\widetilde{T}_{min}$, $\widetilde{T}_{min} \in \{5, 7 \ldots, 35\}$, with $\Delta = 5$ under the (a) binomial and (b) multinomial distributions. The distribution of $T_a$, $a \in A$, is uniform.}
	\label{fig: binomial and multinomial, Tmin}
\end{figure}


\begin{figure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}
		\includegraphics[width = \linewidth]{graphs_clear/exp2/graph_costs.jpg}
		\caption{Binomial distribution.}	
		\label{fig: binomial, costs}
	\end{subfigure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}	
		\includegraphics[width = \linewidth]{graphs_clear/exp3/graph_costs.jpg}
		\caption{Multinomial distribution.}
		\label{fig: multinomial, costs}
	\end{subfigure}
	\caption{Average nominal expected costs (in an increasing order) and associated worst-case expected costs under the (a) binomial and (b) multinomial distributions with  $\widetilde{T}_{min} = 25$ and $\Delta = 5$. The distribution of $T_a$, $a \in A$, is uniform.}
	\label{fig: binomial and multinomial, costs}
\end{figure}
 
\subsubsection{Dependence on the sample size} \label{subsec: dependence on the sample size}
Henceforth, the data-generating distribution is supposed to be binomial. We explore how the distribution of $T_a$, $a \in A$, affects the out-of-sample performance of the baseline DRO approach and Hoeffding bounds. Specifically, we set $\Delta = 10$ and apply the Binomial$_1$ and Binomial$_2$ distributions from Table \ref{tab: choice of T}. The corresponding plots of the nominal relative loss (\ref{eq: nominal relative loss}) as a function of $\widetilde{T}_{min}$ are depicted in Figures \ref{fig: binomial T} and \ref{fig: binomial T reversed}. 

We conclude that the DRO approach outperforms Hoeffding bounds only in the case of Binomial$_2$ distribution. 
This observation can be justified in the context of Figure \ref{fig: binomial, costs}. That is, under the Binomial$_2$ distribution the smaller costs are observed more often and thus, the worst-case expected costs induced by the DRO approach are improved for sufficiently small expected costs. The opposite holds for Hoeffding bounds and larger expected costs under the Binomial$_1$ distribution.  

\begin{figure}
	\begin{center}
		\includegraphics[width = 0.5\linewidth]{graphs_clear/exp4/graph_loss.jpg}
		\caption{Average relative loss (\ref{eq: nominal relative loss}) and MADs as a function of $\sigma$, $\sigma \in \{1,3, \ldots, 49\}$, under the discretized normal distribution with $\widetilde{T}_{min} = 25$ and $\Delta = 5$. The distribution of $T_a$, $a \in A$, is uniform.}
		\label{fig: normal distribution, variance}
	\end{center}
\end{figure}

\begin{figure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}
		\includegraphics[width = \linewidth]{graphs_clear/exp5a/graph_loss.jpg}	
		\caption{Binomial distribution with Binomial$_1$ $T_a$, $a \in A$.}	
		\label{fig: binomial T}
	\end{subfigure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}
		\includegraphics[width = \linewidth]{graphs_clear/exp5b/graph_loss.jpg}
		\caption{Binomial distribution with Binomial$_2$ $T_a$, $a \in A$.}
		\label{fig: binomial T reversed}
	\end{subfigure}
	\caption{Average relative loss (\ref{eq: nominal relative loss}) and MADs as a function of $\widetilde{T}_{min}$, $\widetilde{T}_{min} \in \{5, 7 \ldots, 35\}$, with $\Delta = 10$ under the binomial distribution. The distributions of $T_a$, $a \in A$, are assumed to be (a) Binomial$_1$ and (b) Binomial$_2$.}
	\label{fig: binomial with binomial T}
\end{figure}

Importantly, the construction of Binomial$_2$ distribution can be also motivated from an application perspective. That is, the decision-maker may collect the data by traversing between the source and destination nodes in the network. In this setting the arcs with smaller expected costs are more preferable, if the decision-maker aims to minimize its expected loss based on the currently observed data. This process falls into the category of online learning problems and is not in the scope of the current paper. However, we note that the described data collection process should be consistent with the exploration-exploitation trade-off~\cite{Sutton2018}. 

Next, we assume that the distribution of $T_a$, $a \in A$, is uniform and customize the gap $\Delta = \widetilde{T}_{max} - \widetilde{T}_{min}$. To this end, in Figure \ref{fig: binomial distribution, variance} we plot the nominal relative loss (\ref{eq: nominal relative loss}) as a function of $\Delta$, $\Delta \in \{0, 2, \ldots, 40\}$, with $\widetilde{T}_{min} = 5$. 

\begin{figure}[htp] 
	\begin{center}
		\includegraphics[width = 0.5\linewidth]{graphs_clear/exp6/graph_loss.jpg}
			\caption{Average relative loss (\ref{eq: nominal relative loss}) and MADs as a function of $\Delta$, $\Delta \in \{0, 2, \ldots, 40\}$, under the binomial distribution with $\widetilde{T}_{min} = 5$.}
		\label{fig: binomial distribution, variance}	
\end{center}
\end{figure}

As a result, the baseline DRO approach outperforms Hoeffding bounds whenever the gap $\widetilde{T}_{max} - \widetilde{T}_{min}$ is sufficiently small. In fact, if the gap is large, then smaller expected costs can be observed less frequently than larger expected costs, which provides an advantage to the approach based on Hoeffding's inequality; recall Figure \ref{fig: binomial, costs}. One may argue that a similar advantage is given to the baseline DRO approach (if an opposite situation holds), but we note that smaller expected costs are more likely to be contained in an optimal solution.   

\subsubsection{Truncation of the data}
In this section we assume that the nominal distribution is discretized normal. Furthermore, we note that the complexity of the truncated DRO$_1$ approach substantially increases with the increase of the network's size. Hence, we focus on comparatively small instances of the shortest path problem with $h = w = 3$. 

The distribution of $T_a$, $a \in A$, is assumed to be uniform and the dependence on $\widetilde{T}_{max}$ is of interest since the truncated methods do not take into account a part of the data. In Figure \ref{fig: normal distribution, sigma = d/4, T max} we depict the nominal relative loss (\ref{eq: nominal relative loss}) as a function of $\Delta$, $\Delta \in \{0, 2, \ldots, 40\}$, with $\widetilde{T}_{min} = 10$ and $\sigma_a = \sigma = \frac{d}{4}$ for each $a \in A$. 

\begin{figure}[htp]
	\begin{center}
		\includegraphics[width = 0.5\linewidth]{graphs_clear/exp7/graph_loss.jpg}
			\caption{Average relative loss (\ref{eq: nominal relative loss}) and MADs as a function of $\Delta$, $\Delta \in \{0, 2, \ldots, 40\}$, under the discretized normal distribution with $\sigma_a = \frac{d}{4}$, $a \in A$. The distribution of $T_a$, $a \in A$, is uniform; $\widetilde{T}_{min} = 10$ and $w = h = 3$.}
		\label{fig: normal distribution, sigma = d/4, T max}
		\end{center}
	\end{figure}
	
The numerical results indicate that both methods with the truncated data set substantially outperform the baseline DRO approach across all values of $\widetilde{T}_{max}$ even though $\widetilde{T}_{min}$ is rather small. This observation appears counterintuitive at first sight, but it can be rationalized by analyzing the in-sample performance of both methods. 

In this regard, we consider a particular instance of the shortest path problem where the truncated DRO$_2$ approach outperforms the baseline DRO approach and $\widetilde{T}_{max} = 30$. We depict the nominal expected costs (sorted in an increasing order) and associated worst-case expected costs for the outlined solution approaches; see Figure \ref{fig: normal distribution, sigma = d/4, costs}. 

\begin{figure}[htp]
	\begin{center}
		\includegraphics[width = 0.5\textwidth]{graphs_clear/exp8/graph_costs.jpg}
		\caption{Nominal expected costs (in an increasing order) and associated worst-case expected costs (in one particular test instance) under the discretized normal distribution with $\sigma = \frac{d}{4}$. The distribution of $T_a$, $a \in A$, is uniform; $\widetilde{T}_{min} = 10$, $\widetilde{T}_{max} = 30$ and $w = h = 3$.}
		\label{fig: normal distribution, sigma = d/4, costs}
	\end{center}
\end{figure}

Note that the baseline DRO method demonstrates a better in-sample performance, which is rather intuitive due to the truncation of the data set.
At the same time, the worst-case expected costs obtained from the truncated data set in some sense \textit{repeat the form} of the nominal expected costs, which results in a better out-of-sample performance. More specifically, in this particular test instance the number of ``pairwise errors'' (i.e., situations where the worst-case expected costs for a pair of actions reverse the order of nominal costs) is given by $44$ and $29$ for the DRO and DRO$_2$ approaches, respectively. These implications provide some practical insights behind a rather poor out-of-sample performance of the baseline DRO approach. Furthermore, we make the following remark.

\begin{remark} \upshape Prediction and prescription problems (\ref{prediction problem}) and (\ref{prescription problem}) do not minimize the out-of-sample performance explicitly. That is, we attempt to minimize the in-sample performance, but ensure that the out-of-sample performance is smaller than the in-sample performance with increasingly high probability as the sample size grows.  \vspace{-9mm}\flushright$\square$
\end{remark}

Another important point is that a truncation of the data is reasonable due to a very specific form of the relative entropy distance. In fact, one may verify that the described procedure in the case of Hoeffding bounds deteriorates the out-of-sample performance (these results are omitted in order to simplify the structure of computational study.) 

Finally, we note that the baseline DRO approach may outperform DRO$_1$ and DRO$_2$ methods, e.g., if the variance $\sigma$ is relatively large or the parameters $T_a$, $a \in A$, are governed by the Binomial$_2$ distribution; see Figures 8a and 8b, respectively. 

 The intuition behind the latter observation is discussed in Section \ref{subsec: dependence on the sample size}; recall Figure \ref{fig: binomial T reversed}.  Meanwhile, the former observation can be explained as follows. The number of random samples in the truncated data set is relatively small. Hence, a large variance may result in biased estimates of the nominal expected costs that, in turn, results in a poor out-of-sample performance.

\begin{figure}[htp]
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}
		\includegraphics[width = \linewidth]{graphs_clear/exp9a/graph_loss.jpg}	
		\caption{Normal distribution with uniform $T_a$, $a \in A$, and $\sigma=\frac{3d}{4}$.}	
		\label{fig: normal, sigma = 3d/4}
	\end{subfigure}
	\begin{subfigure}[b]{0.5\textwidth}
		\centering
		\captionsetup{justification=centering}
		\includegraphics[width = \linewidth]{graphs_clear/exp9b/graph_loss.jpg}
		\caption{Normal distribution with Binomial$_2$ $T_a$, $a \in A$, and $\sigma=\frac{d}{4}$.}
		\label{fig: normal, reversed T}
	\end{subfigure}
	\caption{Average relative loss (\ref{eq: nominal relative loss}) and MADs as a function of $\Delta$, $\Delta \in \{0, 2, \ldots, 40\}$, under the discretized normal distribution with (a) uniform $T_a$, $a \in A$, and $\sigma=\frac{3d}{4}$ (b) Binomial$_2$ $T_a$, $a \in A$, and $\sigma=\frac{d}{4}$; $\widetilde{T}_{min} = 10$ and $w = h = 3$.}
	\label{fig: normal distribution, sigma = 3d/4 and reversed T}
\end{figure}


Summarizing our numerical results, advantages of the baseline DRO approach induced by (\ref{distributionally robust prescriptor}) can be described as follows:

\begin{itemize}
	\item The proposed approach outperforms Hoeffding bounds induced by (\ref{hoeffding prescriptor}) whenever relatively small expected costs can be observed sufficiently often or the gap between the sample sizes for individual actions is comparatively small;
	\item The out-of-sample performance of the baseline DRO approach improves with the decrease of variance of the data-generating distribution;
	\item If the variance is sufficiently small, then the truncation of the data set may substantially improve the out-of-sample performance of the baseline DRO approach. Surprisingly, this effect can be achieved ``for free'', i.e., without affecting complexity of the baseline model.
\end{itemize}
\section{Conclusion} \label{sec: conclusion}
In this paper we address linear combinatorial optimization problems, where the cost vector is governed by some unknown probability distribution $\mathbb{Q}^*$. The decision-maker attempts to minimize its expected loss under $\mathbb{Q}^*$ by leveraging a finite set of independent samples from this distribution. In contrast to the related study in \cite{VanParys2020}, the data set is not complete, i.e., the number of random observations for particular components of the cost vector may vary. Since the expected value under $\mathbb{Q}^*$ cannot be optimized directly we focus on prediction rules that approximate the expected loss for a fixed decision and prescription rules that seek optimal decisions for a given prediction rule.

The overall goal is to find the least conservative prediction and prescription rules, which also hedge against underestimated losses whenever the sample size tends to infinity. We demonstrate that the associated prediction and prescription problems admit a weakly optimal solution that can be obtained by solving a specified distributionally robust optimization problem (DRO). More precisely, the decision-maker may optimize the worst-case expected loss across all probability distributions with given component-wise relative entropy distances from the empirical marginal distributions. Importantly, this problem can be solved rather effectively whenever there exist effective algorithms for the deterministic version of a combinatorial optimization problem. 

In addition to asymptotic performance guarantees, we provide finite sample guarantees based on very recent results in the large deviation theory \cite{Agrawal2020, Mardia2020}. We conduct numerical experiments, which certify that under some parameter settings the proposed approach can be effective against several benchmark approaches. More specifically, our approach provides better out-of-sample tests whenever the actions with low expected costs can be observed sufficiently often or the maximal gap between the sample sizes for individual actions is comparatively small. Furthermore, for data-generating distributions with a sufficiently small variance a naive truncation of the data set may substantially improve the out-of-sample performance of the baseline DRO model without affecting its complexity.       

Admittedly, our theoretical results exploit the structure of the expected loss criterion. Therefore, it would be interesting to consider some other risk measures, especially those that do not account for any covariance  information among the cost coefficients (but may account for some variance information.) In fact, we shrink the probability space to the set of univariate distributions mostly due incompleteness of the data, i.e., standard estimates of joint probabilities are not applicable directly to partially observable data. 

In addition, a rather natural step forward is to extend the current results to the case of a continuous support. We refer the reader to Section 5 in \cite{VanParys2020} for the related discussion in the case of complete data. Finally, it is still an open question whether the proposed distributionally robust optimization approach provides strongly optimal prediction (prescription) rules or not.  

\textbf{Acknowledgments.} The article was prepared within the framework of the Basic Research Program at the National Research University Higher School of Economics (Sections 1-2). The research of Sergey S. Ketkov is funded by RFBR according to the research project №20-37-90060 (Section 3-5). 

\textbf{Conflict of interest:} The authors declare that they have no conflict of interest.

\bibliographystyle{ieeetr}
