Bounds for Probabilities of Extreme Events Defined by Two Random Variables

1. Introduction

Modeling extreme events like the U.S. stock market crash in 1929 and catastrophe insurance losses is certainly a major task for risk managers. In such events, typically, complete information of the underlying distribution is not available. Instead, one has only partial information such as estimates of mean, variance, covariance, or range based on a relatively small sample. Moreover, the sample may contain no extreme observations. While it is almost impossible to obtain accurate tail risk measures based on incomplete information, one can use limited information to compute bounds on the tail risk measures. This paper presents a method to obtain bounds on tail probabilities using only moment information for two kinds of extreme events based on two random variables. These problems are called semiparametric bound problems or generalized Tchebyshev bound problems.

Scarf (1958) applies these ideas in inventory management and Lo (1987) applies them in mathematical finance. Other applications in finance focus on option pricing in the well-known Black and Scholes (1973) setting (Merton 1973; Levy 1985; Ritchken 1985; Boyle and Lin 1997; Bruckner 2007; Schepper and Heijnen 2007) and other asset pricing and portfolio problems (Gallant, Hansen, and Tauchen 1990; Hansen and Jagannathan 1991; Ferson and Siegel 2001, 2003).

The purpose of this paper is to apply the semiparametric bounds approach to estimate the tail probability of joint events which, in many cases, cannot be reliably estimated by traditional statistical methods (e.g., the parametric approach). In particular, our approach is useful in the situations when it is very difficult or inappropriate to make distributional assumptions about random variables of interest, among others, due to scarcity and/or very high volatility of available data. Our approach also aims at tackling the problem of estimating the likelihood of extreme (tail) events for which we have very few observations of outliers. Traditional methods do not work for such tasks because these approaches typically produce a good fit in those regions in which most of the data reside but at the expense of a good fit in the tails (Hsieh 2004).

To address this problem, instead of assuming full knowledge of the distributions of the random variables of interest, we show how to numerically compute upper and lower bounds on the probabilities Pr(X₁ ≤ t₁ and X₂ ≤ t₂) and Pr(w₁X₁ + w₂X₂ ≤ a) for some appropriate values of t₁, t₂, w₁, w₂, a ∈ ℝ, when only second order moment information (means, variances, and covariance) and the support of random variables X₁ and X₂ are known. Our approach explicitly considers correlations between variables when estimating the bounds. Incorporating variable correlations is important because many models (e.g., models of risk-based capital and enterprise risk management) often involve several random variables, most of which are correlated. For example, let X₁ and X₂ stand for a random discount factor and a random future insurance payment. If the insurance payment X₂ is subject to economic inflation, it will be correlated with the interest rate which determines the discount factor X₁. As another example, the variables X₁ and X₂ can be the returns of two stocks, both of which respond to security market forces.

Following the work of Smith (1990), Cox (1991), Brochett et al. (1995), Zuluaga (2004), Popescu (2005), and Bertsimas and Popescu (2005), we obtain a range of possible values for each of our tail risk measures, corresponding to every distribution that has given moments on a given support. This range can be considered as a 100% confidence interval for the tail risk measure. Generally, semiparametric bounds are robust bounds that any reasonable model must satisfy. It is worth pointing out that the bounds provide “bestcase” estimates and “worst-case” estimates of the probabilities of extreme events. They would be useful for very risk-loving and very risk-averse investors. Moreover, in situations when distributional assumptions can be made, they provide a mechanism for checking the consistency of such assumptions, as well as an initial estimate for cumulative probabilities regardless of any model specifications.

The remainder of the paper is organized as follows. In Section 2, we formally state the semiparametric bound problems considered here and explain the methodology for solving them. Sections 3 and 4 show how the desired semiparametric bounds can be numerically computed with readily available optimization solvers. We present relevant numerical experiments to illustrate the application of our results. In Section 5, we discuss the possible extension of our methodology to obtain bounds when only confidence intervals on moments are given. Section 6 concludes the paper.

2. Preliminaries and notation

For a function φ(x₁,x₂) of two random variables with joint cumulative distribution function F(X₁, X₂), its expected value is

\[ \mathbb{E}_{F}\left[\phi\left(X_{1}, X_{2}\right)\right]=\int_{\mathcal{D}} \phi\left(X_{1}, X_{2}\right) d F\left(x_{1}, x_{2}\right), \]

where the set 𝒟 ⊆ ℝ² is the support of random variables X₁ and X₂ and ∫_𝒟 dF(x₁,x₂) = 1.

The semiparametric upper bound of \(\mathbb{E}_F\left[\phi\left(X_1\right.\right.\), \(\left.\left.X_2\right)\right]\) given up to second order moment information can be expressed as follows:

\[ \begin{array}{l} \begin{array}{ll} \bar{p}=\max & \mathbb{E}_{F}\left[\phi\left(X_{1}, X_{2}\right)\right] \\ \text { subject to } & \mathbb{E}_{F}\left(X_{i}\right)=\mu_{i}, \quad i=1,2, \\ & \mathbb{E}_{F}\left(X_{i}^{2}\right)=\mu_{i}^{(2)}, \quad i=1,2, \\ & \mathbb{E}_{F}\left(X_{1} X_{2}\right)=\mu_{12}, \end{array}\\ F\left(x_{1}, x_{2}\right) \text { a probability distribution on } \mathcal{D} \text {, } \end{array} \tag{2.1} \]

where \(\mu_i, \mu_i^{(2)}\) and \(\mu_{12}\) are the given first and second order non-central moments of \(X_i(i=1,2), \mathcal{D}\) is the given support of the distribution, and \(\bar{p}\) denotes the upper bound value. In order to simplify our presentation, let’s assume for the moment that the “point estimates” of the moments of the interested random variables are known. Later, in Section 5, we will show how our results can be adapted in a straightforward fashion to take into account the situation in which confidence intervals rather than point estimates of the moments are known.

The corresponding semiparametric lower bound problem is analogous, except that the objective function is

\[ \underline{p}=\min \mathbb{E}_{F}\left[\phi\left(X_{1}, X_{2}\right)\right], \tag{2.2} \]

with the same constraints as (2.1).

Notice that from the definitions of \(\bar{p}\) and \(\underline{p}\) in problems (2.1) and (2.2), the interval \([\underline{p}, \bar{p}]\) is a sharp (or tight) \(100 \%\) confidence interval on the expected value of \(\phi\left(X_1, X_2\right)\) for all joint distributions of \(X_1\) and \(X_2\) with the given moments and support. It follows that for any \(\bar{p}^{\prime} \geq \bar{p}\) and \(\underline{p}^{\prime} \leq \underline{p}\), the interval \(\left[\underline{p}^{\prime}, \bar{p}^{\prime}\right]\) is also a \(100 \%\) confidence interval, although not necessarily sharp. Our aim is to numerically compute useful \(100 \%\) confidence intervals for relevant choices of the function \(\phi\left(X_1, X_2\right)\), balancing computational effort and tightness of the confidence interval, using recent advances in optimization.

