Practical algorithms for value-at-risk portfolio optimization problems

Abstract

This article compares algorithms for solving portfolio optimization problems involving value-at-risk (VaR). These problems can be formulated as mixed integer programs (MIPs) or as chance-constrained mathematical programs (CCMPs). We propose improvements to their state-of-the-art MIP formulations. We also specialize an algorithm for solving general CCMPs, featuring practical interpretations. We present numerical experiments on practical-scale VaR problems using various algorithms and provide practical advice for solving these problems.

Keywords:

• value-at-risk
• portfolio optimization
• mixed integer programming
• mixed integer linear programming
• chance-constrained mathematical programs

JEL Classification:

• C02
• C61
• G11

1. Introduction

Value-at-risk (VaR) and conditional value-at-risk (CVaR) are both risk measures that describe tail risks of financial portfolios. VaR can be defined as a percentile of the loss distribution and CVaR as the average loss for losses exceeding the corresponding VaR. Readers are encouraged to refer to Rockafellar and Uryasev (2002) for rigorous discussions on VaR and CVaR for general distributions. The study of coherent risk measures (Artzner et al. 1999) and the development of a convex CVaR optimization framework (see Krokhmal et al. 2002, Rockafellar and Uryasev 2002, and references therein) show advantages of CVaR over VaR. In particular, scenario-based CVaR problems can be reformulated as linear programs (LPs). Even when the number of assets and the number of scenarios are large, there are special formulations and algorithms that can solve CVaR problems more efficiently than the originally proposed LP reformulation (see Hong and Liu 2009, Lim et al. 2010, Ogryczak and Śliwiński, 2011, for example). Nevertheless, since its introduction as a risk measure VaR has played an important role in both practical applications and academic research (see Jorion 1997, Pritsker 1997, for example). In particular, VaR has been written into industry regulations such as Basel Accord II (BCBS 2004) as a capital adequacy requirement. Despite the aforementioned favorable developments for CVaR, VaR retains a similar status in Basel Accord III (BCBS 2010), which is scheduled to be implemented from 2013 to 2018. In addition to its practical importance, VaR also attracts significant academic research attention. For example, Kou and Peng (2012) provide counterarguments against common criticisms, and Kou et al. (2013) recommend a new risk measure for regulatory purpose, the tail conditional median, which is essentially VaR at a higher confidence level. As a consequence, VaR optimization both is a valuable tool in practice and an interesting area of research.

VaR problems can be solved heuristically or optimally. Readers are encouraged to refer to Gilli and Këllezi (2002), Larsen et al. (2002), and Hong et al. (2014) for various heuristics for VaR problems. In this article we explore ways to solve VaR problems to global optimality. Since VaR is defined as a percentile of the portfolio loss distribution, optimization problems with a VaR constraint or objective can be naturally formulated as chance-constrained mathematical programs (CCMPs). A common approach to solve CCMPs is to reformulate them into mixed integer programs (MIPs). Two MIP formulations for VaR portfolio optimization problems are proposed by Benati and Rizzi (2007) and are commented on by Lin (2009). Moreover, there have been advances in recent years for solving general CCMPs directly (see Qiu et al. 2012, Song et al. 2012, Luedtke 2013, for example). We examine both the traditional approach and the new advances and compare their effectiveness in solving practical VaR problems.

This paper is organized as follows. Section 2 reviews the definition of VaR and presents common MIP formulations for VaR portfolio optimization problems. Section 3 proposes some improvements for the MIP formulations. The application of a decomposition algorithm for general CCMPs to VaR problems is shown in Section 4. Section 5 summarizes results of numerical experiments. Section 6 concludes with practical advice for solving large-scale VaR problems.

2. VaR problems and MIP formulations

Given a universe of n assets and a set of m scenarios, denote the loss (or the negative return) incurred in scenario i by per unit of asset j by . We denote the loss matrix by for , . Each row of the loss matrix, , represents the assets' losses in scenario i, . We denote the probability of realizing scenario i by , and the decision variable by , where is the set of all feasible portfolios. We assume that the loss of asset j in any scenario is proportional to , . Then the portfolio loss random variable is represented by the vector , where represents the portfolio loss in the ith scenario.

