Regularized Stochastic Dual Dynamic Programming for convex nonlinear optimization problems
Vincent Guigues, Miguel Lejeune, Wajdi Tekaya
TL;DR
This paper introduces regularized variants of dual dynamic programming algorithms, REDDP and SDDP-REG, for nonlinear stochastic optimization, demonstrating significant computational speedups in portfolio management problems with market impact costs.
Contribution
It extends dual dynamic programming methods with regularization to nonlinear and stochastic problems, providing convergence proofs and demonstrating substantial efficiency improvements.
Findings
REDDP is up to 184 times faster than traditional DDP.
SDDP-REG outperforms standard SDDP in portfolio problems with market impact costs.
Numerical results confirm the efficiency and convergence of the proposed algorithms.
Abstract
We define a regularized variant of the Dual Dynamic Programming algorithm called REDDP (REgularized Dual Dynamic Programming) to solve nonlinear dynamic programming equations. We extend the algorithm to solve nonlinear stochastic dynamic programming equations. The corresponding algorithm, called SDDP-REG, can be seen as an extension of a regularization of the Stochastic Dual Dynamic Programming (SDDP) algorithm recently introduced which was studied for linear problems only and with less general prox-centers. We show the convergence of REDDP and SDDP-REG. We assess the performance of REDDP and SDDP-REG on portfolio models with direct transaction and market impact costs. In particular, we propose a risk-neutral portfolio selection model which can be cast as a multistage stochastic second-order cone program. The formulation is motivated by the impact of market impact costs on large…
|
|
|
|||||
|---|---|---|---|---|---|---|---|
|
with | ||||||
|
|||||||
|
with | ||||||
|
| 10 | 50 | 100 | 150 | 200 | 250 | 300 | 350 | |
|---|---|---|---|---|---|---|---|---|
| DDP | 3 | 69 | 268 | 780 | 1304 | 2400 | 4289 | 5348 |
| DDP-REG-PREV-REG2 | 1 | 4 | 8 | 13 | 17 | 30 | 25 | 29 |
| DDP-REG-PREV-REG1-0.2 | 1 | 4 | 12 | 21 | 28 | 23 | 60 | 56 |
| DDP-REG-PREV-REG1-0.9 | 1 | 4 | 8 | 12 | 17 | 21 | 25 | 29 |
| DDP-REG-AVG-REG2 | 1 | 4 | 8 | 13 | 17 | 30 | 25 | 29 |
| DDP-REG-AVG-REG1-0.2 | 1 | 5 | 8 | 12 | 17 | 21 | 47 | 55 |
| DDP-REG-AVG-REG1-0.9 | 1 | 4 | 9 | 13 | 17 | 22 | 26 | 30 |
| 10 | 50 | 100 | 150 | 200 | 250 | 300 | 350 | |
|---|---|---|---|---|---|---|---|---|
| DDP-REG-PREV-REG2 | 3.0 | 17.3 | 33.5 | 60.0 | 76.7 | 80.0 | 171.6 | 184.4 |
| DDP-REG-PREV-REG1-0.2 | 3.0 | 17.3 | 22.3 | 37.1 | 46.6 | 104.4 | 71.5 | 95.5 |
| DDP-REG-PREV-REG1-0.9 | 3.0 | 17.3 | 33.5 | 65.0 | 76.7 | 114.3 | 171.6 | 184.4 |
| DDP-REG-AVG-REG2 | 3.0 | 17.3 | 33.5 | 60.0 | 76.7 | 80.0 | 171.6 | 184.4 |
| DDP-REG-AVG-REG1-0.2 | 3.0 | 13.8 | 33.5 | 65.0 | 76.7 | 114.3 | 91.3 | 97.2 |
| DDP-REG-AVG-REG1-0.9 | 3.0 | 17.3 | 29.8 | 60.0 | 76.7 | 109.1 | 165.0 | 178.3 |
| 10 | 50 | 100 | 150 | 200 | 250 | 300 | 350 | |
|---|---|---|---|---|---|---|---|---|
| DDP | 10 | 26 | 39 | 58 | 66 | 83 | 100 | 104 |
| DDP-REG-PREV-REG2 | 3 | 3 | 3 | 3 | 3 | 4 | 3 | 3 |
| DDP-REG-PREV-REG1-0.2 | 3 | 3 | 3 | 3 | 3 | 3 | 6 | 5 |
| DDP-REG-PREV-REG1-0.9 | 3 | 3 | 3 | 3 | 3 | 3 | 3 | 3 |
| DDP-REG-AVG-REG2 | 3 | 3 | 3 | 3 | 3 | 4 | 3 | 3 |
| DDP-REG-AVG-REG1-0.2 | 3 | 3 | 3 | 3 | 3 | 3 | 5 | 5 |
| DDP-REG-AVG-REG1-0.9 | 3 | 3 | 3 | 3 | 3 | 3 | 3 | 3 |
| T | Variant | CPU time (s) | Number of iterations |
|---|---|---|---|
| SDDP | 14 | 4 | |
| SDDP-REG-PREV-REG2 | 6 | 2 | |
| SDDP | 29 | 5 | |
| SDDP-REG-PREV-REG2 | 21 | 4 | |
| SDDP | 40 | 6 | |
| SDDP-REG-PREV-REG2 | 18 | 3 |
| Variant | CPU time (s) | Number of iterations |
|---|---|---|
| SDDP | 3895 | 50 |
| SDDP-REG-PREV-REG2 | 3921 | 50 |
| Variant | CPU time (s) | Number of iterations | ||
|---|---|---|---|---|
| 12 | SDDP | 20 | 8 | |
| 12 | SDDP-REG-PREV-REG2 | 7 | 3 | |
| 12 | SDDP | 6 | 3 | |
| 12 | SDDP-REG-PREV-REG2 | 7 | 3 | |
| 20 | SDDP | 43 | 10 | |
| 20 | SDDP-REG-PREV-REG2 | 19 | 5 | |
| 20 | SDDP | 11 | 3 | |
| 20 | SDDP-REG-PREV-REG2 | 11 | 3 | |
| 24 | SDDP | 55 | 10 | |
| 24 | SDDP-REG-PREV-REG2 | 13 | 3 | |
| 24 | SDDP | 57 | 11 | |
| 24 | SDDP-REG-PREV-REG2 | 13 | 3 |
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Regularized Stochastic Dual Dynamic Programming for convex nonlinear optimization problems
Abstract.
We define a regularized variant of the Dual Dynamic Programming algorithm called DDP-REG to solve nonlinear dynamic programming equations. We extend the algorithm to solve nonlinear stochastic dynamic programming equations. The corresponding algorithm, called SDDP-REG, can be seen as an extension of a regularization of the Stochastic Dual Dynamic Programming (SDDP) algorithm recently introduced which was studied for linear problems only and with less general prox-centers. We show the convergence of DDP-REG and SDDP-REG. We assess the performance of DDP-REG and SDDP-REG on portfolio models with direct transaction and market impact costs. In particular, we propose a risk-neutral portfolio selection model which can be cast as a multistage stochastic second-order cone program. The formulation is motivated by the impact of market impact costs on large portfolio rebalancing operations. Numerical simulations show that DDP-REG is much quicker than DDP on all problem instances considered (up to 184 times quicker than DDP) and that SDDP-REG is quicker on the instances of portfolio selection problems with market impact costs tested and much faster on the instance of risk-neutral multistage stochastic linear program implemented (8.2 times faster).
[TABLE]
[TABLE]
Keywords. Stochastic Optimization, Stochastic Dual Dynamic Programming, Regularization, Portfolio Selection, Market Impact Costs.
AMS subject classifications. 90C15, 90C90.
Multistage stochastic optimization problems are used to model many real-life applications where a sequence of decisions has to be made, subject to random costs and constraints arising from the observations of a stochastic process, see [36, 41] for an overview on multistage stochastic programs. Solving such problems is challenging and often requires some assumptions on the underlying stochastic process, on the problem structure, and some sort of decomposition. In this paper, we are interested in problems for which deterministic or stochastic dynamic programming equations can be written. In this latter case, we will focus on situations where the underlying stochastic process is discrete interstage independent, the number of stages is moderate to large, and the state vector is of small size.
Two popular solution methods to solve stochastic dynamic programming equations are Approximate Dynamic Programming (ADP) and Stochastic Dual Dynamic Programming (SDDP) [32]. Several enhancements of SDDP have been proposed such as the extension to interstage dependent stochastic processes [22, 16], the introduction and analysis of risk-averse variants [20, 21, 23, 34, 40, 42], cut selection strategies [33, 8, 18, 5], and convergence proofs of the algorithm and variants in [35, 13, 17, 19, 5]. However, a known drawback of the method is its slow convergence rate. To cope with this difficulty, a regularized variant of SDDP was recently proposed in [4] for Multistage Stochastic Linear Programs (MSLPs). This variant consists in computing in the forward pass of SDDP the trial points using a Tikhonov regularization [43]. More precisely, the objective is penalized with a quadratic term depending on a prox-center updated at each iteration. On the tests reported in [4], the regularized method converges faster than the classical SDDP method on risk-neutral instances of MSLPs. On the basis of these encouraging numerical results, several natural questions arise:
- a)
When specialized to deterministic problems, how does the regularized method behave? How to extend the method when nonlinear objective and constraints are present and under which assumptions to ensure the convergence of the method?
- b)
How can the regularized algorithm be extended to solve Multistage Stochastic NonLinear Problems (MSNLPs) and under which assumptions to ensure the convergence of the algorithm to an optimal policy?
- c)
What other prox-centers and penalization schemes can be proposed? Find a MSLP for testing the new prox-centers and penalization schemes. Can we observe on this application a faster convergence of the regularized method, as for the application considered in [4]?
- d)
Find a relevant application, modeled by a MSNLP, to test the regularized variant of SDDP.
The objective of this paper is to study items a)-d) above. Our findings on these topics are as follows:
- a)
Regularized Dual Dynamic Programming. We propose a regularized variant of Dual Dynamic Programming (DDP, the deterministic counterpart of SDDP) called DDP-REG, for nonlinear optimization problems. For DDP-REG, in Theorem 1.4, we show the convergence of the sequence of approximate first step optimal values to the optimal value of the problem and that any accumulation point of the sequence of trial points is an optimal solution of the problem. The same proof, with weaker assumptions can be used to show the convergence of this regularized variant of DDP applied to linear problems.
We then consider instances of a portfolio problem with direct transaction costs with a large number of stages and compare the computational time required to solve these instances with DDP and DDP-REG. In all experiments, the computational time was drastically reduced using DDP-REG. More precisely, we tested 6 different implementations of DDP-REG and for problems with , and time periods, the range (for these 6 implementations) of the reduction factor of the overall computational time with DDP-REG was respectively , , , , , , , and .
- b)
SDDP-REG: Regularized SDDP. We define a Regularized SDDP method for MSNLPs which samples in the backward pass to compute cuts at trial points computed, as in [4], in a forward pass, penalizing the objective with a quadratic term depending on a prox-center. In Theorem 3.2, we show the convergence of this algorithm. More precisely, we show (i) the convergence of the sequence of the optimal values of the approximate first stage problems and that (ii) any accumulation point of the sequence of decisions can be used to define an optimal solution of the problem. It will turn out that (ii) improves already known results for SDDP.
- c)
On prox-centers, penalization parameters, and on the performance of the regularization for MSLPs. We propose new prox-centers and penalization schemes and test them on risk-neutral and risk-averse instances of portfolio selection problems.
- d)
Portfolio Selection with Direct Transaction and Market Impact Costs. The multistage optimization models studied in this paper are directly applicable in finance and in particular for the rebalancing of portfolios that incur transaction costs. Transaction costs can have a major impact on the performance of an investment strategy (see, e.g., the survey [7]). Two main types of transaction costs, implicit and explicit, can be distinguished.
Explicit or direct transaction costs are directly observable, are directly charged to the investor, and are generally modelled as linear or piecewise linear. In reality, it is however not possible to trade arbitrary large quantities of securities at their current theoretical market price.
Implicit or indirect costs, often called market impact costs, result from imperfect markets due for example to market or liquidity restrictions (e.g., bid-ask spreads), depend on the order-book situation when the order is executed, and are not itemized explicitly, thereby making it difficult for investors to recognize them. Yet, for large orders, they are typically much larger than the direct transaction costs. Market impact costs are equal to the difference between the transaction price and the (unperturbed) market price that would have prevailed if the trade had not occurred [44]. Market impact costs are typically nonlinear (see, e.g., [1, 3, 12, 14, 44]), and much more challenging to model than direct transaction costs. Market impact costs are particularly important for large institutional investors, for which they can represent a major proportion of the total transaction costs [27, 44]. They can be viewed as an additional price for the immediate execution of large trades.
There is a widespread interest in the modeling and analysis of market impact costs as they are (one of) the main reducible parts of the transaction costs [27]. In this study, we propose a series of dynamic - deterministic and stochastic (risk-neutral and risk-averse) - optimization models for portfolio optimization with direct transaction and market impact costs.
We compare the computational time required to solve with SDDP-REG and SDDP instances of risk-neutral and risk-averse portfolio problems with direct transaction costs. We also compare the computational time required to solve with SDDP-REG and SDDP risk-neutral instances of portfolio problems with market impact costs using real data and stages. To our knowledge, no dynamic optimization problem for portfolio optimization with conic market impact costs has been proposed so far. Also, we are not aware of other published numerical tests on the application of SDDP to a real-life application modelled by a multistage stochastic second-order cone program with a large number of stages.
The paper is organized as follows. In section 1, we present a class of convex deterministic nonlinear optimization problems for which dynamic programming equations can be written. We propose the variant DDP-REG of DDP to solve these problems and show the convergence of the method in Theorem 1.4. In section 2, we introduce the type of stochastic nonlinear problems we are interested in and propose SDDP-REG, a regularized decomposition algorithm to solve these problems. In section 3, we show in Theorem 3.2 the convergence of SDDP-REG. The portfolio selection models described in item d) above are discussed in section 4. Finally, the last section 5 presents the results of numerical simulations that illustrate our results. We show that DDP-REG is much quicker than DDP on all problem instances considered (up to 184 times quicker than DDP) and that SDDP-REG is quicker on the instances of nonlinear stochastic programs tested and much faster on the instance of risk-neutral multistage stochastic linear program implemented (8.2 times faster).
We use the following notation and terminology:
-
The usual scalar product in is denoted by for . The corresponding norm is .
-
is the relative interior of set .
-
.
-
dom() is the domain of function .
-
is the normal cone to set at point .
-
is the Average Value-at-Risk with confidence level , [37].
-
is the diameter of set .
-
The notation represents the matrix \left(\begin{array}[]{l}A\\ B\end{array}\right).
1. Regularized dual dynamic programming: Algorithm and convergence
1.1. Problem formulation and assumptions
Consider the problem
[TABLE]
where is given by
[TABLE]
is a convex function, , and is given.
For this problem, we can write dynamic programming equations defining recursively the functions on as
[TABLE]
with the convention that . Clearly, is the optimal value of (1.1). More generally, we have
[TABLE]
We make the following assumptions: setting and
[TABLE]
(H0) there exists such that for
- (a)
is nonempty, convex, and compact;
- (b)
is proper, convex, and lower semicontinuous;
- (c)
Setting , for , the -th component function is a convex lower semicontinuous function;
- (d)
and for every , there exists such that and ;
- (e)
if , there exists
[TABLE]
such that , and .
The DDP algorithm solves (1.1) exploiting the convexity of recourse functions :
Lemma 1.1**.**
Consider recourse functions , given by (1.2). Let Assumptions (H0)-(a), (H0)-(b), (H0)-(c), and (H0)-(d) hold. Then for , is convex, finite on , and Lipschitz continuous on .
Proof: See the proof of Proposition 3.1 in [17].
The description of the subdifferential of given in the following proposition will be useful for DDP, DDP-REG, and SDDP-REG:
Proposition 1.2**.**
*Lemma 2.1 in [17]. Let Asssumptions (H0) hold. Let be an optimal solution of (1.2). Then for every for every , if and only if *
[TABLE]
where I(x_{t-1},x_{t}(x_{t-1}))=\Big{\{}i\in\{1,\ldots,p\}\;:\;g_{t,i}(x_{t-1},x_{t}(x_{t-1}))=0\Big{\}}.
Proof: See [17].
1.2. Dual Dynamic Programming
We first recall DDP method to solve (1.2). It uses approximations of . At iteration , let functions such that
[TABLE]
be given and define for the function as
[TABLE]
Clearly, (1.4) implies that: It is assumed that the functions can be evaluated at any point . The DDP algorithm works as follows:
DDP (Dual Dynamic Programming).
Step 1)
Initialization. Let , satisfying (1.4) be given. Set .
Step 2)
Forward pass. Setting , for , compute
(1.5)
Step 3)
Backward pass. Define . For , solve
(1.6)
using Proposition 1.2 take a subgradient of at , and store the new cut
for , making up the new approximation .
Step 4)
Do and go to Step 2).
1.3. Regularized Dual Dynamic Programming
For the regularized DDP to be presented in this section, we still define
[TABLE]
[TABLE]
However, since the function computed by regularized DDP is different from the function computed by DDP, the functions obtained with respectively regularized DDP and DDP are different. The regularized DDP algorithm is given below:
Regularized DDP (DDP-REG).
Step 1)
Initialization. Let , satisfying (1.4) be given. Set .
Step 2)
Forward pass. Setting , for , compute
(1.8)
where the prox-center is any point in and where is
for some exogenous nonnegative penalization with if or .
Step 3)
Backward pass. Define . For , solve
(1.9)
using Proposition 1.2 take a subgradient of at , and store the new cut
for , making up the new approximation .
Step 4)
Do and go to Step 2).
Observe that the backward passes of the regularized and non-regularized DDP are the same. The algorithms differ from the way the trial points are computed: for regularized DDP a proximal term is added to the objective function of each period to avoid moving too far from the prox-center.
1.4. Convergence analysis
The following lemma will be useful to analyze the convergence of regularized DDP:
Lemma 1.3**.**
Let Assumptions (H0) hold. Then the functions generated by DDP-REG are Lipschitz continuous on , satisfy , and and are bounded for all .
Proof: It suffices to follow the proof of Lemma 3.2 in [17].
We have for DDP-REG the following convergence theorem which is a special case of Theorem 3.2 shown in section 3 (obtained considering deterministic processes ).
Theorem 1.4**.**
Consider the sequences of decisions and approximate recourse functions generated by DDP-REG. Let Assumptions (H0) hold and assume that for , we have and for every . Then we have ,
[TABLE]
and for ,
[TABLE]
Also, (i) , the optimal value of (1.1), and (ii) any accumulation point of the sequence is an optimal solution of (1.1).
If convergence of DDP-REG holds for any sequence of prox-centers in and of penalty parameters converging to zero for every , the performance of the method depends on how these sequences are chosen. DDP is obtained taking for every .
For all numerical experiments of section 5.2, DDP-REG was much faster than DDP. Some natural candidates for and , used in our numerical tests, are the following:
- •
Weighted average of previous values: with nonnegative weights and . Note that because all are in the convex set . Special cases include the average of previous values and the last trial point for , .
- •
where or for , .
If for a given stage , is a polytope and we do not have the nonlinear constraints given by constraint functions (i.e., the constraints for this stage are linear), then the conclusions of Lemmas 1.1, 1.3, and Theorem 1.4 hold under weaker assumptions. More precisely, for such stages , we assume (H0)-a), (H0)-(b), and instead of (H0)-(d), (H0)-(e), the weaker assumption (H0)-(c’):
(H0)-(c’) There exists such that:
- (c’).1)
- (c’).2)
for every , the set is nonempty.
2. Regularized Stochastic Dual Dynamic programming
2.1. Problem formulation and assumptions
Consider a stochastic process where is a discrete random vector with finite support containing in particular as components the entries in in a given order where are random vectors and are random matrices.
Let denote the sigma-algebra , be the set of -measurable functions, and be the conditional expectation at .
With this notation, we are interested in solving problems of form
[TABLE]
for some functions taking values in , where is given and where
[TABLE]
for some vector-valued function and some nonempty compact convex set .
We make the following assumption on :
- (H1)
is interstage independent and for , is a random vector taking values in with discrete distribution and finite support while is deterministic.
To alleviate notation and without loss of generality, we have assumed that the number of possible realizations of , the size of , and dimension of do not depend on .
Under Assumption (H1), coincides with its unconditional counterpart where is the expectation computed with respect to the distribution of . To ease notation, we will drop the index in . As a result, for problem (2.11), we can write the following dynamic programming equations: we set and for , define
[TABLE]
with
[TABLE]
Problem (2.11) can then be written
[TABLE]
with optimal value denoted by .
Recalling definition (1.3) of the -fattening of a set, we make the following Assumption (H2): setting , there exists such that for :
is nonempty, convex, and compact.
- 2)
For every , the function is proper, convex, and lower semicontinuous.
- 3)
For every , each component of the function is a convex lower semicontinuous function.
- 4)
we have
- 4.1)
for every ,
- 4.2)
for every , for every , the set is nonempty.
- 5)
If , for every , there exists
[TABLE]
such that .
The following proposition, proved in [17], shows that Assumption (H2) guarantees that for , recourse function is convex and Lipschitz continuous on the set for every . SDDP-REG and its convergence analysis are based on this proposition.
Proposition 2.1**.**
Under Assumption (H2), for , for every , the recourse function is convex, finite and Lipschitz continuous on .
Proof: We refer to the proof of Proposition 3.1 in [17] where similar value functions are considered.
Assumption (H2) will also be used to derive explicit formulas for the cuts to be built for recourse functions in SDDP-REG applied to the nonlinear problems we are interested in.
2.2. Algorithm
Recalling Assumption (H1), the distribution of is discrete and the possible realizations of can be organized in a finite scenario tree with the root node associated to a stage [math] (with decision taken at that node) having one child node associated to the first stage (with deterministic). In this section, we describe SDDP-REG algorithm, a regularization of SDDP which can be seen as an extension of the regularization proposed in [4] (which applied to linear problems) to nonlinear problems.
To describe this algorithm, we need some notation: is the set of nodes, Nodes is the set of nodes of the scenario tree for stage and for a node of the tree, we denote by:
- •
the set of its children nodes (the empty set for the leaves);
- •
a decision taken at that node;
- •
the transition probability from the parent node of to ;
- •
the realization of process at node 111Note that to simplify notation, the same notation is used to denote the realization of the process at node Index of the scenario tree and the value of the process for stage Index. The context will allow us to know which concept is being referred to. In particular, letters and will only be used to refer to nodes while will be used to refer to stages.: for a node of stage , this realization contains in particular the realizations of , of , and of ;
- •
is the history of the realizations of process from the first stage node to node : for a node of stage , the -th component of is for , where is the function associating to a node its parent node (the empty set for the root node).
At iteration of the algorithm, trial points are computed for a set of sampled nodes of the scenario tree replacing recourse functions by the approximations available at the beginning of this iteration and penalizing the objective with a quadratic term with prox-center for all the nodes of stage . The nodes selected at iteration are denoted (with , and for , a node of stage , child of node ) and correspond to a sample of . For , a cut
[TABLE]
is computed for at (see the algorithm below for details). To alleviate notation, we will write .
Gathering the cuts computed until iteration , we get at the end of this iteration for the polyhedral lower approximations , given by
[TABLE]
To describe and analyze the algorithm, it is convenient to introduce the function given by
[TABLE]
where
[TABLE]
The SDDP-REG algorithm is given below:
**SDDP-REG (Regularized SDDP).
**
- Step 1)
Initialization. Let , satisfying be given and . Set , .
- Step 2)
Forward pass.
Sample a scenario from the distribution of corresponding to a set of nodes .
For ,
Find an optimal solution of
[TABLE]
where , is any point in and where
is the function given by
[TABLE]
with if or .
End For
- Step 3)
**Backward pass.
**Set and .
For ,
**For **every child node of solve
[TABLE]
and compute, using Proposition 1.2, a subgradient
at .222Note that the proposition can be applied because Assumption (H2) holds and thus the assumptions of the proposition are satisfied for value function .
**End For
** The new cut is obtained computing
[TABLE]
End For
- Step 4)
Do and go to Step 2).
2.3. On the prox-centers and penalizations
Though are now random variables, the remarks of section 1.4 on the choice of the prox-centers for DDP-REG still apply for SDDP-REG. Indeed, convergence of SDDP-REG holds for any sequence of prox-centers in and of penalty parameters converging to zero for every , but the performance of the method depends on how these sequences are chosen. The following choices for and will be used in our numerical tests of SDDP-REG:
- •
Weighted average of previous values: with nonnegative weights and . As a special case, was used in [4].
- •
where (used in [4] with constant) or for , .
3. Convergence analysis of SDDP-REG
The approximate recourse functions available at the end of iteration of SDDP-REG define a policy allowing us to compute decisions for every node of the scenario tree with the following loops: **Simulation of SDDP-REG in the end of iteration .
** For ,
**For **every node of stage ,
**For **every child node of node , compute an optimal solution of
[TABLE]
where .
**End For
** **End For
End For
**
We will assume that the sampling procedure in SDDP-REG satisfies the following property:
- (H3)
for every , for every , and for every , For every , and ,
[TABLE]
The following lemma will be useful in the sequel:
Lemma 3.1**.**
Consider the sequences , and generated by SDDP-REG. Under Assumptions (H1), (H2), then almost surely, for , the following holds:
- (a)
* is convex with on for all ;*
- (b)
the sequences and are bounded;
- (c)
for , is Lipschitz continuous on .
Proof: The proof is similar to the proof of Lemma 3.2 in [17].333In [17] a forward, instead of a forward-backward algorithm, is considered. In this setting, finiteness of coefficients and is not guaranteed for the first iterations (for instance are ) but the proof is similar. We give the main steps of the proof which is by backward induction on starting with where the statement holds by definition of . Assuming for that is Lipschitz continuous on with , then setting , for every we have which gives
[TABLE]
where for the last inequality, we have used Proposition 1.2 which can be applied since Assumption (H2)-5) holds. Therefore, defines a valid cut for and . Assumptions (H2)-1)-4) and finiteness of on imply that and are bounded for every , and allow us to obtain a uniform upper bound on , i.e., a Lipschitz constant valid for all functions .
Theorem 3.2 shows the convergence of the sequence to the optimal value of (2.11) and that any accumulation point of the sequence can be used to define an optimal solution of (2.11).
Theorem 3.2** (Convergence analysis of SDDP-REG).**
Consider the sequences of stochastic decisions and of recourse functions generated by SDDP-REG to solve dynamic programming equations (2.12)-(2.13)-(2.14). Let Assumptions (H1), (H2), and (H3) hold and assume that and that for every , we have . Then
- (i)
almost surely, for , the following holds:
[TABLE]
- (ii)
Almost surely, the limit of the sequence of the approximate first stage optimal values and of the sequence is the optimal value of (2.11). Also, let be any accumulation point of the sequence . Now define with by where is given by . Then is an optimal solution to (2.11).
Proof: See the Appendix.
Remark 3.3** (Extension of SDDP-REG to risk-averse nonlinear problems.).**
Using [17], SDDP-REG algorithm can be extended to nested risk-averse formulations of risk-averse multistage stochastic nonlinear programs of form
[TABLE]
where is a coherent and law invariant conditional risk measure. The convergence proof of this variant of SDDP-REG can be easily obtained combining the convergence proof of risk-averse decomposition methods from [17] with the convergence proof of Theorem 3.2.
Similary to DDP-REG, if for a given stage , is a polytope and we do not have the nonlinear constraints given by constraint functions (i.e., the constraints for this stage are linear), then the conclusions of Proposition 2.1, Lemma 3.1, and Theorem 3.2 hold under weaker assumptions. More precisely, for such stages , we assume (H2)-1), (H2)-2), and instead of (H2)-4), (H2)-5), the weaker assumption (H2)-3’): there exists such that 3.1’) for every , \mathcal{X}_{t-1}^{\varepsilon}{\small{\times}}\mathcal{X}_{t}\subset\mbox{dom}\;f_{t}\Big{(}\cdot,\cdot,\xi_{t,j}\Big{)}; and 3.2’) for every , for every , the set is nonempty.
4. Multistage portfolio optimization models with direct transaction and market impact costs
4.1. Multistage portfolio selection models with direct transaction costs
This section presents risk-neutral and risk-averse multistage portfolio optimization models with direct transaction costs over a discretized horizon of stages. We model the direct transaction costs incurred by selling and purchasing securities as being proportional to the amount of the transaction [9].
Let be the number of risky assets and asset be cash. Next is the dollar value of asset at the end of stage , is the return of asset at , is the amount of asset sold at the end of , is the amount of asset bought at the end of , and are respectively the proportional selling and purchasing transaction costs. Each component , of is a known parameter. The expression is the budget available at the start of the investment planning horizon. The notation is a parameter defining the maximal amount that can be invested in each financial security . To allow for a direct application of SDDP-REG to solve the problem, we will assume that returns are interstage independent. Simple modifications of SDDP-REG could be used to deal with returns following generalized linear models with finite memory (as in [16]) and with Markov chain approximations of the returns (as in [28]).
For , given a portfolio and , we define the set as the set of satisfying
[TABLE]
and for ,
[TABLE]
Constraints (4.24a) define the amount of security held at each stage and take into account the proportional transaction costs. Constraints (4.23) are the cash flow balance constraints and define how much cash is available at each stage. Constraints (4.24b) preclude selling an amount larger than the one held. Constraints (4.24c) do not allow the position in security at time to exceed a specified limit , while (4.24d) prevents short-selling and enforces the non-negativity of the amounts purchased and sold.
Risk-neutral model. With this notation, the dynamic programming equations of a risk-neutral portfolio model of form (2.12), (2.13), (2.14) can be written444It is indeed immediately seen that (4.25)-(4.26) is of form (2.12), (2.13), (2.14), writing the maximization problems as minimization problems and introducing the extended state .: for , setting we solve the problem
[TABLE]
while at stage , we solve
[TABLE]
where for , . With this model, we maximize the expected return of the portfolio taking into account the transaction costs, non-negativity constraints, and bounds imposed on the different securities.
Risk-averse model. As we recall from the previous section, SDDP-REG can be easily extended to solve risk-averse problems of form (3.22). We can therefore define a nested risk-averse counterpart of the risk-neutral portfolio problem we have just introduced and solve it with SDDP-REG. This model is obtained replacing the expectation in the risk-neutral portfolio problem above by the (unconditional, due to Assumption (H1)) risk measure given by
[TABLE]
where , is the confidence level of the Average Value-at-Risk, and is computed with respect to the distribution of . Therefore, a risk-averse portfolio problem with direct transaction costs is written as follows: at stage , setting , we solve
[TABLE]
while at stage , we solve
[TABLE]
where for , .
4.2. Conic quadratic models for multistage portfolio selection with market impact costs
Due to market imperfections, securities can seldom be traded at their current theoretical market price, which leads to additional costs, called market impact costs. If the trade is very large and involves the purchase (resp., selling) of a security, the price of the share may rise (resp., drop) between the placement of the trade and the completion of its execution [30]. As more of a security is bought or sold, the proportional cost increases due to the scarcity effect. Market impact costs are particularly important for large institutional investors, for which they represent a major proportion of the total transaction costs [27, 44]. Often, large trading orders are not executed at once, but are instead split into a sequence of smaller orders executed within a given time window. Taken individually, these small orders exert little or no pressure on the market [45], which can curb market impact costs. The downside is that the execution of the entire trade order is postponed, which may lead to a loss in opportunities caused by (unfavorable) changes in market prices.
The change in a security price is impacted by the size of the transaction and is often modelled as a concave monotonically increasing function of the trade size [2]. In that vein, Lillo et al. [25] and Gabaix et al. [11] model market impact costs as a concave power law function of the transaction size. Bouchaud et al. [6] use a logarithmic function of the transaction size and assert that the market impact is temporary and decays as a power law. Moazeni et al. [29] propose linear market impact costs and evaluate the sensitivity of optimal execution strategies with respect to errors in the estimation of the parameters. Mitchell and Braun [27] study the standard portfolio selection problem in which they incorporate convex transaction costs, including market impact costs, incurred when rebalancing the portfolio. They rescale the budget available after paying transaction costs, which results into a fractional problem that can be reformulated as a convex one. Frino et al. [10] approximate impact costs with a linear regression based on quantized transaction sizes, while Zagst and Kalin [45] use a piecewise linear function. Loeb [26] shows that market impact costs are a function of the square root of the amount traded. Similarly, Torre [44] models the price change as proportional to the square root of the order size. This led to the so-called square-formula which defines the market impact costs as proportional to the square root of the ratio of the number of shares traded to the average daily trading volume of the security [12]. The square-root formula is widely used in the financial practice [12] to provide a pre-trade estimate of market impact costs and is preconized by Andersen et al. [3] as well as by Grinold and Kahn [15]. The latter observe that this approach is consistent with the trading rule-of-thumb according to which it costs roughly one day’s volatility to trade one day’s volume. In Barra’s Market Impact Model Handbook [44], it is showed that the square-root formula fits transaction cost data remarkably well. An empirical study conducted by Almgren et al. [2] advocates to set the price change as proportional to a 3/5 power law of block size.
In this study, the modeling of the market impact costs is based on the square-formula. More precisely, we follow the approach proposed by Grinold and Kahn [15] and Andersen [3], and model the market impact costs as proportional to a 3/2 power law of the transacted amount (see (4.30)).
Let be the volume of security in the considered transaction and be the overall market volume for security at . Additionally, is the monetary value of asset transacted at . The market impact costs for asset are defined as:
[TABLE]
where and are non-negative parameters that must be estimated. The market impact costs capture the fact that the price of an asset increases or decreases if one buys or sells very many shares of this asset. The total market impact costs depend on both the cost per unit and the square root of the amount traded (which is aligned with the empirical tests reported in [26]):
[TABLE]
For , given a portfolio and , we define
[TABLE]
Constraints (4.31)-(a) define how much cash is held at each period and take into account the market impact costs. Constraints (4.31)-(b) define the total amount of security traded at time . The nonlinear constraints (4.31)-(c) follow from (4.30) and permit to define the total market impact costs incurred for security at time . Same as in the previous section, we assume that returns are interstage independent.
The risk-neutral multistage portfolio optimization problem with market impact costs writes as follows: for , setting we solve the problem
[TABLE]
while at stage , we solve
[TABLE]
where for , .
It is easy to see that the left-hand side of the constraints (4.31)-(c) are convex functions, which implies that Assumptions (H2)-3) are satisfied and SDDP-REG can be applied to solve the portfolio problem under consideration. For implementation purposes, it is convenient to rewrite constraints (4.31)-(c) as a conic quadratic constraint:
Theorem 4.1**.**
For , the convex feasible sets
[TABLE]
can be equivalently represented with the rotated second-order constraints (4.34a) and the linear constraints (4.34b)-(4.34d):
[TABLE]
Proof: This representation is proved in [3].
For , given and , denoting
[TABLE]
and using Theorem 4.1, our portfolio optimization problem with market impact costs (4.32)-(4.33) can be rewritten substituting in (4.32)-(4.33) for , the constraints by
[TABLE]
This formulation of the portfolio problem can be solved using SDDP-REG with all subproblems of the forward and backward passes being conic quadratic optimization problems.
5. Numerical experiments
In this section, we evaluate the computational efficiency of the DDP-REG and SDDP-REG algorithms presented in sections 1 and 2, and benchmark them with standard, non regularized versions of the deterministic and stochastic DDP algorithms. The analysis starts (section 5.2) with the deterministic setting and the DDP-REG algorithm tested on a portfolio optimization problem with direct transaction costs, and continues in section 5.3 with the stochastic case and the SDDP-REG algorithm tested on risk-neutral and risk-averse formulations involving either direct transaction or market impact costs. In practice, portfolio selection problem parameters (the returns) are not known in advance and stochastic optimization models are used for these applications. We use such models in section 5.3. However, to compare DDP and DDP-REG, we assume that the parameters of the portfolio problems, namely the returns, are known over the optimization period. This allows us to easily generate feasible problem instances that can be solved with DDP and DDP-REG and to know what would have been the best return for these instances.
5.1. Data and parameter settings
The problem instances and the algorithms are modelled in Python and the problems are solved with MOSEK 8.0.0.50 solver [31]. The experiments are carried out using a single thread of an Intel(R) Core(TM) i5-4200M CPU @ 2.50GHz machine.
The following settings are used for the parameters of the portfolio optimization problems described in section 4. The budget available is n=6\eta=\nuare set to 1%. The return data of six securities were collected from WRDS [[39](#bib.bib39)] for the period ranging from July 2005 to May 2016. The monthly fixed cash return is equal to0.2%u_{i}=20%T\Theta_{t}(using the notation of section [2](#S2)), is set toM=60$.
As we recall from section 1.4 for DDP-REG and from subsection 2.3 for SDDP-REG, we need to define sequences of prox-centers and of penalization parameters to define instances of DDP-REG and SDDP-REG. In our study, we will use the prox-centers and penalization parameters given in Table 1 (we recall that no penalization is used for and for , i.e., for all ). This table also contains the names used for the corresponding DDP-REG and SDDP-REG variants.
We recall that in [4], only the variant SDDP-REG-PREV-REG1- was tested for linear programs. In this section, we test all deterministic variants from Table 1 for linear programs and variant SDDP-REG-PREV-REG2 for multistage stochastic linear and nonlinear programs.
5.2. Deterministic instances
In this section, we consider the deterministic counterpart of the portfolio optimization problem with direct transaction costs presented in section 4.1 using the parameters given in the previous section and 8 different values for the number of time periods: , and . We solve these problems using DDP and the following 6 variants of DDP-REG (using the notation of Table 1): DDP-REG-PREV-REG1-0.2 (DDP-REG-PREV-REG1- with ), DDP-REG-PREV-REG1-0.9 (DDP-REG-PREV-REG1- with ), DDP-REG-PREV-REG2, DDP-REG-AVG-REG1-0.2 (DDP-REG-AVG-REG1- with ), DDP-REG-AVG-REG1-0.9 (DDP-REG-AVG-REG1- with ), and DDP-REG-AVG-REG2.
Stopping criterion. When studying the convergence of DDP-REG in section 1, we have not discussed the stopping criterion. At each iteration, this algorithm can compute an approximate lower bound on the optimal value of the problem which is given at iteration by (using the notation of section 1), the optimal value of the approximate problem for the first time period. Observe that we can make of an exact lower bound if we take (such strategy was used in our tests). DDP-REG can also compute at iteration the upper bound on the optimal value. Given a tolerance (taken equal to in our experiments), the algorithm stops when the difference between the upper and lower bound is less than (in this case, we have computed an -optimal solution to the problem). Note, however, that since our portfolio problems are maximization problems, the approximate first stage problem provides an upper bound on the optimal value and provides a lower bound.
We have checked that on all instances, all algorithms correctly compute the same optimal value and that the upper and lower bounds were converging to this optimal value. For illustration, Figure 1 displays the evolution of the upper and lower bounds and of the optimality gap across the iterative process with DDP for the instance with .
The CPU time needed to solve the different instances with DDP and our 6 variants of DDP-REG is given in Table 2 and the corresponding reduction factor in CPU time for these DDP-REG variants is given in Table 3. The number of iterations of the algorithms is given in Table 4. We observe that on all instances DDP-REG variants are much faster and need many fewer iterations than DDP. Most importantly, the benefits of regularization increase as the problem gets larger and the number of stages raises. When is large there is a drastic improvement in CPU time with DDP-REG variants. For instance, for , and , the reduction factor in CPU time varies (among the 6 DDP-REG variants) respectively in the interval , , and . Remarkably, the solution time with the regularized algorithm DDP-REG is not monotonically increasing with the number of stages, which points out to the scalability of the algorithm and the possibility to use it for even larger problems. As an illustration, the difference in time and number of iterations between DDP and DDP-REG-PREV-REG2 is shown in Figure 2, which highlights that the time and iteration differential increase with the number of stages.
5.3. Stochastic instances
In this section, we evaluate the computational efficiency of the SDDP-REG algorithm presented in section 2, and benchmark it with the standard, non regularized version of the SDDP algorithm. We have implemented the regularization scheme SDDP-REG-PREV-REG2 given that penalization scheme REG2 performed best for the deterministic instances (see section 5.2). The algorithms are tested on three types of problem instances with periods: risk-neutral portfolio models of subsection 4.1, risk-averse portfolio models of subsection 4.1, and risk-neutral porfolio model with market impact costs from subsection 4.2.
Stopping criterion. For risk-neutral SDDP, we used the following stopping criterion. The algorithm stops if the gap is . The gap is defined as where and correspond to upper and lower bounds, respectively. The upper bound corresponds to the optimal value of the first stage problem (recall that we have a maximization problem), obtained taking, as for DDP-REG, (if , we get a sequence of approximate upper bounds, which, as we have seen, converges almost surely to the optimal value of the problem). The lower bound corresponds to the lower end of a 95%-one-sided confidence interval on the optimal value for policy realizations, see [40] for a detailed discussion on this stopping criterion. Risk averse SDDP was terminated after a fixed number of iterations ().
5.3.1. Risk-neutral multistage linear problem with direct transaction costs (4.25)-(4.26)
We report in Table 5 the computational time and number of iterations required for SDDP and SDDP-REG-PREV-REG2 to solve the instance of portfolio problem (4.27)-(4.28) obtained taking and the problem parameters given in subsection 5.1. We observe that as in the deterministic case, the regularized decomposition method converges much faster (it is about twice as fast for ) and requires many fewer iterations. We also refer to Figure 3 where the evolution of the upper and lower bounds and the gap (in % of the upper bound) are represented for SDDP and SDDP-REG-PREV-REG2 for . We see that the gap decreases much faster with SDDP-REG-PREV-REG2.
5.3.2. Risk-averse multistage linear problem with direct transaction costs (4.27)-(4.28)
We implemented risk-averse models (4.27)-(4.28) taking and , running the algorithms for 50 iterations. The CPU time is reported in Table 6. Since both problems are run for the same number of iterations and since the regularized variant requires solving quadratic problems instead of just linear programs in the forward passes, it was expected to have a larger computational time with the regularized variant. However, the difference is small. We also report in Figure 4 the evolution of the upper bounds for SDDP and SDDP-REG-PREV-REG2. We see again that the upper bound decreases much faster with SDDP-REG-PREV-REG2.
5.3.3. Conic risk-neutral multistage stochastic problem with market impact costs from section 4.2
We consider two variants of the portfolio problem with market impact costs given in section 4.2 in which we set the market impact unit cost , to respectively 3 basis points (we recall that a basis point is ) for the first model and for the second. The CPU time and number of iterations to solve these problems with SDDP and SDDP-REG-PREV-REG2 are given in Table 7. The evolution of the upper and lower bounds and of the gap along the iterations of the algorithms are reported in Figures 5 and 6 for . We observe that when are small the regularized variant is much quicker and the gap decreases much faster. When increases, in particular for the value , more money is invested in cash and the computational time and gap evolution with the non-regularized and regularized variants of SDDP are similar.
6. Conclusion
We presented and studied regularized variants of DDP and SDDP which are extensions of [4] to nonlinear problems and tested several prox-centers. On the one hand, for deterministic problems, the important reduction in CPU time when passing from DDP to DDP-REG, in the vein of [24], was expected. In the stochastic case, we would a priori need different prox-centers for all nodes of the scenario tree (see also [38]). However, such regularized variant is not computationally tractable. Therefore, the proposed SDDP-REG offers a tractable regularized variant of SDDP whose convergence can be shown for vanishing penalties and which can converge quicker than SDDP on some problem instances and for some choices of prox-centers as shown in our experiments. We also observe that it is possible, as in [4], to partition the decision for stage into state and control variables and to take as prox-centers for stage and iteration the state component of . The convergence of both this variant of SDDP-REG as well as the variant that uses prox-centers attached to nodes of the scenario tree mentioned above can be shown following the steps of the convergence proof of SDDP-REG given in section 3.
An interesting topic on the regularization of SDDP is to define a regularization that can be proved to have better complexity than standard SDDP.
Acknowledgments
The first author’s research was partially supported by an FGV grant, CNPq grant 307287/2013-0, FAPERJ grants E-26/110.313/2014, and E-26/201.599/2014. This research was initiated during the visit of the second author at FGV. The second author thanks Dr. Guigues for the support provided through FAPERJ grant E-26/201.599/2014.
Appendix
Proof of Theorem 3.2. In this proof, all equalities and inequalities hold almost surely. We show , , , by induction backwards in time. follows from the fact that . Now assume that holds for some . We want to show that holds. Take a node . Let be the set of iterations such that the sampled scenario passes through node . Due to Assumption (H3), the set is infinite. We first show that
[TABLE]
Take . We have and recalling (2.15), we have . Using definition (2.20) of , it follows that
[TABLE]
Now let such that where is defined by (2.16) with replaced by . Using (6.37) and the definition of , we get
[TABLE]
Now using the definitions of and we obtain
[TABLE]
Observing that for every the decision , we obtain, using definition (2.13) of , that
[TABLE]
Combining this relation with (6.39) gives for
[TABLE]
Next,
[TABLE]
where the above inequality comes from the fact , i.e., is feasible for optimization problem (3.21) with objective function and optimal solution . We get
[TABLE]
where is the diameter of (finite, since is compact). Plugging (6.42) and (6.40) into (6.38) yields for any
[TABLE]
Using the induction hypothesis , we have for every child node of node :
[TABLE]
Now recall that is convex on the compact set (Proposition 2.1), for every child node of node , and the functions , are Lipschitz continuous with on (Lemma 3.1). It follows that we can use Lemma A.1 in [13] to deduce from (6.44) that for every
[TABLE]
Combining this relation with (6.43) and using the fact that , we obtain
[TABLE]
To show , it remains to show that
[TABLE]
Relation (6.46) can be shown following the end of the proof of Theorem 4.1 in [17], by contradiction and using the Strong Law of Large Numbers (the same arguments were first used in a similar context in Theorem 3.1 of [13] but for a different problem formulation and sampling scheme). The key to the proof being the fact that the sampled nodes for iteration are independent on the decisions computed at the nodes of the scenario tree for that iteration and on recourse functions . This achieves the proof of (i).
(ii) Recalling that the root node with decision taken at that node has a single child node with corresponding decision at iteration , the computations in (i) show that for every 555Though when deriving these relations in (i) we had fixed , the inequalities we now re-use for (ii) are valid for any ., we have
[TABLE]
We have shown in (i) that . Plugging this relation into (6.47) shows that
[TABLE]
Now take an accumulation point of the sequence and let be an infinite set of iterations such that for every , .666The existence of an accumulation point comes from the fact that the decisions belong almost surely to a compact set. Using once again computations from (i), we get for any , , , ,
[TABLE]
which can be written
[TABLE]
Since (due to (i)), the above relation shows that
[TABLE]
We will now use the continuity of which follows from (H2) (see Lemma 3.2 in [17] for a proof). We have
[TABLE]
where for the last inequality we have used the continuity of . To prove (ii) it suffices to observe that the sequence belongs to the set
[TABLE]
and this set is closed since is lower semicontinuous and is closed. Thus, , which, together with (6.49), shows that is an optimal solution of and achieves the proof of (ii).
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1[1] R. Almgren , Optimal execution with nonlinear impact functions and trading-enhanced risk , Applied Mathematical Finance, 10 (2003), pp. 1–18.
- 2[2] R. Almgren, C. Thum, and H. Li , Equity market impact , Risk, 18 (2005), pp. 57–62.
- 3[3] E. Andersen, J. Dahl, and H. Friberg , Markowitz portfolio optimization using mosek , 2012.
- 4[4] T. Asamov and W. Powell , Regularized decomposition of high-dimensional multistage stochastic programs with markov uncertainty , Accepted in SIAM Journal on Optimization, (2015).
- 5[5] M. Bandarra and V. Guigues , Single cut and multicut SDDP with cut selection for multistage stochastic linear programs: convergence proof and numerical experiments , Available at https://arxiv.org/abs/1902.06757 , (2019).
- 6[6] J. Bouchaud, Y. Gefen, M. Potters, and M. Wyart , Fluctuations and response in financial markets: The subtle nature of’random’ price changes , Quantitative Finance, 4 (2004), pp. 176–190.
- 7[7] A. Cadenillas , Consumption-investment problems with transaction costs: Survey and open problems , Mathematical Methods of Operations Research, 51 (2000), pp. 43–68.
- 8[8] V. de Matos, A. Philpott, and E. Finardi , Improving the performance of stochastic dual dynamic programming , Journal of Computational and Applied Mathematics, 290 (2015), pp. 196–208.