In particular, given \(t_1, t_2 \in \mathbb{R}^{+}\)and non-negative random variables \(X_1\) and \(X_2\), we compute \(100 \%\) confidence intervals on the probability of the extreme events \(X_1 \leq t_1\) and \(X_2 \leq t_2\), by setting \(\phi\left(X_1, X_2\right)=\mathbb{I}_{\left\{X_1 \leq t_1\right.}\) and \(\left.X_2 \leq t_2\right\}\) and \(\mathcal{D}=\mathbb{R}^{+2}\), where \(\mathbb{I}_{\mathcal{S}}\) is the indicator function of the set \(\mathcal{S}\). Similarly, given \(w_1, w_2, a \in \mathbb{R}\), we compute \(100 \%\) confidence intervals on the probability \(\operatorname{Pr}\left(w_1 X_1+\right. w_2 X_2 \leq a\) ), by setting \(\phi\left(X_1, X_2\right)=\mathbb{I}_{\left\{w_1 X_1+w_2 X_2 \leq a\right\}}\) and \(\mathcal{D}=\mathbb{R}^2\). In the second case, we strengthen the bounds in problems (2.1) and (2.2) by adding an additional moment constraint, \(\mathbb{E}_F\left[\left(X_1-X_2\right)^{+}\right] =\gamma\) where \(x^{+}=\max \{x, 0\}\). That is, we strengthen the bounds by only considering distributions of \(X_1, X_2\) that can replicate the expected payoff \(\gamma\) of an exchange option on \(X_1\) and \(X_2\). This illustrates how additional information can improve the semiparametric bounds. More details are shown in Sections 3 and 4.

The following is the dual of the upper bound problem (2.1) (see, e.g., Karlin and Studden 1966; Bertsimas and Popescu 2002; and Zuluaga and Peña 2005):

\[ \begin{aligned} \bar{d}=\min & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \\ & \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right) \geq \phi\left(x_{1}, x_{2}\right) \\ \text { for all } & \left(x_{1}, x_{2}\right) \in \mathcal{D}. \end{aligned} \tag{2.3} \]

The dual of the lower bound problem (2.2) is

\[ \begin{aligned} \underline{d}=\max & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \\ & \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right) \leq \phi\left(x_{1}, x_{2}\right) \\ \text { for all } & \left(x_{1}, x_{2}\right) \in \mathcal{D}, \end{aligned} \tag{2.4} \]

where the quadratic polynomial p(x₁,x₂) is defined as

\[ \begin{aligned} p\left(x_{1}, x_{2}\right)= & y_{00}+y_{10} x_{1}+y_{01} x_{2}+y_{20} x_{1}^{2} \\ & +y_{02} x_{2}^{2}+y_{11} x_{1} x_{2} . \end{aligned} \]

It is not difficult to see that weak duality holds between (2.1) and (2.3), or between (2.2) and (2.4); that is, p ≤ d (or p ≥ d) (Bertsimas and Popescu 2005, Theorem 2.1, p. 785). Furthermore, strong duality holds, i.e., p = d (or p = d), if the following conditions are satisfied (Zuluaga and Peña 2005, Proposition 4.1(ii)):

If problem (2.1) is feasible and there exist y₀₀, y₀₁, y₁₀, y₂₀, y₀₂, y₁₁such that

\[ p\left(x_{1}, x_{2}\right)>\phi\left(x_{1}, x_{2}\right), \quad \text { for all } \quad\left(x_{1}, x_{2}\right) \in \mathcal{D}, \]

then p = d. Similarly, if problem (2.2) is feasible and there exist y₀₀, y₀₁, y₁₀, y₂₀, y₀₂, y₁₁ such that

\[ p\left(x_{1}, x_{2}\right)<\phi\left(x_{1}, x_{2}\right), \quad \text { for all } \quad\left(x_{1}, x_{2}\right) \in \mathcal{D}, \]

then p = d.

Notice that for the two problems to be solved in Sections 3 and 4, \(\phi\left(x_1, x_2\right)\) is an indicator function bounded on \([0,1]\). Therefore, the inequality \(p\left(x_1, x_2\right)>\phi\left(x_1, x_2\right)\) for problem (2.1) holds if we set \(y_{00}>1\) and \(y_{i j}=0\) for all \((i, j) \neq(0,0)\). Similarly, setting \(y_{00}<0\) and \(y_{i j}=0\) for all \((i, j) \neq (0,0)\), the inequality \(p\left(x_1, x_2\right)<\phi\left(x_1, x_2\right)\) holds for the lower bound problem (2.2). Thus, as long as problem (2.1) or problem (2.2) is feasible when \(\phi\left(x_1, x_2\right)\) is an indicator function, strong duality \(\bar{p}=d\) (or \(\underline{p}=\underline{d}\) ) holds. Therefore one can solve (2.3) (or (2.4)) to obtain the desired semiparametric bounds. Before explaining how to solve (2.3) and (2.4), we introduce the following wellknown definition and theorems relevant to the discussion to follow.

Definition 1 (SOS polynomials).

A polynomial

\[ p(x)=p\left(x_{1}, \ldots, x_{n}\right)=\sum_{i_{1}, \ldots, i_{n} \in \mathbb{N}} a_{\left(i_{1}, \ldots, i_{n}\right)} x_{1}^{i_{1}} \cdots x_{n}^{i_{n}} \]

is said to be a sum of squares (SOS) if

\[ p(x)=\sum_{i}\left[q_{i}(x)\right]^{2} \]

for some polynomials \(q_i(x)=q_i\left(x_1, \ldots, x_n\right)\).

Theorem 1 (Diananda 1962). Let \(p\left(x_1, \ldots, x_n\right)\) be a quadratic polynomial. If \(n \leq 3\), then \(p\left(x_1, \ldots, x_n\right) \geq 0\), for all \(x_1, \ldots, x_n \geq 0\) if and only if \(p\left(x_1^2, \ldots, x_n^2\right)\) is an SOS polynomial.

Theorem 1 states that to check if

\[ \begin{aligned} p\left(x_{1}, x_{2}\right)= & y_{00}+y_{10} x_{1}+y_{01} x_{2}+y_{20} x_{1}^{2} \\ & +y_{02} x_{2}^{2}+y_{11} x_{1} x_{2} \end{aligned} \]

is positive for all x₁,x₂ ≥ 0, one can check whether

\[ \begin{aligned} p\left(x_{1}^{2}, x_{2}^{2}\right)= & y_{00}+y_{10} x_{1}^{2}+y_{01} x_{2}^{2}+y_{20} x_{1}^{4} \\ & +y_{02} x_{2}^{4}+y_{11} x_{1}^{2} x_{2}^{2} \end{aligned} \]

is an SOS. Here we present Diananda’s Theorem in a form (shown as Theorem 1) that will be suitable for our purposes, instead of presenting it in its original form. Parrilo (2000) and Zuluaga (2004) discuss the equivalence of the original version of Diananda’s Theorem and Theorem 1.

Loosely speaking, in order to solve (2.3), or (2.4), we break the constraint

\[ \begin{array}{c} p\left(x_{1}, x_{2}\right) \geq(\text { or } \leq) \phi\left(x_{1}, x_{2}\right) \\ \text { for all }\left(x_{1}, x_{2}\right) \in \mathcal{D} \end{array} \]

into a number of constraints of the form

\[ \begin{array}{l} p_{i}\left(x_{1}, x_{2}\right) \geq 0, \quad \text { for all } \quad\left(x_{1}, x_{2}\right) \in \mathbb{R}^{+2},\\ i=1, \ldots, m, \end{array} \tag{2.5} \]