Given , VaR at confidence level α, or α-VaR, of the portfolio loss random variable l is defined as its α-percentile. Rigorous definitions of α-VaR for general distributions can be found in Rockafellar and Uryasev (2002). In our settings the α-VaR is computed as follows.

Proposition 2.1.

Assume that for each portfolio the distribution of the portfolio loss is concentrated in finitely many points, that is, the vector has finite dimension. Let those losses be ordered and let the corresponding probabilities be . The α-VaR of the loss is given by where is the unique index such that .

To simplify exposition we assume that all scenarios are equally probable, that is, for . Then the unique index in Proposition 2.1 is given by . The extension to general discrete distributions is straightforward. We are interested in finding a portfolio with minimum α-VaR of portfolio losses among all feasible portfolios, that is, (1) The set of feasible portfolios is described by a given set of portfolio constraints. Examples of portfolio constraints considered in practice are budget constraint (), no short-selling constraints , , etc. One can formulate (1(1) ) into a mathematical program using implication constraints, which requires additional binary variables and an auxiliary variable. In particular, by introducing binary variables and an auxiliary variable we can reformulate (1(1) ) into the following mathematical program with implication constraints. (2a) (2b) (2c) (2d)

Constraints (2b(2a) ) and (2c(2a) ) ensure that , so for any feasible solution the value of z is an upper bound for . Moreover, since and (2) is a minimization problem, if is an optimal solution, then must be the least upper bound for , which is itself. In particular, equals the th largest portfolio loss, that is, for some scenario . If the probability of realizing each scenario is different, that is, , , the only change needed is replacing (2c(2a) ) by .

A common way to solve (1(1) ) is to reformulate (2) into a MIP using the so-called ‘big-M’ constraints given by (3) If one replaces (2b(2a) ) with (3(2d) ), then the MIP formulation for (1(1) ) is obtained, provided that the coefficient vector is sufficiently large. It is clear that for any , if . If , then (2b(2a) ) does not impose any restriction on the variables . Thus if is large enough to ensure the constraint is not active when , then (3(2d) ) is a valid reformulation of (2b(2a) ). However, big-Ms that are too large result in poor continuous relaxations in branch-and-bound methods and hence hinder computational performance. Suitable choices of , for constraints (3(2d) ) is the main focus in the next section.

3. Big-Ms tailored for VaR problems

In general, finding suitable big-Ms for MIP formulations can be a challenging task. For portfolio optimization problems such as (1(1) ) there is often a natural way of bounding the portfolio losses using bounds of individual assets' losses, which can help calculating big-Ms. For example, to deactivate constraint (3(2d) ) it requires for some scenario . Without knowing the optimal solution an upper bound on the difference of scenario losses is needed. Knowing the bounds of portfolio losses for some , then , is valid for (3(2d) ).

Due to limited liability, common stocks have a maximum loss of . The big-Ms considered by Benati and Rizzi (2007) are based on this bound. Similarly, one can derive natural bounds for an insurance portfolio or any long-only credit portfolio because individual losses in these cases are non-negative. Prior to the developments in chance-constrained mathematical programming in the last few years (see Qiu et al. 2012, Luedtke 2013, among others), this natural way of calculating big-Ms was the state-of-the-art approach for solving VaR problems as MIPs. We present a different choice of big-Ms that is tailored to VaR problems. This choice entails additional computational costs but we expect it to improve the overall computational efficiency. The following discussions about big-Ms not only reveal a choice of big-Ms for VaR problems but also provide meaningful interpretations for a specialized decomposition algorithm, which is presented in the next section.

Definition 3.1.

For any two scenarios define the relative excessive loss REL of i relative to j as (4)

Consider the following pivoting procedure: For a scenario , called the pivot scenario, we calculate for all , and then sort them to obtain a permutation of such that (5)

