Anton Friedl, Jiří J. Klemeš, Stefan Radl, Petar S. Varbanov, Thomas Wallek (Eds.) Proceedings of the 28th European Symposium on Computer Aided Process Engineering June 10th to 13th, 2018, Graz, Austria. © 2018 Elsevier B.V. All rights reserved. https://doi.org/10.1016/B978-0-444-64235-6.50092-9
A Relaxed Knapsack-Problem Based Decomposition Heuristic for Large-Scale Multistage Stochastic Programs Zuo Zeng, Selen Cremaschi* Department of Chemical Engineering, Auburn University, Auburn, AL 36849, USA
[email protected]
Abstract Most optimization problems with endogenous uncertainty can be modelled as multistage stochastic programs (MSSPs). The MSSPs grow exponentially with the increase in the numbers of scenarios and time periods, and quickly become computationally intractable for real-world sized problems. This paper presents a Relaxed Knapsack-problem based decomposition Algorithm (RKDA) to efficiently generate tight dual bounds for large-scale MSSPs. The algorithm builds upon the Knapsack-problem based Decomposition Algorithm (KDA) (Christian and Cremaschi, 2015), which generates feasible solutions for MSSPs with endogenous uncertainty by decomposing the original multi-period MSSP into a series of knapsack problems and by solving these problems at appropriate decision points of the planning horizon. The new approach, RKDA, employs KDA to solve a relaxed version of the original MSSP. The relaxed MSSP is obtained by removing its resource constraints. We applied KDA and RKDA to generate initial primal and dual bounds for the pharmaceutical clinical trial planning problem with two-, three-, four-, five-, six-, seven- and ten-products and three clinical trials with different planning horizon lengths. The results revealed that the initial relative gaps ranged from 2.10% to 4.66%, and KDA and RKDA yielded these gaps in 0.1 to 2377.0 CPU seconds. Keywords: stochastic programming, endogenous uncertainty, knapsack-problem based decomposition heuristic, primal and dual bounds
1. Introduction Optimization problems with decision-dependent, i.e., endogenous, uncertainty are commonly observed in process industry, e.g., R&D pipeline management (Colvin and Maravelias, 2010), synthesis of process networks with uncertain process yields (Tarhan and Grossman, 2008), and artificial lift infrastructure planning (Zeng and Cremaschi, 2017). In these problems, decisions may impact the resolution time (Type II) and/or the probability distributions of the uncertain parameters (Type I) (Apap and Grossmann, 2017). This paper restricts its scope to optimization problems with Type II uncertainty. One approach to address these problems is multistage stochastic programming, which is a scenario-based method that considers decisions and recourse actions in multiple stages as uncertainty is realized. Each scenario defines a unique set of realizations for all uncertain parameters. In MSSP formulations of problems with endogenous uncertainty, decision variables are defined independently for each scenario. To prevent decision variable values that anticipate future outcomes of uncertain parameters, non-anticipativity constraints (NACs) are explicitly incorporated to MSSP formulations. Size of decision variables and
520
Z. Zeng and S. Cremaschi
NACs grow exponentially as the number of uncertain parameters and the number of outcomes for each uncertain parameter increase. This growth results in both space and time complexities for solving MSSPs for real-world sized problems. In general, the solution approaches for large-scale MSSPs with endogenous uncertainty rely on heuristic, approximation, and decomposition methods, such as sample average approximation algorithm (Solak et al., 2010), improved Lagrangean decomposition framework (Gupta and Grossmann, 2014), sequential scenario decomposition approach (Apap and Grossmann, 2017), a branch and bound algorithm (Christian and Cremaschi, 2017) and a rolling-horizon heuristic approach (Colvin and Maravelias, 2009). It has been shown that moderate-size problems can be solved to optimality. Here, we presents a Relaxed Knapsack-problem based Decomposition Algorithm (RKDA) to generate tight dual bounds for large-scale MSSPs with endogenous uncertainty. The algorithm builds upon the KDA (Christian and Cremaschi, 2014). Instead of enumerating all uncertain parameter space as scenarios, the KDA decomposes the original multi-period MSSP into a series of knapsack problems and solves these problems at appropriate decision points of the planning horizon. It provides tight primal bounds for the original problem, and generates these bounds several orders of magnitude faster than solving the original MSSP. The RKDA generates a tight dual bound by employing a relaxed version of the KDA to solve a relaxed version of the original MSSP. The relaxed MSSP is obtained by removing its resource constraints. The application of RKDA is demonstrated by solving instances of clinical trial planning problem, which is introduced in Section 2. Then, we describe KDA and RKDA (Section 3). In Section 4, we apply KDA and RKDA to different instances of the clinical trial planning problem with scenarios ranging from 64 up to 1,048,576. Finally, the concluding remarks are summarized in Section 5.
2. Clinical Trial Planning Problem and its MSSP formulation The goal of the clinical trial planning problem is to identify the schedule of clinical trials for a set of candidate drugs, i∈I={1,2,…,|I|}, that maximizes the expected net present value (ENPV) of the clinical trial pipeline. The outcomes of the clinical trials for each drug are uncertain and decision-dependent, because the outcome of a clinical trial can only be realized once the trial is completed. A drug has to complete three clinical trials in order (j∈J={PI,PII,PIII}), and it can either pass (P) a clinical trial or fail (F) it. A drug that has failed a trial cannot continue to the subsequent clinical trial(s). Let θi be the uncertain parameter associated with this uncertainty for drug i. We then define the realizable values of this uncertain parameter with set Θi={PI(F),PII(F),PIII(F),PIII(P)}, where PIII(P) corresponds to the outcome that drug i passes trials PI,PII, and PIII. The scenarios for the MSSP model are constructed as Cartesian product of uncertain parameter outcomes. Therefore, the total number of scenarios for a clinical trial planning problem with |I| drugs each of which should complete all three clinical trials is |S| = 4|I|, where S is the scenario set. The planning horizon is defined as t∈T={1,2,…,|T|} (period t starts at time t-1 and ends at time t). Each drug has a known cost Cij, required resources ρijr (where r∈R={1,2,…,|R|}), and fixed duration of τij for each clinical trial. After successfully completing all clinical trials, the revenue from drug i is realized. Potential revenue is defined using the maximum revenue, revimax, associated with drug i. Patent life of a drug is assumed to starts shrinking
A Relaxed Knapsack-Prob Based Decomposition Heuristic
521
once the drug enters the pipeline, and associated losses are represented by two penalty terms: γiD (loss of market share) and γiL (loss of patent life). A slightly revised version of the MSSP formulation by Colvin and Maravelias (2010) is given in Figure 1. The decision variable Xi,j,t,s is equal to one if clinical trial j of drug i is started at time t in scenario s. Equations 2-5 calculate the ENPV. Equations 6-8 are sequencing constraints. Per Colvin and Maravelias (2010), a clinical trial j’ of drug i cannot be started at time t in scenario s, (1) if drug i fails a prerequisite trial (j
1. For Eq. 11, we define the subset B of S×S as scenarios s and s’, which are distinguishable in the outcome of one (drug, clinical trial) pair (is,s’,js,s’).
𝑖𝑖 𝑠𝑠,𝑠𝑠′ ,𝑗𝑗 𝑠𝑠,𝑠𝑠′
′,𝑠𝑠
𝑋𝑋𝑖𝑖,𝑗𝑗 ,
,𝑠𝑠
= 𝑋𝑋𝑖𝑖,𝑗𝑗 ,
,𝑠𝑠 ′
≤ 𝜌𝜌𝑟𝑟
�∨ �
𝑚𝑚𝑚𝑚𝑚𝑚
𝑡𝑡
𝑌𝑌𝑖𝑖 𝑠𝑠,𝑠𝑠′ ,𝑗𝑗 𝑠𝑠,𝑠𝑠′ ,
𝑡𝑡
�
′ ≤𝑝𝑝−𝜏𝜏
𝑡𝑡
¬
𝑡𝑡
�
𝑡𝑡
′ ,𝑠𝑠
∀𝑖𝑖, 𝑗𝑗, 𝑠𝑠 𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝜌𝜌𝑖𝑖,𝑗𝑗 ,𝑟𝑟 𝑋𝑋𝑖𝑖,𝑗𝑗 ,
𝑡𝑡
𝑡𝑡
𝑋𝑋𝑖𝑖,𝑗𝑗 −1,
𝑡𝑡
𝑡𝑡
𝑡𝑡 𝑡𝑡
𝑡𝑡
𝑡𝑡 𝑡𝑡
′ −𝜏𝜏 𝑡𝑡
𝑡𝑡 𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡 𝑡𝑡
−𝜏𝜏 𝑖𝑖,𝑗𝑗
𝑋𝑋𝑖𝑖,𝑗𝑗 ,1,1 = 𝑋𝑋𝑖𝑖,𝑗𝑗 ,1,𝑠𝑠
𝑖𝑖,𝑗𝑗 −1
𝑖𝑖,𝑗𝑗 −1 ,𝑠𝑠
+
−
(3)
(4)
(5) (6)
∀𝑖𝑖, 𝑗𝑗 > 1, ′, 𝑠𝑠
𝑡𝑡
∑𝑖𝑖 ∑𝑗𝑗 ∑
′≤ ′>
′ >𝜏𝜏
∀𝑠𝑠
(7)
′≤
(8)
∀𝑟𝑟, , 𝑠𝑠
(9)
�
−𝜏𝜏 𝑖𝑖 𝑠𝑠,𝑠𝑠′ ,𝑗𝑗 𝑠𝑠,𝑠𝑠′
(10)
𝑌𝑌𝑖𝑖 𝑠𝑠,𝑠𝑠′ ,𝑗𝑗 𝑠𝑠,𝑠𝑠′ ,
′,𝑠𝑠 �
∀𝑖𝑖, 𝑗𝑗, ≥ 1, (𝑠𝑠, 𝑠𝑠 ′ ) ∈ 𝐁𝐁 (11) 𝑡𝑡
−𝜏𝜏
′ ≤|𝑇𝑇|
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
∀𝑠𝑠
∑ ′ 𝑋𝑋𝑖𝑖,𝑗𝑗 , ′ ,𝑠𝑠 ≤ ∑ ′ 𝑖𝑖,𝑗𝑗 −1 𝑋𝑋𝑖𝑖,𝑗𝑗 −1, ′ ,𝑠𝑠 𝑋𝑋𝑖𝑖,𝑗𝑗 , ,𝑠𝑠 = 0 ∀(𝑖𝑖, 𝑗𝑗, , 𝑠𝑠) ∈ 𝐍𝐍
) − 𝛾𝛾𝑖𝑖 𝐿𝐿 � + 𝜏𝜏𝑖𝑖,𝑃𝑃𝑃𝑃𝑃𝑃𝑃𝑃 �𝑋𝑋𝑖𝑖,𝑃𝑃𝑃𝑃𝑃𝑃𝑃𝑃, ,𝑠𝑠 �
′ 𝑜𝑜𝑜𝑜𝑜𝑜𝑜𝑜 ∑ ′ ≤|𝑇𝑇| 𝑋𝑋𝑖𝑖,𝑃𝑃𝑃𝑃, ′ ,𝑠𝑠 � 𝑖𝑖,𝑃𝑃𝑃𝑃 𝑓𝑓𝑖𝑖,𝑃𝑃𝑃𝑃 �1 − 𝑟𝑟𝑟𝑟 𝑟𝑟 𝑖𝑖,𝑗𝑗 ∀𝑠𝑠 , 𝑓𝑓𝑖𝑖,𝑗𝑗 +1 𝑋𝑋𝑖𝑖,𝑗𝑗 , ,𝑠𝑠
>|𝑇𝑇|−𝜏𝜏 𝑖𝑖,𝑗𝑗
,𝑠𝑠
′ ,𝑠𝑠
𝑡𝑡
∑𝑖𝑖 ∑𝑗𝑗 ∈{𝑃𝑃𝑃𝑃,𝑃𝑃𝑃𝑃𝑃𝑃} ∑
𝑡𝑡
�−𝑋𝑋𝑖𝑖,𝑗𝑗 ,1,𝑠𝑠 + ∑
� + ∑𝑖𝑖 𝑟𝑟
𝐶𝐶𝐶𝐶𝐶𝐶𝑠𝑠 = ∑𝑖𝑖,𝑗𝑗 , 𝑐𝑐𝑐𝑐 𝐶𝐶𝑖𝑖,𝑗𝑗 𝑋𝑋𝑖𝑖,𝑗𝑗 , ∑ 𝑋𝑋𝑖𝑖,𝑗𝑗 , ,𝑠𝑠 ≤ 1 ∀𝑖𝑖, 𝑗𝑗, 𝑠𝑠
𝑋𝑋𝑖𝑖,𝑗𝑗 ,
𝑡𝑡
′ ,𝑠𝑠
′≤ ′
𝑡𝑡
𝑋𝑋𝑖𝑖,𝑗𝑗 ,
𝑜𝑜𝑜𝑜𝑜𝑜𝑜𝑜 𝑓𝑓𝑖𝑖,𝑗𝑗 𝑖𝑖,𝑗𝑗
−∑
𝑟𝑟
′
𝑖𝑖,𝑗𝑗 −1 ,𝑠𝑠
− 𝛾𝛾𝑖𝑖 𝐷𝐷 ∑𝑗𝑗 =𝑃𝑃𝑃𝑃𝑃𝑃,𝑃𝑃𝑃𝑃𝑃𝑃𝑃𝑃 (−𝑋𝑋𝑖𝑖,𝑗𝑗 ,1,𝑠𝑠 +
𝑡𝑡
𝐹𝐹𝐹𝐹𝐹𝐹
𝑡𝑡
𝑡𝑡
∑
′ ≤|𝑇𝑇|
′ −𝜏𝜏
,𝑠𝑠
𝑟𝑟𝑟𝑟
= ∑𝑖𝑖 ∑𝑗𝑗 ≠𝑃𝑃𝑃𝑃 𝑟𝑟
𝑠𝑠
𝑋𝑋𝑖𝑖,𝑃𝑃𝑃𝑃𝑃𝑃𝑃𝑃,
𝑋𝑋𝑖𝑖,𝑗𝑗 −1, 𝑟𝑟𝑟𝑟
𝐹𝐹
𝑖𝑖,𝑗𝑗 −1
𝑚𝑚𝑚𝑚𝑚𝑚
𝑡𝑡
𝑡𝑡
𝑡𝑡
𝑡𝑡
∑
𝑖𝑖
𝑡𝑡
′≤ ′ >𝜏𝜏
𝑟𝑟𝑟𝑟
𝑡𝑡
𝑟𝑟𝑟𝑟
𝑅𝑅𝑅𝑅𝑅𝑅𝑠𝑠 = ∑𝑖𝑖 ∑ �𝑟𝑟
(1) (2)
𝑠𝑠 − 𝐶𝐶𝐶𝐶𝐶𝐶𝑠𝑠 )
𝑡𝑡
𝐹𝐹𝐹𝐹𝐹𝐹
max 𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸 s.t. 𝐸𝐸𝐸𝐸𝐸𝐸𝐸𝐸 = ∑𝑠𝑠 𝑝𝑝𝑠𝑠 (𝑅𝑅𝑅𝑅𝑅𝑅𝑠𝑠 + 𝐹𝐹
Figure 1. An MSSP formulation for the clinical trial planning problem
3. Proposed Heuristic Algorithms 3.1. Knapsack-problem based decomposition algorithm (KDA) The KDA generates and solves a series of |R|-dimensional 0-1 knapsack problems based on the realizations of uncertainty. The KDA uses the concept of eligible item set K, where items k∈K are created by enumerating all allowable decisions associated with endogenous uncertain parameters at time period t according to the logic of the sequencing constraints. Each item has a value, Vk, and weight(s), Wkr, associated with it. Item values and weight(s) are estimated using expected gains and resource requirements, respectively.
Z. Zeng and S. Cremaschi
522
Suppose that a clinical trial planning problem only includes two drugs (A and B). The KDA starts by generating the eligible item set at t=1, the root node. The eligible items are (A, PI) and (B, PI), and they are determined using the logic of sequencing constraints (Eqns. 6-8). Next, item values and weights are calculated. The value for (A, PI) is the expected value of its revenue assuming all remaining clinical trials are completed continuously, and it is calculated via Eqns. 12 and 13. Note that t=1 for Eq. 12. The item values approximate the impact of starting the drug, clinical trial pair at the current time period would have on objective function value of the original MSSP (Eqns. 2-5). Equation 12 reduces the potential revenue, revA, due to loss of active patent life for the total duration of the current and remaining clinical trials. The potential revenue is calculated by deducting the linearly depreciated clinical trial costs CA,j’ from the maximum revenue revAmax (Eq. 13). Item weights for (A, PI) are equal to ρA,PI,r. The |R|dimensional capacity vector {Wrmax} contains the available resources for investment, ρrmax. 𝑉𝑉(𝐴𝐴,𝑃𝑃𝑃𝑃) = 𝐸𝐸�𝑟𝑟𝑟𝑟𝑟𝑟𝐴𝐴 − 𝛾𝛾𝐴𝐴 𝐿𝐿 �𝑡𝑡 + ∑𝑗𝑗 ′ ≥𝑃𝑃𝑃𝑃 𝜏𝜏𝐴𝐴,𝑗𝑗 ′ + 1�� 𝑟𝑟𝑟𝑟𝑟𝑟𝐴𝐴 =
𝑟𝑟𝑟𝑟𝑟𝑟𝐴𝐴𝑚𝑚𝑚𝑚𝑚𝑚
− ∑𝑗𝑗 ′ ≥𝑃𝑃𝑃𝑃 𝐶𝐶𝐴𝐴,𝑗𝑗 ′ �1 −
𝑗𝑗 >𝑗𝑗 ′′ 0.025 ∑𝑗𝑗 ′′ >𝑃𝑃𝑃𝑃 𝜏𝜏𝐴𝐴,𝑗𝑗 ′′ −1 �
(12) (13)
The knapsack problems that are generated and solved are given in Eqns. 14-16. The objective is to maximize the total value (Eq. 14), and the weight constraints (Eq. 15) are equivalent to the resource constraints of the original MSSP (Eq. 9). max ∑𝑘𝑘∈𝐸𝐸 𝑉𝑉𝑘𝑘 𝑦𝑦𝑘𝑘 s.t. 𝑊𝑊𝑘𝑘,𝑟𝑟 𝑦𝑦𝑘𝑘 ≤ 𝑊𝑊𝑟𝑟𝑚𝑚𝑚𝑚𝑚𝑚 ∀𝑟𝑟 ∈ 𝐑𝐑 𝑦𝑦𝑘𝑘 ∈ {0, 1} ∀𝑘𝑘 ∈ 𝐊𝐊
(14) (15)
(16) Once the root node knapsack problem is solved, the number of new knapsack problems generated and their generation times are based on the outcomes of the root node knapsack problem solution. For our example, if the solution selects (A, PI) and (B, PI), time period is set to max(τA,PI,τB,PI) and 22=4 additional knapsack problems are generated. These problems correspond to the four possible outcomes: (1) drugs A and B pass clinical trial PI, (2) drug A fails PI while drug B passes it, (3) drug A passes trial PI while drug B fails it, and (4) both drugs fail trial PI. For each knapsack problem, the algorithm updates set K and calculates corresponding item values and weights at the corresponding time period. The item values are updated based on the realizations of uncertainties that correspond to each individual knapsack problem. The algorithm continues to solve and generate new child knapsack problems until the end of the planning horizon or no eligible items remain. The KDA yields a feasible solution, and hence, generates a primal bound for the MSSP. Detailed discussions of KDA can be found in Christian and Cremaschi (2015, 2017).
3.2. Relaxed knapsack-problem based decomposition approach (RKDA) The RKDA generates a valid dual bound for the original MSSP by solving a relaxed version of the original MSSP to optimality. The original MSSP is relaxed by removing the resource constraints (Eq. 9). We refer to this relaxed MSSP formulation as RMSSP. The optimal solution of RMSSP provides a dual bound for the original MSSP. The weight constraints of the knapsack problems in the KDA are equivalent to the resource constraints of the MSSP. Therefore, to solve the RMSSP, the weight
A Relaxed Knapsack-Prob Based Decomposition Heuristic
523
constraints of the knapsack problems are removed yielding the relaxed knapsack problems given in Eqns. 17-18. The original KDA is modified to generate and solve the relaxed knapsack problems every time period at which the outcome of an uncertain parameter is realized. We refer to this modified KDA as relaxed KDA (RKDA). The solution obtained by RKDA to the RMSSP satisfies the sequencing constraints (Eqns. 68) and the NACs (Eqns. 10-11) of the RMSSP, and recommends starting all eligible clinical trials at the earliest time allowable. Therefore, the RKDA yields the optimal solution of the RMSSP, which is a tight dual bound for the original MSSP. (17) (18)
𝑉𝑉
max ∑𝑘𝑘∈𝐸𝐸 𝑘𝑘 𝑦𝑦𝑘𝑘 s.t. 𝑦𝑦𝑘𝑘 ∈ {0, 1}
∀𝑘𝑘 ∈ 𝐊𝐊
4. Results and Discussion We applied KDA and RKDA to generate initial primal and dual bounds for the clinical trial planning problem with two-, three-, four-, five-, six-, seven- and ten-products and three clinical trials with different planning horizons. We also solved the original MSSPs of these problems to 0.1% optimality gap when possible. The models and algorithms were implemented in Pyomo and solved using CPLEX 12.6.3 on a standard node of Auburn University Hopper Cluster. Problem data are available upon request. Number of trials, time periods, and scenarios for the problems are summarized in Table 1 along with the optimum ENPVs and the MSSP solution times for the problems that were solved to optimality. As expected, the MSSP solution times grow exponentially with the number scenarios. For example, CPLEX 12.6.3 required more than 30 hours to solve the six-product problem. The MSSP enumerates all scenarios, and quickly becomes computationally intractable due to its space and time complexity. Therefore, the MSSP formulations of seven- and ten-product problems cannot be generated in RAM and hence cannot be solved using CPLEX 12.6.3 due to their space complexities. Table 1. Main characteristics and the MSSP solutions of the problems
Case Name 2-prod 3-prod 4-prod 5-prod 6-prod 7-prod 10-prod
Number of Trials 2 3 3 3 3 3 3
Number of Time Periods 5 12 6 6 6 8 10
Number of Scenarios 9 64 256 1,024 4,096 16,384 1,048,576
ENPV ($M) 1104 1189 1696 2082 2450 ---
MSSP Solution Time (CPUs) 0.1 4 7 185 29,357 ---
Primal and dual bounds obtained for the problems using KDA and RKDA are compiled in Table 2. An initial relative gap for each problem is calculated using the primal and dual bounds. The gaps ranged between 2.10% and 4.66% (Table 2). The KDA and RKDA yielded these relative gaps very quickly, in 0.1 to 2377 CPU seconds. The Pearson correlation coefficient between the total solution time and number of scenarios is 1.00, revealing the strong linear relationship between them. It is worth noting that KDA and RKDA are able to generate a feasible solution and quantify its quality for even the ten-product problem under one CPU hour.
Z. Zeng and S. Cremaschi
524
5. Conclusion This paper discussed decomposition based heuristic algorithms, KDA and RKDA, to quickly obtain a feasible solution, a primal and a dual bound for large scale MSSP problems with endogenous uncertainty. Unlike the deterministic equivalent of MSSP formulations, the KDA and RKDA do not suffer from space complexity, and hence, can be used to generate a feasible solution and to obtain its solution quality for problems that are currently computationally intractable. Furthermore, the proposed algorithms is an initial step towards optimally solving large-scale MSSP for real-world sized problems. Table 2. Primal and dual bounds, relative gap, and solution times for the KDA and RKDA
Case Name 2-prod 3-prod 4-prod 5-prod 6-prod 7-prod 10-prod
KDA Primal Bound ($M) 1097 1178 1677 2052 2407 2874 4078
RKDA Dual Bound ($M) 1120 1221 1721 2127 2516 3008 4237
Relative Gap 2.10% 3.65% 2.62% 3.65% 4.53% 4.66% 3.90%
Solution Time (CPUs) 0.1 2 3 7 18 69 2377
References Apap, R. M., & Grossmann, I. E. (2017). Models and computational strategies for multistage stochastic programming under endogenous and exogenous uncertainties. Computers & Chemical Engineering, 103, 233-274. Christian, B., & Cremaschi, S. (2014). A Quick Knapsack Heuristic Solution for Pharmaceutical R&D Pipeline Management Problems. Computer-Aided Chem Engineering, 33, 1285-1290. Christian, B., & Cremaschi, S. (2015). Heuristic solution approaches to the pharmaceutical R&D pipeline management problem. Computers & Chemical Engineering, 74, 34-47. Christian, B., & Cremaschi, S. (2017). Variants to a knapsack decomposition heuristic for solving R&D pipeline management problems. Computers & Chemical Engineering, 96, 18-32. Christian, B., & Cremaschi, S. (2017). A branch and bound algorithm to solve large-scale multistage stochastic programs with endogenous uncertainty. AIChE Journal, DOI 10.1002/aic.16019. Colvin, M., & Maravelias, C. T. (2009). Scheduling of testing tasks and resource planning in new product development using stochastic programming. Computers & Chemical Engineering, 33(5), 964-976. Colvin, M., & Maravelias, C. T. (2010). Modeling methods and a branch and cut algorithm for pharmaceutical clinical trial planning using stochastic programming. European Journal of Operational Research, 203, 205-215. Gupta, V., & Grossmann, I. E. (2014). A new decomposition algorithm for multistage stochastic programs with endogenous uncertainties. Computers & Chemical Engineering, 62, 62-79. Solak, S., Clarke, J. P. B., Johnson, E. L., & Barnes, E. R. (2010). Optimization of R&D project portfolios under endogenous uncertainty. European Journal of Operational Research, 207(1), 420-433. Tarhan, B., & Grossmann, I. E. (2008). A multistage stochastic programming approach with strategies for uncertainty reduction in the synthesis of process networks with uncertain yields. Computers & Chemical Engineering, 32(4), 766-788. Z. Zeng, & S.Cremaschi (2017) Artificial lift infrastructure planning for shale gas producing horizontal wells. FOCAPO/CPC, Tuscan, AZ.