where \(p_i, i=1, \ldots, m\) are suitable quadratic polynomials whose coefficients are linear functions of the coefficients of \(p\left(x_1, x_2\right)\). Theorem 1 implies that (2.5) is equivalent to

\[ p_{i}\left(x_{1}^{2}, x_{2}^{2}\right) \text { is an SOS polynomial, } \quad i=1, \ldots, m . \]

As we will show in detail in Sections 3 and 4, this allows us to reformulate problems (2.3) and (2.4) as SOS programs; that is, as an optimization problem, the variables are coefficients of polynomials, the objective is a linear combination of the polynomial coefficients, and the constraints are given by the polynomials being SOS. A detailed discussion about SOS programming is beyond the scope of this article, but the key fact is that these SOS programs can be readily solved by recently developed SOS programming solvers such as SOSTOOLS (Prajna, Papachristodoulou, and Parrilo 2002), GloptiPoly (Henrion and Laserre 2003), or YALMIP (Löfberg 2004). These SOS programming solvers enable us to find the desired bounds on problems (2.1) and (2.2). This approach has been widely used to solve semiparametric bound problems in other areas (see, e.g., Bertsimas and Popescu 2002; Boyle and Lin 1997; and Laserre 2002).

Parrilo (2000) and Todd (2001) show that any SOS program can be reformulated as a semi-definite program (SDP). Specifically, SOS programming solvers work by reformulating an SOS program as an SDP, and then applying SDP solvers such as SeDuMi (Sturm 1999). However, the SDP formulations of SOS programs can be fairly involved. To make it easy to present and reproduce our results, throughout the article we use SOS programming formulations instead of directly reformulating problems (2.1) and (2.2) as SDPs.

3. Extreme probability bounds

In this section, we consider the problem of finding upper and lower bounds on the probability Pr(X₁ ≤ t₁ and X₂ ≤ t₂) of two non-negative random variables X₁ and X₂, attaining values lower than or equal to t₁,t₂ ∈ ℝ⁺ respectively, without making any assumption on the distribution of X₁ and X₂, other than the knowledge of the first and second order moments of their joint distribution (means, variances, and covariance).

3.1. SOS programming formulations

The upper semiparametric bounds for this problem come from problem (2.1) with \(\phi\left(X_1, X_2\right) =\mathbb{I}_{\left\{X_1 \leq t_1\right.}\) and \(\left.X_2 \leq t_2\right\}\) and \(\mathcal{D}=\mathbb{R}^{+2}\) (Section 2):

\[ \begin{aligned} \bar{p}_{\text {Extreme }}=\max & \mathbb{E}_{F}\left[\mathbb{I}_{\left\{X_{1} \leq t_{1} \text { and } X_{2} \leq t_{2}\right\}}\right] \\ \text { subject to } & \mathbb{E}_{F}\left(X_{i}\right)=\mu_{i}, \quad i=1,2, \\ & \mathbb{E}_{F}\left(X_{i}^{2}\right)=\mu_{i}^{(2)}, \quad i=1,2, \\ & \mathbb{E}_{F}\left(X_{1} X_{2}\right)=\mu_{12}, \\ & F\left(x_{1}, x_{2}\right) \text { a probability } \\ & \text { distribution on } \mathbb{R}^{+2} . \end{aligned} \tag{3.1} \]

Similarly, the lower semiparametric bounds for this problem can be obtained by setting the objective function of problem (2.2) as follows:

\[ \underline{p}_{\text {Extreme }}=\min \mathbb{E}_{F}\left[\mathbb{I}_{\left\{X_{1} \leq t_{1} \text { and } X_{2} \leq t_{2}\right\}}\right], \tag{3.2} \]

with the same constraints as (3.1).

Before obtaining the SOS programming formulation of these problems, let us first examine their feasibility in terms of the moment information. Using Theorem 1 and convex duality (Rockafellar 1970), one can show that problems (3.1) and (3.2) are feasible (i.e., they have solutions), if the moment matrix Σ is a positive definite matrix (i.e., all eigenvalues are greater than zero) and all elements of Σ are greater than zero, where the moment matrix Σ is

\[ \Sigma=\left[\begin{array}{ccc} 1 & \mu_{1} & \mu_{2} \\ \mu_{1} & \mu_{1}^{(2)} & \mu_{12} \\ \mu_{2} & \mu_{12} & \mu_{2}^{(2)} \end{array}\right]. \]

Now we derive SOS programs to numerically approximate p_Extreme and p _Extreme with SOS programming solvers.

3.1.1. Upper bound

To derive an SOS program for problem (3.1), we begin by stating its dual explicitly:

\[ \begin{aligned} \bar{d}_{\text {Extreme }}=\min & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \\ & \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right) \geq \mathbb{I}_{\left\{x_{1} \leq t_{1} \text { and } x_{2} \leq t_{2}\right\}}, \\ \text { for all } & x_{1}, x_{2} \geq 0 . \end{aligned} \tag{3.3} \]

To formulate problem (3.3) as an SOS program, we proceed as follows. First notice that the constraint in (3.3) is equivalent to

\[\begin{align} p\left(x_1, x_2\right) \geq 1, \quad &\text{for all} \quad 0 \leq x_1 \leq t_1, 0 \leq x_2 \leq t_2 \\ p\left(x_1, x_2\right) \geq 0, \quad &\text{for all} \quad x_1, x_2 \geq 0. \end{align} \tag{3.4}\]

While the second constraint of (3.4) can be directly reformulated as an SOS constraint using Theorem 1, the first constraint is difficult to reformulate as an SOS constraint. That is, there is no linear transformation from 0 ≤ x₁ ≤ t₁, 0 ≤ x₂ ≤ t₂ to ℝ⁺² (that would allow us to use Theorem 1). Thus, we change the problem to obtain an SOS program that either exactly or approximately solves problem (3.4). Specifically, consider the following problem related to (3.4):

\[ \begin{aligned} \bar{d}_{\text {Extreme }}^{\prime}=\min & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \\ & \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right) \geq 1 \\ \text { for all } & x_{1} \leq t_{1}, x_{2} \leq t_{2} \\ & p\left(x_{1}, x_{2}\right) \geq 0 \\ \text { for all } & x_{1} \geq 0, x_{2} \geq 0. \end{aligned} \tag{3.5} \]

We relaxed the requirement that x₁ and x₂ are non-negative in the first constraint. Notice that the constraints in (3.5) are stricter than those in (3.4) since the first constraint of (3.5) includes more values of x₁ and x₂. Thus, d′_Extreme is a (not necessarily sharp) upper bound on d′_Extreme; that is, d′_Extreme ≥ d′_Extreme.

After we apply the substitution x₁ → t₁ − x₁,x₂ → t₂ − x₂ to the first constraint of (3.5), the constraints of (3.5) can be rewritten as

\[\begin{align} p\left(t_{1}-x_{1}, t_{2}-x_{2}\right)-1 \geq 0, \quad &\text{for all} \quad x_{1}, x_{2} \geq 0 \\p\left(x_{1}, x_{2}\right) \geq 0, \quad &\text{for all} \quad x_{1}, x_{2} \geq 0. \end{align}\tag{3.6} \]

To finish, we apply Theorem 1 to the constraints (3.6) and conclude that (3.5) is equivalent to the following SOS program:

\[ \begin{align} \bar{d}_{\text {Extreme }}^{\prime}=\min &\left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } \quad &p\left(t_{1}-x_{1}^{2}, t_{2}-x_{2}^{2}\right)-1 \\ &\qquad \text{is an SOS polynomial}\\ &p\left(x_{1}^{2}, x_{2}^{2}\right) \\ &\qquad \text{is an SOS polynomial}. \end{align}\tag{3.7} \]

The SOS program (3.7) can be readily solved with an SOS programming solver. Thus, if problem (3.1) is feasible, we can numerically obtain a (not necessarily sharp) semiparametric upper bound on the extreme probability, Pr(X₁ ≤ t₁, X₂ ≤ t₂) ≤ d′_Extreme, by solving problem (3.7) with an SOS solver.

3.1.2. Lower bound

The dual of the lower bound problem (3.2) can be expressed as

\[ \begin{aligned} \underline{d}_{\text {Extreme }}=\max & \left(y_{00}+y_{10} \mu_1+y_{01} \mu_2+y_{20} \mu_1^{(2)}\right. \\ & \left.+y_{02} \mu_2^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & p\left(x_1, x_2\right)<\pi_{\left\{x_1<I_1 \text { and } x_2<I_2\right\}} \\ \text { for all } & x_1, x_2>0 \end{aligned} \tag{3.8} \]

The constraint in problem (3.8) is equivalent to

\[ \begin{array}{lll} p\left(x_1, x_2\right) \leq 1, & \text { for all } & 0 \leq x_1 \leq t_1, 0 \leq x_2 \leq t_2 \\ p\left(x_1, x_2\right) \leq 0, & \text { for all } & x_1 \geq t_1, x_2 \geq 0 \\ p\left(x_1, x_2\right) \leq 0, & \text { for all } & x_1 \geq 0, x_2 \geq t_2 . \end{array}\tag{3.9} \]

Proceeding in the same way as for the upper bound problem, we now change the problem to obtain an SOS program that either exactly or approximately solves problem (3.8). Specifically, consider the following problem related to (3.8):

\[ \begin{aligned} \underline{d}_{\text {Extreme }}^{\prime}=\max & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right) \leq 1, \quad \text { for all } \quad x_{1} \leq t_{1}, x_{2} \leq t_{2} \\ & p\left(x_{1}, x_{2}\right) \leq 0, \quad \text { for all } \quad x_{1} \geq t_{1}, x_{2} \geq 0 \\ & p\left(x_{1}, x_{2}\right) \leq 0, \quad \text { for all } \quad x_{1} \geq 0, x_{2} \geq t_{2} . \end{aligned} \tag{3.10} \]

Notice that the constraints in (3.10) are stricter than those in (3.8). Thus, d′_Extreme is a (not necessarily sharp) lower bound on d_Extreme; that is, d′_Extreme ≤ d′_Extreme

Applying the substitutions x₁ → t₁ − x₁,x₂ → t₂ − x₂ to the first constraint of (3.10) and x1 → t₁ + x₁,x₂ → t₂ + x₂ to the second and third constraints respectively, it follows that problem (3.10) is equivalent to the following SOS program when Theorem 1 is applied:

\[ \begin{array}{rll} \underline{d}_{\text {Extreme }}^{\prime}=\max & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}\right) \\ \text { subject to } & 1-p\left(t_{1}-x_{1}^{2}, t_{2}-x_{2}^{2}\right) \qquad \qquad \text { is an } \operatorname{SOS} \text { polynomial } \\ & -p\left(t_{1}+x_{1}^{2}, x_{2}^{2}\right) \qquad \qquad \qquad \quad \text { is an } \operatorname{SOS} \text { polynomial } \\ & -p\left(x_{1}^{2}, t_{2}+x_{2}^{2}\right) \qquad \qquad \qquad \quad \text { is an } \operatorname{SOS} \text { polynomial. } \end{array} \tag{3.11} \]

Following the same route, we can also derive the upper and lower bounds on the joint survival probability Pr(X₁ ≥ t₁ and X₂ ≥ t₂) of two nonnegative random variables X₁ and X₂. The details are shown in Appendix A.

3.2. Example of extreme probability bounds

We select from the NAIC database a major property/casualty insurance company, which we call insurer A. Suppose the insurer faces the problem of managing its risk of unexpectedly high claims and simultaneously unanticipated poor asset returns. This leads insurer A to calculate the bounds on Pr(R ≤ t₁, M ≤ t₂) given moment information, where ℝ is the company’s return on its invested assets and M is the margin on its insurance business.

The return \(R_i\) of asset \(i\) in insurer A’s portfolio is equal to \(P_{i, t} / P_{i, t-1}-1\) where \(P_{i, t-1}\) and \(P_{i, t}\) denote the prices of asset \(i\) at the beginning and the end of the period. Insurer A’s asset portfolio return \(R\) is the weighted average return of six asset classes: stocks, government bonds, corporate bonds, real estates, mortgages, and short-term investments; that is

\[ \begin{aligned} R & =\sum_{i=1}^{6} w_{i} R_{i}=\sum_{i=1}^{6} w_{i}\left(\frac{P_{i, t}}{P_{i, t-1}}-1\right) \\ & =\sum_{i=1}^{6} w_{i} \frac{P_{i, t}}{P_{i, t-1}}-1=X_{1}-1, \end{aligned} \]

where \(w_i\) is the weight of asset class \(i(i=1,2\), \(\ldots, 6)\) in the portfolio and \(X_1=\sum_{i=1}^6 w_i P_{i, t} / P_{i, t-1}\). The following inequalities are equivalent:

\[ R \leq t_{1} \Longleftrightarrow X_{1} \leq t_{1}+1. \tag{3.12} \]

We make this shift from asset returns to price ratios to apply our SOS results because we need non-negative random variables.

The margin on the insurance business is defined as

\[ M=1-\mathrm{CR}=1-\mathrm{LR}-\mathrm{ER}, \]

where CR is the combined ratio, LR is the loss ratio, and ER is the expense ratio.^[1] So M is the profit from the underwriting business.

In order to reformulate the condition M ≤ t₂ so that the condition fits our SOS results, we replace M ≤ t₂ with X₂ ≤ t₂ + 1 where X₂ = M + 1. Using this with (3.12) we get the following:

\[ \operatorname{Pr}\left(R \leq t_{1}, M \leq t_{2}\right)=\operatorname{Pr}\left(X_{1} \leq t_{1}+1, X_{2} \leq t_{2}+1\right) . \]

The weights w_i of different asset categories were calculated from the quarterly data of the National Association of Insurance Commissioners (NAIC). We used the quarterly returns of the Standard & Poor’s 500 (S&P500), the Lehman Brothers intermediate term total return, the domestic high-yield corporate bond total return, the National Association of Real Estate Investment Trusts (NAREIT) total return, the Merrill Lynch mortgage backed securities total return, and the U.S. 30-Day T-Bill as proxies for insurer A’s stock returns, government bond returns, corporate bond returns, real estate returns, mortgage returns and short-term investment returns, respectively. Based on insurer A’s quarterly losses, expenses, and premiums, we calculate the moments of X₁ and X₂ as follows:

\[ \begin{aligned} \mathrm{E}\left(X_{1}\right) & =1.0442, & & \mathrm{E}\left(X_{1}^{2}\right)=1.0967 \\ \mathrm{E}\left(X_{2}\right) & =1.1555, & & \mathrm{E}\left(X_{2}^{2}\right)=1.3715 \\ \mathrm{E}\left(X_{1} X_{2}\right) & =1.2086, & & \operatorname{Cov}\left(X_{1}, X_{2}\right)=0.0021 \\ \operatorname{Var}\left(X_{1}\right) & =0.0063, & & \operatorname{Var}\left(X_{2}\right)=0.0364 \\ \rho & =0.1387 . & & \end{aligned} \]

Insurer A’s average margin on its insurance business (E(M) = 0.1555) is higher than its average asset return (E(R) = 0.0442), while the margin is more volatile (Var(M) > Var(R)). Moreover, the asset return and insurance margin are somewhat positively correlated (0.1387). This implies that occasionally insurer A’s insurance business and investment performances move in the same direction.

Next we compute bounds on the tail probability Pr(R ≤ t₁, M ≤ t₂) using SOS programming. Then we compare it to the bivariate normal cumulative joint probability with the same moments. The upper left plot in Figure 1 shows the upper bounds of the joint probability Pr(R ≤ t₁, M ≤ t₂) for different values of t₁ and t₂, and the upper right one is the corresponding bivariate normal cumulative joint probabilities. Since we are looking at low values of t₁ and t₂ corresponding to joint extreme events, it is not surprising that our calculated lower bound is zero over this range of their values. The ratios of the upper bounds to the bivariate normal cumulative joint probabilities are shown in the lower graphs.

Figure 1.The upper left plot shows the upper bound of the joint probability \(\operatorname{Pr}\left(R \leq t_1, M \leq t_2\right)\) where \(R\) is the invested asset return and \(M\) is the insurance business margin of insurer \(\mathbf{A}\). The upper right one is the bivariate normal cumulative probabilities with the same moments. The ratio of the upper bound to the bivariate normal cumulative joint probabilities is shown in the lower left graph. The lower right one is a zoom-in plot of the ratio, illustrating a special case of \(\operatorname{Pr}(R \leq 0, M \leq 0)\). The vertical axis of the upper graphs is the probability. It is the ratio for the lower graphs. The two axes at the bottom in all graphs represent the value of return \(R\) and the value of insurance margin \(M\).

The ratio is large when t₁ and t₂ are low. For example, consider the event that insurer A has negative investment earnings and simultaneously it has an aggregate loss on its insurance business. In addition to the case with zero investment and insurance returns, this is stated as ℝ ≤ 0, M ≤ 0. From the lower right graph of Figure 1, we see that for t₁ = 0 and t₂ = 0, the upper bound is about 7.2 times higher than the cumulative joint normal probability. This means that the actual joint distribution may have a much fatter tail than the joint normal distribution. In other words, an extreme event may be more likely to occur than the normal distribution suggests.

4. Value-at-risk probability bounds

Here we find upper and lower bounds on the probability that a portfolio w₁X₁ + w₂X₂ (w₁, w₂ ∈ ℝ⁺) attains values lower than or equal to a ∈ ℝ, given up to the second order moment information (means, variances, and covariance) on the random variables X₁, X₂ (X₁, X₂ ∈ ℝ).

4.1. SOS programming formulations

Finding the sharp upper and lower semiparametric bounds for this problem can be formulated by setting \(\phi\left(X_1, X_2\right)=\mathbb{I}_{\left\{w_1 X_1+w_2 X_2 \leq a\right\}}\), and \(\mathcal{D}=\mathbb{R}^2\). To obtain tighter bounds (see numerical results in Section 4.2), we include the information of the expected payoff \(\gamma\) of an exchange option on the assets; that is, we add the moment constraint \(\mathbb{E}_F\left[\left(X_1-X_2\right)^{+}\right]=\gamma\left(\right.\) where \(x^{+}= \max \{0, x\}\) ) to illustrate how to incorporate additional information. This is the resulting semiparametric upper bound problem:

\[ \begin{aligned} \bar{p}_{\mathrm{VaR}}=\max & \mathbb{E}_{F}\left[\mathbb{I}_{\left\{w_{1} X_{1}+w_{2} X_{2} \leq a\right\}}\right] \\ \text { subject to } & \mathbb{E}_{F}\left(X_{i}\right)=\mu_{i}, \quad i=1,2, \\ & \mathbb{E}_{F}\left(X_{i}^{2}\right)=\mu_{i}^{(2)}, \quad i=1,2, \\ & \mathbb{E}_{F}\left(X_{1} X_{2}\right)=\mu_{12}, \\ & \mathbb{E}_{F}\left[\left(X_{1}-X_{2}\right)^{+}\right]=\gamma, \\ & F\left(x_{1}, x_{2}\right) \text { a probability distribution on } \mathbb{R}^{2} . \end{aligned} \tag{4.1} \]

The corresponding lower bound problem has the same constraints as (4.1) and its objective function is

\[ \underline{p}_{\mathrm{VaR}}=\min \mathbb{E}_{F}\left[\mathbb{I}_{\left\{w_{1} X_{1}+w_{2} X_{2} \leq a\right\}}\right] . \tag{4.2} \]

Problems (4.1) and (4.2) have solutions if and only if the moment matrix Σ is a positive semi-definite matrix and p_Exch < γ < p_Exch, where p_Exch and p_Exch are the upper and lower bounds of problems (2.1) and (2.2) with φ(x₁,x₂) = (x₁ − x₂)⁺, which can be readily computed using SOS techniques (Zuluaga and Peña 2005).

The dual of problem (4.1) is

\[ \begin{aligned} \bar{d}_{\mathrm{VaR}}=\min & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \\ & \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}+y_{0} \gamma\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right)+y_{0}\left(x_{1}-x_{2}\right)^{+} \\ & \geq \mathbb{I}_{\left\{w_{1} x_{1}+w_{2} x_{2} \leq a\right\}}, \\ \text { for all } & x_{1}, x_{2} \in \mathbb{R} . \end{aligned} \tag{4.3} \]

Similarly, the dual of problem (4.2) is:

\[ \begin{aligned} \underline{d}_{\mathrm{VaR}}=\max & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \mu_{1}^{(2)}\right. \\ & \left.+y_{02} \mu_{2}^{(2)}+y_{11} \mu_{12}+y_{0} \gamma\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right)+y_{0}\left(x_{1}-x_{2}\right)^{+} \\ & \leq \mathbb{I}_{\left\{w_{1} x_{1}+w_{2} x_{2} \leq a\right\}} \\ \text { for all } & x_{1}, x_{2} \in \mathbb{R}. \end{aligned} \tag{4.4} \]

A straightforward generalization of Proposition 1, and the discussion after Proposition 1, shows that if problems (4.1) and (4.2) are feasible, then p_VaR = d_VaR and p_VaR = d_VaR. Thus, if problems (4.1) and (4.2) are feasible, we can solve (4.3) and (4.4) to obtain the desired bounds. Finally, notice that setting y₀ = 0 in (4.3) and (4.4) is equivalent to solving the semiparametric bounds (4.1) and (4.2) without using information about the exchange option expected payoff.

4.1.1. Upper bound

The upper bound problem (4.3) is equivalent to