The following proposition shows how the RELs and the pivoting procedure relate to the feasible set of (2). A valid choice of big-Ms for constraints (3(2d) ) is presented based on these relationships.

Proposition 3.2.

Given any pivot scenario the following constraints are valid for any feasible solution to problem (2). (6a) (6b) (6c)

A formal proof for Proposition 3.2 is given in Appendix 1. We present here the intuitive interpretations for constraints (6). By definition of α-VaR in our setting there can be at most scenarios with positive excessive losses over α-VaR. However, for any and any pivot scenario i, by definition of RELs and the pivoting procedure there are at least scenarios whose losses are greater than or equal to . Therefore α-VaR must be at least , which validates constraint (6a). If , the th scenario loss is less than or equal to α-VaR. Therefore the difference between the ith scenario loss and the α-VaR cannot exceed that between the ith scenario's and the th scenario's loss, whose maximum is . The above interpretations remain valid for any upper bound of α-VaR instead of α-VaR. For any feasible solution to (1(1) ), z is an upper bound of α-VaR for portfolio x thus constraints (6a) and (6b) are valid. Based on Equation (6a) we propose the following big-Ms.

Proposition 3.2.

With the following big-Ms, Equation (3(2d) ) are valid constraints for Equation (2): (7) where the second equality holds because .

Proof.

For any scenario i, if then Equation (3(2d) ) requires that . If , then Equation (3(2d) ) becomes Equation (6a). Because Equation (6a) is valid for any feasible solution to (2), the big-M constraint (3(2d) ) is not active. Therefore, Equation (7(6c) ) is valid because this choice of big-M makes Equation (3(2d) ) a valid reformulation of Equation (2b(2a) ).

In general, computing Equation (7(6c) ) requires solving optimization problems. Efficient algorithms for solving these problems, if available, should always be considered. For example, if has one linear constraint and bound constraints then Equation (4(3) ) is a continuous knapsack problem with one constraint, which can be solved efficiently with an algorithm. Suppose is described by a small number of N constraints. In general, Equation (4(3) ) has n variables and N constraints, which is usually a small optimization problem. Qiu et al. (2012) propose an iterative algorithm to tighten the big-Ms. In that proposal each iteration of the algorithm requires solving m optimization problems with variables and constraints, which may be computationally expensive in financial applications where and .

A decomposition algorithm is presented in the next section. We will see the similarities between constraints (6b) and (6c) and the valid inequalities used in the decomposition algorithm in the following discussions.

4. A branch-and-cut decomposition algorithm for VaR problems

Note that Equation (1(1) ) can be viewed as a special case of CCMP because it can be cast as (8) In this section, we apply the decomposition algorithm proposed by Luedtke (2013) for general CCMPs to solve Equation (8(7) ). Some components of the algorithm can be simplified due to the VaR problems' special structures.

The outer loop of decomposition algorithm proposed by Luedtke (2013) is similar to a branch-and-bound method. A generic branch-and-bound method for solving MIPs can be summarized as follows. At the beginning, a continuous relaxation of a master problem is solved where all integrality constraints on variables are relaxed. This initial master problem is commonly referred to as the root node. Two child nodes are then created, which inherit the master problem and one of the constraints or where is an integral variable but has a fractional optimal value at the current node and (respectively, ) is the largest (respectively, smallest) integer that is less than (respectively, greater than) . Based on the optimal solution of a child node, new child nodes are created. In general, a node is fathomed when it is infeasible, when all discrete variables happen to have optimal integer values, or when the optimal objective is worse than the best optimal objective from all previously found integral solutions. A fathomed node will not be processed any further. A node is open if it has not been fathomed yet. The branch-and-bound method terminates when there are no open nodes.

As stated by Luedtke (2013), the only important difference between the decomposition algorithm and a standard branch-and-bound algorithm is how nodes are processed. The decomposition algorithm generates valid inequalities based on optimal solutions to current node. The current node is solved repeatedly if a valid inequality is violated and added to the current node. The initial master problem for the decomposition algorithm may not be a complete description of the original chance-constrained problem so it is possible to generate valid inequalities even when the solution to the current node is integral. For ease of discussion, if the solution to the current node is integral, then the valid inequalities generated with respect to this solution are referred to as lazy constraints. Otherwise when the solution to current node is fractional then the valid inequalities generated are called user cuts. This terminology is consistent with those in the CPLEX documentation (IBM ILOG 2013).

Applying the decomposition algorithm proposed by Luedtke (2013) to solve Equation (8(7) ) requires solving a master problem of the following form: (9a) (9b) (9c) (9d) where is a polyhedron that contains the feasible region of Equation (2), and are such that . Note that the description of C may not contain all implication constraints at the beginning, but it will be augmented when processing the branch-and-bound tree so that correctness of the algorithm is guaranteed.

The following components are needed to fully specify our decomposition algorithm:

1. Given a solution to the current node, identify pivot scenarios from which violated valid inequalities may be obtained.

2. Solve subproblems that are required to generate valid inequalities.

3. Identify a good violated inequality from the set of inequalities generated from step (ii).

Component (i) can be specified in a straightforward manner. As shown by Luedtke (2013), given a solution to the current node, any scenario whose implication constraint is violated can be considered to generate valid inequalities. For our problem specifically, we consider a scenario i as a pivot scenario if . There may be multiple scenarios that satisfy this condition. Luedtke (2013) shows that if the candidate solution is not feasible for Equation (2), then for any of these scenarios at least one valid inequality, which will be presented below, is violated and can be added to C.

To specify component (ii), suppose scenario i is chosen to be a pivot scenario. If one follows the steps of the general decomposition algorithm in the original proposal (Luedtke 2013) then the required single-scenario subproblems for scenario i are for all . We suggest the following simplifications to these subproblems: (10) because and appears in only one constraint . Our suggestion helps to reduce one variable and one constraint in every subproblem. This reduction not only helps to solve the subproblems more efficiently but also lends the subproblems practical interpretations in the context of VaR-optimization, as discussed in Section 3.

To condense notation, the dependence on pivot scenario i for and is suppressed in the following discussion. After pivoting on scenario i (i.e. calculating its RELs and sorting them) the following inequalities (11) are valid for the feasible set of Equation (2). Inequalities (11(10) ) can be viewed as an application of Lemma 1 of Luedtke (2013) to the feasible set of Equation (2). Alternatively, they can also be seen as a way to formulate both Equations (6a) and (6b) in one inequality for scenarios . In particular, if , then Equation (11(10) ) resembles Equation (6a). If , then Equation (11(10) ) resembles Equation (6b). The following theorem presents a family of valid inequalities by mixing Equation (11(10) ). The valid inequalities generated in component (ii) are presented in the following theorem.

Theorem 4.1.

For any pivot scenario i, let be any non-empty index set such that for where . Then the inequalities (12) are valid for the feasible set of Equation (2).

Theorem 4.1 is a special case of Theorem 3 of Luedtke (2013) and readers are encouraged to refer to his discussion on the strength of these inequalities. It is beneficial to have in Equation (12(11) ) because in that case the inequalities are facet-defining in a certain sense (Luedtke 2013). Moreover, setting for some yields the inequality (13) which dominates inequality (11(10) ) for this particular k. We observe that by setting inequalities (12(11) ) resemble the form of big-M constraints (3(2d) ) with big-Ms defined in Equation (7(6c) ). Therefore the specialized decomposition algorithm is able to generate big-M constraints when needed.

It is shown by Luedtke (2013) that constraints (12(11) ) are both a necessary and sufficient description for the implication constraints (2b(2a) ). The discussion on big-Ms in Section 2 sheds some light on the interpretations of Equation (12(11) ). They are a compact way of expressing the necessary conditions (6c). Specifically, any that satisfies Equation (12(11) ) also satisfies Equation (6c).