\[\small{ \begin{aligned} \bar{d}_{\mathrm{VaR}}=\min & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \sigma_{1}^{2}+y_{02} \sigma_{2}^{2}+y_{11} \sigma_{12}+y_{0} \gamma\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right)+y_{0}\left(x_{1}-x_{2}\right) \geq 1, \quad \text { for all } x_{1}, x_{2} \quad \text { with } w_{1} x_{1}+w_{2} x_{2} \leq a, x_{1} \geq x_{2} \\ & p\left(x_{1}, x_{2}\right) \geq 1, \qquad \qquad \qquad \quad \ \ \text { for all } x_{1}, x_{2} \quad \text { with } w_{1} x_{1}+w_{2} x_{2} \leq a, x_{1} \leq x_{2} \\ & p\left(x_{1}, x_{2}\right)+y_{0}\left(x_{1}-x_{2}\right) \geq 0, \quad \text { for all } x_{1}, x_{2} \quad \text { with } x_{1} \geq x_{2} \\ & p\left(x_{1}, x_{2}\right) \geq 0, \quad \qquad \qquad \qquad \ \ \text { for all } x_{1}, x_{2} \quad \text { with } x_{1} \leq x_{2} . \end{aligned} \tag{4.5}} \]

In order to use Theorem 1, we will use the following transformations:

Applying the upper left transformation in (4.6) to the first and third constraints of problem (4.5) and applying the upper right transformation in (4.6) to the second and fourth constraints of problem (4.5), the constraints in (4.5) are equivalent to

\[ \begin{array}{rlll} p\left(z_{1}+z_{2}, z_{2}\right)+y_{0} z_{1} \geq 1, & \text { for all } \quad z_{1}, z_{2} \quad \text {with} \quad w_1\left(z_1+z_2\right)+w_2 z_2 \leq a, z_1 \geq 0\\ p\left(z_{1}, z_{1}+z_{2}\right) & \text { for all } \quad z_{1}, z_{2} \quad \text{with} \quad w_1 z_1+w_2\left(z_1+z_2\right) \leq a, z_2 \geq 0 \\ p\left(z_{1}+z_{2}, z_{2}\right)+y_{0} z_{1} \geq 0, & \text { for all } \quad z_{1}, z_{2} \quad \text {with} \quad z_1 \geq 0\\ p\left(z_1, z_1+z_2\right) \geq 0 & \text { for all } \quad z_{1}, z_{2} \quad \text {with} \quad z_2 \geq 0. \end{array} \tag{4.7} \]

Now applying the lower left and right transformations in (4.6) to the first two constraints of (4.7) respectively, these two constraints are equivalent to

\[ \begin{array}{l} p\left(t_{1}+\frac{a-w_{1} t_{1}}{w_{1}+w_{2}}-t_{2}, \frac{a-w_{1} t_{1}}{w_{1}+w_{2}}-t_{2}\right)+y_{0} t_{1} \geq 1, \\ \quad \text { for all } t_{1} \geq 0, t_{2} \geq 0 \\ p\left(\frac{a-w_{2} t_{2}}{w_{1}+w_{2}}-t_{1}, \frac{a-w_{2} t_{2}}{w_{1}+w_{2}}-t_{1}+t_{2}\right) \geq 1, \\ \quad \text { for all } t_{1} \geq 0, t_{2} \geq 0 . \end{array} \tag{4.8} \]

Finally, the last two constraints in (4.7) are equivalent to

\[ \begin{array}{rll} p\left(z_{1}+z_{2}, z_{2}\right)+y_{0} z_{1} \geq 0, & \text { for all } & z_{1} \geq 0, z_{2} \geq 0 \\ p\left(z_{1}-z_{2},-z_{2}\right)+y_{0} z_{1} \geq 0, & \text { for all } & z_{1} \geq 0, z_{2} \geq 0 \\ p\left(z_{1}, z_{1}+z_{2}\right) \geq 0, & \text { for all } & z_{1} \geq 0, z_{2} \geq 0 \\ p\left(-z_{1},-z_{1}+z_{2}\right) \geq 0, & \text { for all } & z_{1} \geq 0, z_{2} \geq 0 . \end{array} \]

After applying Theorem 1, we obtain the SOS formulation for the upper bound of problem (4.3):

\[ \begin{array}{c} \bar{d}_{\mathrm{VaR}}=\min \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \sigma_{1}^{2}\right. \\ \left.+y_{02} \sigma_{2}^{2}+y_{11} \sigma_{12}+y_{0} \gamma\right) \end{array} \tag{4.9} \]

subject to the following being SOS polynomials:

\[ \begin{array}{l} p\left(t_{1}^{2}+\frac{a-w_{1} t_{1}^{2}}{w_{1}+w_{2}}-t_{2}^{2}, \frac{a-w_{1} t_{1}^{2}}{w_{1}+w_{2}}-t_{2}^{2}\right)+y_{0} t_{1}^{2}-1 \\ p\left(\frac{a-w_{2} t_{2}^{2}}{w_{1}+w_{2}}-t_{1}^{2}, \frac{a-w_{2} t_{2}^{2}}{w_{1}+w_{2}}-t_{1}^{2}+t_{2}^{2}\right)-1 \\ p\left(z_{1}^{2}+z_{2}^{2}, z_{2}^{2}\right)+y_{0} z_{1}^{2} \\ p\left(z_{1}^{2}-z_{2}^{2},-z_{2}^{2}\right)+y_{0} z_{1}^{2} \\ p\left(z_{1}^{2}, z_{1}^{2}+z_{2}^{2}\right) \\ p\left(-z_{1}^{2},-z_{1}^{2}+z_{2}^{2}\right) . \end{array} \]

4.1.2. Lower bound

The lower bound problem (4.4) is equivalent to

\[\small{ \begin{aligned} \underline{d}_{\mathrm{VaR}}=\max & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \sigma_{1}^{2}+y_{02} \sigma_{2}^{2}+y_{11} \sigma_{12}+y_{0} \gamma\right) \\ \text { subject to } & p\left(x_{1}, x_{2}\right)+y_{0}\left(x_{1}-x_{2}\right) \leq 1, \quad \text { for all } \quad x_{1}, x_{2} \quad \text { with } \quad x_{1} \geq x_{2} \\ & p\left(x_{1}, x_{2}\right) \leq 1, \quad \qquad \qquad \quad \ \ \text { for all } \quad x_{1}, x_{2} \quad \text { with } \quad x_{1} \leq x_{2} \\ & p\left(x_{1}, x_{2}\right)+y_{0}\left(x_{1}-x_{2}\right) \leq 0, \quad \text { for all } \quad x_{1}, x_{2} \quad \text { with } \quad w_{1} x_{1}+w_{2} x_{2} \leq a, x_{1} \geq x_{2} \\ & p\left(x_{1}, x_{2}\right) \leq 0, \qquad \qquad \quad \ \ \quad \text { for all } \quad x_{1}, x_{2} \quad \text { with } \quad w_{1} x_{1}+w_{2} x_{2} \leq a, x_{1} \leq x_{2} \end{aligned} \tag{4.10}} \]

Here we will use the following extra transformations:

Following steps analogous to those taken in Section 4.1.1 for problem (4.5), we obtain that problem (4.4) is equivalent to