Lastly, we specify component (iii). For each pivot scenario i there are exponentially many ( to be precise) valid inequalities of the form (12(11) ). Separating a ‘good’ valid inequality from this exponential family is the next task. An separation algorithm is proposed by Günlük and Pochet (2001) and is employed in Luedtke (2013). We show that the complexity of the separation algorithm proposed by Günlük and Pochet (2001) improves from to when it is implemented as a subroutine of the decomposition algorithm. This is not only true in our VaR problem but also true in general when implementing the decomposition algorithm in Luedtke (2013) to solve CCMPs.

We first formulate the separation problem mathematically. Adding to both sides of Equation (12(11) ) we obtain . Then we have by replacing the RHS with a telescopic sum. Finally, we group the like terms and arrive at the following inequality: (14) Although Equation (14(13) ) is algebraically equivalent to Equation (12(11) ), its LHS is independent of the choice of the index set R. Given a candidate solution to the current node, that is, may be fractional, one can minimize the RHS of Equation (14(13) ) to identify a most violated inequality in Equation (12(11) ), if any, or determine that none can be found. The minimization problem to be solved is as follows. (15a) (15b) (15c)

Next we examine optimal solutions for Equation (15). Observe that the index set for any feasible R is also feasible and has objective value no greater than that of R. Therefore there is at least one optimal solution in which is a member. Following the arguments in Günlük and Pochet (2001) we can find an optimal index set which satisfies the conditions and . An algorithm for identifying such an index set is shown in lines 6 to 9 in Algorithm 2 in Appendix 2.

We state a basic version of the specialized algorithm for solving VaR problems in Algorithm 1 in the appendix. The algorithm is a specialization of the one proposed by Luedtke (2013) with a detailed separation subroutine Algorithm 2. The general mechanism of the algorithm can be found in its original proposal and is omitted to avoid repetition. The user cuts added (Line 17 of Algorithm 1) do not affect the correctness of the algorithm but may be used to improve computational efficiency. For example, one can include user cuts every other 50 nodes in the branch-and-bound tree or include only user cuts that have significant violations. In addition, this decomposition algorithm can be used in conjunction with the special branching scheme proposed by Qiu et al. (2012). The special branching scheme adds a cut of the form to nodes that have and all subsequent child nodes. Such local cuts enforce the implication which does not affect the correctness of the decomposition algorithm. Readers are encouraged to refer to Qiu et al. (2012) for details of the special branching scheme.

5. Numerical experiments

Numerical experiments are conducted using CPLEX 12.5.1 on a desktop running Ubuntu 12.04.4 LTS, with 8 cores 3.4 GHz Intel Core i7, 32 GB RAM and 32 GB swap space. We set the number of threads to one for fair comparison. Some advanced CPLEX callbacks are necessary to implement the decomposition algorithm and the special branching scheme described at the end of Section 4. However, some advanced CPLEX features, such as dynamic search, are disabled when advanced callbacks are used. Unless stated otherwise, advanced CPLEX features are disabled for fair comparison.

We considered VaR-minimization problems of the form for reinsurance portfolios. An anonymous reinsurance company generously provided us with simulated reinsurance loss and premium data. The data was hypothetically simulated to represent realistic risk and returns of real business. The characteristics of these loss matrices include non-negativity, sparsity (e.g. to non-zero entries), positive skewness, and high kurtosis. It is shown by Qiu et al. (2012) that sparsity is favorable for computation efficiency. The premiums , are the market premiums of the risk contracts. The decision variable , represents the fraction of risk contract j that the reinsurer is willing to take. The portfolio premium is then cx. A reward constraint impose a minimum level of portfolio premium. In our experiments we set , , and .

In our first set of experiments we compare the following algorithms:

1. NaturalM: MIP formulation with a natural choice of big-Ms. Since reinsurance losses and the decision variables , the portfolio loss in scenario i is bounded by , . The discussion at the beginning of Section 2 explains that , is a natural choice of big-Ms.

2. TightM: MIP formulation with big-Ms defined in Equation (7(6c) ). In our experiments contains only one linear constraint and bound constraints. Therefore, the REL subproblems can be solved with a algorithm.

3. Decomp: Decomposition algorithm with lazy constraints only. The initial master problem has an empty description in the cut set C. This is a minimal form of the decomposition algorithm because no user cut is added. When a strong valid inequality of the form (12(11) ) is violated by a integer solution to the current node, we add not only this inequality but also any inequality of the form (13(12) ) that is violated. As observed by Luedtke (2013), including inequalities (13(12) ) is a simple way to add additional sparse valid inequalities in one round.

4. Decomp+: Decomposition algorithm with lazy constraints as well as user cuts at the root node. At the root node we continue generating valid inequalities with violation greater than 0.01 until none can be found.

5. TightMB: TightM with special branching scheme implemented, as discussed at the end of Section 4.

6. DecompB: Decomp with special branching scheme implemented.

We minimize the -VaR for 10 instances of reinsurance problems with 5000 scenarios and 25 risk contracts. A time limit of 3 hours for each instance is set, which is a realistic computational budget. We are interested in both the CPU time if an instance is solved and the final optimality gap if an instance is halted due to the time limit.

Table 1 summarizes the CPU time of the 10 instances using different algorithms. We define the optimality gap as where is the objective value of best incumbent found up to the current node and is the best lower bound. Table 2 summarizes the optimality gaps after 3 hours if an instance is not solved within the time limit. NaturalM serves as a benchmark for the comparison because it is the state-of-the-art method. We see from tables 1 and 2 that TightM, TightMB, Decomp, and Decomp+ all have better performance than NaturalM: More instances are solved within the time limit, it takes less CPU time to solve an instance, and there is a smaller optimality gap when a problem is halted. It takes only about 20 seconds to solve subproblems to calculate the tight big-Ms. Recognizing the knapsack problem structure enables us the calculate the tight big-Ms so efficiently. Our experiments show that TightM and TightMB are the fastest algorithms to solve a given instance within time limit. On average, TightM and TightMB can solve an instance more than 20 times faster than the NaturalM does. However, since this observation takes into account only problems that are solved within an arbitrary time limit, it may be an overstatement about the improvement of TightM and TightMB over NaturalM. For instances which are halted, either TightM, TightMB, or Decomp+ has the smallest optimality gap. In 1 out of 10 instances Decomp takes less CPU time than TightM does and Decomp+ has a smaller optimality gap than TightM does.

Table 1. CPU time (in seconds) for solving -VaR minimization problems with 5000 scenarios and 25 risk contracts (including time for calculating big-Ms).

Instance ID	NaturalM	TightM	TightMB	Decomp	Decomp+	DecompB
1	(–)	1135.0	1080.6	(–)	3480.8	(–)
2	(–)	(–)	(–)	(–)	(–)	(–)
3	(–)	(–)	(–)	(–)	(–)	(–)
4	7035.6	29.1	29.0	168.8	97.7	2589.7
5	1699.6	45.1	45.5	150.8	175.0	(–)
6	(–)	473.5	465.0	1649.5	1439.2	(–)
7	(–)	3232.9	2775.9	2686.7	3506.8	6320.7
8	4757.7	470.7	475.7	934.2	797.3	(–)
9	(–)	(–)	(–)	(–)	(–)	(–)
10	(–)	(–)	(–)	(–)	(–)	(–)

Table 2. Optimality gap after 3 hours for solving -VaR minimization problems with 5000 scenarios and 25 risk contracts.

Instance ID	NaturalM	TightM	TightMB	Decomp	Decomp+	DecompB
1	23.69%	(–)	(–)	0.79%	(–)	16.91%
2	24.16%	4.97%	5.39%	6.86%	5.33%	7.55%
3	43.40%	16.02%	13.92%	13.34%	10.78%	33.61%
4	(–)	(–)	(–)	(–)	(–)	(–)
5	(–)	(–)	(–)	(–)	(–)	0.50%
6	11.56%	(–)	(–)	(–)	(–)	7.84%
7	16.45%	(–)	(–)	(–)	(–)	(–)
8	(–)	(–)	(–)	(–)	(–)	5.72%
9	23.52%	7.28%	6.68%	10.90%	8.60%	24.53%
10	18.35%	0.87%	1.63%	5.72%	2.29%	25.06%