\[ \begin{aligned} \underline{d}_{\mathrm{VaR}}=\max & \left(y_{00}+y_{10} \mu_{1}+y_{01} \mu_{2}+y_{20} \sigma_{1}^{2}\right. \\ & \left.+y_{02} \sigma_{2}^{2}+y_{11} \sigma_{12}+y_{0} \gamma\right) \end{aligned} \tag{4.12} \]

subject to the following being SOS polynomials:

\[ \begin{array}{l} 1-p\left(z_{1}^{2}+z_{2}^{2}, z_{2}^{2}\right)-y_{0} z_{1}^{2} \\ 1-p\left(z_{1}^{2}-z_{2}^{2},-z_{2}^{2}\right)-y_{0} z_{1}^{2} \\ 1-p\left(z_{1}^{2}, z_{1}^{2}+z_{2}^{2}\right) \\ 1-p\left(-z_{1}^{2},-z_{1}^{2}+z_{2}^{2}\right) \\ -p\left(t_{1}^{2}+\frac{a-w_{1} t_{1}^{2}}{w_{1}+w_{2}}+t_{2}^{2}, \frac{a-w_{1} t_{1}^{2}}{w_{1}+w_{2}}+t_{2}^{2}\right)-y_{0} t_{1}^{2} \\ -p\left(\frac{a-w_{2} t_{2}^{2}}{w_{1}+w_{2}}+t_{1}^{2}, \frac{a-w_{2} t_{2}^{2}}{w_{1}+w_{2}}+t_{1}^{2}+t_{2}^{2}\right) . \end{array} \]

4.2. Example of value-at-risk probability bounds

Given a specified tail probability β, the weights w₁ and w₂, as well as the moment information on X₁ and X₂, the value-at-risk (VaR) bound problem finds the upper and lower bounds on a where Pr(w₁X₁ + w₂X₂ ≤ a) = β. To solve this problem, we find bounds on Pr(w₁X₁ + w₂X₂ ≤ a) for different values of a and then solve the inverse problem for bounds on a given β.

We first calculate the semiparametric VaR probability bounds given only the mean, variance, and covariance of the two components of the portfolio. Then we add one more constraint that the expected value of the exchange option between \(X_1\) and \(X_2, \mathbb{E}_F\left[\left(X_1-X_2\right)^{+}\right]\), must equal \(\gamma\). The following example shows that the VaR bounds with exchange option information are tighter since we add more constraints to the optimization.

We analyze a portfolio investing in the S&P 500 Index and the Dow Jones U.S. Small-Cap Index. Suppose we invest 1/3 of our assets in the S&P 500 Index, 1/3 in the Dow Jones U.S. Small-Cap Index, and 1/3 in a risk-free fund paying a flat 0.01 percent per day. Thus, our portfolio daily return is (1/3)X₁ + (1/3)X₂ + (1/3)0.01.

The moments are based on the daily historical log-returns from February 24, 2000 to October 24, 2007. There are 1,923 observations in our sample. Let \(X_1\) and \(X_2\) be the log-return of the S&P 500 Index and Dow Jones U.S. SmallCap Index in percentage per day \(\left(X_i(t)=\right. 100 \log \left(S_i(t+1) / S_i(t)\right)\) for day \(t\) ). Their moments are as follows:

\[ \begin{aligned} \mathrm{E}\left(X_{1}\right) & =0.0059, & & \mathrm{E}\left(X_{1}^{2}\right)=1.2158 \\ \mathrm{E}\left(X_{2}\right) & =-0.2117, & & \mathrm{E}\left(X_{2}^{2}\right)=112.8609 \\ \mathrm{E}\left(X_{1} X_{2}\right) & =1.4161, & & \operatorname{Cov}\left(X_{1}, X_{2}\right)=1.41736 \\ \operatorname{Var}\left(X_{1}\right) & =1.2158, & & \operatorname{Var}\left(X_{2}\right)=112.8160 \\ \rho & =0.1210, & & \mathrm{E}\left(\left(X_{1}-X_{2}\right)^{+}\right)=0.4464 . \end{aligned} \tag{4.13} \]

Apparently, the Dow Jones U.S. Small-Cap Index is much more volatile than the S&P 500 (112.8160 percent vs. 1.2158 percent).

We now calculate the upper and lower bounds for the probability when the portfolio return falls below a given level a, i.e.,

\[ \operatorname{Pr}\left((1 / 3) X_{1}+(1 / 3) X_{2}+(1 / 3) 0.01 \leq a\right) . \]

The corresponding bounds are shown in Figure 2. The lines with –o-- represent the upper and lower bounds on the VaR probability, without using the exchange option information. These bounds are obtained by setting y₀ = 0 in Equations (4.9) and (4.12). The lines with \(-*-\) represent the upper and lower bounds on the VaR probability, using the exchange option information. Obviously the exchange option information tightens the VaR probability bounds significantly. These semiparametric upper and lower bounds apply to all possible joint probabilities, including the bivariate normal joint probability. The VaR probability corresponding to a normal distribution with the same first and second order moments is drawn with the broken line in the middle. Interestingly, the normal VaR probability lies outside the tighter bounds using the exchange option information. This means that the normal model does not satisfy the constraint \(\mathbb{E}_F\left[\left(X_1-X_2\right)^{+}\right]=\gamma\).

Figure 2.Comparison of VaR probability bounds with and without exchange option information.

Here we use the VaR probability bounds in Figure 2 to obtain the upper and lower bounds of the VaR itself, with and without the exchange option price constraint. Figure 2 gives us an idea of how likely the return of this portfolio will be lower than a in 1 day under different conditions. Consider a 5% VaR. We look at the horizontal line through the 0.05 level on the vertical axis—it intersects the –o-- curves at a values of −16 and 1. The best we can say is that

\[ -16 \%<\mathrm{VaR}_{0.05}<1 \% \]

per day. Then reading the curves with the exchange option information (\(-*-\)), we find that

\[ -6.2 \%<\mathrm{VaR}_{0.05}<0 \% \]

per day. Clearly the additional information greatly improves our knowledge of future possible outcomes.

5. Semiparametric bounds given confidence intervals on moments

Thus far in this paper we have assumed that the moments information is given in the form of point estimates of the moments of the random variables of interest. In practice, it is typical to have not only a point estimate but also a confidence interval estimate on a moment, which provides information about the accuracy of the point estimate. Therefore, an important question is how to adapt the results presented so far in order to handle the situation in which the moment information is given in the form of confidence intervals on the moments. To answer this question, consider the following general upper semiparametric bound problem:

\[ \begin{array}{cl} \max & \left\{\mathbb{E}_{F}\left[\phi\left(X_{1}, \ldots, X_{n}\right)\right]\right\} \\ \text { subject to } & \mathbb{E}_{F}\left[f_{j}\left(X_{1}, \ldots, X_{n}\right)\right]=\sigma_{j}, \\ & j=1, \ldots, m, \\ & F\left(x_{1}, \ldots, x_{n}\right) \text { a probability } \\ & \text { distribution on } \mathcal{D} \subseteq \mathbb{R}^{n}, \end{array} \tag{5.1} \]