The addition of a special branching scheme to TightM yields little additional improvement in either CPU time or optimality gap, as evidenced by comparing columns TightM and TightMB in tables 1 and 2. We also observe a significant deterioration in CPU time and optimality gap when the special branching scheme is implemented in conjunction with Decomp when comparing columns Decomp and DecompB in tables 1 and 2. Even when DecompB manages to solve an instance, its CPU time is much longer than that of Decomp. We currently have no explanation why the special branching scheme deteriorates the decomposition algorithm so significantly.

For NaturalM and TightM, it is for the sake of fair comparison that the advanced CPLEX features were disabled. In practice, however, these features should be enabled for NaturalM and TightM. The second to fifth columns of table 3 summarize the CPU time and final optimality gap for the same 10 instances when CPLEX features are enabled (CPLEX default settings). We see that enabling CPLEX features significantly improves performance for both NaturalM and TightM. More problems are solved within the given time limit and problems are solved much faster. Under the CPLEX default settings TightM can be 2–14 times faster than NaturalM. When a problem is halted, TightM always has a smaller optimality gap than that of NaturalM. With TightM, 9 out 10 instances have final optimality gaps less than . Compared to NaturalM with CPLEX features disabled, the computational enhancement from calculating the tight big-Ms is more significant than enabling the CPLEX features. Comparing the second column in table 3 and the third column (TightM) in table 1 we see that TightM with CPLEX features disabled takes less CPU time to solve an instance than NaturalM with the features enabled.

Table 3. CPU time (in seconds, including big-M calculation time) and optimality gap for solving -VaR minimization problems with 25 risk contracts.

	5000 scenarios and 3-hours limit				10 000 scenarios and 12-hours limit
	NaturalM		TightM		NaturalM		TightM
Instance ID	CPU time	Opt. gap	CPU time	Opt. gap	CPU time	Opt. gap	CPU time	Opt. gap
1	2376.8	(–)	149.9	(–)	(–)	9.75%	(–)	2.10%
2	(–)	6.93%	(25 071.3)	2.04%	(–)	14.44%	(–)	7.97%
3	(–)	10.88%	(134 222.2)	6.18%	691.3	(–)	97.0	(–)
4	53.6	(–)	29.9	(–)	(–)	11.21%	10 128.5	(–)
5	111.8	(–)	41.0	(–)	13 568.1	(–)	757.1	(–)
6	1162.2	(–)	119.7	(–)	(–)	9.52%	34 225.9	(–)
7	3993.1	(–)	281.8	(–)	(–)	10.12%	(–)	1.95%
8	609.6	(–)	182.0	(–)	(–)	2.77%	16 144.2	(–)
9	(–)	8.68%	(17 423.2)	1.66%	(1807.9)	15.95%	(–)	8.10%
10	(–)	3.35%	1859.8	(–)	(–)	17.63%	(–)	5.95%

Notes: CPLEX features enabled. Number of scenarios and time limit are as specified in table. CPU time in a parenthesis denotes the eventual solve time of an instance.

The algorithms can find a high-quality solution quickly but devote most of the computation time to proving optimality. We define the incumbent gap as where and are the objective of the current best incumbent and of the optimal solution, respectively. Table 4 shows incumbent gaps for various algorithms after 1 minute and 10 minutes. Although solving VaR-minimization problems of the size we considered may take hours or even days, as shown in the fourth column of table 3, the algorithms can find a high-quality feasible solution in a few minutes. For example, using TightM with CPLEX features enabled (the fifth column in table 4), incumbents within of optimal solutions are found in 1 minute for all 10 instances. Moreover, after 10 minutes TightM with CPLEX features enabled found incumbents within of the optimal solutions for 4 instances, yet 3 of them took more than 12 hours to eventually be solved to optimality.

Table 4. Percentage gaps between best incumbent and optimal solution for solving -VaR minimization problems with 5000 scenarios and 25 risk contracts.

	Incumbent gap after 1 minute				Incumbent gap after 10 minutes
	Features disabled		Features enabled		Features disabled		Features enabled
Instance ID	NaturalM	TightM	NaturalM	TightM	NaturalM	TightM	NaturalM	TightM
1	0.28%	0.00%	0.00%	0.00%	0.00%	0.00%	0.00%	Solved
2	1.76%	1.10%	0.70%	0.78%	0.21%	0.00%	0.02%	0.01%
3	1.67%	1.81%	0.44%	0.65%	1.67%	0.00%	0.00%	0.00%
4	0.00%	Solved	Solved	Solved	0.00%	Solved	Solved	Solved
5	8.65%	0.00%	Solved	Solved	0.00%	Solved	Solved	Solved
6	2.61%	0.63%	0.00%	0.00%	0.38%	0.00%	0.00%	Solved
7	1.07%	0.70%	0.39%	0.21%	0.00%	0.00%	0.00%	Solved
8	0.94%	0.76%	0.00%	0.75%	0.00%	Solved	Solved	Solved
9	0.37%	0.37%	0.37%	0.63%	0.37%	0.37%	0.00%	0.00%
10	0.89%	0.44%	0.79%	0.00%	0.44%	0.00%	0.44%	0.00%

In the second set of experiments we aim to provide practical advice for solving large scale problems. All CPLEX advanced features are enabled. We considered 10 instances of reinsurance problems with 10 000 scenarios and 25 risk contracts. A -VaR minimization problem contains 50 tail scenarios, that is, . The number of scenarios and number of tail scenarios of this size are similar to problems considered in practice. In this case we set a time limit of 12 hours to replicate real-life situations where an over-night computational budget is allowed. The results are summarized in the sixth to the ninth columns of table 3. Out of 10 large scale problems, TightM solved 5 instances within the time limit and had a final optimality gap less than for 7 instances. NaturalM, however, solved only 2 instances within the time limit and had 3 instances whose optimality gaps were less than . For these large instances TightM is 7–18 times faster than NaturalM when an instance is solved within the time limit. When an instance is halted, TightM always has a smaller optimality gap than NaturalM.

6. Conclusion

We examined several algorithms to solve VaR-minimization problems. Our experiments suggest that calculating tight big-Ms improves overall performance. When CPLEX's advanced features are disabled, variations of the specialized decomposition algorithm perform comparably to the MIP formulation with tight big-Ms. We hope to see significant improvements in computational performance if we were able to enable some of the CPLEX features for the decomposition algorithm. The efficacy of these methods when VaR-constrained problems remains to be investigated. At the current stage, we recommend TightM as a practical method to solve VaR problems.

Acknowledgments

The authors thank Dr James Luedtke and Dr Simge Küçükyavuz for offering insights into this work. Mingbin Feng and Andreas Wächter were supported partially by the NSF grant DMS-1216920.

Appendix 1. Proof for Proposition 3.2

Proof.

Let i be a pivot scenario. Inequality (6a) is backed by (2b(2a) ), (2c(2a) ) and the definition of RELs. Since Equation (2b(2a) ) implies and Equation (2c(2a) ) implies , thus . By definition of RELs and the ordering (5(4) ) we know that . Thus . Assume to the contrary that Equation (6a) is violated by any solution , that is, . Then we have , which is a contradiction. Therefore Equation (6a) is valid for all feasible solution to Equation (2).

The implication constraint (6b) results from the definition of REL and the implication constraint (2b(2a) ). By Equation (2b(2a) ) we know that for all . So implies that , which validates the implication constraint (6b).

Constraints (6c) is a consequence of (6a), (6b) and the ordering (5(4) ). This completes the proof for Proposition 3.2.

Appendix 2. Specialized decomposition algorithm for VaR problems