where \(\sigma_j, j=1, \ldots, m\),, represent the moments. Assume now that instead of knowing the point estimates of the moments, we have an estimate of the moments in the form of confidence intervals. Loosely speaking, assume that the estimates are given in the form \(\hat{\sigma}_j \in\left[\hat{\sigma}_j^{-}, \hat{\sigma}_j^{+}\right]= \left[\hat{\sigma}_j-\delta_j, \hat{\sigma}_j+\delta_j\right]\), where \(\hat{\sigma}_j, j=1, \ldots, m\), is the point estimate. With these estimates, one can compute a \(100 \%\) confidence interval on the expected value of \(\phi\left(X_1, \ldots, X_n\right)\) over all distributions of the random variables with moments within the confidence intervals by solving the problem

\[ \begin{array}{cl} \max & \left\{\mathbb{E}_{F}\left[\phi\left(X_{1}, \ldots, X_{n}\right)\right]\right\} \\ \text { subject to } & \hat{\sigma}_{j}^{-} \leq \mathbb{E}_{F}\left[f_{j}\left(X_{1}, \ldots, X_{n}\right)\right] \leq \hat{\sigma}_{j}^{+}, \\ & j=1, \ldots, m, \\ & F\left(x_{1}, \ldots, x_{n}\right) \text { a probability } \\ & \text { distribution on } \mathcal{D} \subseteq \mathbb{R}^{n} . \end{array} \tag{5.2} \]

Using the same duality arguments discussed in Section 2, under suitable conditions (similar to Proposition 1), the objective value of the problem above can be found by solving the following dual problem:

\[ \begin{aligned} \min & \left\{y_{0}+\sum_{j=1}^{m}\left(y_{j}^{+} \hat{\sigma}_{j}^{+}-y_{j}^{-} \hat{\sigma}_{j}^{-}\right)\right\} \\ \text {subject to } & y_{0}+\sum_{j=1}^{m}\left(y_{j}^{+}-y_{j}^{-}\right) f_{j}\left(x_{1}, \ldots, x_{n}\right) \\ & \geq \phi\left(x_{1}, \ldots, x_{n}\right) \\ \text { for all } & \left(x_{1}, \ldots, x_{n}\right) \in \mathcal{D}, \\ & y_{j}^{+}, y_{j}^{-} \in \mathbb{R}^{+}, \quad j=1, \ldots, m, \\ & y_{0} \in \mathbb{R} . \end{aligned} \tag{5.3} \]

Note that the dual of the original problem (5.1) is:

\[ \begin{array}{cl} \min & \left\{y_{0}+\sum_{j=1}^{m} y_{j} \sigma_{j}\right\} \\ \text { subject to } & y_{0}+\sum_{j=1}^{m} y_{j} f_{j}\left(x_{1}, \ldots, x_{n}\right) \\ & \geq \phi\left(x_{1}, \ldots, x_{n}\right) \\ \text { for all } & \left(x_{1}, \ldots, x_{n}\right) \in \mathcal{D} \\ & y_{j} \in \mathbb{R}, \quad j=0, \ldots, m. \end{array} \tag{5.4} \]

In both problems (5.3) and (5.4), the objective function is linear, and all the constraints are linear except for the first constraint of the problem. In fact, the difficulty of solving any of these problems comes from the first constraint, which is the “same” for both problems. For the particular semiparametric bound problems we have considered here, we have shown how to address this first constraint using SOS techniques. Because the difference between these two problems is the addition of extra variables that are linear in the objective and in the extra constraints, having an SOS formulation for the original problem given point estimates of moments means that an SOS formulation for the problem with confidence intervals on moments can be obtained in straightforward fashion by accordingly changing the objective of the SOS formulation of the original problem, and adding the adequate linear constraints. The resulting SOS formulation can then be efficiently solved using SOS optimization softwares such as SOSTOOLS. As an example, the SOS formulation (3.5) to obtain a semiparametric upper bound for the extreme probability Pr(X₁ ≤ t₁ and X₂ ≤ t₂) can be modified as follows, in order to consider the situation in which only the confidence intervals on the moments are available:

\[ \begin{aligned} \min & \left(y_{00}+y_{10}^{+} \mu_1^{+}+y_{01}^{+} \mu_2^{+}+y_{20}^{+} \mu_1^{(2)^{+}}+y_{02}^{+} \mu_2^{(2)^{+}}+y_{11}^{+} \mu_{12}^{+}\right. \\ & \left.-y_{10}^{-} \mu_1^{-}-y_{01}^{-} \mu_2^{-}-y_{20}^{-} \mu_1^{(2)^{-}}-y_{02}^{-} \mu_2^{(2)^{-}}-y_{11}^{-} \mu_{12}^{-}\right) \\ \text {subject to } & p\left(t_1-x_1^2, t_2-x_2^2\right)-1 \quad \text { is an SOS polynomial, } \\ & p\left(x_1^2, x_2^2\right) \qquad \qquad \quad \ \ \ \quad \text { is an SOS polynomial, } \\ & y_{10}^{+}, y_{10}^{-}, y_{01}^{+}, y_{01}^{-}, y_{20}^{+}, y_{20}^{-}, y_{02}^{+}, y_{02}^{-}, y_{11}^{+}, y_{11}^{-} \in \mathbb{R}^{+}, \\ & y_{00} \in \mathbb{R} . \end{aligned} \tag{5.5} \]

Above, the {·}⁺ and {·}⁻ in the moments, represent the upper and lower bounds of the confidence interval used to estimate the moments. Similar straightforward modifications can be made for all the SOS formulations of semiparametric bound problems considered in the paper.

It is important to note that in our discussion above we have assumed that the confidence intervals have been obtained independently for each moment. Developing SOS formulations for relevant semiparametric bound problems considering more complex (dependent) confidence intervals will be the topic of future work.

6. Conclusions

In this paper, we have illustrated a new optimization technique known as sum of squares (SOS) programming to find optimal bounds for the probability of extreme events involving two random variables, given only the first and second order moment information. An interesting aspect is that we work solely under the physical measure. This avoids the difficulty of estimating moments of the risk-neutral distribution.

We extend the application of classical moment problems (or semiparametric methods) to finance, insurance, and actuarial science by examining two extreme probability problems, both taking into account correlations between random variables. The first problem allows us to put “100% confidence intervals” on the probability of joint extreme events. The second finds VaR probability bounds on the sum of two variables, given up-to-the-second moment information. In each case the moment information is given by point estimates, which are based on historical observations or judgments from scenario analysis. We provide the examples to illustrate the potential usefulness of moment methods in assessing probability of rare events. We also show that the proposed method can be modified in a straight-forward fashion to obtain semiparametric bounds based on confidence intervals rather than point estimates of the moments.

There are other applications where our approach could be useful. For example, this approach can be used to estimate the default probability of fixed-income securities where incomplete knowledge on the enterprise and economic factors drives the credit risk. In other areas such as inventory and supply chain management, this approach can be applied to find inventory policies that will be applicable to different (unknown) demand distributions in the future. Even when the distributions of the random variables are assumed to be known, this approach can be implemented to measure sensitivity of a joint probability, VaR, or other variables to model misspecification as in Lo (1987) and Hobson, Laurence, and Wang (2005).

Some important issues for future research clearly deserve more investigation. For example, it will be interesting to analyze the bounds on tail distribution given the moments of extreme values instead of the moments of the whole distribution. A further question involves to what extent our results would change if we incorporate distribution class information (e.g., continuous, symmetric, unimodal, etc.) in our bound problems. We leave these questions for future research.