Integrated size and price optimization for a fashion retailer

Universit¨ at Bayreuth Fakult¨ at f¨ ur Mathematik, Physik und Informatik Lehrstuhl f¨ ur Wirtschaftsmathematik Prof. Dr. J¨ org Rambau Dissertation zur Erlangung des Grades ”Doktor der Naturwissenschaften“ an der Universit¨ at Bayreuth Integrated Size and Price Optimization for a fashion retailer vorgelegt von: Dipl. Math. Miriam Kießling Lehrstuhl f¨ ur Wirtschaftsmathematik Universit¨ at Bayreuth 95440 Bayreuth Tel.: 0921/55-7354 [email protected] vorgelegt am: 25. September 2012 Betreuer: Herr Prof. Dr. J¨ org Rambau Zusammenfassung Die vorliegende Arbeit ist das Ergebnis einer mehrj¨ ahrigen Zusammenarbeit mit einem deutschen Textildiscounter. Das Ziel war die Entwicklung eines entscheidungsunterst¨ utzenden Systems f¨ ur die Belieferung der etwa 1300 Filialen in Deutschland. Diese weist einige Besonderheiten auf: Die Filialen werden mit vorverpackten Kombinationen von Gr¨ oßen eines Artikels, sogenannten Lot-Typen, beliefert. Durch die Zusammenstellung dieser Lot-Typen, die bereits in dem Niedriglohnland erfolgt, in dem die Ware produziert wird, sollen die Handgriffe im Zentrallager und damit die Kosten in Deutschland reduziert werden. Um auch den weiteren Aufwand im Zentrallager m¨ oglichst klein zu halten, werden nur maximal vier bis f¨ unf verschiedene Lot-Typen f¨ ur einen Artikel verwendet. Außerdem wird pro Artikel jede Filiale nur mit einem Lot-Typ in einer Vielfachheit beliefert. Da es sich um Modeartikel handelt, die in der Regel nicht nachbestellt werden k¨ onnen, ist die Popularit¨ at des jeweiligen Produkts von vornherein unbekannt. Bedarfe k¨ onnen nur sehr grob, das heißt durchschnittlich und auf Warengruppenebene, gesch¨ atzt werden. ¨ Uberund Unterbelieferungen lassen sich nicht vermeiden. Eine Einflussnahme auf den Verkaufsprozess ist durch Preisreduzierungen m¨ oglich. Um ¨ Uberbelieferungen zu kompensieren, kann der Preis eines Artikels w¨ ochentlich auf vordefinierte vom Startpreis abh¨ angige Preisstufen herabgesetzt werden. Preisreduzierungen erfolgen f¨ ur einen Artikel in allen Filialen und f¨ ur alle Gr¨ oßen simultan. In Rahmen der Kooperation wurden mathematische Problemformulierungen entwickelt, mit dem Zweck, Kosten f¨ ur die Abweichung von Belieferung und gesch¨ atztem Bedarf zu minimieren. Der eigentliche Verkaufsprozess wurde bei der Ermittlung dieser Kosten nicht oder nur sehr grob betrachtet. Wir beziehen nun die M¨ oglichkeit von Preisreduzierungen bereits bei der Entscheidung ¨ uber die Belieferung ein. Das Ergebnis ist das zweistufige stochastische Programm ISPO: Die sogenannte Erststufenentscheidung ist die Festlegung einer Belieferungsstrategie, die Zweitstufenentscheidung oder der Rekurs, die Entscheidung ¨ uber Preisreduzierungen im Verkaufsverlauf. ISPO liefert eine ertragsmaximierende Belieferungsstrategie sowie sich darauf beziehende optimale Preisreduzierungsstrategien f¨ ur betrachtete Szenarien. ISPO ist zu komplex um es mit Standardverfahren zu l¨ osen. Die Entwicklung von speziellen L¨ osern war notwendig. Zum einen pr¨ asentieren wir einen exakten L¨ oser zum Benchmarking und zum anderen eine schnelle Heuristik f¨ ur den praktischen Einsatz beim Industriepartner. Der exakte L¨ oser basiert auf der Idee m¨ ogliche Preisreduzierungsstrategien zu enumerieren. Damit kann ISPO auf eine fr¨ uhere Problemformulierung zur Optimierung der Belieferungsstrategie, die mit Standardverfahren gel¨ ost werden kann, zur¨ uckgef¨ uhrt werden. i ii In der Praxis ist eine L¨ osung von ISPO nur durch Enumeration aller m¨ oglichen Preisreduzierungsstrategien zeitlich unm¨ oglich. Daher wird die Idee zu einem problembezogenen Branch&Bound Verfahren erweitert. In diesem Zusammenhang entwickeln wir duale Schranken f¨ ur allgemeine zweistufige stochastische Optimierungsprobleme, die auf der sogenannten wait-andsee solution aus der stochastischen Optimierung basieren. Wir zeigen, dass unsere Schranken im Allgemeinen sch¨ arfer sind. Die Heuristik sucht ausgehend von einer zul¨ assigen Zweitstufenentscheidung eine dazu optimale Erststufenentscheidung und alterniert dann bis zur Konvergenz zwischen zweiter und erster Stufe. Die Optimalit¨ atsl¨ ucke ist klein genug um einen praktischen Einsatz zu rechtfertigen. In der Praxis werden die bez¨ uglich ISPO optimalen Preisreduzierungsstrategien nicht verwendet. Stattdessen werden aktuelle Verkaufszahlen ausgenutzt. Gem¨ aß diesen und einer angepassten Bedarfssch¨ atzung wird w¨ ochentlich eine neue optimale Preisreduzierungsstrategie f¨ ur den verbleibenden Verkaufszeitraum ermittelt. Daf¨ ur pr¨ asentieren wir einen Algorithmus, der auf dynamischer Programmierung beruht und nicht optimale L¨ osungen durch sogenannte Dominanztests von vornherein auszuschliessen versucht. ISPO, genauer gesagt unsere Heuristik, zusammen mit der w¨ ochentlichen Aktualisierung der Preisreduzierungsstrategien bildet unser entscheidungsunterst¨ utzendes System zur integrierten Gr¨ oßenund Preisoptimierung DISPO. Wir testeten DISPO in einem f¨ unfmonatigen Feldversuch, durchgef¨ uhrt als statistisches Experiment, beim Praxispartner. Hierbei wurden Paare ¨ ahnlicher Filialen miteinander verglichen: In einer Filiale, der Testfiliale, wurde die von ISPO vorgeschlagene Belieferungsstrategie umgesetzt und w¨ ochentlich, wie oben beschrieben, die Preisreduzierungsstrategie aktualisiert. In der anderen Filiale, der Kontrollfiliale, wurde ein fr¨ uheres Modell zur Festlegung der Belieferungsstrategie eingesetzt, in dem der Verkaufsprozess nicht integriert ist. Preisreduzierungen in den Kontrollfilialen wurden vom Projektpartner angeordnet. In den Testfilialen, f¨ ur die DISPO eingesetzt wurde, erzielten wir einen um mehr als 1,5Prozentpunkte h¨ oheren realisierten Rohertrag als in den Vergleichsfilialen. Abstract This thesis is the result of a collaboration with a German fashion retailer which lasted for several years. The aim was the development of a decision-support system for the supply of the about 1300 branches in Germany. There are some specialties about the situation at our industrial partner: The branches are supplied by prepackaged size-assortments of a product which we call lot-types. With the objective to economize handling cost, these lot-types are already composed at the respective low-wage country where the article is also produced. The expense at the German central warehouse is further reduced by allowing only four or five different lot-types for the delivery of one product. Moreover, each branch is supplied by a certain quantity of a single lot-type. For the most fashion articles replenishment is not possible. The sales success of a product is a priori unknown. Historical sales data can only be used on a higher aggregation level, e.g., the average historical demand on the commodity group level. Demand estimation is therefore very vague. Underand oversupplies are unavoidable. Influence over the sales process is possible by marking down prices. To compensate for an oversupply of a product, weekly the price can be reduced to predefined price steps which depend on the starting price of the product. Mark-downs for an article are performed simultaneously for all branches and sizes. Within the cooperation mathematical problem formulations with the aim to minimize measures for the deviation of supply from estimated demand had been developed. In these measures the selling process is not or only very vaguely regarded. Now we include the possibility of marking down prices during the selling time already when deciding on the supply. The result is the two-stage stochastic program ISPO: The so-called first stage decision is the determination of a supply policy. The second stage decision, or recourse, is the decision on mark-downs during the selling time. ISPO yields an expected revenue maximizing supply strategy and corresponding optimal mark-down strategies for the considered scenarios. ISPO it too complex to solve it via standard approaches. Customized methods had to be devised to solve ISPO. On the one side we present an exact solver for benchmarking. On the other side a fast heuristic was developed for practical use at our partner. The basic idea of our exact solver is to enumerate all possible mark-down strategies. With this it is possible to reduce ISPO to a former formulation for the optimization of supply, which can be solved via standard approaches. In practice enumeration of all valid mark-down strategies for the purpose of solving ISPO is for reasons of time impossible. Therefore the idea is extended to a customized Branch&Bound approach. In this context we derived dual bounds for general twostage stochastic programs which are based on the so-called wait-and-see solution from stochastic programming. We show that in general our bounds are tighter. The heuristic, beginning with a valid second stage decision, determines an optimal iii iv first stage decision and alternates between solving the first stage and the second stage until convergence is reached. The optimality gap is small enough to justify a practical use at the industrial partner. In practice the by ISPO proposed mark-down strategies are not applied; instead latest sales figures are exploited. According to these and an updated demand estimation weekly a new optimal mark-down strategy for the remaining selling time of the product is determined. For this purpose we propose an algorithm which relies on dynamic programming and tries to exclude non-optimal solutions a priori by dominance checks. ISPO, more precisely our heuristic approach, together with the weekly adaption of the mark-down strategy forms our decision support system for integrated size and price optimization DISPO. We tested DISPO in a five-month field study, performed as a statistical experiment, at our partner where pairs of similar branches were compared. At one branch of each pair, the test branch, supply and mark-down decisions came from ISPO. With respect to latest sales figures the mark-down decisions were weekly updated via our dynamic programming approach. At the other branch, the control branch, these decisions were not integrated: Supply was determined according to a strategy resulting from a former model that disregarded the selling process and mark-downs were handled manually by our partner. For the branches at which the decisions of ISPO were implemented an average raise of 1.5 percentage points of relative revenue was observed. Contents List of Symbols x 1 Introduction 1 1.1 Related Work .............................. 2 1.2 Our contribution ............................. 3 1.3 Outline of the thesis ........................... 4 1.4 Preliminary remarks ........................... 4 1.4.1 Basics from mixed-integer linear programming ........ 4 1.4.2 Labelling of own results .................... 5 1.4.3 Computational results ..................... 5 2 Collaboration with the industrial partner – historical progress 6 2.1 Lot-types and lots ............................ 6 2.2 The Lot-type Design Problem LDP ................... 7 2.2.1 Problem specification ...................... 7 2.2.2 Problem formulation ...................... 7 2.2.3 The Score-fix-adjust heuristic ................. 8 2.2.4 Implementation at the industrial partner ............ 9 2.3 The Stochastic Lot-type Design Problem SLDP ............ 10 2.3.1 Problem specification ...................... 10 2.3.2 Modelling the SLDP ...................... 10 2.3.3 Solving the SLDP by the LDP ................. 12 2.3.4 A column generation approach ................. 13 2.4 Reasons for integrating price optimization ............... 14 2.5 Price optimization ............................ 15 2.5.1 Problem specification ...................... 15 2.5.2 Problem formulation ...................... 15 2.5.3 Justificating the problem formulation ............. 17 2.5.4 Price optimization with receding horizon – POP-RH ..... 17 2.6 Integrated size and price optimization ................. 18 3 Demand estimation 19 3.1 Literature review ............................ 19 3.2 Empirical estimation .......................... 21 3.2.1 Relative demand estimation .................. 22 3.2.2 Regarding different scenarios ................. 24 3.2.3 Splitting up the demand to sales periods ............ 25 3.2.4 Price-dependent demand .................... 26 v CONTENTS vi 3.2.5 Combining the estimated factors ................ 26 3.2.6 Updating the scenario ..................... 26 3.3 Logistic regression ........................... 27 3.3.1 Maximum likelihood estimation ................ 27 3.3.2 Binary logistic Regression ................... 28 3.3.3 Ordinal logistic regression ................... 29 3.4 Applying ordinal logistic regression .................. 29 3.4.1 Data sample ........................... 29 3.4.2 Choosing the model ...................... 30 3.4.3 Result .............................. 32 3.5 Comparison of different estimation methods .............. 34 3.5.1 Methodology .......................... 34 3.5.2 Results ............................. 35 4 Price Optimization 37 4.1 Extending POP by mark-down costs – a mixed-integer nonlinear program 38 4.1.1 Problem formulation ...................... 38 4.1.2 Nonlinearity by mark-down costs ............... 39 4.2 Enumerating price trajectories ..................... 39 4.3 Excursus: Dynamic programming ................... 42 4.3.1 General dynamic program ................... 42 4.3.2 The dynamic programming algorithm ............. 43 4.3.3 Deterministic Systems ..................... 43 4.3.4 Solving shortest path problems ................. 44 4.3.5 Resource constraint shortest path problems and dominance . . 45 4.4 Dynamic generation of mark-down strategies ............. 46 4.5 Pruning the enumeration tree – dominating partial mark-down strategies 48 4.6 Implementation ............................. 52 4.7 An accompanying example ....................... 53 4.8 POP-DYN applied on the accompanying example ........... 54 4.9 Computational results .......................... 58 4.10 Conclusion of the chapter ........................ 59 5 Stochastic Optimization 60 5.1 Two-stage stochastic programs ..................... 61 5.2 Solving stochastic programs ...................... 62 5.2.1 The L-shaped method for two-stage linear stochastic programs ............................ 62 5.2.2 Solving two-stage mixed-integer stochastic programs ..... 62 5.3 Common bounds for two-stage stochastic programs .......... 63 5.3.1 Dual bounds .......................... 63 5.3.2 Primal bounds ......................... 67 5.4 Multi-stage stochastic programs .................... 68 6 The Integrated Size and Price Optimization Problem (ISPO) 69 6.1 Problem specification .......................... 69 6.2 ISPO as a two-stage stochastic mixed-integer program (SMIP) in its extensive form .............................. 70 6.3 Complexity of ISPO .......................... 73 6.4 Solving ISPO with standard approaches ................ 73 CONTENTS vii 7 Reducing ISPO to the SLDP 75 7.1 Single supply revenues ......................... 75 7.1.1 Runtime of Algorithm 6 .................... 80 7.1.2 An Example .......................... 84 7.1.3 Computational results ..................... 85 7.2 Establishing lot-type revenues ..................... 87 7.3 Fixing price trajectories in the ISPO – an SLDP ............ 87 7.4 Solving the ISPO by enumerating SLDPs ............... 88 8 Dual Bounds 90 8.1 Dual bounds from wait-and-see solutions ............... 90 8.1.1 The wait-and-see solution ................... 90 8.1.2 Extending wait-and-see solutions ............... 91 8.1.3 Relaxations of extended wait-and-see solutions ........ 94 8.2 Application to ISPO ........................... 94 8.2.1 Relaxing the lot-type constraint – single supply relaxations . . 94 8.2.2 Extended wait-and-wee solutions for ISPO .......... 98 8.2.3 Computational results ..................... 100 8.3 Conclusion of the chapter ........................ 101 9 Solving the Integrated Size and Price Optimization Problem 102 9.1 An exact Branch&Bound approach ................... 102 9.1.1 The algorithm .......................... 103 9.1.2 Some implementational aspects ................ 104 9.1.3 Computational results ..................... 107 9.1.4 ISPO-BAB applied to the accompanying example ....... 109 9.2 A heuristic approach – ISPO-PingPong ................ 112 9.2.1 Reversible recourse ....................... 112 9.2.2 The main algorithm ....................... 112 9.2.3 Fixing price trajectories .................... 112 9.2.4 Solving the SLDP(WE).................... 113 9.2.5 Solving the POP ........................ 114 9.2.6 Computational results ..................... 114 9.2.7 Similarities to familiar approaches ............... 116 9.3 Computational results for real-world instances ............. 118 9.4 General goodness of ISPO-PingPong .................. 121 9.5 Conclusion of the chapter ........................ 122 10 DISPO in practical application – real-world experiments 124 10.1 Performing statistical experiments ................... 124 10.1.1 Blind experiments ....................... 125 10.1.2 Statistical significance ..................... 125 10.1.3 Statistical tests in general .................... 125 10.1.4 Wilcoxon signed-rank test ................... 126 10.2 Performing our field-studies as statistical experiments ........................ 128 10.3 POP-RH in real-world studies ..................... 130 10.3.1 Performing price optimization with receding horizon – POP-RH 130 10.3.2 Sales increase by mark-downs ................. 130 10.3.3 Earnings increase by mark-downs ............... 131 CONTENTS viii 10.4 Potential of ISPO ............................ 134 10.5 DISPO – the field study ......................... 136 10.5.1 Preparation ........................... 136 10.5.2 Setup of the field study ..................... 136 10.5.3 Evaluation ........................... 138 10.5.4 Results of the field study .................... 139 11 Conclusion 147 A ISPO-PingPong – further results 149 B Sales increase by mark-downs – further results 153 C Single supply revenues for the accompanying example 164 D Demand estimation via logistic regression 167 E Instances 169 CHAPTER 1. INTRODUCTION 5 1.4.2 Labelling of own results This thesis is a result of more than a five years long cooperation (2006 to 2011) at our industrial partner. We joined the project at 2009. Therefore not all topics in this thesis concerning the collaboration are a result of our own or our complete own work. We will use three terms to differentiate between us and the other colleagues: •former DISPO-team: involved persons (in alphabetical order) were Konstantin Gaul, Tobias Kreisel, PD Dr. Sascha Kurz, Alexander Lawall, Prof. Dr. J¨ org Rambau, •DISPO-team: involved persons were Miriam Kießling, Tobias Kreisel, PD Dr. Sascha Kurz, Alexander Lawall, Prof. Dr. J¨ org Rambau, •we: Miriam Kießling. 1.4.3 Computational results In this thesis we will state several computational results. If not otherwise specified they were provided by a machine with Intel(R) Xeon(R) processor with 2.33 GHz and 62 GB of RAM. We implemented all stated algorithms in C++. Whenever we will use a state-of-the-art solver for mixed-integer linear programs in our results we use IBM ILOG CPLEX, version 12.2 (as alternative also SCIP could be chosen; we tested the programs for SCIP in version 2.0.1 combined with SOPLEX-1.5.0 how it is included in the ZIBOPTSUITE-2.0.1). Chapter 2 Collaboration with the industrial partner – historical progress The Integrated Size and Price Optimization Problem ISPO is an enhancement of former models for the optimization of supply or size optimization that were developed during more than a five years long cooperation with our industrial partner. In this chapter we outline the main results related to the time before we developed DISPO. We show the historical progress from deterministic size optimization to integrated size and price optimization and outline the basic ideas of the implemented approaches. In Section 2.1 we treat the terms lot-type and lot in detail. The first model which assumes deterministic demand – the Lot-type Design Problem LDP presented by Gaul, Kurz and Rambau [GKR09], which is currently as a standard implemented at our industrial partner, is outlined at first in Section 2.2. We describe the so-called SFA heuristic which was introduced in [GKR10] as a solving method for the LDP. By additionally regarding stochastic demand and lot-opening costs we arrive at the SLDP – the Stochastic Lot-type Design Problem in Section 2.3. We show how to reduce the SLDP to the LDP which makes it possible to apply the SFA heuristic on it. Because both models do not contain all relevant properties of the sales process at our industrial partner as we outline in Section 2.4, we extend the SLDP by regarding the possibility of mark-downs, i.e. integrating price optimization. Price optimization as it was implemented by the former DISPO-team is treated in Section 2.5. We give a short outlook on integrated size and price optimization in Section 2.6. 2.1 Lot-types and lots Our industrial partner already before the cooperation supplied its branches with lottypes. This is done to economize handling costs. A lot-type describes a prepackage which contains items of one product in different sizes and numbers. Mathematically we are given a lot-type by a n-tuple where nequals the number of sizes. The entries describe the number of items per size where we assume that the sizes are ordered increasingly. For example if we want to specify a prepackage containing 1item of size S, 3items of size M and 2items of size L, we do this by the lot-type (1,3,2). Before 6 CHAPTER 2. HISTORICAL PROGRESS 7 the collaboration all branches were supplied by so-called standard prepackages. These are prepackages, i.e lot-types, which contain always just 1item of the extreme sizes and 2items for the middle sizes. An example would be the lot-type (1,2,2,2,1). A branch can be supplied by a specific number of prepackages, provided all prepackages are specified by the same lot-type. (To avoid handling costs it is not allowed to mix differing prepackages for supplying a branch.) For example supplying 2times Lot-type (1,3,2) means supplying 2 items of Size S, 6 items of Size M and 4 items of Size L. 2.2 The Lot-type Design Problem LDP The Lot-type Design Problem LDP was the first formulation for optimization of supply at our industrial partner and was first presented in [GKR09]. 2.2.1 Problem specification We consider an article with a given set of sizes S. We want to supply each branch of a set of branches Bwith one lot-type from a set of lot-types Lin a multiplicity from a set of multiplicities M={1, . . . , mmax}. The lot-types are given by four parameters: the minimum supply per lot-type and size vmin, the maximum supply per lot-type and size vmax, the minimum supply per lot-type vlmin and the maximum supply per lot-type vlmax. At the maximum κdifferent lot-types can be used for supplying the branches. The overall supply must lie in between a lower bound Iand an upper bound ¯ I. The supply shall meet the dependent demand db,s, s ∈S, b ∈Bfor each branch band size sas good as possible. 2.2.2 Problem formulation The LDP is formulated as follows. Problem 1 (LDP [GKR09]). min X b∈BX `∈LX m∈M distLDP b,l,m ·xb,l,m (2.1) subject to X `∈LX m∈M xb,l,m = 1 ∀b∈B, (2.2) X `∈L y`≤κ, (2.3) X m∈M xb,l,m ≤yl∀b∈B, ` ∈L, (2.4) Ib,s =X `∈LX m∈M m·`s·xb,`,m ∀b∈B, s ∈S, (2.5) I=X b∈BX s∈S Ib,s,(2.6) I∈[I,I],(2.7) xb,`,m ∈ {0,1} ∀b∈B, ` ∈L, m ∈M, (2.8) y`∈ {0,1} ∀`∈L. (2.9) CHAPTER 2. HISTORICAL PROGRESS 8 Binary variables xb,`,m indicate if Lot-type `is delivered to Branch bin Multiplicity m. If this is answered by “yes” the variable takes value one, otherwise zero. The binary variable y`takes value one if at least one branch is supplied by Lot-type `, otherwise it takes value zero. With Constraint (2.2) it is ensured that every branch is supplied by exactly one lottype in one multiplicity. Constraint (2.4) connects the variables xb,`,m and ylin such a way that y`can take value one only if lot-type `is delivered to at least one branch. The adherence of the upper bound κfor the number of different lot-types is enforced by Constraint (2.3). With the constraints (2.5) and (2.6) it is ensured that the overall supply adheres to the lower and the upper bound; the variables Ib,s here describe the supply for Branch band Size s,Ithe overall supply and lsthe number of items of Size sin Lot-type `. The objective coefficients distLDP b,l,m measure the deviation between supply and demand. In our case we restrict ourselves to the L1-Norm, which is also implemented at the industrial partner. For other measurements see [GKR09]. With a demand db,s for Size sin Branch b, see Chapter 3for the estimation method, the objective coefficients are given by distLDP b,l,m := X s∈S |db,s −m·ls|.(2.10) Remark 1 (Lower and upper bounds for the overall supply [GKR09]).For each product our partner first decides on an overall capacity Dbefore the items are distributed to the particular branches. There are two reason why this amount is softened to the interval I≤D≤Iin Constraint (2.7): If for example the overall capacity Dfor an article was prime and there were two or more sizes for the considered product than the LDP would be infeasible. Such from a theoretical point of view the soft bound is needed to guarantee feasibility of the problem. The other reason is practical. Our partner does not always obtain the ordered supply from the supplier. As a rule there is a deviation between the ordered and the actual delivered amount. Deviations up to 5% from the ordered volume may occur. Remark 2 (Complexity of the LDP [GKR09]).The LDP is NP-hard. This is shown by reducing the p-median problem on it after restricting the set of multiplicities to the case M={1}and adapting the lower and upper bound in Constraint (2.7)in such a way that it is not a real restriction. Depending on the number of sizes, allowed lot-types and multiplicities the solving process for real-world instances by using state-of-the-art solvers for mixed-integer linear programs (MIPs) as SCIP or CPLEX can be very time-consuming (more than 5 hours) and therefore is not suitable for practical purposes. For that reason Gaul, Kurz and Rambau implemented the so-called SFA heuristic. 2.2.3 The Score-fix-adjust heuristic In [GKR10] the Score-fix-adjust (SFA) heuristic for the LDP was proposed. The name of the heuristic stems from the three basic steps the heuristic consists of. 1. Score: The lot-types get scores in terms of how good they meet the demands of the branches. CHAPTER 2. HISTORICAL PROGRESS 9 2. Fix: For a given time period κ-subsets of the set of lot-types are traversed according to the scoring from the previous step. For each branch the best fitting lot-type from the considered subset and the related best fitting multiplicity is fixed. 3. Adjust: The multiplicities are adjusted to adhere to the bounds Iand ¯ Ifor the overall supply. In the Score-step each lot-type `gets points according to the objective coefficient distLDP b,l,m. This is done in the following way. For every branch band lot-type `first the best fitting multiplicity m(b, `)is determined. This is the multiplicity m∈Mfor which Ps∈S|db,s −m·ls|is minimal. According to Ps∈S|db,s −m(b, `)·ls|the lot-types `∈Lare ordered decreasingly. This yields an ordering of the lot-types in terms of how good they meet the demand of Branch b. Starting from this for each branch the three locally best fitting lot-types can be determined. A score of 100 to the best fitting lot-type, a score of 10 to the second best fitting lot-type and a score of 1 to the third best fitting lot-type is added. (Of course this can be generalized to the first tbest fitting lot-types and different scoring schemes.) In the Fix-step the best κ-subsets of lot-types – best in terms of the highest sums of scores over all branches – are traversed. This is done for a predefined time period trusting that the most promising selections of lot-types were checked. For the considered κ-subset L0⊆Lfor a branch bthe lot-type `0∈L0is fixed which minimizes Ps∈S|db,s −m(b, `)·l0 s|. Such, a preliminary supply policy with a corresponding overall supply I0is specified. If I0∈[I,I]the supply policy is valid. Otherwise the Adjust-step, see below, has to be performed to establish feasibility. If the supply policy yields a smaller objective value of the LDP than the already considered ones or if it is the first considered one, we update our best found solution correspondingly. The Adjust-step assures the adherence of the bounds Iand ¯ Ifor the overall supply. Fixing the best fitting lot-type from the considered κ-subset with best matching multiplicity for each branch in the previous step might violate Constraint (2.7). There are two cases of infeasibility: 1. I0< I 2. I0>¯ I In the first case supply is increased until the lower bound is met. This is done in a greedy way: The branch for which increasing the currently fixed multiplicity by one is valid and leads to the smallest additional costs in terms of the objective function is determined. The multiplicity is increased by one and fixed. This procedure is iterated unless the lower bound for the overall supply Iis met. The proceeding in the second case is similarly. Supply iteratively is reduced until the upper bound is met. From all branches that are at least supplied with a lot-type in multiplicity 2 we choose the branch for which decreasing the multiplicity by one leads to the smallest additional costs. For 36 real instances the authors performed the SFA heuristics with a computation time of one second. This led to a mean optimality gap of 0.327% while the highest gap amounts to 2.114%. 2.2.4 Implementation at the industrial partner Gaul, Kurz and Rambau [GKR10] performed a preliminary study at our industrial partner to evaluate if the LDP performs better than manual planning of supply. Previously, CHAPTER 2. HISTORICAL PROGRESS 10 the branches were all supplied by the same standard lot-type as we described in Section 2.1. In [KRSW08] the authors could show that the size dependent demand among different sizes at our industrial partner actually varies and that the LDP together with their presented demand estimation method, see Chapter 3, could increase the gross yield about 0.85 percentage points. Since 2006 the LDP, more precisely, the SFA heuristic is implemented at the partner and used for nearly all fashion articles which are supplied in terms of lots. 2.3 The Stochastic Lot-type Design Problem SLDP The former DISPO-team enhanced the LDP to the Stochastic Lot-type Design Problem, the SLDP. The SLDP can be understood as an intermediate model between the deterministic LDP and the final stochastic model ISPO integrating price optimization. In the SLDP different scenarios with scenario probabilities and scenario dependent demands estimated from historical data are treated. While in the LDP abstract costs in form of the L1-norm are considered and overand undersupply are treated the same, now monetary asymmetric costs are imposed for the deviation between supply and demand. With these costs the SLDP can be seen as a first step in integrating the sales process in the size optimization. Monetary measurement now allows also to take other costs into account. On the one hand pick costs which arise from arranging the lot-types to lots and on the other hand lot-opening costs which arise from the fact that each additional supplied lot-type leads to higher logistic effort. 2.3.1 Problem specification We consider an article with a given set Sof sizes. We want to deliver each branch from a set Bof branches with one lot-type from a set Lof lot-types in a multiplicity from the set M={1, . . . , mmax}of multiplicities. At the maximum κdifferent lot-types are allowed to use for supply. For the ith supplied new lot-type from the set Llot-opening cost δiarise . For every handgrip needed for putting together the lot-types to lots pick cost pcost arise. The overall supply must lie in between a lower bound Iand an upper bound ¯ I. Now we consider a set Eof different scenarios with scenario probabilites Prob(e),∀e∈E. With given demands de b,s for each size s, branch band scenario ean oversupply is penalized by acquisiton price ap minus salvage value πpmax , an undersupply by starting price π0minus ap. This means that it is assumed, that each undersupply would lead to a loss of the full starting price while each oversupplied item can just be sold for the salvage value. The aim is to minimize the expected overall costs, i.e. the sum of the handling costs – lot-opening and pick cost – together with the expected costs for oversupply and undersupply. In terms of demand estimation and the estimation of the probabilities Prob(e), e ∈Esee Chapter 3. 2.3.2 Modelling the SLDP Before we introduce the entire model we first focus on the coefficients in the objective: The expected dependent demand db,s for Branch band Size sis given by db,s =Pe∈EProb(e)·de b,s, where db,s equals the dependent demand in the LDP. Now asymmetric costs for overand undersupply are introduced. An oversupply is penalized by acquisition price minus salvage value ap −πpmax , an undersupply by starting CHAPTER 2. HISTORICAL PROGRESS 11 price minus acquisition price π0−ap. Then the cost arising by suppling Branch bwith Lot-type `in Multiplicity mare given by distSLDP b,`,m which is defined as distSLDP b,`,m := X s∈S (max{m·ls−db,s,0}·(ap−πpmax )+max{db,s−m·ls,0}·(π0−ap)). (2.11) The Stochastic Lot-type Design Problem SLDP is modeled as follows: Problem 2 (SLDP). min X b∈BX `∈LX m∈MdistSLDP b,`,m +m·pcostxb,`,m + κ X i=1 δi·zi(2.12) subject to X `∈LX m∈M xb,`,m = 1 ∀b∈B, (2.13) X m∈M xb,`,m ≤y`∀b∈B, ` ∈L, (2.14) X `∈L y`≤ κ X i=1 zi,(2.15) zi≤zi−1i= 1 . . . , κ, (2.16) Ib,s =X `∈LX m∈M m·`s·xb,`,m ∀b∈B, s ∈S, (2.17) I=X b∈BX s∈S Ib,s,(2.18) I∈[I,I],(2.19) xb,`,m ∈ {0,1} ∀b∈B, ` ∈L, m ∈M, (2.20) y`∈ {0,1} ∀`∈L, (2.21) zi∈ {0,1}i= 1, . . . , κ. (2.22) Most constraints are similar to the ones of the LDP. At this point we explain only the differences and refer the reader to Problem 1. To take handling costs into account we introduce the binary variables zi, i = 1 ...,κ which indicate if at least idifferent lot-types are opened. Constraint (2.15) links the variables ziwith yl. By Constraint (2.16) it is ensured that zican take value one only if zi−1also does. The additional costs for opening new lot-types are added in the objective function and for every delivered lot pick costs m·pcost arise in addition to the costs for overand undersupply (2.12). Corollary 1 (Complexity of the SLDP).The SLDP is NP-hard. Proof. If we set δito zero for i= 1, . . . , κ in the SLDP, we obtain an LDP with changed objective coefficients because the constraints (2.15) and (2.16) in this case are equivalent to Constraint (2.3). That means, we can reduce the LDP in polynomial time to the SLDP. Because the LDP is NP-hard – as stated in Remark 2– the SLDP is, too. CHAPTER 2. HISTORICAL PROGRESS 12 2.3.3 Solving the SLDP by the LDP The SLDP simplifies to an LDP if we set all lot-opening costs to zero. We now show that in similar way we are able to determine the optimal solution of the SLDP as the best solution resulting from solving κLDPs. For i= 1, . . . ,κ we consider the following formulation of the LDP. Problem 3 (LDP(i)). min X b∈BX `∈LX m∈MdistSLDP b,`,m +m·pcostxb,`,m (2.23) subject to X `∈LX m∈M xb,l,m = 1 ∀b∈B, (2.24) X `∈L y`≤i, (2.25) X m∈M xb,l,m ≤yl∀b∈B, ` ∈L, (2.26) Ib,s =X `∈LX m∈M m·`s·xb,`,m ∀b∈B, s ∈S, (2.27) I=X b∈BX s∈S Ib,s,(2.28) I∈[I,I],(2.29) xb,`,m ∈ {0,1} ∀b∈B, ` ∈L, m ∈M, (2.30) y`∈ {0,1} ∀`∈L. (2.31) The LDP(i) is an LDP with the restriction that at most iinstead of κdifferent lottypes are allowed for supply, Constraint (2.25). The coefficients of the variables xb,`,m in the objective function are these from the SLDP. Opening-costs for new lot-types are not regarded. Having solved the LPD(i) with the optimal solution (x∗(i), y∗(i)) we can compute the corresponding overall opening costs by adding P`∈Ly∗ `(i) X j=1 δj to the optimal objective value z∗ LDP(i). Thus, by solving the LDP(i) for each 1≤ i≤κseparately we obtain the optimal supply for each possible allowed number of different lot-types. Adding the opening costs to the related objective value yields the optimal objective value of the SLDP. The LDP(i) for which the objective value plus the corresponding opening costs is minimal among 1≤i≤κthen yields the optimal solution of the SLDP. Theorem 1 (Deducing the optimal solution of the SLDP from the LDP).With z∗ LDP(i) we denote the optimal objective value of the LDP(i). The corresponding optimal solutions are denoted by x∗(i)and y∗(i). With z∗ SLDP we denote the optimal objective value of the SLDP and with x∗,y∗and z∗the related values of the variables. We define i∗:= arg min i=1,...,κ   z∗ LDP(i)+P`∈Ly∗ `(i) X j=1 δj  .(2.32) CHAPTER 2. HISTORICAL PROGRESS 13 Then the optimal objective function value of the SLDP is given by z∗ SLDP =z∗ LDP(i∗)+P`∈Ly∗ `(i∗) X j=1 δj.(2.33) It is x∗=x∗(i∗),y∗=y∗(i∗)and z∗ j= 1 for j= 1,...,P`∈Ly∗ `(i∗)and z∗ j= 0 for P`∈Ly∗ `(i∗)< j ≤κ. Proof. It is i+the number of used lot-type according to the optimal solution of the SLDP. If we would set κ=i+in the SLDP this would yield the same optimal solution. We call the SLDP restricted to maximal κ=i+different lot-types SLDP(i+). The corresponding optimal objective value is denoted by z∗ SLDP(i+). The LDP(i+) yields a supply x∗(i+)that minimizes Pb∈BP`∈LPm∈MdistSLDP b,`,m +m·pcostx∗(i+)b,`,m for maximal i+different lot-types not regarding lot-opening costs. Because the ziand the lot-opening costs are independent from the selected lot-types and depend only on the number of them the LDP(i+) yields the same optimal solutions in terms of the supply than the SLDP(i+). To obtain the same solution we could set the lot-opening costs δiin the SLDP(i+) to zero, compute the optimal supply and later on add the costs δifor the i+used lot-types. This is the same as solving the LDP(i+) and adding the corresponding lot-opening costs. Overall that means z∗ SLDP =z∗ SLDP(i+)=z∗ LDP(i+)+P`∈Ly∗ l(i+) X j=1 δj. It is z∗ SLDP ≥z∗ LDP(i) + PP`∈Ly`(i∗) j=1 δjfor all i= 1, . . . , κ, i 6=i+. Otherwise the SLDP would yield an optimal solution with less or more than i+lot-types. By setting i∗=i+the claim follows. Remark 3. In order to compute the optimal solution of the SLDP we propose to solve the LDP(i)s in ordering i=κ, . . . , 1. If the LDP(i)yielded a supply policy with just i−< i lot-types we would not have to solve the LDP(j)for i−≤j < i. The LDP(j)s would yield the same optimal solution as the LDP(i). So traversing the LDP(i)s in order i=κ, . . . , 1may reduce the computational effort. By reducing the SLDP to the LDP now it is possible to apply solving methods for the LDP – as the described SFA heuristic – to the SLDP. Later on, in Chapter 9, we will mention how this property can be exploited when solving the Integrated Size and Price Optimization Problem ISPO which is discussed in Chapter 6. 2.3.4 A column generation approach In [KKR11a] an exact column generation approach for the LDP is presented by the DISPO-team. The approach is guided by two main ideas1. •Considering the restricted master problem (RMP) with only a subset L0⊂Lof lot-types 1For further information about column generation we refer the reader to [LD05],[LD11] or [L¨ ub10] CHAPTER 2. HISTORICAL PROGRESS 14 •Solving the LDP for the most promising subset of lot-types ¯ L⊆L0exactly We will sketch the main parts of the approach. For further details we refer the reader to [KKR11a]. A restricted master problem RMP of the LP relaxation of the LDP is considered. The only difference to the LP relaxation is that only a subset of the lot-types are considered. Thus, because the optimal solution may not be contained, the restricted master problem yields an upper bound for the LP relaxation of the original problem. 1. At first, a starting solution (x∗, y∗)of the LDP is determined. This can be done via an adapted version of the SFA heuristic. The used lot-types, i.e. lot-types with y∗ `= 1 then are added to the initial subset L0of lot-types. Additionally the three best fitting lot-types for each branch – as they result from the score-step of the SFA-heuristic, see 2.2.3 – are added to L0. The lot-types from the set L0are the only lot-types that are considered in the RMP at the beginning. 2. With (xRMP, yRMP)we denote the optimal solution of the RMP. The set of most promising lot-types ¯ Lis the subset of all lot-types from the set L0with yRMP l≥ε where εis a small constant, for example ε= 0.15. If the optimal objective value of the RMP is smaller than the objective value of the current best integer solution (x∗, y∗)the LDP restricted to ¯ Lis solved exactly and possibly the currently best integer solution (x∗, y∗)is updated. (If the RMP yields an optimal value higher than the to (x∗, y∗)corresponding objective value the set L0cannot contain the optimal subset of lot-types. Because the RMP is a relaxation of the LDP that contains only the lot-types L0the optimal objective value is a lower bound for the LDP restricted to the set L0of lot-types. Such, we are not able to obtain a better integer solution than (x∗, y∗)by only regarding the lot-types from the set L0.) Cover cuts are added to the RMP to forbid that the optimal solution of the RMP yields ¯ Las the set of most promising lot-types again. This implies a branching on the set ¯ Land the rest of lot-types L0currently considered in the RMP. 3. Whenever the optimal function value of the RMP is higher than or equals the objective value of the current best solution (x∗, y∗)the pricing step is performed in which – if possible – new lot-types are added to the RMP – i.e. L0is updated and the RMP is solved again and so on. Whenever the optimal objective function value of the RMP is smaller than the to (x∗, y∗)related objective function value of the LDP, then we update the subset of most promising lot-types L0and branch on this subset, i.e. perform Step 2. If the optimal objective function value of the RMP exceeds or equals the objective value of the LDP corresponding to (x∗, y∗) and no more lot-types are/can be added to the RMP than we end up at this point and return (x∗, y∗)as optimal solution. The results in [GKR09] show that for real-world instances the maximum amount of time for solving can be reduced from 36 minutes to 4 seconds. Even very large instances – for which state-of-the-art MIP solvers fail – can be solved in less than 16 minutes. 2.4 Reasons for integrating price optimization The introduction of monetary costs in Problem 2is a first step in integrating the sales process in the size optimization. But by penalizing oversupply with acquisition price CHAPTER 3. DEMAND ESTIMATION 21 Applicability To omit right-censored data in our case is not possible. Because the most sizes per branch are supplied by maximal one or two items and nearly all observations are rightcensored this would lead to a tiny sample which could not serve as basis for demand estimation whatsoever. Thus, we have to deal with censored observations. For our right-censored and ordinal data linear regression as mentioned above is not suitable. Ordinary linear regression may lead to negative values for the dependent variable, i.e. in our case the demand. Because the supply and consequently the censoring-point for each size and branch differs the tobit model is not an alternative, rather a general censored regression model. But also this is not conceived for integer numbers of outcomes. The EM approach by Vulcano et al. described above in principle appears promising. However, the approach needs a set of comparable products as input. One could consider all items of the same commodity group or sub commodity group as comparable but this is not the case. The articles differ in color, fashion, etc. and most important price. But even if our industrial partner could commit us lists of comparable products there would be another difficulty: Because there is no reordering of products and the products have different sales starts we have to regard that the list of comparable articles may change over time. Moreover, incomparability caused by mark-downs would have to be regarded. The Kaplan-Meier approach as applied in [HLRO11] can not directly be adopted. In our situation the items are not reordered. So far, we did not see how the method could be adapted to our situation where we also have to include the possibility of mark-downs during the selling time. From all mentioned methods from literature the for us most promising approach is the ordinal logistic regression model. This can explicitly deal with small integer outcomes. Moreover we can include price dependencies and dependencies in terms of the popularity of the observed product. By including the current stock as independent variable we estimate sales not demand and have not to deal with right-censored observations. 3.2 Empirical estimation All models and results in this thesis are based on an empirical demand estimation developed by the former DISPO-team. Parts of this method were already implemented to estimate the demand in terms of the LDP the SLDP and the Price Optimization Problem POP. Because there are no publications including a detailed description we outline the method at this point. The overall demand Dis considered as an exogenous quantity.1 The estimator circumvents the difficulties arising from varying popularity among different articles by considering the conditional probability that – if an item is sold – this happens in Branch band Size s. The conditional probability is also called the mean amount of demand δb,s for Branch band Size s. Additionally we consider the mean amount of demand δbfor Branch b. It is the conditional probability that if an item is sold this happens in Branch b. 1This coincides with the proceeding at our industrial partner. Before the supply according to lots takes place, the purchaser decides on the overall number Dof supplied items. CHAPTER 3. DEMAND ESTIMATION 22 With these quantities the fractional mean demand db,s per branch band size sfor the sets Bof branches and Sof sizes is computed. We outline the approach in Subsection 3.2.1. We show how to include the set Eof scenarios in Subsection 3.2.2 before we split up db,s to the sales periods – in our case weeks – K\ {kmax}and include price dependency for the price indices P\ {pmax}in the subsections 3.2.3 and 3.2.4. In Subsection 3.2.5 we combine the single estimates to arrive at the estimate we use for DISPO. Concluding, in Subsection 3.2.6, we describe how to adapt the demand estimation according to a scenario in effect, as it is the case when we perform price optimization with receding horizon. Demand is always estimated by pooling observations from products from the same commodity group. A smaller division for example to sub commodity groups in our case would lead to an insufficient amount of data. Commodity groups are for example “women overgarments classic” or “women overgarments fashion” or “men trousers”. 3.2.1 Relative demand estimation To exclude the influence of the popularity of the different articles and also to reduce the influence of lost-sales we consider sales just until the day when 50% of the observed product’s overall supply (over all branches and sizes) is sold. The advantage of doing so is that we obtain a measurement of the sales speed for different sizes and branches. If we considered the complete selling time per article varying behavior among branches and sizes – because due to mark-downs nearly all items would be sold out – would not longer be recognizable. For article a∈ A we denote the amount of sales (over all branches and sizes) until the point in time when 50% of the overall supply for aare sold by sal50 a. The amount of sales for Branch bfor the same time frame is denoted by sal50a band for Branch b and Size sby sal50a b,s. Whenever in this subsection we talk about sales we refer always to the point in time until 50% of the supplied items of the considered article are sold. We define with ˆ δa b:= sal50a b sal50a· |B|(3.1) the scaled relative demand for branch bfor Article a– scaled in such a way that the mean of ˆ δa bfor Article aover all branches from the set Btakes value one. Because not every observed product may be delivered to all branches this scaling is necessary to guarantee comparability of the observations. With ˜ δb:= Pa∈A ˆ δa b |A| (3.2) we define the mean relative demand for Branch bin terms of the set of articles A. The mean amount of demand δbfor Branch bis given by the relative frequency δb:= ˜ δb Pb0∈B˜ δb0 .(3.3) Analogously we define the scaled relative demand for Size sin Branch bfor Article aby ˆ δa b,s := sal50a b,s sal50a b · |S|(3.4) CHAPTER 3. DEMAND ESTIMATION 23 and the mean relative demand for Size sin Branch bfor the set of articles Aby ˜ δb,s := Pa∈A ˆ δa b,s |A| .(3.5) The mean amount of demand δb,s for Branch band Size sis given by δb,s := ˜ δb,s Ps0∈S˜ δb,s0 .(3.6) Finally, the mean demand db,s for Branch band Size sby db,s =D·δb·δb,s.(3.7) That means we split up the estimated overall supply Dto the particular sizes and branches in terms of the observed relative frequencies of sales. Example 1 (relative demand estimation).We want to illustrate the approach on a small example. We assume observations of sales according to the third column of the following table. We consider the case of three different observed articles and branches B={b1, b2, b3}. The scaled relative demands per Article aand Branch bare stated in the last column. a b sal50a bˆ δa b a1b16 1.50 a1b22 0.50 a1b34 1.00 a2b13 1.00 a2b24 1.33 a2b32 0.66 a3b12 0.86 a3b22 0.86 a3b33 1.29 In the next table in the second column we stated the mean relative demands per branch in terms of the set Aof observed articles. In the third columns the corresponding mean amounts of demand are stated. b˜ δbδb b11.12 0.37 b20.90 0.30 b30.97 0.33 In the following table in the second column the sales per article for the particular branches and sizes s1, s2, s3, s4are stated. The corresponding scaled relative demands are stated in the fourth column. a b (sal50a b,s1,sal50a b,s2,sal50a b,s3,sal50a b,s4) (ˆ δa b,s1,ˆ δa b,s2,ˆ δa b,s3,ˆ δa b,s4) a1b1(2,2,2,0) (1.33,1.33,1.33,0.00) a1b2(1,1,0,0) (2.00,2.00,0.00,0.00) a1b3(0,2,2,0) (0.00,2.00,2.00,0.00) a2b1(1,2,0,0) (1.33,2.67,0.00,0.00) a2b2(1,1,1,1) (1.00,1.00,1.00,1.00) a2b3(1,0,1,0) (2.00,0.00,2.00,0.00) a3b1(1,1,0,0) (2.00,2.00,0.00,0.00) a3b2(0,1,1,0) (0.00,2.00,2.00,0.00) a3b3(0,2,1,0) (0.00,2.67,1.33,0.00) The mean relative demands per branch and size are stated in the next table in the second column. In the third column the mean amount of demand per branch and size is stated. CHAPTER 3. DEMAND ESTIMATION 24 b(˜ δb,s1,˜ δb,s2,˜ δb,s3,˜ δb,s4) (δb,s1, δb,s2, δb,s3, δb,s4) b1(1.55,2.00,0.44,0.00) (0.39,0.50,0.11,0.00) b2(1.00,1.67,1.00,0.33) (0.25,0.42,0.25,0.08) b3(0.67,1.56,1.78,0.00) (0.17,0.39,0.44,0.00) On the basis of these computations for a given estimated overall supply of D= 20 we compute the mean demand db,s per branch and size: For example, the mean demand for Size s1in Branch b1is given by db,s = 20 ·0.37 ·0.39 = 2.89. The mean demands for all branches and sizes are stated in the following table. b(db,s1, db,s2, db,s3, db,s4) b1(2.89,3.70,0.81,0.00) b2(1.50,2.52,1.50,0.48) b3(1.12,2.57,2.90,0.00) 3.2.2 Regarding different scenarios As a next step we include the consideration of different scenarios. The resulting mean demand de b,s per branch b, size sand scenario eis applied in the SLDP, see Section 2.3. Additionally we estimate the scenario probabilities Prob(e)which are applied in the SLDP and DISPO. At first we determine the realized scenario for each observed article a∈ A. We observe sales until two weeks after sales start. For article a∈ A the number of these sales over all branches and sizes is given by sal2 aand the supply by sup2a. The relation rel2 a=sal2 a sup2a∈[0,1] then indicates the popularity of the article. It is rel2a= 0 if no item is sold in the first two weeks and rel2 a= 1 if all items are already sold out after two weeks. We categorize three different scenarios as they are stated in the following table. rel2 ae <0.33 low seller ≥0.33,≤0.66 normal seller >0.66 high seller The scenario probabilities Prob(e)are given by the relative frequencies of the observed scenario in the historical data. Now we determine how the demand for the low and the high scenario behaves against the normal scenario. This is done the following way: We categorize the set of articles Aby their scenarios. For every scenario e∈ {low seller, normal seller, high seller}thus we obtain a set Ae. The sales until 3 months after sales start for each article over all supplied branches and sizes are given by sala, the supply by supa. The mean relative sales releover all articles for one particular scenario eare given by rele:= Pa∈Ae sala supa |Ae|.(3.8) We compute the change of demand df efor scenario eas df e:= rele relnormal seller (3.9) Thus, df = 1 for the normal seller scenario. The mean demand for Scenario e, Branch band Size sis given as CHAPTER 3. DEMAND ESTIMATION 25 de b,s := df e·db,s.(3.10) 3.2.3 Splitting up the demand to sales periods The next step is to split up the demand to the sales periods. This is done by estimating sales rates from the historical data for all periods – in our case weeks. With sala kwe denote the number of sales for Article aover all branches and sizes in Period k. We denote with stoa kthe overall stock for Article aat the beginning of Period k. For Period kthe relative sales per period rsa kfor Article aare given by rsa k:= sala k stoa k .(3.11) The mean relative sales per period rsa kfor Period kare given by rsk:= Pa∈A rsa k |A| .(3.12) The value rskequals the mean relative amount of sold items depending on the stock at the beginning of Period k. We now convert the rsk, k = 1, . . . , kmax −1to a factor which describes the amount of sold pieces per size and branch depending on the supply. We compute this factor, we call it the sales rate per period srkfor Period k, by Algorithm 1. Algorithm 1 Sales rates per period Require: mean relative sales rsk, k ∈K\ {kmax} Ensure: sales rate srk, k ∈K\ {kmax} 1: init sum = 0 2: init stock = 1 3: for all k= 0, . . . , kmax −1do 4: nr =rsk·stock 5: ˜srk=nr 6: stock =stock −nr 7: end for 8: for all k∈K\ {kmax}do 9: srk=˜srk Pj∈K˜srk 10: end for In the for-loop in Step 3of Algorithm 1the percentage remaining stock (beginning with a stock of one item or 100%) is computed according to the mean relative sales. From the current stock and the mean relative sales the values ˜srkarises. In Step 9of the algorithm the values srkare computed by scaling the values ˜srkin such a way that the sum over all resulting srktakes value one. Example 2 (sales rates).In this example we assume kmax = 4 that means 4real sales periods. We assume that the historical data yields mean relative sales as they are stated below. k0123 rsk0.5 0.7 0.2 0.4 CHAPTER 3. DEMAND ESTIMATION 26 Now we apply Algorithm 1. At first stock is set to value one. In Period 0according to rs0we sell 50% of the current stock or 0.5items. The updated stock is set to 0.5 and ˜sr0= 0.5. In Period 1we start with a stock of 0.5. With rs1= 0.7, we assume that 70% of the current stock that means 0.35 items are sold. This yields a new stock of 0.15. It is ˜sr1= 0.35. We proceed analogously for the last two periods and obtain ˜sr2= 0.03 and ˜sr3= 0.048. We scale the values ˜srkaccording to Step 9of Algorithm 1and obtain the results stated in the following table. k0123 srk0.539 0.377 0.032 0.052 3.2.4 Price-dependent demand An important factor in terms of our demand estimation is the influence of the sales price. To estimate the impact of a mark-down in week kfrom price πp1to πp2with p2> p1the relative sales per week, as defined in the last subsection, for an observed article ain the week before the mark-down rsa k−1and in the week of the mark-down rsa k are compared. The observed increase of sales is given by the factor ˜ elasa πp1→πp2=rsa k rsa k−1. With noπp1→πp2beeing the number of articles for which a mark-down from πp1to πp2 was observed, the mean elasticity elasπp1→πp2is given by Pa∈A ˜ elasa πp1→πp2 noπp1→πp2. 3.2.5 Combining the estimated factors The scenario-, timeand price-dependent demand for a product with starting price π0 for branch band size sgiven by de k,b,s,p =db,s ·df e·srk·Y p0∈P:p0≤p elasπ0→πp0.(3.13) 3.2.6 Updating the scenario In the approach POP-RH, see Subsection 2.5.4, we use latest sales figures to determine the scenario in effect. Before we perform POP-RH to adapt our mark-down policy for the subsequent periods we update demand estimation by adapting the factor for the change of demand from Subsection 3.2.2. This is done by comparing predicted overall sales with realized overall sales. It is salobs kthe amount of sales over all branches and sizes for the last period k. The amount of predicted sales over all branches and sizes for Period kis given by salreal k. Then our updated change of demand df ekis given by df ek=salobs k salreal k . We compute the dependent demands for the next period k+ 1 by dek k+1,b,s,p =db,s ·df ek·srk+1 ·Y p0∈P:p0≤p elasπp0→πp0.(3.14) CHAPTER 3. DEMAND ESTIMATION 27 3.3 Logistic regression With the aim to compare the empirical estimation method outlined in the last chapter with a more common parametric approach we performed ordinal logistic regression. In Subsection 3.3.1 we outline the basic concepts of maximum likelihood estimation – this is the common approach to estimate logistic regression models. We introduce binary logistic regression in Subsection 3.3.2 before we extend it to ordinary logistic regression in Subsection 3.3.3. We are mainly guided by [Har10] and [Rya08]. 3.3.1 Maximum likelihood estimation Maximum likelihood estimation (MLE) is a common statistical method to estimate parameters in generalized linear models. The so-called likelihood function is defined as the joint probability function of the random variables. Assumed we are given oobservations Yi, i = 1, . . . , o and we want to estimate unknown parameters β= (β1, . . . , βn). We denote with fi(y, β)the density function of the random variable yfor the i-th observation. The likelihood for the i-th observation is given by Li(β) := fi(Yi, β). With the assumption that the observations are independent from each other the likelihood function is given by L(β) = o Y i=1 Li(β).(3.15) The log likelihood function is the logarithm of the likelihood function: ln(L(β)) = n X i=1 ln(Li(β)) (3.16) The maximum likelihood is the value of βthat maximizes ln(L(β)) as a function of β. It can be computed by considering the first derivates of Li(β), the so-called score-vector, i.e. the gradient, and the matrix of the second derivates, also known as observed information matrix i.e. the Hessian-matrix. For a maximizing βall first partial derivates have to take value zero while the observed information matrix has to be negative definite. The log likelihood function is applied because in much cases the derivate of the log likelihood function – because of the properties of the logarithm, for example changing from product to sum – is easier to compute as the likelihood itself. And if βmaximizes the log-likelihood function then it also maximizes the likelihood function. If it is not possible to determine the maximum likelihood exactly, so-called iterative trial-and-error methods are used. The most common method is the so-called NewtonRaphson or simply Newton method. This method is a standard approach in nonlinear optimization. The score vector U(β)is locally approximated by a linear function of βin a small region. In the following we denote with I(β)the observed information matrix. With a starting estimate of β(0) of the maximum likelihood βthe linear approximation is given by U(β) = U(β(0))−I(β(0))(β−β(0)).(3.17) Equating to 0and solving by βyields β=β(0) +I−1(β(0))U(β(0)).(3.18) CHAPTER 3. DEMAND ESTIMATION 28 Thus, we determined the null of the linear approximation of the maximum likelihood function in β. At the i-th step we obtain the next estimate by β(i+1) =β(i)+I−1(β(i))U(β(i)).(3.19) If it is the case that the log likelihood worsened at β(i+1) then “step-halfing” is applied, that mean β(i+1) is replaced by β(i)+β(i+1) 2. This is done until the likelihood still is worse than the likelihood at β(i). The method is iterated until convergence. For more details about maximum likelihood estimation see for example [Har10] or [Rya08]. For further reading in terms of the Newton method, we refer to [Rus06]. 3.3.2 Binary logistic Regression The assumption for binary logistic regression is, that the outcome – the dependent yvariable only takes value zero or one. Independent variables xi, i = 1, . . . , k may be binary, integer, continuous or multinomial. We denote the probability that the outcome ytakes value one under the observations xi, i = 1, . . . , k by P(y= 1|x1, x2, . . . , xk) =: P(y= 1). Analogously the probability for ytaking value zero is denoted by P(y= 0|x1, x2, . . . , xk) =: P(y= 0). Binary logistic regression tries to estimate these probabilities. For this purpose the so called odds is defined as odds := P(y= 1) 1−P(y= 1) =P(y= 1) P(y= 0) (3.20) where both probabilities P(y= 1) and P(y= 0) should not take value zero. The odds which takes values ∈[0,∞[is transformed to the so-called logit by the logarithmic function. logit := ln(odds) = ln P(y= 1) 1−P(y= 1)(3.21) In contrast to the odds the logit can take all real values. The logit is estimated by a linear function ln P(y= 1) 1−P(y= 1)=α+β1x1+β2x2+. . . +βkxk.(3.22) Rearranging Equation (3.22) yields P(y= 1) = 1 1 + exp(−(α+Pk i=1 βixi)).(3.23) where the right side is known as the logistic function. This is a nonlinear function. The coefficients αand βi,= 1, . . . , n are estimated by maximum likelihood estimation. Because there are only two possibilities for the outcomes – zero or one – each observation can be seen as a Bernoulli experiment. For a sample of size nand observations o1, . . . , onthe likelihood function therefore is given by L(β) = n Y i=1 Poi i(1 −Pi)1−oi.(3.24) where it is Pi=P(oi= 1) = 1 1+exp(Pk i=1 βixi). CHAPTER 3. DEMAND ESTIMATION 29 3.3.3 Ordinal logistic regression If the dependent variable ycan also take other values than zero or one the binary logistic regression model may be extended to an ordinal logistic regression model. The restriction is that the possible outcomes stand together in an ordinal relationship. If this is not the case one would favor a multinomial logistic regression model. The proceeding in that case is different. For further reading about multinomial regression see for example [Har10]. We assume levels j= 0,1, . . . , jmax for the dependent variable. Now the odds for j= 0,1, . . . , jmax is defined as odds =P(y≥j) 1−P(y≥j)(3.25) and the logit as logit = ln P(y≥j) 1−P(y≥j)=αj+β1x1+β2x2+. . . +βkxk.(3.26) This yields the following coherence: P(y≥j) = 1 1 + exp(−(αj+Pk i=1 βixi)) (3.27) In this model the regression coefficients βiare independent of j. Only the αjwhich describe the intercept depend on j. For a specific jthe formulation is consistent with the binary regression model, see (3.23). The ordinal model can also be estimated via maximum likelihood estimation. For detailed information we refer the reader to [McC80]. 3.4 Applying ordinal logistic regression Now we apply ordinal logistic regression to estimate the sales per periodˆ=week, size, branch and price at our industrial partner. For this purpose we use the environment GNU R for statistical computing2. 3.4.1 Data sample In our first trials GNU R in terms of memory could not deal with the transaction data of the hole tested commodity group. Thus, we will perform all tests on a smaller data set. We consider historical data from the commodity group “women overgarments classic” containing transaction data for about one year starting in September 2009. We restrict the set of about 1400 branches by randomly choosing 30 branches and the set of sizes by randomly choosing 3sizes. (Mostly there are 6or 7different sizes per article in this commodity group.) Moreover we will only regard articles which can be observed at least 13 weeks beginning from the sales start. We divided the remaining 116 articles randomly aiming at ratio 70 : 30. This yields a set of 81 articles for estimation and a set of 35 articles for validation. We categorize our observations in observations per article, week, branch and size. For our test set of 81 articles this results in 81 ·13 ·30 ·3 = 94 770 data points. 2GNU R version 2.10.2 CHAPTER 3. DEMAND ESTIMATION 30 Frequencies of Responses 0 1 2 3 4 90326 4118 299 17 10 Table 3.1: Estimation – Responses Analogously we get 40 950 observations for the validation set. In our data set there are five types of responses. The dependent variable, the number of sales, takes value 0,1,2,3or 4. The frequencies are stated in Table 3.4.1. 3.4.2 Choosing the model Analogously to the empirical estimation we include popularity of the article, week, price, size and branch in the estimation and interpret the branches and sizes as categorical incomes. In contrast to the empirical method now we also consider current stock as input. At this point we want to estimate realized sales and not demand. We have to decide on how the price is adopted in the estimation. We are given the starting price and the realized price per week. We tried three different models in terms of the adoption of the price which include 1. only realized price (Model 1), 2. only ratio realized price divided by starting price (Model 2), 3. ratio and starting price (Model 3). as independent variables. We will present the complete Model 3 before we comment on statistical tests in terms of the model fitness for all models. We denote by sa the number of sold items and by at the popularity of the corresponding article. Analogously to the empirical estimation the popularity is computed by dividing the sales of the first two weeks by the stock at the beginning of the first week. With we we denote the corresponding week and with st the stock at beginning of the week and by pr the sales price and by sp the starting price. The additional binary variables {cb},{cs},{cb,s}take value one if and only if the observation concerns Branch band/or Size s. The probability P(sa ≥j)with j > 0after Model 3 is given as follows (P(sa ≥0) by definition takes value one): P(sa ≥j) = 1 1 + exp(−αj−βX),(3.28) where βX =βat at +βwe we +βst st +βsp sp +βpr pr sp +X b∈B βb{cb}+X s∈S βs{cs}+X β∈B X β∈S {cb,s} We estimate the coefficients of the model under use of GNU R. For maximum likelihood estimation we used the function lrm (logistic regression model) with default settings from the Design-package which was implemented by Frank Harrell, see also [Har10]. To evaluate the fitness of the model we state results from some statistical tests. These tests also were performed under use of the Design-package. Chapter 4 Price Optimization In this chapter we elaborate on the Price Optimization Problem, how it takes part in DISPO and present algorithms to solve it. We restrict us to the case that the set of scenarios contains only one scenario. This is the case when we perform price optimization with receding horizon POP-RH, Subsection 2.5.4: The set of scenarios consists only of the current scenario in effect. In Section 4.1 we introduce a mixed-integer programming formulation for the Price Optimization Problem for one fixed scenario ˆe, the POPˆe. To describe the situation at our industrial partner as correctly as possible, we extend price optimization as it is performed by the former DISPO-team by adding mark-down costs depending on the current stock. This leads to nonlinearity of the underlying mixed-integer program. We adapt the mark-down policy every period/week, see Subsection 2.5.4 or Figure 1.1. The last sales for the week can be observed on Saturday evening. On Monday morning already the decision for mark-downs has to be made. At our industrial partner weekly more than 4000 articles have to be considered. Therefore the solving process of price optimization must not last too long. Even without regarding stock-depending costs for mark-downs solving price optimization with state-of-the-art MIP solvers – because of the long computation times – is not suitable in terms of these real-world requirements. In the MIP formulation of the POPˆethe variables are finely grained: For every period and price index there exists a binary variable that indicates if the related price is assigned to the period or not. One idea would be to enumerate all possible combinations of these variables – the so-called price trajectories, Section 4.2. In Section 4.3 we introduce some basics of dynamic programming which we apply in Section 4.4 on the POP-RH: We outline how to generate price trajectories for a given supply dynamically. For this purpose we develop dominance rules to exclude truncated price trajectories which will not lead to an optimal solution from further consideration, Section 4.5. The result is a so-called label setting algorithm for POPˆe. The detailed implementation is outlined in Section 4.6. We illustrate the algorithm on a small example, Section 4.8, that we introduce in Section 4.7 and that will accompany us in the remainder of the thesis. We state computational results in Section 4.9 and conclude the chapter in Section 4.10. We will apply the enumeration of price trajectories in our exact Branch&Bound solver for ISPO which is presented in Section 9.1. The dynamic generation of price trajectories – besides price optimization with receding horizon, POP-RH – takes place in our heuristic approach for ISPO, see Section 9.2. 37 CHAPTER 4. PRICE OPTIMIZATION 38 4.1 Extending POP by mark-down costs – a mixed-integer nonlinear program In this section we state price optimization for one fixed scenario ˆe, Subsection 4.1.1. We will not go into details in terms of the problem specification at this point and refer the reader to Problem 5. Because we fix scenario ˆein the POPˆewe drop the scenario index for the corresponding variables at this point. In Subsection 4.1.2 we go into the case that mark-down costs depending on the corresponding stock arise as it is the case at our industrial partner. Inclusion of these costs leads to nonlinearity of the introduced model. 4.1.1 Problem formulation We formulate the Price Optimization Problem for a fixed scenario ˆewith the inclusion of mark-down costs as follows: Problem 5 (POPˆe). max X k∈K exp(−ρk)X b∈BX s∈S rk,b,s −µkβk(4.1) subject to X p∈P uk,p = 1 ∀k∈K, (4.2) uk,0= 1 ∀k∈K:k < kobs,(4.3) ukmax,pmax = 1,(4.4) uk−1,p1+uk,p2≤1∀k∈K:k > 0, p1, p2∈P:p2< p1,(4.5) βk≥uk−1,p1+uk,p2−1∀k∈K:k > 0,∀p1, p2∈P:p26=p1, (4.6) v0,b,s =Ib,s ∀b∈B, s ∈S, (4.7) vk−1,b,s −vk,b,s =X p∈P wk−1,b,s,p ∀k∈K:k > 0, b ∈B, s ∈S, (4.8) X p∈P wk,b,s,p ≤vk,b,s ∀k∈K, b ∈B, s ∈S, (4.9) wk,b,s,p ≤uk,p ·dk,p,b,s ∀k∈K:k < kmax, b ∈B, s ∈S, p ∈P:p<pmax,(4.10) wkmax,b,s,pmax =vkmax,b,s ∀b∈B, s ∈S, (4.11) rk,b,s =X p∈P πp·wk,b,s,p ∀k∈K, b ∈B, s ∈S, (4.12) uk,p ∈ {0,1} ∀k∈K, p ∈P, (4.13) βk∈ {0,1} ∀k∈K, p ∈P, (4.14) wk,b,s,p ≥0∀k∈K, b ∈B, s ∈S, p ∈P, (4.15) vk,b,s ≥0∀k∈K, b ∈B, s ∈S, (4.16) rk,b,s ≥0∀k∈K, b ∈B, s ∈S. (4.17) CHAPTER 4. PRICE OPTIMIZATION 39 A mark-down in period kis indicated by the dependent binary variable βk, which is forced to one by Inequality (6.19) if the price compared to the previous period has changed. In the objective the mark-down costs µkfor every period are subtracted from the revenue. 4.1.2 Nonlinearity by mark-down costs We want to put a finer point to the mark-down costs µkfor Period k. Mark-down costs divide into two parts. On the one side there are fixed mark-down cost. If there is a mark-down in a period always cost of µfoccur independent from the number of items that have to be marked down. On the other side there are variable mark-down cost: Every single item that has to be marked down causes cost of µv. At our partner fixed mark-down costs arise from all actions that are necessary to inform the branches about mark-downs, the variable mark-down costs accrue from pricing the items in the branches by the sales personnel. There is an exception for the last/sellout period kmax. We assume that in the sellout process still qkmax mark-downs are necessary. Only variable mark-down costs arise for the sellout period. Altogether the mark-down costs µkin the real sales period kare given by µk=µf+µvX b∈BX s∈S vk,b,s,∀k∈K\ {kmax}.(4.18) For the sellout period kmax the mark-down costs are given by µkmax =qkmax µvX b∈BX s∈S vkmax,b,s.(4.19) Extending the formulation of Problem 5by Constraints (4.18) and (4.19) leads to a mixed-integer nonlinear program because in the objective function we have to multiply the binary variables βkvia µkwith the real variables vk,b,s. In the remainder we will regard the mark-down costs as formulated by the constraints (4.18) and (4.19). Henceforth, whenever we will mention the Price Optimization Problem, we refer to the POPˆe 4.2 Enumerating price trajectories We can enumerate the possible price trajectories by branching on the decisions of the price optimization stage. The variables for the price optimization stage are finely grained – for every period kand price index pthere is a binary variable uk,p which indicates whether the price with index pis assigned to k. Now we consider more widescale decisions. A natural idea is to condense the mark-down decisions in each time period to an entire price trajectory for the complete selling time. Definition 1 (price trajectory, revenue of a price trajectory).We define a price trajectory t= (t0, . . . , tkmax )as a kmax +1-tuple where each entry is a price index p∈P. It is tk= 0 for k < kobs and tkmax =pmax. Moreover it is tk≤tk+1,∀k∈K\ {kmax}. That is a price trajectory equals a valid assignment of the uk,p variables in Problem 5 where tk=pif and only if uk,p = 1. The revenue aresulting from a price trajectory CHAPTER 4. PRICE OPTIMIZATION 40 which equals the related objective value of Problem 5is given by a:= kmax−1 X k=0 exp(−ρk) πtkX b∈BX s∈S min   max   Ib,s − k−1 X j=0 dj,b,s,tj,0  , dk,b,s,tk   −βk µf+µvX b∈BX s∈S max{Ib,s − k−1 X j=0 dj,b,s,tj,0} ! + exp(−ρkmax)(πkmax −qkmax µv)·X b∈BX s∈S max (Ib,s − kmax−1 X k=0 dk,b,s,tk,0) (4.20) According to the constraints (4.7), (4.8), (4.9) and (4.10)max{Ib,s −Pk−1 j=0 dj,b,s,tj,0} equals the stock vk,b,s for Size sin Branch bat the beginning of Sales period k. The number of sold items for the real sales periods is the minimum of the current stock and demand, (4.9) and (4.10) together with the objective function (4.1). At the sellout period kmax all remaining items are sold. Subtracting mark-down costs from the related yield results in the last line of Equation (4.20). Now we want to deduce the general number of all valid price trajectories. For this purpose we consider a small example. Example 3 (number of price trajectories).It is kmax = 4, i.e. |K|= 5, and pmax = 3, i.e. |P|= 4, and kobs = 2. We depict the sales periods by their indices and encode a mark-down after Period kin Period k+ 1 to the next price index tk+ 1 by the symbol ?, a mark-down to the after next price index tk+ 2 by ?? and so on. Then these are all valid price trajectories together with their encoding: price trajectory encoding 0 0 0 0 3 0 1 2 3 ? ? ? 4 0 0 0 1 3 0 1 2 ?3? ? 4 0 0 0 2 3 0 1 2 ? ? 3?4 00113 01?2 3 ? ? 4 00123 01?2?3?4 00223 01? ? 2 3 ?4 Using the encoding scheme from Example 3we can establish the number of all valid price trajectories. Theorem 2 (number of price trajectories [KKR11b]).The number of all valid price trajectories for Problem 5is given by kmax −kobs +pmax −1 pmax −1.(4.21) Proof. We encode the feasible price trajectories by inserting pmax symbols for markdowns, like e.g. ?as in Example 3. In the first kobs observation periods no mark-down is possible, so there is no symbol ?between the related places in the encoding. The last mark-down after Period kmax −1to the salvage value is determined. So we have to distribute our remaining pmax −1mark-downs/symbols ?among the remaining kmax −kobs +pmax −1places. This yields the claim. CHAPTER 4. PRICE OPTIMIZATION 41 period/depth-1 0 1 2 3 4 0 0 1 02 01 2 1 2 2 3 3 3 3 3 3 Figure 4.1: POP – enumeration tree The valid price trajectories can be established by walking through an enumeration tree. The nodes of the tree in depth k+1 with 0≤k < kmax −1correspond with fixed prices/price indices for the first k+1 periods. The leaves at depth kmax +1 correspond with valid price trajectories. In every depth k+ 1 with 0≤k < kmax −1we consider all extensions with price indices p=tk−1, .. . , pmax −1. At depth kmax + 1 – the sellout period – we have to fix the salvage value, i.e. price index pmax. We consider the enumeration tree for a small example. Example 4 (enumeration tree for POP).It is kmax = 4, kobs = 2 and pmax = 3. In Figure 4.1 the corresponding enumeration tree is depicted. The number inside a node at depth k+ 1 is equivalent to the entry tkin the price trajectory. Now we state some corollaries in terms of the enumeration tree which follow from Theorem 2. Corollary 2 (number of nodes per depth).The number of nodes in depth k+ 1 with kobs ≤k < kmax in the enumeration tree amounts to k−kobs +pmax pmax −1.(4.22) Proof. This Corollary follows easily from Theorem 2because it is equivalent to consider the enumeration tree for the same data, but kmax =k+ 1. Corollary 3 (number of nodes).The overall number of nodes in the enumeration tree amounts to kobs + kmax−1 X k=kobs k−kobs +pmax pmax −1+kmax −kobs +pmax −1 pmax −1.(4.23) Proof. This claim – more precisely the second term – follows from Corollary 2. The first term describes the number of nodes for the observation time kobs. Because for periods kwith k < kobs the starting price has to be maintained there is always one node for these periods. The third term equals the number of nodes at depth kmax + 1 which is the number of all valid price trajectories. CHAPTER 4. PRICE OPTIMIZATION 42 Corollary 4 (number of induced price trajectories).The number of induced price trajectories by a node at depth k+1 with kobs ≤k < kmax and tk=pin the enumeration tree is given by kmax −k+pmax −p−2 pmax −p−1.(4.24) Proof. We consider truncated price trajectories ending up with tk=p. The number of induced price trajectories is the number of all valid extensions. So it is equivalent to determine the number of price trajectories for kmax =kmax −k,kobs = 1 and pmax =pmax −p. 4.3 Excursus: Dynamic programming As already mentioned in the introduction a general approach for inventory and pricing problems is dynamic programming. Dynamic programming is based on the Bellman’s optimality principle which roughly says that for a dynamic system (Section 4.3.1) every optimal solution consists of optimal partial solutions. This leads to a backwards dynamic programming algorithm which we outline in Section 4.3.2. While this algorithm computes the optimal partial solutions backwards in time for the special case of deterministic problems an algorithm performing forwards in time can be stated. This is outlined in Section 4.3.3. Deterministic dynamic problems can be reduced to shortest path problems and common algorithms for solving shortest path problems can be applied, Section 4.3.4. Sometimes the state-space for dynamic programs is restricted by resource constraints. Thus, in Section 4.3.5 we consider the case of a resource constraint shortest path problem. Because the in Section 4.3.4 proposed methods only regard the length of the partial path they are not suitable for this problem formulation. The explained approach is extended to a label setting algorithm. Now each partial path gets a label which includes the length of the path and the still available amounts of the resources. The state space is reduced by comparing the labels. If one label can not lead to a better solution than the other one it is said to be dominated. Dominated labels can be excluded from further consideration. In this section we are mainly guided by [Ber05]. 4.3.1 General dynamic program We consider a system of the form xk+1 =fk(xk, uk, wk), k = 0,1, . . . , N −1(4.25) where kis a discrete time index, xkis a state of the system for stage k,ukis the decision variable or control which is selected at time kand wkis a random parameter. The number of stages is stated by Nwhich is also denoted as horizon. With the function fkthe dynamic of the system is described. Additionally we are given a cost function gk(xk, uk, wk). The total costs are given by gN(xN) + N−1 X k=0 gk(xk, uk, wk).(4.26) gN(xN)is also called terminal cost. CHAPTER 4. PRICE OPTIMIZATION 43 With Skwe denote the state space of xk. It is xk∈Skand analogously we consider a space Ckwhere uk∈Ck. The disturbance wkis an element from a space Dk. A control is called admissible if uk∈U(xk)where U(xk)⊂Ck. That means the admissibility of a control at stage kdepends on the state xkin this stage. The control ukis selected with the knowledge of the current state xk. A policy or control law is a sequence of functions π={µ0, . . . , µN−1}(4.27) where µkmaps the state xkinto controls uk=µk(xk). The goal is to minimize the expected cost Jπ(x0)of πstarting at stage x0which is given by Jπ(x0) = E(gN(xN) + N−1 X k=0 gk(xk, µk(xk), wk)).(4.28) We only consider admissible policies, that means policies with µk(xk)∈Uk(xk) ∀xk∈Sk. The set of all admissible policies is denoted by Π. An optimal policy π∗is a policy πthat minimizes the costs, that means Jπ∗(x0) = min π∈ΠJπ(x0).(4.29) 4.3.2 The dynamic programming algorithm The techniques to solve dynamic programs are based on the principle of optimality stated first by Richard Bellman [Bel10]. Definition 2 (principle of optimality).It is π∗= (µ∗ 0, µ∗ 1, . . . , µ∗ N−1)an optimal policy. The truncated policy (µ∗ i, µ∗ i+1, . . . , µ∗ N−1)is also optimal for the subproblem to minimize the expected cost E(gN(xN) + N−1 X k=i gk(xk, µk(xk), wk))(4.30) from Stage ito Stage N. We now denote with Jk(xk)the optimal expected cost for starting at Stage k. With the above principle for every initial state x0the optimal cost J∗(x0)equals J0(x0)and is given by the last step of the following algorithm. The algorithm proceeds backwards in time from Stage N−1to Stage 0: JN(xN) = gN(xN),(4.31) Jk(xk) = min u∈Uk(xk) Ewk{gk(xk, uk, wk)+Jk+1(fk(xk, uk, wk))}, k = 0,1, . . . , N−1. (4.32) 4.3.3 Deterministic Systems In this Section we focus on deterministic problems. These are problems where the disturbance wktakes only one value. This may result from the approximation of a CHAPTER 4. PRICE OPTIMIZATION 44 stochastic problem. For deterministic problems for a given policy (µ0, . . . , µN−1)and the initial state x0the future states are predictable by xk+1 =fk(xk, µk(xk)), k = 0,1, . . . , N −1(4.33) and the corresponding controls are given by uk=µk(xk), k = 0,1, . . . , N. (4.34) A deterministic dynamic program can be seen as a shortest path problem in a directed graph with nodes corresponding to stages. The source scorresponds to State x0 while the sink is an artificial terminal node tthat describes the state after adding the terminal costs. The inner nodes correspond to the stages 1,2, . . . , N. There are only arcs between nodes corresponding to state xkand xk+1, k = 0, . . . , N −1. These arcs describe a transition of the form xk+1 =fk(xk, uk). The length of an arc is given by the transition cost gk(xk, uk). Moreover every node related to state xNis connected with the sink t. The corresponding length of the arc is the terminal cost gN(xN). With this reduction solving a dynamic program to optimality is the same as finding the shortest path in the corresponding graph. This leads to an forward algorithm for the dynamic program what means that we compute optimal partial solutions beginning from Stage 0and ending up at Stage N. With ak ij we denote the cost of transition from Stage kand State i∈Skto State j∈Sk+1. The terminal cost of State i∈SNare denoted by aN ij . It is ˜ JN(j) = a0 sj, j ∈S1(4.35) and ˜ Jk(j) = min i∈SN−k [aN−k ij +˜ Jk+1(i)].(4.36) The optimal cost are given by ˜ J0(j) = min i∈SN [aN ij +˜ J1(i)].(4.37) 4.3.4 Solving shortest path problems In the previous section we stated the context of deterministic dynamic programming and shortest path problems and a forward algorithm which can be seen as a general approach to solve shortest path problems. We consider a graph where we want to find the shortest path from a source node sto a sink node t. The length of the path results as the sum of the lengths of the traversed arcs. The problem can be solved to optimality by a so-called label correcting algorithm. The idea is to discover shorter paths from the source sto every other node jand to maintain the length of the shortest path found so far in a variable djwhich is called the label of j. We start from the source s, Step 2of Algorithm 2, and extend our partial step-bystep to a path ending up at the sink t. For this purpose we consider all possible arcs starting at the end node iof our partial path, Step 5. Whenever a shorter path from the sink to a node jis found, the label is corrected in Step 7, i.e. we always consider only the shortest partial path from the source to node j. Because according to the Bellman’s optimality principle each path consists of optimal partial paths we will end up with a shortest path from the source sto the sink t. CHAPTER 4. PRICE OPTIMIZATION 45 Algorithm 2 Label correcting 1: init dj=∞for all nodes j 2: init OPEN={s} 3: while OPEN6=∅do 4: choose node ifrom OPEN 5: for all childs jof ido 6: if di+aij <min{dj,UPPER}then 7: dj=di+aij 8: if j /∈OPEN and j6=tthen 9: place jin OPEN 10: else 11: if j=tthen 12: UPPER=di+aij 13: end if 14: end if 15: end if 16: end for 17: remove ifrom OPEN 18: end while There are different ways to perform the label correcting algorithm. For example, one could traverse the nodes from the set OPEN in a breadth-first search, also known as Bellman-Ford method with complexity O(nm)where nis the number of nodes and m the number of arcs. Or one can perform a depth-first search with the same complexity but with the advantage that the amount of needed memory is less. By a best-first search, also denoted by Dijkstra’s method the complexity only amounts to O(nlogn+m). But Dijkstra’s algorithm in general works only correctly if there are no negative arc lengths. (Otherwise Bellman’s optimality principle may be violated, because with negative arc lengths an optimal path has not necessarily to consist of shortest partial paths.) If the graph contained negative cycles then the label correcting approach would not terminate. Traversing negative cycles would always reduce the length. The Bellman-Ford method can detect negative circles. (In the case of dynamic programming where there are only forwards arcs – from stage kto stage k+ 1 – no cycles can occur.) The label correcting method can also be extended to a Branch&Bound method where in comparison with the bound UPPER solutions are discarded that have no chance to be optimal. For further reading about shortest path problems we refer to [CGR93]. 4.3.5 Resource constraint shortest path problems and dominance We now deal with shortest path problems with one or more additional resource constraints. For each resource jan initial stock R(j) init is given. Each arc in the graph, see Section 4.3.3, consumes always an amount of the given resources. Now, a path is only valid if the totally amounts of each resources does not violate the resource restriction – i.e. the sum of the consumed amounts of resource jover all arcs in the path must not exceed R(j) init . For resource constraint shortest path problems the algorithms stated in the last section are not convenient. In a simple shortest path problem the shortest path is always the best. For resource constraint shortest path problems this path might violate the resource constraints. Thus, we can not exclude longer paths being optimal as it is implied by Step 7of Algorithm 2. An idea would be to save all possible paths. But according to the problem size this might be inefficient in terms of time and impossible in terms of memory. Handler and Zang [HZ80] showed that the resource constraint shortest path prob- CHAPTER 4. PRICE OPTIMIZATION 46 lem in general – also in our case where no cycles in the graph appear – is NP hard by reducing the knapsack problem to it. The same is shown by Garey and Johnson [GJ79]. But they reduced the partition problem on it. Irnich and Desaulniers [ID05] among others covered so-called dominance rules for resource constraint shortest path problems. The label correcting algorithm above is adapted to a so-called label setting algorithm: Not only the shortest path to a node is regarded but also longer paths which cannot be excluded from being optimal. In our case we are given a resource constraint shortest path problem with nresources. The initial stock for each resource j= 1, . . . , n is given by R(j) init . For each resource j= 1, . . . , n and each arc from Node k1to Node k2a consumption or weight c(j) k1k2is given. With L= (d, R, k)(4.38) we define a label for a node. The first element dof the triple is the length of the path from the sink sto the related node k.R= (r(1), r(2), . . . , r(n))is an n-tuple where r(j)is the still available amount of resource j. We start with the label Ls= (ds, Rs, s)with ds= 0,(4.39) r(j) s=R(j) init , j = 1, . . . , n. (4.40) A label Li1= (di1, Ri1, k1)is extended to a label Li2by setting Li2= (di2, Ri2, k2)(4.41) where di2=di1+ak1k2,(4.42) r(j) i2=r(j) i1−c(j) k1k2, j = 1, . . . , n. (4.43) A label Ln1– which stands for a partial path starting at the sink sand ending at node k1– is said to be dominating over a label Ln2ending at the same node if dn1< dn2 and r(j) n1≥r(j) n2for all j= 1, . . . , n. Because of the higher amount of the resource and both labels ending up at the same node we can extend the path described by label Ln1 in each way we can extend Ln2. And, because dn1< dn2each path containing the to Ln2related partial path cannot be shorter than each path containing the to Ln1related partial path. Thus, we can exclude the label Ln2from further consideration. We say, Ln2is dominated by Ln1. Many further references regarding dominance and dominance rules can be found in literature. We just mention some few of them exemplary. A definition of dominance was given by Manne [Man58] already in the year 1958. In [JC11] dominance rules in combinatorial optimization and their characteristics are generally defined and studied. Fischetti and Salvagnin presented a dominance procedure for general mixed-integer linear programs in [FT88]. General results for applying dominance in Branch&Bound algorithms are presented in [Iba77]. 4.4 Dynamic generation of mark-down strategies In the POPˆewe deal with fixed supply for each branch and size, Constraint (4.7). We apply a dynamic programming approach, see the last section, to solve Problem 5. In CHAPTER 4. PRICE OPTIMIZATION 53 4.7 An accompanying example In order to illustrate basic ideas and algorithms we introduce a small (manageable) instance for ISPO on which we will draw on in the remaining course of the thesis. We consider just two branches and sizes. The selling time amounts to five periods. There are four different prices, salvage value included. This implies at most two possible markdowns during the sales process. Concisely, our data is as follows: • |B|= 2, B ={1,2} • |S|= 2, S ={S,L} •kmax = 4 •kobs = 2 •pmax = 4,π0= 10.99, π1= 5.99, π2= 1.99, π3= 0.99 •µf= 1, µv= 0.10 •ρ= 0.01 •qkmax = 2 •ap = 0.5 •E={low seller,normal seller,high seller}with Prob(low seller)=0.2,Prob(normal seller) = 0.5,Prob(high seller) = 0.3. •κ= 1 •δ1= 1.5 •pcost = 0.01 •I= 5, I = 10 •vmin = 1, vmax = 1,vlmin = 2,vlmax = 3 In the subsequent example we will assume a fixed supply per branch b∈Band size s∈Swhich is given by the following table: b/sS L 1 5 5 2 3 8 In the following tables we state the demands per price index pand period kfor every pair of branches and sizes (b, s)for the “normal seller” scenario. (1,S) k/p0 1 2 0 2 0 0 1 1.5 0 0 2 1 1.5 2 3 0.5 1 1.5 (1,L) k/p0 1 2 0 3 0 0 1 2 0 0 2 1 2 3 3 0.7 0.8 0.9 CHAPTER 4. PRICE OPTIMIZATION 54 (2,S) k/p0 1 2 0 1 0 0 1 0.9 0 0 2 0.8 0.9 1 3 0.7 0.8 0.9 (2,L) k/p0 1 2 0 2 0 0 1 1 0 0 2 0.9 1.2 1.5 3 0.5 1 1.4 The demands for the other scenarios are given by dlow seller k,b,s,p = 0.5·dnormal seller k,b,s,p ∀k∈ {0,1,2,3}, s ∈ {S,L}, b ∈ {1,2}, p ∈ {0,1,2}, (4.56) dhigh seller k,b,s,p = 1.3·dnormal seller k,b,s,p ∀k∈ {0,1,2,3}, s ∈ {S,L}, b ∈ {1,2}, p ∈ {0,1,2}. (4.57) 4.8 POP-DYN applied on the accompanying example We depict the basic ideas of our Algorithm 3, POP-DYN, with help of the example from the last section. Example 5 (solving POPˆeby Algorithm 3).All valid (partial) mark-down strategies for the example from Section 4.7 are pictured in the enumeration tree in Figure 4.2. Revenue and stock can be found beside the nodes. Nodes which can be pruned by dominance are colored gray. We just consider the scenario “normal seller”. We start by extending the partial mark-down strategy ˜ P= (˜a, ˜ t, ˜r)with ˜a= 0,˜ t= () and ˜rb,s =Ib,s,∀s∈S, ∀b∈B. Because kobs = 2 we have to maintain the starting price with Price index 0 in the two first periods. For Branch 1and Size S the demand in Period 0 amounts to 2, the current stock is 5. That means all the demand for this branch and size is met and we obtain a revenue of exp(−0·0.10) ·2·10.99 = 21.98 for this branch and size. The revenue for Branch 1and Size L is given by exp(−0·0.01) ·3·10.99 = 32.97. For Branch 2we earn 10.99 and for Size L 21.98. There is no mark-down in this period. Hence, the revenue for all branches and sizes is the sum of the single revenues, namely 87.92. We get to Period 1. For Branch 1 and Size S there are 3 items left, for Size L 2. In Branch 2 there are still 2items of Size S and 6 items of Size L available. Thus, all demand of this period can be met and our updated revenue at the end of Period 1 is 87.92 + exp(−1·0.01) ·10.99 ·(1.5+2+0.9 + 1) = 146.68. Because we traverse the enumeration tree by depth-first-search we first consider the extension of the current partial mark-down strategy with Price Index 0. In Branch 1 the demand for Size S can be met, Size L is sold out. The demands for both sizes in Branch 2 can be fully met. So the revenue for the current partial mark-down strategy amounts to 146.68 + exp(−2·0.01) ·10.99 ·(1 + 0 + 0.8+0.9) = 175.76. Now we perform the dominance check. Because we stand in the first branch of the enumeration tree pruning is not yet possible. But we are able to update our bounds abound 2,p , p = 0,1,2which still take value −∞. For the current period k= 2 and price index pthere could be only one additional mark-down in the extension (see 4.55) and with the maximum mark-down costs that can arise in Period kmax we get aµ= exp(−3·0.01) ·(0.5+0+0.3+4.1) ·0.1 + exp(−4·0.01) ·((0.5+0+0.3+4.1) ·2·0.1+2·1) = 3.34. CHAPTER 4. PRICE OPTIMIZATION 55 So the updated bound amounts to abound 2,0= 175.76 −3.34 = 172.42. For p= 1 we have to regard that in comparison to a mark-down strategy of the same length with the current last price index 0 an additional mark-down may be possible, that means we update our bound to abound 2,1= 175.76 −exp(−3·0.01) ·(0.1·(0.5+0+0.3+4.1) + 1) + exp(−4·0.01) ·((0.5+0+0.3+4.1) ·2·0.1+2·1) = 171.45. All new bounds for Period 2 are outlined in the following table: p0 1 2 abound 2,p172.42 171.45 171.45 At Period 3 we extend the current partial mark-down strategy again with Price Index 0. In Branch 1 the demand of Size S can be met while Size L is as always sold out. In Branch 2 we meet the demand only partly and sell all remaining items of Size S. The demand for Size L can be fully met. Therefore the revenue is given by 175.76 + exp(−3·0.01) ·10.99 ·(0.5+0+0.3+0.5) = 189.63. Now we come to extend the partial mark-down strategy to a complete one. For every remaining item on the one hand we earn the salvage value on the other hand we have to pay two times variable mark-down cost. We obtain a complete strategy with revenue 189.63 + exp(−4·0.01) ·((0.99 −2·0.1) ·(0+0+0+3.6) −2·1) = 190.44. We set P∗= (190.44,(0,0,0,0,3)). The next step in the depth-first-search is to extend the partial mark-down strategy with partial price trajectory (0,0,0) by Price index 1. This implies a mark-down at the beginning of Period 3 and together with the yield earned by the sales the corresponding revenue amounts to 175.76 + exp(−3·0.01) ·(5.99 ·(0.5+0+0.3 + 1) −2·1−2· 0.1·(0.5+0+0.3+4.1)) = 184.78. Extending the current partial mark-down strategy to a complete one yields a revenue of 185.21. We have to continue the algorithm without the possibility to prune nodes until our partial mark-down strategy will end up in Period 2 with Price Index 1. The revenue for this strategy amounts to 166.09. Because 166.09 <171.45 = abound 2,1the current branch of the tree can be pruned. This is the same with Period 2 and Price Index 2. We end up with the optimal mark-down strategy P∗= (190.44,(0,0,0,0,3)). Algorithm 3 POP-DYN Require: complete data of an instance of the POPˆe Ensure: optimal mark-down strategy P∗= (a∗, t∗) 1: init abound k,p =−∞,∀k < kmax −1,∀p<pmax, init ˜a= 0 ,˜ t= (),˜rb,s =Ib,s ∀b∈B, ∀s∈S, init P∗= (a∗, t∗)with a∗=−∞ and t∗uninitialised 2: set ˜ P= (˜a, ˜ t, ˜r) 3: POP-DFS( ˜ P) 4: return P∗ CHAPTER 4. PRICE OPTIMIZATION 56 Algorithm 4 POP-DFS Require: a partial mark-down strategy ˜ P= (˜a, ˜ t, ˜r)with ˜ t= (t0,...,tk) Ensure: valid non-dominated extensions of ˜ P 1: if k < kmax −1then 2: if POP-DOM( ˜ P)=true then 3: return 4: end if 5: for all valid extensions ˜ Pext of ˜ P(see Definition 5)do 6: POP-DFS( ˜ Pext ) 7: end for 8: else 9: extend ˜ Pto obtain the complete mark-down strategy P= (a, t)(see Definition 6) 10: if a>a∗then 11: P∗=P 12: end if 13: end if Algorithm 5 POP-DOM Require: a partial mark-down strategy ˜ P= (˜a, t, ˜r)with ˜ t= (t0,...,tk), for each price index p:p<pmax a bound abound k,p Ensure: possibly updated bound abound k,p for tk≤p<pmax, ˜ Pdominated? true or false 1: if ˜a≤abound k,tkthen 2: return true 3: else 4: for all p:p≥tk, p < pmax do 5: compute the maximal additional mark-down costs ˜aµfor extending ˜ Pby all valid extensions starting with the price index p(see Lemma 1) 6: if ˜a−˜aµ> abound k,p then 7: abound k,p = ˜a−˜aµ 8: end if 9: end for 10: return false 11: end if CHAPTER 4. PRICE OPTIMIZATION 57 period/depth-1 0 1 2 3 4 0˜a=87.92 ˜rb,s S L 1 3 2 2 2 6 0˜a=146.68 ˜rb,s S L 1 1.5 0 2 1.1 5 1˜a=166.09 ˜rb,s S L 1 0 0 2 0.2 3.8 0˜a=175.76 ˜rb,s S L 1 0.5 0 2 0.3 4.1 2˜a=152.75 ˜rb,s S L 1 0 0 2 0.1 3.5 0˜a=189.63 ˜rb,s S L 1 0 0 2 0 3.6 1˜a=184.78 ˜rb,s S L 1 0 0 2 0 3.1 2˜a=178.56 ˜rb,s S L 1 0 0 2 0 2.7 1˜a=173.06 ˜rb,s S L 1 0 0 2 0 2.8 2˜a=167.82 ˜rb,s S L 1 0 0 2 0 2.4 2˜a=155.65 ˜rb,s S L 1 0 0 2 0 2.1 3a=190.44 3 a=185.21 3 a=178.69 3a=173.27 3 a=167.72 3 a=155.32 Figure 4.2: Example for the dynamic generation of mark-down strategies, Algorithm 3, POP-DYN CHAPTER 4. PRICE OPTIMIZATION 58 Instance scenario low scenario normal scenario high nd d nd d nd d time(s) time(s) %n time(s) time(s) %n time(s) time(s) %n 1 3.47 1.10 34.70 3.29 0.68 21.23 3.45 0.43 13.47 2 3.66 0.71 21.89 4.07 0.37 10.18 3.85 0.24 7.98 3 3.55 0.99 30.97 3.48 0.80 21.60 3.95 0.47 14.79 4 3.55 1.16 22.53 4.07 0.73 22.33 3.37 0.53 14.93 5 4.06 0.79 24.23 3.24 0.52 16.76 3.41 0.36 11.79 6 3.57 1.17 35.80 4.05 0.85 22.84 3.48 0.50 16.32 7 4.57 0.74 21.96 4.46 0.52 15.37 3.50 0.40 11.86 8 2.86 0.90 28.26 2.92 0.51 20.50 3.31 0.32 12.66 9 3.23 0.85 27.60 2.88 0.51 19.62 3.30 0.28 11.64 12 2.96 0.93 30.01 3.35 0.62 22.62 2.80 0.38 24.49 ∅3.55 0.93 27.80 3.58 0.61 19.31 3.44 0.39 13.99 Table 4.1: Dynamic generation of mark-down strategies – computational results 4.9 Computational results In Table 4.1 we compare two implementations of Algorithm 3, one with the usage of dominance checks, Algorithm 5, and one without it. We took articles with known supply from our set I, Appendix E, of real-world instances. We performed the tests for all three different scenarios, see also Chapter 3: good seller, bad seller and normal seller. To compare the number of visited inner nodes we use the results of the corollaries 2,3and 4. We compare runtime in seconds “time” and the percentage of visited inner nodes of the enumeration tree.1With kmax = 13 and pmax = 4 the overall number of nodes in the enumeration tree according to Corollary 3amounts to 1730. Without the leaves there remain 1730 −364 = 1366 b= 100% visited nodes in the complete enumeration tree for applying Algorithm 3without dominance checks “nd”. The percentages of this number for applying the procedure with dominance checks via Algorithm 5“d” are given in the columns with label “%n”. We see that with the application of our dominance rule depending on the considered scenario we can reduce the mean number of visited inner nodes to a number between 12.83 and 27.55 percents of the number of inner nodes we would have to visit if we enumerated all trajectories. The computation times can be reduced by a factor between 3.38 and 8.29. For “higher” scenarios we get better results with the usage of dominance checks. We used the same instances from the set Ito compare Algorithm 3with solving the mixed-integer programming formulation of Problem 5directly via CPLEX. To avoid nonlinearity we set the mark-down costs – both, fixed and variable – to zero. The results for the “low” scenario can be found in Table 4.2. While the solving process for the MIP averagely needs more than 8 hours our dynamic programming approach with dominance checks can yield the optimal solution averagely in less than one second. Compared to the results stated in Table 4.1 we also see that the number of cut off inner nodes is averagely more than 14 percentage points higher as against the case where mark-down costs are regarded. 1We exclude the leaves in the enumeration tree from this consideration because extensions to complete mark-down strategies are easy to compute by including the overall stock in the labels. CHAPTER 4. PRICE OPTIMIZATION 59 Instance nd d MIP time(s) time(s) %nodes time 1 3.56 0.48 14.79 30067.8 2 4.09 0.44 13.69 27423.3 3 4.08 0.44 15.45 28637.7 4 3.59 0.49 15.45 30303.6 5 3.49 0.53 14.79 28752.7 6 3.58 0.41 12.81 25800.5 7 4.56 0.18 4.90 36300.9 8 2.95 0.43 14.13 32898.0 9 2.99 0.36 13.47 30194.3 12 3.39 0.36 14.35 32356.2 ∅3.63 0.41 13.38 30273.5 Table 4.2: POP – dynamic generation versus MIP 4.10 Conclusion of the chapter Every week our industrial partner for more than 4000 products has to decide on markdowns. Hence, a fast approach for solving the Price Optimization Problem is necessary. In praxis there remain only two days to decide for an eventual mark-down. Our results show that therefore applying state-of-the-art solvers on the mixed-integer programming formulation of price optimization are not an alternative. A standard approach for pricing problems is dynamic programming. We applied this idea on the Price Optimization Problem how it takes part in DISPO and extended the approach by dominance rules to a label setting algorithm. Now we can solve the Price Optimization Problem for one article averagely in less than one second. Our label setting algorithm with dominance checks against an enumeration of all mark-down strategies would reduce the runtime for deciding mark-downs for weekly 4000 articles from four hours to one hour. In this form our algorithm POP-DYN with dominance checks could be applied as a standard approach for deciding on mark-downs at our industrial partner. Chapter 5 Stochastic Optimization In this chapter we outline some basics of stochastic programming. We are mainly guided by [BL97]. In contrast to deterministic programs stochastic programs can contain random data. Thus, stochastic programming extends the field of mathematical programming by programming under uncertainty. Uncertain input data are reproduced by random variables with known distribution. Stochastic programs are applied when not all of the input data is known at the time of decision making. The aim is, e.g., to optimize the expected costs over all possible scenarios. Economical problems often have to deal with uncertainties. Demand for products, resources, etc. are not always known a priori. By treating this uncertainties in a stochastic program a more realistic problem formulation is anticipated which shall lead to better decisions at the end. In Section 5.1 we will introduce so-called two-stage stochastic programs. Twostage stochastic programs consist of a so-called first stage decision, a decision which has to be made before the scenario in effect is known. The second stage – or recourse decision – responds to the realization of the scenario. If the random events follow a discrete distribution with a finite number of scenarios the two-stage stochastic program can be formulated as a so-called deterministic equivalent. In Section 5.2 we sketch out approaches from literature to solve the deterministic equivalent. Our focus in this chapter is on dual bounds for general two-stage stochastic programs, Section 5.3 – later on we will apply them to the Integrated Size and Price Optimization Problem in the context of a customized Branch&Bound approach. In Section 5.4 we sketch out the idea of multi-stage stochastic programs. Again a first-stage decision has to be made before anything about future behavior is known. But in contrast to two-stage programs – in which we only deal with one recourse decision – multi-stage programs involve sequences of recourse decisions over time depending on the realizations of the particular outcomes. 60 CHAPTER 5. STOCHASTIC OPTIMIZATION 61 5.1 Two-stage stochastic programs We start with an example of a popular stochastic program, the newsvendor problem, as it is for example stated in [BL97]. Example 6 (newsvendor problem).In the morning a newsvendor buys xnewspapers at a price cper paper from a publisher to sell them at the price of qon the street. The number xof bought newspapers is bounded above by u. The newsvendor sells as many papers as possible for the sales price q. At the end of the day he can return the remaining newspapers to the publisher at a price of rwith r < c . The demand per day is varying and described by a random variable ξ. The newsvendor problem can be formulated as a two-stage stochastic linear program. The first stage is the decision on how many newspapers the newsvendor should buy from the publisher. As second stage or recourse the newsvendor can compensate a wrong first-stage decision by returning overbought newspapers to the publisher. We get to the general formulation of a two-stage stochastic program. Definition 8 (two-stage stochastic program with recourse). min xcTx+EξQ(x, ξ)(5.1) subject to Ax =b, (5.2) x≥0.(5.3) It is Q(x, ξ) = min{qT ξy|Wξy=hξ−Tξx, x, y ≥0}.(5.4) The function Q(x, ξ)is also called recourse function. The recourse is called fixed if the so-called recourse matrix Wξdoes not depend on any uncertainties, then it is Wξ=W. The recourse is called complete if there is a valid second-stage decision for every first-stage decision and relative complete if for every valid first-stage decision for every scenario a valid second-stage decision exists. The matrix T is also denoted as technology matrix. In the case that the random events follow a discrete distribution with a finite number of scenarios it is possible to reduce the two-stage stochastic program to a deterministic program. Then the expected value EξQ(x, ξ)can be computed explicitely: For every variable yof the second stage one introduces a random variable for every particular scenario and obtains an equivalent linear program – the so-called deterministic equivalent. Definition 9 (two-stage stochastic problem in its extensive form – deterministic equivalent). min xcTx+X ξ∈Ξ pξqT ξyξ(5.5) subject to Ax =b, (5.6) Tξx+Wξyξ=hξ,∀ξ∈Ξ,(5.7) x≥0, y ≥0.(5.8) The set Ξcontains all scenarios which follow from the discrete distribution. The probability for the occurrence of scenario ξis given by pξ. CHAPTER 5. STOCHASTIC OPTIMIZATION 62 The number of the second stage variables in terms of the classical formulation from Definition 8in this models multiplies with the number of scenarios. This may lead to a large deterministic program. It may be that this deterministic equivalent can not be handled with standard approaches from linear programming. 5.2 Solving stochastic programs In this section we sketch out common solving methods for two-stage stochastic programs. One of the most prominent method for linear programs (without integer variables) is the so-called L-shaped method. We outline the basic idea in Subsection 5.2.1. For general mixed-integer programs an important property the L-shaped method is based on is violated: The convexity of the recourse function. In literature applications of general approaches from mixed-integer programming like Branch&Bound are often proposed to handle two-stage mixed-integer linear programs. We state some references in Subsection 5.2.2. 5.2.1 The L-shaped method for two-stage linear stochastic programs The L-shaped method for two-stage linear stochastic programs – in general also Benders’ decomposition [Ben62] – exploits the special structure of the deterministic equivalent.1 The Bender’s decomposition method is closely linked to the Dantzig-Wolfe decomposition [DW60]: It equals the Dantzig-Wolfe decomposition on the dual linear program: While in the case of a Dantzig-Wolfe decomposition in a column-generation algorithm [LD11] sequentially “promising” variables are added, in the L-shaped method step-by-step cuts are added. For this purpose the L-shaped method exploits the property that the recourse function Q(x, ξ)is piecewise linear and convex in xfor a fixed ξ. Positive linear combinations of convex functions are convex and so the expected value Q(x) := Eξ(Q(x, ξ)) is as well. The main idea of the L-shaped method is to approximate the term Q(x)in the objective function by piecewise linear functions. Two types of constraints are sequentially added: On the on hand optimality cuts which are linear approximations of Q(x)and on the other hand feasibility cuts which restrict the set {x|Ax =b, x ≥0}. For details and the exact algorithm we refer the reader to [BL97], [KW94] or [Pr´ e95]. 5.2.2 Solving two-stage mixed-integer stochastic programs In the case of (mixed-)integer stochastic programs the convexity of the function Q(x) – which is exploited by the L-shaped method – is no longer guaranteed. In literature some extensions of the Benders’ decomposition which deal with integrality of variables can be found. There are approaches for special problem structures. Wollmer [Wol80] extended Benders’ decomposition for stochastic problems with binary first and continuous second stage variables. This algorithm was extended by Laporte and Louveaux [LL93] for binary first and binary or continuous recourse variables. 1The name L-shaped method stems from the block structure of the deterministic equivalent. Chapter 6 The Integrated Size and Price Optimization Problem (ISPO) In this chapter we present our model ISPO for integrated size and price optimization. The model is an advancement of the model SLDP, presented in Chapter 2. Because demand is a priori unknown and depends on the scenario in effect, ISPO is formulated as a stochastic program. More precisely, a two-stage stochastic program with fixed recourse where the so-called first-stage decision is the supply in terms in lots and the recourse the price optimization where oversupply is compensated by marking down prices. The target is to find a supply policy that maximizes the expected profit over all scenarios – “bad seller”, “normal seller” and “good seller”. We will specify the problem in Section 6.1 and outline ISPO in its extensive form in Section 6.2. In Section 6.3 we show that ISPO is NP-hard by reducing to it the SLDP. In Chapter 4we presented a dynamic programming approach for solving the Price Optimization Problem where the supply is fixed. One could also think about applying dynamic programming for solving ISPO for all possible supply strategies. In Section 6.4 we will see that because of the special situation at our industrial partner – the supply in terms of lots – the state space is too large to apply dynamic programming. Moreover state-of-the-art MIP solvers cannot deal with the problem size of ISPO. 6.1 Problem specification An instance of ISPO among others consists of a set Bof branches, a set Lof lot-types and a set M={1, . . . , mmax}of multiplicities for the given lot-types. Also a set of sizes Sis specified. With lswe denote the number of items of Size sin Lot-type l. The set of lot-types is given by four parameters. The minimum number of items per size vmin with vmin ≥1, the maximum number of items per size vmax, the minimum number of items per lot-type vlmin and the maximum number of items per lot-type vlmax. An upper bound for the maximum number of supplied different lot-types κis given. Lot opening costs δifor using an additional i-th lot-type i= 1, . . . , κ are also specified. Moreover pick costs pcost arise for each supplied lot-type. An acquisition price ap has to be paid for each supplied item. The lower and upper bounds for overall supply are given by Iand I. 69 CHAPTER 6. THE INTEGRATED SIZE AND PRICE OPTIMIZATION PROBLEM (ISPO)70 In terms of the sales success of the considered article a set of scenarios Ewith scenario probabilites Prob(e),∀e∈Eis given. A set of sales periods K={0, . . . , kmax} is specified. We call 0, . . . , kmax −1the real sales periods where kmax is the sellout period. For a duration of kobs observation periods the article will not be marked down. We are given a set of price indices P={0, . . . , pmax}where the related price steps πp, p ∈Pare in descending order. The first price step π0is the starting price, the last price step πpmax is the salvage value. Moreover a factor ρfor weekly discounting is given. For every mark-down uniquely fixed mark-down cost µfand mark-down cost µvper marked down item arise. In the sellout period kmax, where all remaining items are sold, we assume that every left over item has to be marked down qkmax times. Here only variable mark-down costs µvper item accrue. For every scenario e, period k, branch b, size sand price index pthe dependent demand is given by de k,b,s,p. Our goal is to maximize the expected revenue for supplying each branch with a lot-type in a multiplicity by taking into account that mark-downs may occur during the selling time. 6.2 ISPO as a two-stage stochastic mixed-integer program (SMIP) in its extensive form Now we present one of the main aspects of this work – our formulation of the Integrated Size and Price Optimization Problem ISPO. We start with a look at the objective function coefficients. For the objective function we introduce handling costs cb,`,m. They are computed as described in the following definition. Definition 13 (handling costs cb,`,m).For a given acquisition price ap for one item and given pick cost pcost for one lot-type the handling costs for supplying Branch bwith Lot-type `in Multiplicity mare given by cb,`,m := m· X s∈S lsap + pcost!.(6.1) For the first stage – the size optimization stage (SOP) – we use binary assignment variables xb,`,m to encode the independent assignment of Branch bto Lot-type `in Multiplicity m. In the second stage – the price optimization stage (POP) – we introduce binary assignment variables for the independent second stage assignment decision ue k,pof Price index pto Period kin Scenario e. In order to account for the profit and the cost, we need additional dependent variables. We list the complete model before we comment on the details. Problem 6 (ISPO). max −X b∈BX `∈LX m∈M xb,`,m ·cb,`,m − κ X i=1 δi·zi(6.2) +X e∈E Prob(e)X k∈K exp(−ρk)X b∈BX s∈S re k,b,s −µkβe k(6.3) CHAPTER 6. THE INTEGRATED SIZE AND PRICE OPTIMIZATION PROBLEM (ISPO)71 Size Optimization Stage (SOP): X `∈LX m∈M xb,`,m = 1 ∀b∈B, (6.4) X m∈M xb,`,m ≤y`∀b∈B, ` ∈L, (6.5) X `∈L y`≤ κ X i=1 zi,(6.6) zi≤zi−1i= 1, . . . , κ, (6.7) Ib,s =X `∈LX m∈M m·`s·xb,`,m ∀b∈B, s ∈S, (6.8) I=X b∈BX s∈S Ib,s,(6.9) I∈[I,I],(6.10) xb,`,m ∈ {0,1} ∀b∈B, ` ∈L, m ∈M, (6.11) y`∈ {0,1} ∀`∈L, (6.12) zi∈ {0,1} ∀i= 1, . . . , κ, (6.13) Coupling via initial inventory: Ib,s −ve 0,b,s = 0,∀b∈B, s ∈S, e ∈E, (6.14) Price Optimization Stage (POP): X p∈P ue k,p = 1 ∀k∈K, e ∈E, (6.15) ue k,0= 1 ∀k∈K:k < kobs, e ∈E, (6.16) ue kmax,pmax = 1 ∀e∈E, (6.17) ue k−1,p1+ue k,p2≤1∀k∈K, e ∈E, p1, p2∈P:p2< p1,(6.18) βe k≥ue k−1,p1+ue k,p2−1∀k∈K, e ∈E, p1, p2∈P:p26=p1, (6.19) ve k−1,b,s −ve k,b,s =X p∈P we k−1,b,s,p ∀k∈K, b ∈B, s ∈S, e ∈E, (6.20) X p∈P we k,b,s,p ≤ve k,b,s ∀k∈K, b ∈B, s ∈S, e ∈E, (6.21) we k,b,s,p ≤ue k,p ·de k,p,b,s ∀k∈K, b ∈B, s ∈S, p ∈P, e ∈E, (6.22) re k,b,s =X p∈P πp·we k,b,s,p ∀k∈K, b ∈B, s ∈S, e ∈E, (6.23) ue k,p ∈ {0,1} ∀k∈K, p ∈P, e ∈E, (6.24) βe k∈ {0,1} ∀k∈K, e ∈E, (6.25) we k,b,s,p ≥0∀k∈K, b ∈B, s ∈S, p ∈P, e ∈E, (6.26) ve k,b,s ≥0∀k∈K, b ∈B, s ∈S, e ∈E, (6.27) re k,b,s ≥0∀k∈K, b ∈B, s ∈S, e ∈E, (6.28) CHAPTER 6. THE INTEGRATED SIZE AND PRICE OPTIMIZATION PROBLEM (ISPO)72 µk=µf+µvX b∈BX s∈S ve k,b,s,∀k∈K\ {kmax}, e ∈E, (6.29) µkmax =qkmax µvX b∈BX s∈S ve kmax,b,s,∀e∈E. (6.30) We start our explanations with the SOP stage: The constraints are the same as in the formulation of the SLDP (Problem 2). For the sake of completeness we will go into them anyway. The binary variables xb,`,m indicate if Lot-type `is delivered to Branch bin Multiplicity m. If this is the case they take value one, and zero otherwise. We force an assignment of a lot-type and a multiplicity to each branch by Equation (6.4). In order to account for the opening costs of the supplied lot-types, we introduce binary variables y`indicating whether or not lot-type `is used at all and binary variables zithat take value one if and only if at least idifferent lot-types are used for supply. Equation (6.5) guarantees that y`= 1 whenever `is assigned to at least one branch b. Inequality (6.6) implies that no more than κlot-types are used. Inequality (6.7) enforces that zi= 1 implies that the number of used lot-types is at least i. We use another dependent variable Ib,s for the inventory in branch band size s, and Equation (6.8) links this variable to the assignment decisions. The total inventory is then given by yet another dependent variable I, computed by Equation (6.9) and enforced to lie in between given bounds by Inequality (6.10). All independent variables have to be binary, see (6.11) through (6.13). Next, let us have a look at the POP stage model that is linked via the start inventories Ib,s to the SOP stage by Equation (6.14). The binary variables uk,p indicate if Price index pis allocated to Period k. If this is the case, they take value one, and zero otherwise. Equation (6.15) enforces the assignment of exactly one price to each period for each particular scenario. For a given number of periods kobs from the beginning of the sales process the starting price is enforced, Equation (6.16). In the last period – the sellout – the salvage value is fixed by Equation (6.17). We forbid increasing prices by Equation (6.18). A mark-down for Scenario ein Period kis indicated by the dependent binary variable βe k, which is forced to one by Inequality (6.19) if the price has changed compared to the previous period. The mark-down costs for the real sales periods are given by Equation (6.29), the mark-down costs for the sellout period by Equation (6.30). The following restrictions model the dynamics of the sales process using dependent variables. The fractional variable ve k,b,s specifies the stock level in Period kin Branch bfor Size sin Scenario e. The fractional variable we k,b,s,p measures the sales in Period kin Branch band Size sfor the price with index pin Scenario e. We capture the yield re k,b,s in Period kfor Size sin Branch bin Scenario eby Equation (6.23). Equation (6.20) describes the change of stock levels from one period to another. Inequality (6.21) models that sales may not exceed stock. In Inequality (6.22) we require that, only if the price with related price index pis chosen, there can be sales at the price index pof at most the demand at the price index p. Because the objective favors larger sales, the sales variables at a price in an optimal solution will attain exactly the minimum of stock and demand at that price. In this POP stage, only the independent price assignment variables need to be binary (6.24). The dependent variables capturing the dynamics of mean stocks, sales, and yields are required to be nonnegative in (6.26) through (6.28). The objective function subtracts the costs for the handling of mlots of Type `in Branch band the lot-type opening costs for using the first, second, . . . , i-th new lottype (6.2) from the expected discounted yields minus the expected discounted costs for CHAPTER 6. THE INTEGRATED SIZE AND PRICE OPTIMIZATION PROBLEM (ISPO)73 mark-downs, (6.3). 6.3 Complexity of ISPO The complexity of ISPO can be derived by reducing the SLDP on it. Corollary 5 (Complexity of ISPO).ISPO is NP-hard. Proof. If we set all prices πpfor all p= 0, . . . , pmax and the fixed and variable markdown costs in ISPO to value zero then the term X e∈E Prob(e)X k∈K exp(−ρk) X b∈BX s∈S re k,b,s −µkβe k! in the objective will take value zero. (This would be the same as completely ignoring the price optimization stage in ISPO). The constraints of the size optimization stage still have to be fulfilled. Because these constraints are the same as for the SLDP we can solve each instance of the SLDP by setting cb,`,m = distSLDP b,`,m +m·pcost. In Corollary 1we mentioned that the SLDP is NP-hard and therefore it follows that ISPO is too. 6.4 Solving ISPO with standard approaches In Chapter 4we presented a dynamic programming approach for the Price Optimization Problem. With a given supply per branch and size as initial state and applying our dominance rules we can solve the problem for all tested instances in less than four seconds to optimality. One could also think about applying dynamic programming to ISPO. In ISPO the initial supply for each branch and size is unknown from the beginning. The optimal supply is the supply that yields – together with the corresponding locally optimal mark-down strategies – the highest expected revenue. For each particular supply strategy theoretically a dynamic programming approach analogous to Chapter 4could be applied. Finally we would choose the supply which maximizes the expected revenue in terms of the locally optimal price trajectories. But the supply has to be in lot-types and the single branches are connected by the overall supply Iand the maximum number κof used lot-types, Constraints (6.6) and (6.10). Moreover mark-downs are applied simultaneously over all branches and sizes. Therefore the state space would be huge. To get an idea how large the state space is, let us consider a small example. As in practice we are given |B|= 1300 branches, but consider only lot-types (1,2,1),(2,1,1) and (1,1,2), set M={1}and restrict the overall supply by I= 5200 and I= 5200. We set κ= 3.1Because every of our 1300 branches can be supplied by each lot-type this yields 31300 >10620 different supply policies. In fact, for real instances with about 729 different lot-types, at least 3multiplicites and larger ranges in terms of the overall supply our state space is much bigger. So there is no possibility to solve ISPO just by dynamic programming to optimality. However dynamic programming is part of our heuristic solver for ISPO which we will present in Chapter 4. Here we restrict ourselves on small subsets of “most promising” lot-types. 1In fact these are no real restrictions. Consider the cardinality of the lot-types, the fact that the sole multiplicity takes value one and the number of branches. CHAPTER 6. THE INTEGRATED SIZE AND PRICE OPTIMIZATION PROBLEM (ISPO)74 Without the fixings (6.16), (6.17) and obviously redundant constraints (6.8), (6.9), (6.23) the number of constraints of ISPO amounts to |E|(|B||S|+|K|(1 + 1.5|P|(|P| − 1) + |B||S|(1 + |P|))) + |B|(1 + |L|) + 2 + κ. Without the obviously redundant variables (I,Ib,s,re k,b,s) there are |L|(|B||M|+ 1) + κ+|E||K|(1 + |B||S|)(|P|+ 1) variables. If we considered practical relevant instances with 3scenarios, 1300 branches, 6 sizes, 14 periods, 5prices, 729 lot-types, 3multiplicities and maximal κ= 4 allowed different lot-types this would lead to 1 990 308 constraints and 4 809 685 variables. For typical computers (8 GB of RAM) the size of these instances exceeds memory capacity and CPLEX fails already at the initialization. An instance of this size was tested on a machine with 128 GB of RAM.2Even after four weeks CPLEX was not able to solve the root relaxation. For an instance with 12 periods, 1000 branches, 7lot-types, 5multiplicities, maximal κ= 5 allowed used lottypes and the rest as above CPLEX needed more than 4 weeks to get to a solution with an optimality gap of 11.33%. So we tried a very small instance with only 30 branches, 5periods and 435 lot-types, 3multiplicities and maximal κ= 4 allowed lot-types and the rest as above.3The optimal solution is found after about five hours but still, after more than four weeks, there is no proof for optimality. To put this into perspective: Our Branch&Bound solver presented in Chapter 9is able to solve such instances in less than 30 seconds. And, more importantly, enables us to tackle real-world instances. 2Quad-Core AMD Opteron(tm) Processor 2384 CPU with 2GHz and 128 GB of RAM 3It is the first instance of our test set Itest 6, see Appendix E. Chapter 7 Reducing ISPO to the SLDP As we have seen in Chapter 6, state-of-the-art solvers fail for real instances of ISPO. In Chapter 4we deduced the number of possible price trajectories for the POPˆe. With this relatively small number – beside dynamic programming – it is possible to enumerate all valid price trajectories a priori. We exploit this property and in the remainder of this thesis only consider price trajectories instead of single assignments “price index to period”. We will show that fixing a price trajectory to each scenario simplifies ISPO to an SLDP with changed objective coefficients. To obtain these, in Section 7.1 we compute for each number of supplied items per branch and size the expected revenue – we call it single supply revenue. Adding up corresponding single supply revenues leads to socalled lot-type revenues which are objective function coefficients in the resulting SLDP, Section 7.2. We state the SLDP that results by fixing price trajectories to scenarios in ISPO in Section 7.3. In Section 7.4 we will show how ISPO theoretically could be solved by enumerating SLDPs for all possible assignments “scenario to price trajectory”. 7.1 Single supply revenues Now we consider the case that in ISPO a price trajectory to a particular scenario is fixed. For the fixed price trajectory we are able to compute the expected revenue for the related scenario for every branch, size and integer supply – the single supply revenue in advance. Before we describe the algorithm for the computation of the single supply revenues we introduce some notation. Notation 1 (map scenario to price trajectory, set of maps).For a scenario eand a price trajectory twe call e→tamap from scenario eto price trajectory t. A map e→t means that we fix Price trajectory t= (t0, . . . , tkmax )to Scenario ein ISPO. In the formulation of Problem 6this would mean setting the corresponding ue k,p with tk=p to 1, ∀k∈Kand the remaining ue k,p to 0. With WE0we denote a set of maps for a subset E0⊆Eof all scenarios, where each scenario e∈E0is mapped to exactly one valid price trajectory. 75 CHAPTER 7. REDUCING ISPO TO THE SLDP 76 Because the sales among the branches and sizes are independent from each other we can determine the expected revenue for each branch and size separately. Our computation of the revenue does not regard the fixed costs for mark-downs a priori. Because they only depend on the corresponding price trajectory and not on the sold items or the stock after a period we will treat them later on. We define the single supply revenue as follows. Definition 14 (single supply revenue).We consider a map scenario to price trajectory e→t. For Branch b∈Band Size s∈Sthe single supply revenue ¯ae→t b,s,n for a number nof supplied items is defined as ¯ae→t b,s,n := −n·ap + kmax−1 X k=0 exp (−ρk)πtkmin   max   n− k−1 X j=0 de j,b,s,tj,0  , de k,b,s,tk   −βkµvmax   n− k−1 X j=0 de j,b,s,tj,0   + exp (−ρkmax) (πkmax −qkmax µv) max (n− kmax−1 X k=0 de k,b,s,tk,0)! (7.1) For supplying 0 items the single supply revenue ¯ab,s,0takes value 0. For each additional supplied item we have to pay the acquisition price and the markdown costs. The more items are supplied the higher the costs for markdowns and for acquisition are. Because of the non-increasing prices and the discounting the yield we earn for each more supplied item will never exceed the yield we obtained for the last item. The single supply revenue ¯ae→t b,s,n is concave in n. Theorem 5 (concavity of the single supply revenue).¯ae→t b,s,n is concave for n≥0. Proof. We show that the function is concave for all real values ˜nwith ˜n > 0. Then we can transfer the concavity to nwith n > 0and ninteger. As a linear function −˜n·ap is concave. The term exp (−ρk) −βkµvmax   ˜n− k−1 X j=0 de j,b,s,tj,0    is concave: As a linear function ˜n−Pk−1 j=0 de j,b,s,tjis convex. Likewise the constant function 0is. Maxima of convex functions are also convex, such is max   ˜n− k−1 X j=0 de j,b,s,tj,0  . It is exp (−ρk),βkand µv≥0. Multiplying a convex function with positive scalars does not change the convexity. But negative convex functions are concave. This yields the claim. CHAPTER 7. REDUCING ISPO TO THE SLDP 77 For the same reasons exp (−ρkmax) −qkmax µvmax (˜n− kmax−1 X k=0 de k,b,s,tk,0)! is a concave function in ˜n. We now show that kmax−1 X k=0 exp (−ρk)πtkmin   max   ˜n− k−1 X j=0 de j,b,s,tj,0  , de k,b,s,tk   + exp (−ρkmax)πkmax max (˜n− kmax−1 X k=0 de k,b,s,tk,0) is concave. Because positive linear combinations of concave functions are also concave the claim of the theorem follows. To simplify the term we denote exp(−ρk)πtkby ˜πk,∀k∈K. Because we consider only one branch, size and scenario we can simplify de k,b,s,tkby ˜ dk. Moreover we include the last sellout period in the first term. This can be done by assuming that the demand ˜ dkmax in the sellout period for ˜nitems always amounts to ˜n. We have to show that kmax X k=0 ˜πkmin   max   ˜n− k−1 X j=0 ˜ dj,0  ,˜ dk   is concave for ˜n > 0. We show that ˜πf˜n0with f˜n0:= min nk∈K:Pk j=0 ˜ dj≥˜n0ois a supergradient for each ˜n0>0, i.e. we show that for every ˜n, ˜n0>0it is kmax X k=0 ˜πkmin   max   ˜n− k−1 X j=0 ˜ dj,0  ,˜ dk   − kmax X k=0 ˜πkmin   max   ˜n0− k−1 X j=0 ˜ dj,0  ,˜ dk  ≤˜πf˜n0(˜n−˜n0). Then concavity follows.1 It is f˜n0as defined above and analogously f˜n:= min{k∈K:Pk j=0 ˜ dj≥˜n}. Then it is kmax X k=0 ˜πkmin   max   ˜n− k−1 X j=0 ˜ dj,0  ,˜ dk   − kmax X k=0 ˜πkmin   max   ˜n0− k−1 X j=0 ˜ dj,0  ,˜ dk   = f˜n−1 X k=0 ˜πk˜ dk+ ˜πf˜n ˜n− f˜n−1 X k=0 dk − f˜n0−1 X k=0 ˜πk˜ dk−˜πf˜n0 ˜n0− f˜n0−1 X k=0 dk  1For detailed information about supergradients and generalized concavity we refer the reader to [ADSZ10]. CHAPTER 7. REDUCING ISPO TO THE SLDP 78 If ˜n= ˜n0the term above takes value 0and the supergradient inequality follows. So in the following we only consider the two cases ˜n > ˜n0and ˜n < ˜n0. Case 1: ˜n > ˜n0 It is f˜n≥f˜n0and because πk≥πk+1,∀k∈K\{kmax}it is also πf˜n≤πf˜n0. Such it is f˜n−1 X k=0 ˜πk˜ dk+ ˜πf˜n(˜n− f˜n−1 X k=0 ˜ dk)− f˜n0−1 X k=0 ˜πk˜ dk−˜πf˜n0( ˜n0− f˜n0−1 X k=0 ˜ dk) = f˜n−1 X k=f˜n0 ˜πk˜ dk+ ˜πf˜n˜n−˜πf˜n f˜n−1 X k=0 ˜ dk−˜πf˜n0˜n0+ ˜πf˜n0 f˜n0−1 X k=0 ˜ dk ≤˜πf˜n0 f˜n−1 X k=f˜n0 ˜ dk+ ˜πf˜n˜n−˜πf˜n f˜n−1 X k=0 ˜ dk−˜πf˜n0˜n0+ ˜πf˜n0 f˜n0−1 X k=0 ˜ dk =˜πf˜n0 f˜n−1 X k=0 ˜ dk+ ˜πf˜n˜n−˜πf˜n f˜n−1 X k=0 ˜ dk−˜πf˜n0˜n0 =˜πf˜n0−˜πf˜nf˜n−1 X k=0 ˜ dk+ ˜πf˜n˜n−˜πf˜n0˜n0 ≤˜πf˜n0−˜πf˜n˜n+ ˜πf˜n˜n−˜πf˜n0˜n0= ˜πf˜n0(˜n−˜n0). Case 2: ˜n < ˜n0 It is f˜n≤f˜n0and because πk≤πk+1,∀k∈K\{kmax}it is also πf˜n≥πf˜n0. Such it is f˜n−1 X k=0 ˜πk˜ dk+ ˜πf˜n(˜n− f˜n−1 X k=0 ˜ dk)− f˜n0−1 X k=0 ˜πk˜ dk−˜πf˜n0( ˜n0− f˜n0−1 X k=0 ˜ dk) = f˜n−1 X k=0 (πk−πf˜n)dk+πf˜n˜n+ f˜n0−1 X k=0 (πf˜n0−πk)dk−πf˜n0˜n0 = f˜n X k=0 (πf˜n0−πf˜n)dk+ f˜n0−1 X k=f˜n+1 (πf˜n0−πk)dk+πf˜n˜n−πf˜n0˜n0 ≤˜n(πf˜n0−πf˜n) + πf˜n˜n−πf˜n0˜n0=πf˜n0(˜n−˜n0) One possibility to compute the single supply revenues for all possible numbers of supplied items nfor a price trajectory and a given scenario would be to consider each nseparately and then apply the formula from Definition 14. But this may lead to unnecessary iterations: For all m < n and all periods kwith Pk j=1 de j,b,s,tj≤mthe yield in this period as well for mas for namounts to exp(−ρk)πtkde k,b,s,tk. With this consideration iterations can be avoided. Nevertheless additional mark-down costs may arise for a higher supply. To handle these we define the aggregated discount. CHAPTER 7. REDUCING ISPO TO THE SLDP 85 We are in the case that 0.5 = ˜ d < ˜r= 1. We update the current revenue and get ˜a= 41.19 + exp (−0.01 ·2) ·0.5·10.99 = 46.58. The current stock is updated to ˜r= 0.5. We continue with period ˜ k= 3 and set the current demand to the demand of Period 3: ˜ d= 1.5. Because there is a mark-down in the current period ˜ k= 3, we have to regard the variable mark-down costs and get ˜a= 46.58 −exp (−0.01 ·3) ·0.5·0.1 = 46.53. Now the current demand exceeds the current stock: It is 1.5 = ˜ d≥˜r= 0.5. The updated current revenue amounts to ˜a= 46.53 + exp (−0.01 ·3) ·0.5·1.99 = 47.49. We set the revenue for suppling 5 items to ¯ae→t 1,S,5= 47.49. The remaining demand now is ˜ d= 1, the remaining stock is set to ˜r= 1. We go on with suppling n= 6 items and subtract the acquisition price from the current revenue: ˜a= 47.49 −0,5 = 46.99. Additionally we have to regard the mark-down costs for the current stock. It is ˜a= 46.99 −exp (−0.01 ·3) ·0.1 = 46.89 and 1 = ˜ d= ˜r. The updated revenue is ˜a= 46.89 + exp (−0.01 ·3) ·1·1.99 = 48.82. The cost for suppling 6 items is set to ¯ae→t 1,S,6= 48.82. The remaining demand is ˜ d= 0. We continue with n= 7 and set the remaining stock to ˜r= 1. The revenue is updated to ˜a= 48.82 −0.5 = 48.32. The mark-down costs for one additional item are subtracted from the current revenue. This results in ˜a= 48.32 −exp (−0.01 ·3) ·0.1 = 48.22. It is 0 = ˜ d < ˜r= 1. Because we are in the last real sales period, we are in the case of Step 20. We have to add the salvage value and subtract the variable mark-down costs for period kmax and get ¯ae→t 1,S,7= 48.22 + exp (−0.01 ·4) ·1·(0.99 −2·0.1) = 48.98 and end up the computation for Branch 1 and Size S. The additional revenue per supplied item for more than 7supplied items amounts to −0.5 + exp (−0.01 ·3) ·0.1 + exp (−0.01 ·4) ·(0.99 −2·0.1) = 0.16. 7.1.3 Computational results We compare Algorithm 6with the more intuitive Algorithm 7to compute the single supply revenues. In Algorithm 7for each integer supply of items nper branch and size with vmin ≤n≤ dPkmax−1 k=0 de k,tk,b,sewe always walk through the periods k= 0, . . . , kmax −1until all items are sold or we reached the last real sales period. Thus, the alternative algorithm has a runtime of O |B||S|kmax &kmax−1 X k=0 de k,b,s,tk'!.(7.11) That means with Algorithm 6we can reduce the pseudo-polynomial runtime against the more intuitive Algorithm 7at least by one factor. In Table 7.1 the results for comparing Algorithms 6with Algorithm 7for a subset from the set of real instances I(see Appendix E) are stated. We took the supply Ib,s from historical transaction data we obtained from our partner. The algorithms were performed consecutively for all three scenarios “low seller”, “normal seller” and “high seller”. Thus, time and number of iterations refer to all three scenarios. The mean values in the last line refer to the whole set I. The number of iterations “#iter” and the computation time in seconds “time(s)” is stated for Algorithm 6and the described alternative. In Algorithm 6on average the number of iterations is 170 586 less. Let us take a look on the overall computation time. It amounts to 49.00 seconds for our proposed CHAPTER 7. REDUCING ISPO TO THE SLDP 86 Algorithm 7 Single supply revenue 2 Require: map e→tfrom Scenario efrom Price trajectory t= (t0,...,tkmax ), complete data of an instance for ISPO Ensure: for every branch b, every size sand integer numbers n:vmin ≤n≤ dPkmax−1 k=0 de k,b,s,tkesingle supply revenue ¯ae→t b,s,n 1: for all b∈Bdo 2: for all s∈Sdo 3: set n=vmin 4: while true do 5: set ˜r=n 6: set ˜a=−ap ·n 7: set ˜p=t0 8: set ˜ k= 0 9: set ˜ d=de 0,b,s,t0 10: while ˜ k < kmax do 11: if ˜ k > 0and t˜ k−16=t˜ kthen 12: ˜a= ˜a−exp −ρ˜ k˜rµv 13: if ˜ d < ˜rthen 14: ˜a= ˜a+exp −ρ˜ k˜ dπ˜p 15: ˜r= ˜r−˜ d 16: ˜ k=˜ k+ 1 17: ˜ d=de ˜ k,b,s,t˜ k 18: ˜p=t˜ k 19: else 20: ˜a= ˜a+ exp −ρ˜ k˜rπ˜p 21: break 22: end if 23: end if 24: end while 25: ˜a= ˜a+exp (−ρkmax) ˜r(πpmax −qkmax µv) 26: n=n+ 1 27: end while 28: end for 29: end for Instance #iter time(s) Alg.7Alg.6Alg.7Alg.6 1328 308 212 244 0.06 0.05 2328 308 212 244 0.06 0.05 3389 730 224 119 0.06 0.06 4365 680 218 762 0.06 0.05 5365 680 218 762 0.07 0.04 6365 680 218 762 0.05 0.06 7365 680 218 762 0.06 0.05 8355 236 217 681 0.05 0.06 9366 060 218 278 0.05 0.06 10 333 383 212 806 0.05 0.05 11 274 576 203 267 0.05 0.05 12 274 576 203 267 0.06 0.05 13 314 439 209 643 0.05 0.05 14 361 619 235 900 0.08 0.06 15 447 372 284 027 0.09 0.08 . . . . . . . . . . . . . . . ∅496 757 326 171 0.0569 0.0606 P280 507 337 427 210 602 49.00 52.15 Table 7.1: Comparison – computation of single supply revenues CHAPTER 7. REDUCING ISPO TO THE SLDP 87 Algorithm and 52.15 second for the alternative. In terms of time this is only a marginal improvement. Yet, we proposed an algorithm for the computation of the single supply revenues which always needs less iterations than a more intuitive algorithm. 7.2 Establishing lot-type revenues To obtain revenues for a map e→tScenario eto Price trajectory tfor a given lot-type `and multiplicity mone has to add up the according to Section 7.1 computed single supply revenues ¯ae→t b,s,n in the following way. Definition 16 (lot-type revenue).The lot-type revenue ˆab,`,m for Branch b, Lot-type ` and Multiplicity mfor a map e→tis given as ˆae→t b,`,m := X s∈S ¯ae→t b,s,m·ls−m·pcost.(7.12) By definition the single supply revenue contains as well acquisition price as yield resulting from fixing Price trajectory tto Scenario e, so also the lot-type revenue does. We included the corresponding pick-costs by subtracting m·pcost. 7.3 Fixing price trajectories in the ISPO – an SLDP Let us consider a subset E0of the scenarios E. We are given a set of maps WE0where each e∈E0is mapped to a price trajectory te. The part of the expected revenue of the ISPO according to the subset E0of scenarios is given by the optimal objective value of the following SLDP. Problem 7 (SLDP(WE0)). max X e∈E0 Prob (e) X b∈BX `∈LX m∈M ˆae→te b,`,mxb,`,m −µf˜ δte kmax−1− κ X i=1 δi·zi! (7.13) subject to X `∈LX m∈M xb,`,m = 1 ∀b∈B, (7.14) X m∈M xb,`,m ≤y`∀b∈B, ` ∈L, (7.15) X `∈L y`≤ κ X i=1 zi,(7.16) zi≤zi−1i= 1 . . . , κ, (7.17) Ib,s =X `∈LX m∈M m·`s·xb,`,m,∀b∈B, s ∈S, (7.18) I=X b∈BX s∈S Ib,s,(7.19) CHAPTER 7. REDUCING ISPO TO THE SLDP 88 I∈[I,I],(7.20) xb,`,m ∈ {0,1} ∀b∈B, ` ∈L, m ∈M, (7.21) y`∈ {0,1} ∀`∈L, (7.22) zi∈ {0,1}i= 1, . . . , κ. (7.23) The constraints do not differ from the formulation of Problem 2. Now we want to consider revenue not costs and changed the objective from minimization to maximization. The objective coefficient ˆae→te b,`,m is the expected revenue for Scenario efor supplying Branch bwith Lot-type `in Multiplicity mexcluding the fixed mark-down costs. We have to subtract these costs for the periods with mark-downs separately. This is done by using the aggregated discount from Definition 15 for every scenario. Additionaly we subtract the lot-opening costs in the objective. Notation 2 (optimal objective value of SLDPWE0).With z∗ SLDP(WE0)we denote the optimal objective function value of Problem 7. Notation 3 (optimal objective value of the LP relaxation of SLDPWE0).With z∗ SLDP-LP(WE0)we denote the optimal objective value of the LP relaxation of the SLDPWE0. 7.4 Solving the ISPO by enumerating SLDPs With the results of the previous sections, theoretically we can solve ISPO by solving an SLDP for every set of maps “scenario to price trajectory” in which each considered scenario is assigned to a price trajectory. The subsequent Corollary follows directly from Theorem 2. Corollary 8 (solving the ISPO by enumerating SLDPs).We consider ISPO how it is stated in the formulation of Problem 6. We can solve ISPO by enumerating kmax −kobs +pmax −1 pmax −1|E| SLDP(WE)s. We will illustrate the proceeding on a small example. Example 9 (ISPO – enumeration tree).For simplicity we assume in this example kobs = 2,kmax = 3 and pmax = 2. According to Theorem 2this leads to 2 different price trajectories: t0= (0,0,0,2) and t1= (0,0,1,2). It is |E|= 3. In Figure 7.1 we see the corresponding enumeration tree. A depth corresponds to a scenario e∈E. The width corresponds to a price trajectory. At the leaves every scenario is fixed to a price trajectory and solving the corresponding SLDP(WE) yields a lower bound for the optimal solution of the ISPO. The optimal solution of the ISPO is the optimal solution of the SLDP(WE) with the highest optimal solution value among all possible set of maps WE. For example, for the set of ordered scenarios E={e1, e2, e3}and the set of ordered price trajectories {t1, t2}at the leaf of the third branch in the enumeration tree in Figure 7.1 we would solve an SLDP(WE) with WE={e1→t1, e2→t2, e3→t1}. CHAPTER 7. REDUCING ISPO TO THE SLDP 89 scenario e1 e2 e3 t1t2 t1t2t1t2 t1t2t1t2t1t2t1t2 Figure 7.1: ISPO – enumeration tree Enumerating all maps WEand solving the related SLDP(WE) is for us firstly theoretically interesting. With real-world instances with kobs = 2,kmax = 13 and pmax = 4 this would yield 364 price trajectories for every scenario what means 3643= 48 228 544 leaves or to be solved SLDP(WE)s. Even if we assumed that the solving process for each SLDP(WE) would take only 1 second, we would have to wait more than 558 days for an optimal solution of the ISPO. But with this reduction we are able to reduce the size of ISPO to the size of – although – many SLDP(WE)s and theoretically can solve the same instances – in terms of the numbers of branches, lot-types and multiplicities – as we can handle with the SLDP. Because also the SLDP(WE) can be reduced to κLDPs as shown in Chapter 2, we can use all known approaches for the LDP in the solving process. Enumerating the SLDP(WE)s in the outlined way is the basic idea of our exact Branch&Bound solver we will present in Chapter 9. Because solving all SLDP(WE)s is not possible, we have to prune the enumeration tree. This is done by dual bounds which are topic of the following chapter. Chapter 8 Dual Bounds In Chapter 9we will present an exact Branch&Bound approach to solve the Integrated Size and Price Optimization Problem. One of the most important points in terms of efficiency of this algorithm is the computation of possibly tight dual bounds. At the nodes of our Branch&Bound we will fix price-trajectories to scenarios. At the leaves then for every scenario a price trajectory is fixed and the corresponding SLDP(WE) is solved to optimality in order to obtain a primal bound for ISPO. To reduce the number of processed nodes and especially of processed leaves, i.e. to be solved SLDP(WE)s, we will use lower bounds based on the wait-and-see solution from stochastic programming, see also Chapter 5. We developed tighter dual bounds based on the wait-and-see solution by combining subsets of scenarios and call the result the extended wait-and-see solution. We present these bounds in Section 8.1 Computational effort is reduced by combining wait-and-see and extended wait-andsee solutions with relaxations of the SLDP(WE0) for subsets E0⊆Eof scenarios. A relaxation we use neglects the fact that the supply has to be in lot-types – i.e. independent integer supplies per branch and size are allowed. Additionally we use LP relaxations. In Section 8.2 we apply these dual bounds to ISPO. We conclude the chapter in Section 8.3. 8.1 Dual bounds from wait-and-see solutions A well known dual bound from stochastic programming is the wait-and-see solution, see also Chapter 5. In Subsection 8.1.1 we will formulate the wait-and-see solution for a general stochastic mixed-integer program. We show how to extend this bound by combining scenarios, Subsection 8.1.2. Computation time can be reduced by relaxing the underlying problems, Subsection 8.1.3. 8.1.1 The wait-and-see solution In this section we will refer to the following formulation of a general stochastic mixedinteger program: 90 CHAPTER 8. DUAL BOUNDS 91 Problem 8 (general two-stage SMIP). max xz(x(Ξ),Ξ) = max xcTx+Eξ∈ΞQ(x, ξ)(8.1) subject to Ax =b, (8.2) x∈Rk×Zn−k.(8.3) Ξis the set of all possible scenarios and Q(x, ξ)is the objective of the second stage for the particular scenario ξand may be nonlinear. With x∗(Ξ) we denote the optimal solution of the problem, z(x∗(Ξ),Ξ) is the optimal objective value. If we restrict Problem 8to a particular scenario ξwe obtain Problem 9. Problem 9 (general two-stage SMIP associated with Scenario ξ). max xzSIN(x(ξ), ξ) = max xcTx+Q(x, ξ)(8.4) subject to Ax =b, (8.5) x∈Rk×Zn−k.(8.6) With x∗(ξ)we denote the optimal solution of Problem 9. The related optimal objective value is denoted by zSIN(x∗(ξ), ξ). Definition 17 (wait-and-see solution for a general SMIP).The wait-and-see solution WS for Problem 8is given by WS =Eξ∈Ξhmax xzSIN(x(ξ), ξ)i=Eξ∈ΞzSIN(x∗(ξ), ξ).(8.7) The wait-and-see solution is a dual bound of the original problem. It follows from the proof of Theorem 4in Chapter 5. Integrality of variables does not affect the result. 8.1.2 Extending wait-and-see solutions Based on the wait-and-see solution from the previous subsection we developed stronger dual bounds for general two-stage stochastic optimization problems by considering subsets of scenarios. We call them extended wait-and-see solutions. For this purpose we define the so-called partial stochastic program on which the extended wait-and-see solution is based. Problem 10 (partial two-stage SMIP).It is pξthe probability of the occurrence of scenario ξ. The partial two-stage stochastic (mixed-integer) program for the subset of scenarios Ξ0⊆Ξis given by max xzPA(x(Ξ0),Ξ0) = max xX ξ∈Ξ0 pξ(cTx+Q(x, ξ)) (8.8) subject to Ax =b, (8.9) x≥0.(8.10) It is x∗(Ξ0)the optimal solution of the problem and with zPA(x∗(Ξ0),Ξ0)we denote the corresponding objective value. CHAPTER 8. DUAL BOUNDS 92 Definition 18 (extended wait-and-see solution).We define the extended wait-and-see solution ˜ WS(Ξ0)of an SMIP for the subset Ξ0⊆Ξof scenarios by ˜ WS(Ξ0) := zPA(x∗(Ξ0),Ξ0) + X ξ∈Ξ\Ξ0 zPA(x∗({ξ}),{ξ})(8.11) =zPA(x∗(Ξ0),Ξ0) + X ξ∈Ξ\Ξ0 pξzSIN(x∗(ξ), ξ).(8.12) For every subset of scenarios from the set Ξthe extended wait-and-see solution is a tighter dual bound than the classical wait-and-see solution. Theorem 7 (extended wait-and-see solution as dual bound).For the optimal solution of Problem 8and the extended wait-and-wee solution ˜ WS(Ξ0)the following inequalities hold: WS ≥˜ WS(Ξ0)≥z(x∗(Ξ),Ξ).(8.13) Proof. The first inequality equals X ξ∈Ξ pξzSIN(x∗(ξ), ξ)≥zPA(x∗(Ξ0),Ξ0) + X ξ∈Ξ\Ξ0 pξzSIN(x∗(ξ), ξ). Subtracting Pξ∈Ξ\Ξ0zSIN(x∗(ξ), ξ)on both sides yields X ξ∈Ξ0 pξzSIN(x∗(ξ), ξ)≥zPA(x∗(Ξ0),Ξ0). With zPA(x∗(Ξ0), ξ)we denote the objective value that results from applying the optimal solution of Problem 10 to Problem 9. By definition it is zSIN(x∗(ξ), ξ)≥zPA(x∗(Ξ0), ξ) for every scenario ξ∈Ξ0. Multiplying both sides with the corresponding scenario probability and adding up for all scenarios ξ∈Ξ0yields the claim. The proof of the second inequality is similar. We denote with zSIN(x∗(Ξ), ξ)and zPA(x∗(Ξ),Ξ0)the objective values that result from applying x∗(Ξ) to the problems 9 and 10. By definition it is zSIN(x∗(ξ), ξ)≥zSIN(x∗(Ξ), ξ) and also by definition it is zPA(x∗(Ξ0),Ξ0)≥zPA(x∗(Ξ),Ξ0). Multiplying the first inequality for each scenario ξ∈Ξ\Ξ0on both sides with the corresponding scenario probability does not change the relation. Adding up the second inequality and the with the scenario probabilities multiplied first inequalities for all scenarios ξ∈Ξ\Ξ0yields the claim. Partial SMIPs and group subproblems In Chapter 5we introduced the so-called group subproblem as it is defined by Sandikc¸i et.al. [SKS12]. The group subproblem yields a hierarchy of bounds for stochastic CHAPTER 8. DUAL BOUNDS 93 programs by combining all subsets of scenarios with same cardinality and computing the expected optimal value over all group subproblems for these subsets. We recognized similarities between the group subproblem and our extended waitand-wee solution. To outline our ideas we first consider a group subproblem where the reference scenario ξrtakes probability zero, i.e. the reference scenario is not contained in our subset Ξof scenarios. For our formulation of a general SMIP, Problem 8, the group subproblem for a subset Ξ0of scenarios can be stated as follows. Problem 11 (group subproblem for a general SMIP, reference scenario not in Ξ). z∗ g(Ξ0) = max xcTx+X ξ∈Ξ0 pξ ρ(Ξ0)Q(x, ξ)(8.14) subject to Ax =b, (8.15) x∈Rk×Zn−k.(8.16) where ρ(Ξ0) := Pξ∈Ξ0pξ. It is ρ(Ξ0)·z∗ g(Ξ0) = max xρ(Ξ0)cTx+ρ(Ξ0)X ξ∈Ξ0 pξ ρ(Ξ0)Q(x, ξ)(8.17) = max xρ(Ξ0)cTx+X ξ∈Ξ0 pξQ(x, ξ)(8.18) = max xX ξ∈Ξ0 pξ(cTx+Q(x, ξ)).(8.19) This is exactly the objective function of our partial two-stage SMIP, Problem 10. That means we can compute the expected value of the group subproblem for i scenarios, with Pi(Ξ) being the set of all subsets of iscenarios from the set of all scenarios Ξ, in this case originally given by EGSO(i) = 1 |Ξ| − 1 i−1X Ξ0∈Pi(Ξ) ρ(Ξ0)z∗ g(Ξ0)(8.20) by EGSO(i) = 1 |Ξ| − 1 i−1X Ξ0∈Pi(Ξ) zPA(x∗(Ξ0),Ξ0).(8.21) Our partial two-stage SMIP yields the same optimal solution as the group subproblem. The objective differs because in relation to the formulation of the group subproblem we scaled the vector cby the summed up scenario probabilities. Such, our partial SMIP can be seen as a special case of the group subproblem where the reference scenario is not contained in the set Ξ. The wait-and-see solution is with Ξ0=∅a special case of the extended wait-andsee solution. Therefore, in the following we will refer solely to extended wait-and-see solutions. CHAPTER 8. DUAL BOUNDS 94 8.1.3 Relaxations of extended wait-and-see solutions If the original problem restricted to one scenario is still hard to solve, we can use the bounding property of the wait-and-see solution together with the bounding property of relaxations and use the optimal values of the relaxed problems 9and 10 for the computation of the extended wait-and-see solutions. The inequalities from the theorems 4and 7will also hold for the relaxed extended wait-and-see solutions, the bounds at most increase. Definition 19 (relaxed extended wait-and-see solution).It is x∗ R(ξ)the optimal solution of a relaxation of Problem 9. With zSIN R(x∗ R(ξ), ξ)we denote the related optimal objective value. With x∗ R(Ξ0)we denote the optimal solution of a relaxation for Problem 10. The related optimal objective value is zPA R(x∗ R(Ξ0),Ξ0). We define the relaxed extended wait-and-see solution ˜ WSR(Ξ0)for the subset Ξ0of scenarios as follows: ˜ WSR(Ξ0) : = zPA R(x∗ R(Ξ0),Ξ0) + X ξ∈Ξ\Ξ0 zPA R(x∗ R({ξ}),{ξ})(8.22) =zPA R(x∗ R(Ξ0),Ξ0) + X ξ∈Ξ\Ξ0 pξzSIN R(x∗ R(ξ), ξ)(8.23) Corollary 9 (bounding by relaxed extended wait-and-see solutions).It is ˜ WSR(Ξ0)≥˜ WS(Ξ0)≥z(x∗(Ξ),Ξ) (8.24) for all subsets of scenarios Ξ0⊆Ξ. Proof. The theorem follows from Theorem 7by exploiting the dual-bound property of relaxations. 8.2 Application to ISPO Now we apply the outlined bounds to the ISPO. Because the SLDP in general is hard to solve in our solver we solely use relaxed (extended) wait-and-see solutions. On the one side LP relaxations and on the other side combinatorial bounds – based on relaxing the restriction that the items have to be supplied in terms of lot-types – are applied. In Subsection 8.2.1 we present our combinatorial bounds. Then, in Subsection 8.2.2 we apply the different dual bounds to ISPO and in Subsection 8.2.3 we present computational results. 8.2.1 Relaxing the lot-type constraint – single supply relaxations By neglecting the restriction that the supply has to be in terms of lots and allowing independent integer supplies of items among branches and sizes we obtain a relaxation of the SLDP(WE0). First we will present the problem formulation of this relaxation as an integer linear program before we outline a fast greedy algorithm to solve it. CHAPTER 8. DUAL BOUNDS 101 and-see solution ILPB. For ELPB this is even always the case. As already mentioned ELPB for E0=Eis the same as the LP relaxation SLDP-LP(WE) of the SLDP(WE). This means that the SLDP(WE) has a small integrality gap, i.e. the optimal objective value of SLDP-LP(WE) is not far from the optimal objective value of the SLDP(WE). 8.3 Conclusion of the chapter We introduced new bounds for two-stage stochastic programs based on the wait-andsee solution from stochastic programming. The extended wait-and-see solutions are tighter than the classical wait-and-see solution. We applied these bounds on the leaves of the enumeration tree of ISPO – here each leaf corresponds to an SLDP(WE). Because the underlying binary programs are hard to solve we relax them – either by allowing single supply instead of lot-types or by the LP relaxation. The results for a set of small test instances show that the bounds based on the single supply relaxation can be obtained very fast. The relaxed extended wait-and-see solutions are always faster to compute than the classical wait-and-see solution. In some cases they even beat the classical wait-and-see solution in terms of the optimality gap; for the LP-relaxed extended wait-and-see solution this is even always the case. Chapter 9 Solving the Integrated Size and Price Optimization Problem Now we present two of the main results of this thesis: our solvers for the Integrated Size and Price Optimization problem ISPO. Since the MIP formulation of ISPO presented in Chapter 6for real instances cannot be solved directly by state-of-the-art MIP solvers, we developed two approaches: An exact algorithm for benchmarking and a fast heuristic for practical use. In the exact algorithm we exploit the fact that for fixed price trajectories ISPO reduces to an SLDP, Chapter 7. Dual bounds on the base of the wait-and-see solution, see Chapter 8, allow us to prune the enumeration tree from Section 7.4. We outline the resulting Branch&Bound algorithm named ISPO-BAB in Section 9.1. The heuristic solver, ISPO-PingPong, is presented in Section 9.2. It exploits the fact that for every valid second stage decision, i.e. for every set of maps “scenario¸to price trajectory” for all scenarios there exists a valid first stage decision, a valid supply policy. We call this property reversible recourse. We present computational results for both solvers on real-world instances in Section 9.3 and give some remarks about the goodness of our proposed heuristic in Section 9.4, before we conclude this chapter in Section 9.5. 9.1 An exact Branch&Bound approach Now we present our exact solver for ISPO – a customized Branch&Bound algorithm. In Chapter 7, Section 7.4, we already showed how to solve ISPO by enumerating all possible combinations of assignments of price trajectories to scenarios. We adopt this idea and develop it further by applying the dual bounds presented in Chapter 8. The result is our exact Branch&Bound algorithm ISPO-BAB. A node at Depth jcorresponds to all maps “scenario to price trajectory” with the images of the first jscenarios fixed. The leaves are the maps with fixed images for all scenarios. In the branching step we extend a partially defined set of maps at a node by maps to all valid price trajectories for the next scenario. 102 CHAPTER 9. SOLVING ISPO 103 9.1.1 The algorithm We present the detailed implementation of the above concept – Algorithm 9. In a first step we compute for each map “scenario to price trajectory” the single supply revenues for each branch and size and every number nfrom the set Nof integer supplies, see Theorems 3and 4, by Step 4. We apply a simple dominance check, Step 6, in Algorithm 10: If for each number, each size and each branch for one price trajectory tthe single supply revenue is always smaller than the single supply revenue for another price trajectory t0in terms of the same scenario, the price trajectory tis dominated by t0. If this is the case we can exclude tfrom further consideration for this scenario. The set of non-dominated price trajectories for scenario eis denoted by Te. Then, in Step 10 of Algorithm 9, we solve for each scenario e∈Eand each price trajectory t∈Tethe SLDP-CB({e→t}), Problem 12, and store the optimal objective function value – it is a summand for the computation of the relaxed waitand-see solution – as Γ (e→t). We update these addends possibly later on by solving the LP-relaxation SLDP-LP({e→t}). In further progress we add up these values to obtain relaxed wait-and-see solutions for the particular SLDP(WE)s. But from these combinatorial bounds we also benefit in another way: These bounds are used to get an idea how good a price trajectory fits to a scenario. By ordering the trajectories decreasingly according to their combinatorial bounds in the Branch&Bound tree we expect to find the optimal solution for ISPO as early as possible which may lead to the possibility to prune a bigger part of the tree. For the purpose of computing dual bounds also at the inner nodes of the tree we save the maximum lower bound for each scenario as Γmax(e)in Step 13 and in the next step the index of the corresponding price trajectory as argmax(e). Then we start a depth-first-search. We start with an empty set WE0of maps “scenario→price trajectory”, Step 18, and set the current depth to value one, Step 19. The depth-first-search itself is outlined as Algorithm 11: We fix the first price trajectory in the order to the scenario related to the current depth, Step 1, and compute the dual bound in terms of a relaxed wait-and-see solution, Step 2. If we are in the first branch of the tree, we abandon to compute improved lower bounds because we will not be able to prune the tree until the primal bound is updated for the first time. Otherwise – if the dual bound is to weak to prune the tree – we try to improve the current bound by updating the bounds Γ(e→t)by LP-relaxations or by computing extended waitand-see solutions in Step 6, more precisely by the algorithms 12,13 and 14. If the dual bound is smaller than the current primal bound, we are able to prune the current branch at this point. Whenever we reach a leaf of the tree at depth |E|and are not able to prune, we solve the corresponding SLDP(WE)and possibly update the primal bound. In Algorithm 12 the summands for the computation of the LP-relaxed wait-and-see solution are computed. We avoid to solve more LP relaxations of the SLDP(WE0) than necessary. That means we mix combinatorial and LP-relaxed bounds. We compute LP relaxations of the SLDP(WE0) just until the resulting “mixed” relaxed wait-and-see solution is smaller than the current primal bound. We also update the corresponding values Γ (e→t)if the SLDP-LP({e→t}) yields a smaller optimal objective function value than the SLDP-CB({e→t}), i.e. if the LP-relaxation of the SLDP({e→t}) has a smaller optimality gap than the related single supply relaxation, Step 6. If this is the case, we may have to adapt the maximum bound for the current scenario Γ (e)max in Step 7. CHAPTER 9. SOLVING ISPO 104 Algorithm 9 ISPO-BAB Require: complete data of an instance of ISPO, set of price trajectories T Ensure: expected revenue maximizing supply x∗in terms of lots and related price trajectories 1: init global primal bound z∗=−∞ 2: for all e∈Edo 3: for all t∈Tdo 4: compute ¯ae→t b,s,n∀b∈B, s ∈S, n ∈N, Algorithm 6 5: end for 6: Te=ISPO-DOM(e) 7: end for 8: for all e∈Edo 9: for all t∈Tedo 10: compute the optimal solution of SLDP-CB({e→t}) and save the related optimal objective function value as Γ (e→t) 11: end for 12: sort the price trajectories Tein descending order according to Γ (e→t)→set of indexed price trajectories Te= (te 1, te 2,...,te |T|) 13: set Γmax(e) = Γ (e→te 1) 14: set argmax(e) = te 1 15: end for 16: sort the scenarios →set of indexed scenarios E={e1,...,e|E|} 17: for all i= 1,...,|Te1|do 18: set WE0=∅ 19: set j= 1 20: ISPO-DFS(1,i,WE0) 21: end for 22: return z∗and related optimal solution of ISPO 9.1.2 Some implementational aspects At this point we outline some implementational details of our Algorithm 9. MIP solvers As MIP solver or LP solver for the SLDP(WE0)s and the SLDP-LP(WE0)s one can choose between ILOG CPLEX and SCIP. By default we use CPLEX because it leads to shorter computation times. An advantage of CPLEX is that if once an SLDP-LP(WE0) is generated we can keep it in memory and only have to update the objective coefficients the next time we have to solve an SLDP-LP(WE0). Moreover we use the possibility to perform so-called warm starts in the simplex algorithm. Warm starts The constraints of the SLDP-LP(WE0) and SLDP(WE0) are independent from the fixed price trajectories, only the coefficients of the objective differ. So every optimal basis of an SLDP-LP(WE0) is also feasible for the SLDP-LP(W˜ E)s with W˜ E6=WE0 at the other nodes of the tree. When solving an SLDP-LP(WE0) we let CPLEX start the primal simplex with the optimal basis of an already solved SLDP-LP(WE0) – if available. For detailed information about the simplex algorithm we refer the reader to [Van07]. Additional Bounding As soon as a primal bound for ISPO is found, we hand over this primal bound to the MIP solver. The MIP-solver uses this bound in its internal Branch&Bound approach to prune branches which can not lead to a better solution. CHAPTER 9. SOLVING ISPO 105 Algorithm 10 ISPO-DOM Require: single supply revenues ¯ae→t b,s,n ∀b∈B, s ∈S, n ∈N, ∀t∈Tfor the considered scenario e Ensure: set of non-dominated ordered price trajectories Te 1: init Te=∅ 2: init domt=true 3: for all i= 1,...,|T|do 4: set t=ti 5: init tdomt0=true 6: init t0domt =true 7: if domt=true then 8: continue 9: end if 10: for all j=i+ 1,...,|T|do 11: set t0=tj 12: if domt0=true then 13: continue 14: end if 15: for all b∈Bdo 16: for all s∈Sdo 17: for all n∈Ndo 18: if ¯aei→t b,s,n >¯aei→t0 b,s,n then 19: t0domt =false 20: else 21: if ¯aei→t b,s,n <¯aei→t0 b,s,n then 22: tdomt0=false 23: end if 24: end if 25: if tdomt0=t0domt =false then 26: goto 39 27: end if 28: end for 29: end for 30: end for 31: if tdomt0then 32: set domt0=true 33: else 34: if t0domt then 35: set domt=true 36: end if 37: end if 38: end for 39: if domt=false then 40: set Te=Te∪ {t} 41: end if 42: end for CHAPTER 9. SOLVING ISPO 106 Algorithm 11 ISPO-DFS Require: depth/scenario index j width/index iof ordered price trajectories bounds Γ (e→te)for all e∈E,te∈Te set of maps WEj−1for Ej−1:= {ei|i<j} current global upper bound z∗ Ensure: possibly updated bounds Γ (e→te)and/or updated primal bound z∗ 1: WEj=WEj−1∪ {ej→tej i} 2: compute the (mixed) relaxed wait-and-see solution ˆ ΓWEjby ˆ ΓWEj=Pe→t∈WEjΓ (e→t) + P|E| o=j+1 Γ (eo)max 3: if ˆ ΓWEj≤z∗then 4: return 5: end if 6: try to improve the dual bound by 7: ˆ ΓWEj=ISPO-ECBWEj,ˆ ΓWEj, jand or 8: ˆ ΓWEj=ISPO-ELPBWEj,ˆ ΓWEj, jand or 9: ˆ ΓWEj=ISPO-LPBWEj,ˆ ΓWEj,nΓ (e→t), e →t∈WEjo 10: if ˆ ΓWEj≤z∗then 11: return 12: end if 13: if j < |E|then 14: for all i= 1,...,|Tej+1 |do 15: ISPO-DFSj+ 1, i, W Ej 16: end for 17: else 18: solve the SLDP(WEj)→optimal objective value z∗ SLDPWEjand optimal supply policy x∗ SLDPWEj 19: if z∗ SLDPWEj> z∗then 20: z∗=z∗ SLDPWEj 21: x∗=x∗ SLDPWEj 22: end if 23: end if Algorithm 12 ISPO-LPB Require: set of maps WEjof price trajectories fixed to scenarios dual bound ˆ ΓWEjfor the given set of maps dual bounds Γ (e→t)for each single map e→tin WEj Ensure: possibly improved dual bound ˆ ΓWEjbased on the LP-relaxed wait-and-see solution possibly improved dual bounds Γ (eo→teo)for o≤j 1: for all o= 1,...,jdo 2: if SLDP-LP({eo→teo})not solved until now then 3: solve SLDP-LP({eo→teo})→optimal objective function value z∗ SLDP-LP({eo→teo}) 4: if z∗ SLDP-LP({eo→teo})<Γ (eo→teo)then 5: set Γold = Γ (eo→teo) 6: Γ (eo→teo) = z∗ SLDP-LP({eo→teo}) 7: if argmax(eo) = teothen 8: arg(eo)max = arg maxt0∈Teo{Γeo→t0} 9: Γmax(eo) = maxt0∈Teo{Γeo→t0} 10: end if 11: ˆ ΓWEj=ˆ ΓWEj−Γold +z∗ SLDP-LP({eo→teo}) 12: if ˆ ΓWEj< z∗then 13: return 14: end if 15: end if 16: end if 17: end for 18: return ˆ ΓWEj CHAPTER 9. SOLVING ISPO 107 Algorithm 13 ISPO-ECB Require: set of maps WEjof price trajectories fixed to scenarios dual bound ˆ ΓWEjfor the given set of maps current depth/scenario index j Ensure: possibly improved dual bound ˆ ΓWEjbased on the single-supply-relaxed extended wait-and-see solution 1: solve SLDP-CB(WEj) 2: if z∗ SLDP-CBWEj+P|E| o=j+1 Γ (eo)max <ˆ ΓWEjthen 3: ˆ ΓWEj=z∗ SLDP-CBWEj+P|E| o=j+1 Γ(eo)max 4: end if 5: return ˆ ΓWEj Algorithm 14 ISPO-ELPB Require: set of maps WEjof price trajectories fixed to scenarios dual bound ˆ ΓWEjfor the given set of maps current depth/scenario index j Ensure: possibly improved dual bound ˆ ΓWEjbased on the LP-relaxed extended wait-and-see solution 1: solve SLDP-LP WEj 2: if z∗ SLDP-LPWEj+P|E| o=j+1 Γ (eo)max <ˆ ΓWEjthen 3: ˆ ΓWEj=z∗ SLDP-LPWEj+P|E| o=j+1 Γmax(eo) 4: end if 5: return ˆ ΓWEj 9.1.3 Computational results We compare different settings for ISPO-BAB in terms of the applied dual bounds. We performed ISPO-BAB for all 166 instances from our test-set Itest 6, see Appendix E. For these tests we first abandon dominance checks, Algorithm 10. We compare different combinations of combinatorial bounds CB, extended combinatorial bounds ECB, LP bounds LPB and extended LP bounds ELPB. Additionally we applied dominance “DOM” in terms of the price trajectories to reduce the size of our Branch&Bound tree a priori. In Figure 9.1 we see how many leaves of the tree can not be pruned by applying the given combinations of bounds. In many cases the number of remaining leaves can be much reduced by applying the extended combinatorial bounds against the case where only combinatorial bounds are used (from averagely 1.96% to 1.32%). If we use combinatorial bounds together with LP bounds or extended LP bounds there are averagely only 0.89% or 0.80% leaves remaining. Combining combinatorial bounds, extended combinatorial bounds and LP bounds together yields averagely 0,68% remaining leaves. If we replace in this case the LP bounds by extended LP bounds we can reduce the number of remaining leaves averagely to 0.65%. We depict the improvements by additionally applying dominance in terms of using combinatorial, extended combinatorial and extended LP bounds in Figure 9.2. We see how strong we can reduce the number of remaining leaves in the tree “not pruned by dom.” by applying dominance checks for the price trajectories a priori – averagely to 1.68 percents. While the number of remaining leaves in the Branch&Bound tree – or to be solved SLDP(WE)s – only reduces about averagely 0.18 percents “diff ILPs” and CHAPTER 9. SOLVING ISPO 108 0 2 4 6 8 10 12 14 percentage remaining leaves instance CB CB+ECB CB+LPB CB+ELPB CB+ECB+LPB CB+ECB+ELPB Figure 9.1: ISPO-BAB – percentage of non-pruned leaves applying different bounds 0 5 10 15 20 25 percentage instance not pruned by dom. diff ECBs diff ELPBs diff ILPs Figure 9.2: ISPO-BAB – Applying dominance for the price trajectories the number of to be solved extended combinatorial bounds “diff ECBs” to 0.44 percents we get averagely about 4.33 percents extended LP bounds “diff ELPBs” fewer to solve by applying dominance checks for price trajectories. Now we want to concentrate our attention on the computation time of ISPO-BAB in terms of using different dual bounds. For every instance from the set Itest 6for the stated bounds the average computation time at the leaves is depicted in Figure 9.3. If we neglect the scaling Figure 9.3 looks like Figure 9.1. Roughly speaking the figures illustrate a typical property of Branch&Bound algorithms: the worse the dual bounds the higher the computation times. In average we get 54,31 seconds for just using combinatorial bounds, 40.33 seconds for additionally using extended combinatorial bounds, 38.21 seconds for combinatorial and LP bounds, 38.36 for the combination of combinatorial bounds and extended LP bounds, 32.83 seconds by combining combinatorial, extended combinatorial and LP bounds and 34.19 seconds for using combinatorial, extended combinatorial and extended LP bounds. Averagely 0.03% fewer ILPs have to be solved by applying ELPBs instead of LPBs. If we add dominance in the last case we can reduce the runtime averagely to 21.73 seconds. CHAPTER 9. SOLVING ISPO 109 0 50 100 150 200 250 300 computation time (s) instance CB CB+ECB CB+LPB CB+ELPB CB+ECB+LPB CB+ECB+ELPB DOM+CB+ECB+ELPB Figure 9.3: ISPO-BAB – solving time applying different bounds 9.1.4 ISPO-BAB applied to the accompanying example We apply Algorithm 9to our accompanying example from Chapter 4. Example 10. A first step in the algorithm is the computation of the optimal objective function value of the SLDP-CB({e→t}) for each price trajectory tand scenario e. Therefore we use the single supply revenues, Algorithm 6. The results for each scenario, branch and size can be found in Appendix C. Based on these data the single supply relaxation SLDP-CB({e→t}) is solved to optimality via Algorithm 6. The results are stated in the following tables. For each price trajectory we state its index and the optimal supply for each pair (b,s) with b∈Band s∈S. The optimal objective function values of the SLDP-CB({e→t}) are always stated in the last column. low seller price trajectory index (1,S) (1,L) (2,S) (2,L) bound (0,0,0,0,3) 0 3 3 2 2 18.86 (0,0,0,1,3) 1 3 3 2 2 18.31 (0,0,0,2,3) 2 3 3 2 2 17.64 (0,0,1,1,3) 3 3 4 1 2 16.86 (0,0,1,2,3) 4 3 4 1 2 15.96 (0,0,2,2,3) 5 3 4 1 2 14.39 normal seller price trajectory index (1,S) (1,L) (2,S) (2,L) bound (0,0,0,0,3) 0 2 4 1 3 51.58 (0,0,0,1,3) 1 2 4 1 3 51.10 (0,0,0,2,3) 2 2 4 1 3 51.10 (0,0,1,1,3) 3 2 4 1 3 51.09 (0,0,1,2,3) 4 2 4 1 3 50.61 (0,0,2,2,3) 5 2 4 1 3 51.09 high seller price trajectory index (1,S) (1,L) (2,S) (2,L) bound (0,0,0,0,3) 0 2 4 1 3 31.00 (0,0,0,1,3) 1 2 4 1 3 30.71 (0,0,0,2,3) 2 2 4 1 3 30.71 (0,0,1,1,3) 3 2 4 1 3 30.71 (0,0,1,2,3) 4 2 4 1 3 30.41 (0,0,2,2,3) 5 2 4 1 3 30.71 We order the price trajectories according to the revenues for the particular scenario and get: CHAPTER 9. SOLVING ISPO 110 scenario 1st 2nd 3rd 4th 5th 6th low 0 1 2 3 4 5 normal 0 1 2 3 5 4 high 0 1 2 3 5 4 To illustrate our Branch&Bound approach we first abandon excluding price trajectories by dominance. For this example we order the scenarios in the sequence “high-normal-low” and start with the depth-first-search, Algorithm 11. The related Branch&Bound tree is depicted in Figure 9.4. The combinatorial bounds CB for the visited nodes are stated. The indices of the mapped price trajectories are always stated in the nodes. In the depth-first-search the nodes are visited in the order according to the numbers right above the nodes. At first we fix the first price trajectory for the “high seller” scenario. For this scenario we obtain an expected revenue of 31.00. For the remaining scenarios yet we have not fixed a price trajectory: In terms of bounding we have always to consider the price trajectory that yields the highest revenue in terms of a single supply. For both, the “normal” and the “low” scenario, this is the first trajectory which yields 51.58 and 18.86 of revenue. Thus, the combinatorial bound amounts to 31.00 + 51.58 + 18.86 + 101.44. The bound stays the same for the next two nodes. At Depth 3we solve the SLDP(WE) for WE={high seller →0,normal seller →0,low seller →0}. This yields a primal bound of 101.38. We continue with fixing the first price trajectory for the low seller scenario. This yields an upper bound of 100.88. Because it is 100.88 <101.44 we prune the current branch. This is also the same for fixing the remaining price trajectories to the low seller scenario at this point. We continue with fixing Price Trajectory 1 to the normal seller scenario. The related objective function value of the SLDP-CB amounts to 51.10. For the non-set low seller scenario again we have to add the expected revenue of 18.86. The expected revenue for the high seller scenario stays 31.00 – because we still fix Price Trajectory 0 to it. Thus, for the set of maps WE0={high seller →0,normal seller →1}our current dual bound is 100.96. Because 100.96 <101.38 we prune the current branch. This is also the same for the following nodes at Depth 2. Now we fix Price Trajectory 1 to the high seller scenario. This results in an expected revenue of 30.71. Again, we have to choose the maximum combinatorial bound for the remaining scenarios to obtain an upper bound for the wait-and-see solutions of the related childnodes. Altogether this yields a revenue of 101.14 and 101.14 <101.38. We can prune the tree also for the remaining nodes at Depth 1 and at the end obtain the optimal solution value of 101.44. In the related optimal solution we just delivered Lot-type (1,1). Both branches obtain this lot-type in Multiplicity 3. In this example the computation of other bounds than the combinatorial bound CB was not necessary because we were able to prune all branches of the tree except the first one directly. We did not apply dominance rules for the price trajectories. If we regarded dominance, the price trajectory with index 0 would dominate all other price trajectories for every particular scenario. Only the first branch of the Branch&Bound tree would remain. CHAPTER 9. SOLVING ISPO 117 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 gap (percent) instance best start trajectories iteration: 0.5 1 1.5 2 0 1 2 3 4 5 6 7 gap (percent) instance first start trajectories iteration: 0.5 1 1.5 2 0 10 20 30 40 50 60 gap (percent) instance random start trajectories iteration: 0.5 1 1.5 2 2.5 0 10 20 30 40 50 60 70 gap (percent) instance worst start trajectories iteration: 0.5 1 1.5 2 2.5 3 Figure 9.8: ISPO-PingPong – Progress of gaps CHAPTER 9. SOLVING ISPO 118 tion problem. The lower-level problem considers a variable xas a parameter to compute the optimal value of a variable ywhile the upper-level problem obtains the optimal value of xby using the value of ycomputed in the lower-level problem [CCMGB10]. In our case – by virtue of reversible recourse – we can see the size optimization stage and also the price optimization stage as both, as upper-level and as lower-level subproblems. 9.3 Computational results for real-world instances We now compare ISPO-BAB with ISPO-PingPong for real instances from the set I, Appendix E. We applied ISPO-BAB by using combinatorial bounds CB, extended combinatorial bounds ECB and extended LP bounds ELPB; LP bounds LPB are only applied in depth 1of the Branch&Bound tree. Moreover, in a preprocessing step the state space of valid price trajectories is restricted by dominance DOM, see Algorithm 10. We performed ISPO-PingPong for two different settings: 1. nrκ= 1,tb = 20 seconds, 2. nrκ= 100,tb = 60 seconds. In the second column in Table 9.3 “non.dom.(%)” we state the number(percentage) of remaining leaves in the enumeration tree after applying the dominance check for the price trajectories. By checking dominance most of the price trajectories can be excluded a priori , averagely more than 99.995 percent and so the number of leaves in the Branch&Bound tree reduces extremely. In the third column we state the number(percentage) of the solved SLDP(WE)s “#ILP(%)”. The percentages in the fourth and the subsequent columns are related to the number of remaining leaves we state at column “non.dom”. We see that in relation only a small number of all possible SLDP(WE)s has to be solved. If we would relate the number of solved SLDP(WE)s to the number of leaves of the hole enumeration tree without the restriction to non dominated price trajectories, the maximum number of 104 solved SLDP(WE)s for Instance 2 would mean that just 0.0002156% of all possible SLDP(WE)s had to be solved. In the next seven columns for all applied dual bounds overall solving time “t. . . ” in seconds “(s)” or hours “(h)” and the number “#. . . ” of the computed related bounds are stated. Because we compute combinatorial bounds for all scenarios and price trajectories a priori, i.e. 100% of these bounds, we did not state this specific number. Let us take look at the computation times for the different bounds. We therefore take the average computation times in the last line of the table and divide it by the related average number of solved bounds. The computation time for combinatorial and extended combinatorial bounds in relation to (extended) LP bounds is very small: For one ECB averagely it amounts to about 0.14 seconds while the average time needed for computing the LP bound amounts to nearly one hour and for the extended LP bound to about 9minutes. But for LP bounds and extended LP bounds in relation to the solution of the binary program SLDP(WE) the computational effort is yet small. The overall computation time minus the time for all bounds except CB yield the approximate average computation time for the SLDP(WE)s. It is more than 110 hours. (The time for traversing the nodes in the tree is neglectable.) Not every SLDP(WE) is solved until the internal Branch&Bound process of the MIP solver starts. In Table 9.3 we did not distinguish between SLDP(WE)s solved up to optimality and SLDP(WE)s for which only the LP relaxation is solved. Because we CHAPTER 9. SOLVING ISPO 119 hand over the current best bound to the MIP solver it may be that the solving process breaks up earlier. Moreover, because the ELPB at depth |E|equals the LP relaxation of the related SLDP(WE) we do not compute the bound ELPB explicitly at depth |E|. So the columns “#ELP(%)” and “tELP(h)” in our case of three scenarios are only related to depth 2. At the columns with label “ILP∗(%)” and “t∗” we state the number of SLDP(WE)s that have to be solved at the leaves and the solving time until the optimal solution was found. We see that mostly the optimal solution is found very fast after averagely 5.9 solved SLDP(WE)s. On the other hand altogether averagely 64.3SLDP(WE)s have to be considered. This speaks for our preprocessing step in which for the particular scenarios the price trajectories according to their single supply relaxations are ordered. The results of ISPO-PingPong for the two settings are stated in the last six columns. We state the solving time in minutes “t/min”, the number of iterations “#iter” and the optimality gap “gap(%)”. Also with the first setting – i.e. at every iteration of the size optimization stage ISPO-SF only one SLDP(WE0) is solved – we mostly can get proper results. For the second setting where we traverse the 100 “best” κ-combinations we get no gap higher than 0.69%. The average solving time amounts to 12.77 minutes. CHAPTER 9. SOLVING ISPO 120 i ISPO-BAB ISPO-PingPong1 ISPO-PingPong2 t/h non.dom. #ILP(%) tCB(s) #ECB(%) tECB(s) #LP(%) tLP(h) #ELP(%) tELP(h) ILP∗(%) t∗/h t/min #iter gap(%) t/min #iter gap(%) 1 47.75 1331(0.0057) 81(6.09) 56.69 101(6.96) 8.44 5(15.15) 2.68 21(17.36) 2.46 2(0.15) 9.60 1.05 2 0.30 22.41 2 0.08 2 88.64 1331(0.0057) 104(7.91) 56.58 123(8.47) 9.55 4(12.12) 3.94 20(16.53) 3.82 30(2.25) 45.80 1.01 1.5 0.26 12.61 2 0.01 3 130.93 1331(0.0057) 85(6.39) 89.36 113(7.78) 19.00 5(15.15) 5.68 29(23.97) 6.19 23(1.73) 46.21 1.05 1.5 0.65 20.84 1.5 0.09 4 214.72 1331(0.0057) 90(6.76) 87.76 110(7.58) 22.86 5(15.15) 5.66 21(17.36) 4.98 2(0.15) 12.61 1.06 2 0.31 24.67 1.5 0.06 5 32.58 1331(0.0057) 34(2.55) 57.32 48(3.31) 3.55 4(12.12) 2.21 15(12.40) 1.92 2(0.15) 29.54 0.94 1 0.69 6.50 1 0.69 6 1.29 1(0.0000) 1(100) 56.96 0(0.00) 0.00 0(0.00) 0.00 0(0.00) 0.00 1(100) 1.23 0.93 1 0.50 6.47 1 0.50 7 51.27 1(0.0000) 1(100) 69.34 0(0.00) 0.00 0(0.00) 0.00 0(0.00) 0.00 1(100) 51.25 1.08 1 0.90 8.55 1 0.18 8 4.97 1331(0.0057) 28(2.10) 50.65 44(3.03) 2.33 4(12.12) 2.11 17(14.05) 0.36 1(0.08) 0.19 0.88 2 1.32 29.56 1.5 0.20 9 3.91 1331(0.0057) 21(1.58) 50.26 35(2.41) 1.75 4(12.12) 1.98 15(12.40) 0.35 1(0.08) 0.18 0.86 1.5 1.68 11.59 1 0.29 11 12.53 1331(0.0057) 58(4.36) 52.58 83(5.72) 5.28 5(15.15) 1.34 26(21.49) 0.81 3(0.23) 0.94 0.86 1.5 1.78 29.73 1.5 0.23 12 10.84 1331(0.0057) 47(3.53) 52.38 68(4.68) 4.35 5(15.15) 1.54 22(18.18) 0.71 1(0.08) 0.65 0.95 1.5 1.69 12.11 1 0.20 13 10.26 1331(0.0057) 45(3.38) 52.06 64(4.41) 4.00 5(15.15) 1.43 20(16.53) 0.67 1(0.08) 0.63 0.94 1.5 1.83 8.74 1 0.21 44 112.69 1331(0.0057) 92(6.91) 56.92 112(7.71) 9.11 5(15.15) 2.56 21(17.36) 2.45 2(0.15) 9.51 1.02 1 1.89 12.95 1 0.23 45 123.71 1331(0.0057) 81(6.09) 58.07 113(7.78) 9.08 5(15.15) 2.62 33(27.27) 3.60 29(2.18) 62.55 1.11 2 0.26 12.43 2 0.00 46 197.57 1331(0.0057) 80(6.01) 89.76 112(7.71) 18.22 5(15.15) 5.78 33(27.27) 7.00 1(0.08) 6.23 0.99 1 0.26 6.66 1 0.00 48 67.77 1331(0.0057) 53(3.98) 58.03 71(4.89) 5.79 4(12.12) 2.34 19(15.70) 2.26 2(0.15) 7.21 1.09 1.5 0.52 12.37 1.5 0.52 49 4.58 1(0.0000) 1(100) 93.19 0(0.00) 0.00 0(0.00) 0.00 0(0.00) 0.00 1(100) 4.56 0.99 1 0.16 7.71 1 0.16 50 67.59 1331(0.0057) 45(3.38) 88.87 65(4.48) 9.98 4(12.12) 5.39 21(17.36) 4.98 1(0.08) 4.55 1.00 1 4.12 8.05 1 0.22 52 4.58 1(0.0000) 1(100) 90.09 0(0.00) 0.00 0(0.00) 0.00 0(0.00) 0.00 1(100) 4.55 0.99 1 4.07 7.93 1 0.21 53 340.49 1331(0.0057) 129(9.69) 90.78 159(10.95) 25.22 5(15.15) 6.84 31(25.62) 6.96 43(3.23) 125.05 1.09 1.5 0.11 14.16 1.5 0.00 61 237.59 1331(0.0057) 90(6.76) 92.71 119(8.20) 20.24 6(18.18) 13.57 30(24.79) 7.20 1(0.08) 4.72 1.12 1.5 0.36 11.99 1.5 0.05 62 291.86 1331(0.0057) 102(7.66) 91.61 129(8.88) 23.25 5(15.15) 6.21 28(23.14) 6.33 1(0.08) 4.61 1.01 1 0.33 6.80 1 0.14 63 296.24 1331(0.0057) 116(8.72) 91.70 141(9.71) 23.81 5(15.15) 4.59 26(21.49) 4.96 2(0.15) 10.74 1.11 2 0.19 13.77 1.5 0.12 64 328.68 1331(0.0057) 118(8.87) 91.19 143(9.85) 26.45 5(15.15) 5.80 26(21.49) 5.80 1(0.08) 4.73 1.11 2 0.20 7.60 1 0.12 65 1.17 1(0.0000) 1(100) 91.29 0(0.00) 0.00 0(0.00) 0.00 0(0.00) 0.00 1(100) 1.17 0.99 1 0.25 7.65 1 0.18 66 154.40 1331(0.0057) 106(7.96) 91.51 129(8.88) 23.54 5(15.15) 5.83 24(19.83) 5.45 2(0.15) 9.32 1.12 2 0.25 14.15 2 0.17 67 323.95 1331(0.0057) 125(9.39) 90.88 156(10.74) 29.11 6(18.18) 9.86 32(26.45) 5.97 2(0.15) 9.37 1.0 1 0.33 6.77 1 0.15 ∅117.13 1085(0.0046) 64.3(23.34) 74.02 82.90(5.71) 11.29 4.00(11.90) 3.70 20.00(16.22) 3.16 5.9(18.95) 17.32 1.01 1.43 0.93 12.77 1.30 0.19 Table 9.1: ISPO-BAB and ISPO-PingPong applied on real instances CHAPTER 9. SOLVING ISPO 121 9.4 General goodness of ISPO-PingPong In the previous section we showed that the optimality gap for our heuristic for realworld instances averagely amounts to 0.19% and in the worst case to 0.69% depending on the setting. Now we apply the heuristic on a small example to show that we cannot guarantee such small gaps in general. As in our accompanying example we consider only two branches B={1,2}and two sizes S={S,L}. Only one lot-type is allowed for supplying the two branches, i.e. κ= 1. We consider all lot-types with at minimum one item per size, i.e. vmin = 1 and at maximum two items per size, vmax = 2. We only allow lot-types with cardinality 3. This yields the lot-types (1,2) and (2,1). The set of multiplicities is given by M={1,2,3}. We set I= 0 and I= 20. In this example we assume that pick-cost and lot-opening costs take value zero. We consider four sales periods including the sellout period kmax = 3. The discounting factor is also set to zero, ρ= 0. Moreover we set the fixed and variable mark-down costs to zero, i.e. µv=µf= 0. We are given four prices including the salvage value. It is P={0,1,2,3}with π0= 10, π1= 9, π2= 8 and π3= 4. We number the price trajectories taccording to the following table: index t 1(0,0,0,3) 2(0,0,1,3) 3(0,0,2,3) For simplicity’s sake we consider only one scenario with probability one. The demands per period kand price index pare given as stated in the following tables. (1,S) k/p012 0 1.0- - 1 0.0- - 2 0.0 3.0 3.1 (1,L) k/p012 0 1.0- - 1 1.0- - 2 0.0 1.1 3.0 (2,S) k/p0 1 2 0 1.0- - 1 1.0- - 2 0.0 1.1 1.1 (2,L) k/p012 0 1.0- - 1 0.0- - 2 0.0 0.0 0.0 We state the single supply revenues computed by Algorithm 6in the table found below. In the last column we are given the additional revenue for the demand-exceeding numbers of items according to Observation 2. Because we are given no mark-down costs, the additional revenue always amounts to πkmax −ap = 5 −4 = −1. (branch,size) t/n1 2 3 4 5 6 6+ (1,S) 1 5.0 4.0 3.0 2.0 1.0 0.0−1.0 2 5.0 9.0 13.0 17.0 16.0 15.0−1.0 3 5.0 8.0 11.0 14.0 13.4 12.4−1.0 (1,L) 1 5.0 10.0 9.0 8.0 7.0 6.0−1.0 2 5.0 10.0 14.0 13.5 12.5 11.5−1.0 3 5.0 10.0 13.0 16.0 19.0 18.0−1.0 (2,S) 1 5.0 10.0 9.0 8.0 7.0 6.0−1.0 2 5.0 10.0 14.0 13.5 12.5 11.5−1.0 3 5.0 10.0 13.0 12.4 11.4 10.4−1.0 (2,S) 1 5.0 4.0 3.0 2.0 1.0 0.0−1.0 2 5.0 4.0 3.0 2.0 1.0 0.0−1.0 3 5.0 4.0 3.0 2.0 1.0 0.0−1.0 CHAPTER 9. SOLVING ISPO 122 At the beginning of ISPO-PingPong a price trajectory for each scenario is fixed. In this example at first we choose the “first” price trajectory 1. For this trajectory the best supply in terms of lot-types is computed. Because κ= 1 we can only choose one lot-type, either (1,2) or (2,1), for supply. According to the single supply revenues an optimal supply for Price trajectory 1is given by choosing Lot-type (1,2) in Multiplicity 1for Branch 1 and in Multiplicity 2 for Branch 2. This yields a revenue of 27. The next step is to reject the fixed price trajectory and if possible to choose one which is optimal for the given supply. Otherwise the heuristic converges. Adding up the corresponding single supply revenues yields that no other price trajectory is better for the given supply. For the price trajectories 2and 3the revenues for the given supply also amount to 27. The approach converges with objective function value 27. An alternative would be to start with the “best” price trajectory – that means the price trajectory twhich yields for our scenario ethe highest objective function value of the single supply relaxation SLDP-CB({e→t}). This is Price trajectory 3. The best revenue in terms of single supplies for this trajectory is 51 which is given by supply (4,5) for Branch 1 and (3,1) for Branch 2. Price trajectory 1yields just an revenue of 30 with supply (1,2) for Branch 1 and (2,1) for Branch 2. The revenue for Trajectory 2 amounts to 50 and results from a delivery of (4,3) for Branch 1 and (3,1) for Branch 2. That means we start with Price trajectory 3and compute the best supply in terms of κ= 1 chosen lot-type. The optimal supply for this trajectory is given by supplying 3times Lot-type (1,2) to both branches. The resulting revenue amounts to 42. Now we have to reject the current price trajectory again and to compute the best trajectory in terms of the current supply. According to the single supply revenues this is also Price trajectory 3. For Trajectory 1the revenue amounts to 18 and for Trajectory 2 it amounts to 38.5. That means the approach converges with a revenue of 42. Solving ISPO for this example to optimality yields an optimal solution value of 46.50 which results from supplying Lot-type (2,1) three times to Branch 1 and two times to Branch 2. Price trajectory 1is optimal. For this example neither using the “best” trajectory nor the “first” trajectory – which is also the “worst” – as start trajectories lead to a solution with small gaps. The gaps amount 41.94% for the “first” and 9.68% for the “best” trajectory. Thus, in general we can not assume such a good performance as on our instances. We want to state some specialties of our real-world instances: For our instances the range between the lower and upper bound as a rule amounts to about maximal 10% of I+I−I 2, see Remark 1. In this example it is much higher, namely 200%. Moreover in the real-world case the deviation among the demands is with mostly zero or one sold items less. Although at this point we can not evidence a general warranty of goodness for ISPO-PingPong, the bounds in terms of our real-world instances are small enough to justify practical use. 9.5 Conclusion of the chapter We presented the Branch&Bound solver ISPO-BAB which solves the Integrated Size and Price Optimization Problem for all tested instances to optimality. We branch on maps “scenario to price trajectory”. Dual bounds are obtained by extensions of the wait-and-see solution from stochastic programming. CHAPTER 9. SOLVING ISPO 123 For practical use at our industrial partner we propose the heuristic ISPO-PingPong. The principle is to alternate size and price optimization until convergence. This is possible because of the special structure of ISPO – the reversible recourse. We show that ISPO-PingPong for the tested real-world instances yields solutions with an average optimality gap of 0.19% – the maximum gap amounts 0.69% – in averagely 12.77 minutes. Chapter 10 DISPO in practical application – real-world experiments The collaboration with our industrial partner gave us the opportunity to test the practical relevance of DISPO in a real-world field study. With real-world experiments we want to verify that DISPO performs better than the method currently in use at the partner, i.e. the LDP together with a manual determination of mark-downs. For that reason the DISPO-team performed a so-called single-blind experiment where test and control branches compete against each other. We give some basics about statistical experiments in Section 10.1 before we apply them in Section 10.2 to our field studies. During the cooperation also field studies only in terms of price optimization were performed. In Section 10.3 we will outline the main results. We show how heavily mark-downs can directly affect the number of sales and that price optimization can also increase the realized revenue. Because performing a field study is expensive in terms of work and money we estimated the potential of improvement a change from LDP to ISPO-based supply would bring along. The results – which we outline in Section 10.4 – convinced our partner and the DISPO-team to perform a five-month field study. With DISPO we could finally increase the realized revenue about more than 1.5 percentage points. Moreover by regarding the field study as a statistical experiment we can give a statement about the significance of the result. We present the details in Section 10.5. 10.1 Performing statistical experiments With real-world experiments we want to figure out if DISPO – or also price optimization as a part of DISPO – performs better than the methods currently implemented at our partner. Moreover we want the results to be statistically significant. In order to apply statistical methods later on we want to introduce some statistical basics in this section. We start with a classification of blind experiments in Subsection 10.1.1. We outline the term statistical significance in Subsection 10.1.2. To make a point about statistical significance we have to perform a test of significance. We mention the most common tests in Subsection 10.1.3 before we focus on the Wilcoxon signed-rank test in Subsection 10.1.4. 124 CHAPTER 10. REAL-WORLD EXPERIMENTS 125 In this section we are mainly guided by [FPP07], [Raj06] and [Kan06]. 10.1.1 Blind experiments A blind experiment is a statistical experiment where not all people involved are informed about certain aspects to avoid bias. Blind experiments are typically applied in medical tests. The group of probands is divided into a test and a control group where the test group gets the medicament and the control group just a placebo. If one wants to examine the effect of a medicament it is usual that the experimentees are not informed in which group they are. Otherwise the experimentees might be affected by this information. If all other involved persons – except the experimentees – have full information about the categorization we talk about a single-blind experiment. In some cases it is useful that also the researchers do not know about the category of the tested persons. They might treat the probands accordingly. In this case we would talk about a double-blind experiment. 10.1.2 Statistical significance If our new method performs worse or better than the method currently in use we want to state how big the role of chance for this result was. If the probability that the result could be caused by pure chance is not small enough we would not give general statements about a better or worse performance. The so-called null hypothesis says that the method leads to no differences (or also to no better/worse performance). An alternative hypothesis that it does. A test is statistically significant if the probability that its outcome is the result of chance is smaller than a predefined significance level. A common choice for the significance level is 5%. To make a point about statistical significance a so-called test staticstic is computed. A test statistic is defined as a measurement of the difference between the data and the statement of the null hypothesis, [FPP07]. The test statistic follows a test distribution. If the probability to obtain the test statistic under the test distribution – the so-called p-value – is smaller than the significance level – then we call the result statistically significant – we reject the null hypothesis and rely on the alternative hypothesis. Rejecting the null hypothesis does not mean that the alternative hypothesis is true. It is only an evidence that the result is not caused by random. 10.1.3 Statistical tests in general To evidence statistical significance there are several statistical tests. There are tests for related samples and tests for unrelated samples. A sample is called unrelated if groups of different individuals are compared. For related samples we compare groups with individuals related pairwise to an individual from the other group. Which test can be applied also depends on the kind of the data: Are the observations nominal, ordinal or given by a distribution? Parametric tests assume a specific distribution while nonparametric tests do not. We distinguish between two-sided and one-sided tests. A two-sided test considers both, a better and a worse performance of the test sample simultaneously. The nullhypothesis says that both methods perform the same way, the alternative hypothesis that they do not. With a one-sided test we are only interested in a better/worse performance of the “new method”. The null-hypothesis says that it performs not better/not worse, CHAPTER 10. REAL-WORLD EXPERIMENTS 126 the alternative that it does. Because in our case we are interested in a worse or better performance we focus on one-sided tests in the following. A common approach for two unrelated normal distributed samples is Student’s ttest. The t-test compares differences between the means of the particular samples and compares them with the corresponding standard error to determine if the two samples arise from the same distribution. For related normal distributed samples the t-test can also be applied in a similar way. For further information see for example [FPP07]. In the case that no distribution for the observations can be assumed (but also for observations from a specific distribution) non-parametric tests can be applied. The tests for non-parametric ordinal data assign ranks to the observations. As test statistic rank-sums are computed. The test distribution is the distribution of the rank-sums. For unrelated samples the Mann-Whitney test is commonly used. At first all observations are ordered increasingly and ranks are assigned in terms of the ordering. Then, by summing up all ranks of one sample the rank sum – the test statistic – is computed. To determine the role of random one computes the p-value as the probability to get the observed rank-sum (or for one-sided tests the observed rank-sum or a higher/lower one) among all other possible rank-sums. For related samples there is a similar approach named Wilcoxon signed-rank test. 10.1.4 Wilcoxon signed-rank test In order to certify statistical significance for two related ordinal samples Wilcoxon signed-rank test is very common. The test is named after Frank Wilcoxon who presented it together with the rank sum test for non-paired observations also called MannWhitney test in [Wil45]. Wilcoxon signed-rank test is an alternative to the Student’s t-test if no normal distribution can be assumed. It yields a statement about the symmetric distribution of the pair differences around the median. The test can be performed as one-sided or two-sided test. For our purposes only the one-sided test is relevant. Therefore we formulate Wilcoxon signed-rank test as one-sided test. It is checked if the differences of the ordered paired observations (test −control) are distributed symmetrically around or right of the median ˜xor symmetrically around or left of the median ˜x. Thus, the null hypothesis in the first case is H0: ˜x≤0,(10.1) and the alternative hypothesis H1: ˜x > 0.(10.2) In the second case, the null-hypothesis is H0: ˜x≥0,(10.3) and the alternative hypothesis H1: ˜x < 0.(10.4) In the first case the null hypothesis is equivalent to the statement that the distributions of the paired observations for the two samples are identical or that the distribution of the test sample is shift to the left. In the second case that they are identical or that the distribution of the test sample is shift to the right. We now describe the different steps for performing the Wilcoxon signed-rank test. They are illustrated on the following small example. CHAPTER 10. REAL-WORLD EXPERIMENTS 133 sample relative realized objective gross yield sales test 0.4774 0.6524 0.7891 control 0.4698 0.6494 0.7945 Table 10.1: Performance metrics – 2nd field study POP-RH Originally there were 4668 products included in the test. To get possibly nonbiased result our evaluation uses only 2298 from them which were supplied to all test and control branches. We state some performance metrics for the test and control branches in Table 10.1. We will focus on the relative realized objective which we will define in the following. For each test-control pair of branches, the realized revenues over all articles in A were compared. That means, in particular, that expensive articles have a larger influence on the result than cheap articles. This point of view is in line with our partner’s point of view. For reasons of comparability we divide the realized revenue by maximum possible revenue in terms of the objective of the POPˆe. This means for an initial stock Ia b,s for the considered branch b, size sand article aand a starting price πa 0we compute the relative realized objective of the mark-down decision ta= (ta 0, ta 1, . . . , ta kmax )with t=taa∈Afor Branch bas RROPOP(b) = objective of POPˆefor bachieved by t maximal possible objective = −X a∈A X s∈S Ia b,s ·apa+X k∈K exp(−ρk)X a∈A X s∈S ˆra k,b,s −˜µa kˆna k Pa∈A Ps∈SIa b,s ·(πa 0−apa).(10.6) Depending on Article a,apadenotes the acquisition price, πa 0the starting price and Ia b,s the initial stock per branch and size. During the sales process, we observed ˆra k,b,s (the realized yield for Branch band Size sin Period k) for Article aand ˆna k(mark-down in Period k– yes or no). Since we only consider a subset of branches we have to take into account that fixed mark-down costs must be scaled with respect to the number of considered branches. This way, we get mark-down costs ˜µa kfor period k. Because all test and control branches were supplied by our partner and now our focus lies on the price optimization stage we do not regard costs for supply in terms of lots as lot-opening costs and pick costs as they appear in ISPO. We see that the mean relative realized objective in the test branches is about 0.76 percentage points higher than in the control branches. Also in terms of the other performance metrics POP-RH beats the manual price optimization. Yet, with a rank-sum of 285 and P30(X≥285) = 14.47% and P30(X≤285) = 86% Wilcoxon signed-rank test yields no significance. Still, the p-value for a better performance caused by pure chance is with 14.47% comparatively small. CHAPTER 10. REAL-WORLD EXPERIMENTS 134 0 5 10 15 20 25 Change of supply by ISPO vs. sellout date by LDP amount of sizes per branch with sellout date weeks until sellout date supply ISPO - supply LDP 0 1 2 -1 -2 3 -3 Figure 10.2: Change of supply by ISPO 10.4 Potential of ISPO By integrating price optimization in the decision on the supply we hope to increase the realized revenue. To get an idea of how good ISPO meets the realized demand of the different sizes we performed a test on 136 real instances for which we were in possession of transaction data. Because at this time all articles at our industrial partner were supplied by the LDP we can use real sales figures to compare the two models LDP and ISPO. In Figure 10.2 we depicted the change of supply by ISPO against the supply that LDP yields. The size of the bubbles is related to the summed up number of sizes over all branches and articles which are sold out in the week marked at the x-axis. The bigger the bubble the more sizes per branch and article are sold out in that week. The position of the bubbles on the y-axis describes the difference between the supply determined by ISPO and the realized supply computed by the LDP. For example, let us consider an article for which a particular size in a particular branch is sold out 5 weeks after sales start. If for this article, branch and size ISPO would yield a supply of 5 and LDP would yield a supply of 3, then this article,branch and size would increase the relative size of the bubble at (5,2). We see, that the size of bubbles decreases in terms of the sellout date for positive values on the y-axis, while the size increases for negative ones. This means, ISPO would deliver more items per branch and size the earlier the supplied items by the LDP were sold out and ISPO would deliver less items per branch and size the later the supplied items by the LDP were sold out. This is exactly the behavior we would expect to obtain a more size conform supply and an indication that ISPO leads to a more size conform supply than the LDP. To measure the possible improvement also in terms of money we compared the CHAPTER 10. REAL-WORLD EXPERIMENTS 135 Figure 10.3: Comparision of different models for size optimization three models LDP, SLPD and ISPO. For the same 136 instances we solved the LDP and SLDP exactly, ISPO – because of the long runtime – heuristically. To compare the three models with respect to monetary effects we evaluated them on the objective of ISPO. For the optimal solutions of LDP and SLDP we computed the associated acquisition costs, pick costs and lot opening costs and then applied price optimization (Algorithm 3) to get the yield for the optimal price trajectory in terms of the supply. Subtracting all costs for supply from the yield provides the related objective value of ISPO. The gaps between the objective of the heuristic solution of ISPO and the objective in terms of ISPO by fixing the supply computed by the LDP and SLDP are depicted in Figure 10.3. These gaps can be seen as a practical application of the value of the stochastic solution VSS, see Section 5.3.2. We measure the loss by non-regarding the second stage in the sales process or particularly in terms of the LDP by non-regarding random. The vertical lines illustrate the relative improvement in terms of money for the particular instance. In terms of the LDP, ISPO could lead to an improvement about 2.27 percentage points of revenue, in terms of the SLDP about 1.52 percentage points. The theoretical potential of increasing the revenue by 2.27 percentage points was encouraging and therefore our partner and the DISPO-team decided to deploy DISPO – ISPO combined with POP-RH – in a real-world field study. CHAPTER 10. REAL-WORLD EXPERIMENTS 136 10.5 DISPO – the field study Parts of this section are presented in similar fashion in [KKR11b]. We performed our real-world field study as described above as a controlled statistical experiment to compare DISPO with the currently applied LDP together with manual price optimization.2 10.5.1 Preparation Additionally to the estimation of the potential of ISPO we outlined in Section 10.4, there were some other work we/the DISPO-team and our partner were concerned with in the time before the field study began. The DISPO-team had to confer with the partner about the exact test set-up. The in Section 10.2 described partition in test and control branches was chosen. To supply the test branches according to our proposal pre-packs had to be opened and items had to be removed or added by our industrial partner. Therefore the partner decided to use its German online-shop to compose the lots for the test branches after the lots arrived at the headquarters. Moreover the participating commodity groups were chosen. These were three commodity groups with comparatively many sizes, see 10.5.2. The involved persons also arranged an appropriate point in time for the field study . We met with those responsible of the sales department to come to an realistic estimation of the fixed and variable mark-down cost and the salvage value. It was also discussed if advertising campaigns and special offers should be allowed. The results will follow in Subsection 10.5.2. For performing POP-RH we could take the most scripts and programs developed by the former DISPO-team that were used in former field studies. To detect potential weak points and to practice the procedure on both sides, we performed a test field study together with the IT department of our partner. Additionally potential inconsistencies in terms of the used data formats should be eliminated. The test data our partner provided us included relevant data for all planned orders for the last season of the year 2010, historical data for demand estimation and current transaction data for all commodity groups. Overall this was about 6.6 GB of data. For this data – after we performed the empirical estimation – we computed a supply policy by ISPO-PingPong and provided our partner the results for the test branches. We also tested the process for performing POP-RH. That means with respect to latest sales figures we computed an optimal mark-down strategy via the POPˆeand informed our partner about proposed mark-downs for the subsequent two weeks. The test field study among others led to some adjustments of scripts and readin routines. 10.5.2 Setup of the field study It was necessary to select a small set of articles for the field study because the orders had already been placed in terms of lot-types, and the adaption of the supply for the test branches to the results of the new method is a too expensive logistic operation to be carried out for each article. 2A comparision between ISPO and complete manual planning was not possible, because nearly all articles at our industrial partner are supplied by the LDP now. CHAPTER 10. REAL-WORLD EXPERIMENTS 137 The test branches were supplied in terms of the – by ISPO-PingPong (Algorithm 15) computed – solution of ISPO. Since there are global constraints for the overall number of supplied items we actually computed the supply for all branches with the new method – and so did our industrial partner with the LDP. Our proposed supply was then implemented only for the test branches; the supply of the control branches (and also all remaining branches) was implemented as computed by the LDP by our project partner. The field study ran from end of May until end of September 2011 for 81 articles from three different commodity groups – women overgarments fashion (wof), women overgarments classic (woc) and women underwear (wu). The sales process for these articles started between May 2011 and mid of June 2011 so that all articles could be observed for a time period of 15 to 17 weeks. Some further relevant properties of the used test articles are stated in Table 10.2. Our demand estimation, see 3.2, is based on historical data in a time frame from September 2009 to September 2010. commodity group number of articles number of sizes wof 9 6 woc 9 3 wu 5 6 Table 10.2: Properties of the test articles. Table 10.3 shows the parameter setting we used in ISPO for the field study. Our lottype opening costs δiand pick costs pcost were estimated on the basis of a thourough cost accounting. This cost accounting also revealed that more than four lot-types can only be handled if the area for internal stock-turnover is increased substantially. The discount factor ρis derived from an estimation of the capital binding cost. Whenever other reasons than interest rates favor faster stock-outs this can be increased. The fact that we did not account for mark-down costs µkjust reflects the fact that at the time of the design of the experiment our partner could simply not provide a realistic value for this. An adaption of ISPO-PingPong was not possible until the beginning of the experiment due to time constraints. parameter setting κ4 pcost 0.0545 δ1100 δi,i > 1 50 E{low,normal,high}(period-0sales ∈[0,10 %)/[10 %,30 %]/(30 %,100 %]) dnormal k,p,b,s from empirical distribution and interpolation of historic sales in commodity group de k,p,b,s α×dnormal k,p,b,s (αfrom historic sales in scenario ecompared to scenario normal) Prob(e)from empirical distribution of historic sales in commodity group kobs 2(i.e., realization of eand earliest mark-down after 2periods) kmax periods (= weeks) until end of season (article dependent) ρ0.000974868 pmax 4 (five prices including start price and salvage value) µk0 πkmax dependent on commodity group ∈[15%,30%] of the starting price π0 Table 10.3: Parameter setting for the field study. CHAPTER 10. REAL-WORLD EXPERIMENTS 138 Advertising campaigns and special offers From time to time our industrial partner performs different campaigns – triggered by exogeneous reasons like new competing stores and the like – in its branches. Articles from a commodity group are all marked down at the same time or they sell three shirts for the price of two and so on. The question was how to deal with this? Prohibiting these measures would certainly make the result less biased. However, if our partner decided to implement DISPO this effects would occur anyway. DISPO has to deal with them. In order to better assess the practicability of DISPO, we decided not to forbid campaigns in the test and control branches: A method the performance of which vitally relies on laboratory conditions with all exogenous disturbances removed cannot be used in practice anyway. Stock transfers There are two types of stock transfers. On the one hand branches perform stock transfers caused by single customer requirements. If the requested product is not available in a branch the employees have the possibility to obtain it from an other branch. On the other hand there are systematical stock transfers which are caused by an undersupply of a size in a branch. Then the sales department decides to restock items from a oversupplied branch to a branch where there is a shortage. While the transfers caused by customer requirements in our partner’s view for the field study can not be forbidden the DISPO-team and our partner decided to prohibit all systematical stock transfers. Because they compensate wrong supply – maximizing the realized revenue by a size-compliant supply is the aim of our method – they would bias the result in such a way that it might be useless. 10.5.3 Evaluation Our test set of articles is denoted by A. For reasons of comparability we consider for each branch the objective value of ISPO divided per merchandise value over all articles from the set A. We distinguish the corresponding variables and parameters for the different articles a∈Aby a superscript a. Apart from that the parameters name are identical to the formulation of ISPO, Problem 6. For an initial stock Ia b,s for the considered branch b, size sand and a starting price πa 0we compute the relative realized objective for branch bof the independent nonanticipative decisions RROISPO(b) = objective of ISPO achieved for b maximal possible objective for b= −X a∈A X `∈LX m∈M xa b,`,m ·ca b,`,m − κ X i=1 ˜ δi·za i+X k∈K exp(−ρk)X a∈A X s∈S ˆra k,b,s −˜µa kˆna k −Pa∈A P`∈LPm∈Mxa b,`,m ·ca b,`,m −Pκ i=1 ˜ δi·za i+Pa∈A Ps∈SIa b,s ·πa 0 .(10.7) Depending on article a, the entity za iindicates that an i-th lot-type was used. During the sales process, we observed ˆra k,b,s (the realized yield for Branch band Size sin Period k) for Article aand ˆna k(mark-down in Period k– yes or no). Since we only consider a subset of branches we have to take into account that pick costs, costs for additional lot types, and fixed mark-down costs must be scaled with respect to the number of considered branches. This way, we get a marginal cost ˜ δi CHAPTER 10. REAL-WORLD EXPERIMENTS 139 -0,2 0 0,2 RRO_testRRO_control ordered test-control-pairs Figure 10.4: RROtest −RROcontrol for ordered test-control-pairs – all 81 articles for the i-th selected lot-type and mark-down costs ˜µa kfor period k. (For a complete notational reference see the problem formulation of ISPO in Chapter 6). For the evaluation of the field study we used transaction data from our industrial partner, analogous to the historical data we use for demand estimation. From this data we obtained the daily sales per branch and size and the realized sales price. But this data yields not all information we needed. We did not have full access to the supply in terms of lot-types in the control branches because the partner did not store the data for all articles. So we had to reconstruct the lot-types from the transaction data. In the most cases we can exactly determine the lot-type which is used in a branch. But we know from former field studies that not for all branches exactly the computed lot-type is also delivered. Sometimes there is one item less ore more supplied for a size. We therefore counted a lot-type as ”new“ lot-type if at least three branches were supplied with this lot-type. Otherwise we assumed that this lot-type resulted from an ”old“ lottype by wrong delivery. Then for this ”wrong“ lot-type only one time pick costs were counted. We want to emphasize that this kind of adaption – if it affects the results – is a disadvantage for the test branches and thus for DISPO. Because every branch has to be supplied by at least one lot-type pick cost always arise at least once. Allowing campaigns and special offers sometimes leads to the fact that the realized sales price differs from the determined price – either in the test branches or in the control branches. Moreover we cannot exclude that some items are bought at a reduced rate because of material defect or the like. There are two possibilities: Either we take the determined price per week for the computation of the RRO or the realized one. We decided for the realized price. We want to test the real-world behavior of DISPO and as already stated in 10.5.2 mark-downs beyond the determined prices are part of it. 10.5.4 Results of the field study The relative realized revenues per branch are shown in Table 10.4 for each test-control pair in the second and third column. We see that on average, the use of the new method gains almost two percentage points compared to the old method. We apply the Wilcoxon signed-rank test. The differences of the observations, here RROtest −RROcontrol – at the fourth column of Table 10.4 are ordered increasingly according to their absolute values (depicted in Figure 10.4). The signed ranks are CHAPTER 10. REAL-WORLD EXPERIMENTS 140 test-control-pair RROtest RROcontrol RROtest −RROcontrol signed rank 10.6333 0.6214 0.0119 2 20.6764 0.6080 0.0683 19 30.5919 0.6072 −0.0154 −5 40.6056 0.5898 0.0159 6 50.6637 0.5663 0.0974 26 60.6228 0.6031 0.0197 8 70.6377 0.6500 −0.0123 −3 80.5832 0.5845 −0.0013 −1 90.5968 0.5731 0.0237 11 10 0.5372 0.6276 −0.0904 −23 11 0.5651 0.5489 0.0163 7 12 0.5333 0.5904 −0.0571 −18 13 0.5782 0.5570 0.0212 9 14 0.6381 0.4940 0.1441 28 15 0.5054 0.5845 −0.0791 −21 16 0.5927 0.4993 0.0934 25 17 0.5872 0.4943 0.0929 24 18 0.6078 0.5691 0.0388 16 19 0.5762 0.6476 −0.0714−20 20 0.5682 0.5323 0.0359 14 21 0.5133 0.4250 0.0883 22 22 0.5272 0.5547 −0.0275 −12 23 0.4015 0.5942 −0.1926 −30 24 0.4628 0.4860 −0.0232 −10 25 0.5168 0.4646 0.0522 17 26 0.5843 0.4621 0.1222 27 27 0.5658 0.4137 0.1521 29 28 0.4989 0.4608 0.0380 15 29 0.5466 0.5607 −0.0141 −4 30 0.5593 0.5272 0.0320 13 ∅0.5692 0.5499 0.0193 5.7 Table 10.4: RROs for the test-control-pairs – all 81 articles stated in the fifth column of the table. At first glance we can see that the differences are not distributed equally. The test branches perform visibly better. This is also the result of the Wilcoxon signed-rank test. For the data we get a rank-sum of 318. The probability for getting an equal or higher rank-sum is P30(X≥318) ≈4.02%. Thus, with a probability of 4.02% for a better performance of the test branches resulted by chance we get for the 81 test articles a significant result for an improvement of DISPO – against the LDP with manual planning of mark-downs. However, we observed that some operational anomalies like failed price cuts in the control branches. In order to estimate the influence of the new method in the most conservative fashion, we removed all articles which may have been affected by systematic disturbances of operations. This led to a second set of articles A0with only 23 articles remaining. The particular RROs per branch are stated in Table 10.5, the corresponding differences RROtest −RROcontrol are depicted in Figure 10.5. We see that in the case of heavily cleaned-up data the RRO for the test branches is still more than 1.5percentage points higher than in the control branches. We repeated the Wilcoxon signed-rank test for this smaller test set. The test now yields a rank-sum of 271, which leads to a probability of P30(X≥271) = 22% that a better performance of the test-branches was observed by chance. Thus, for the heavily cleaned-up data we still observe a relevant effect (1.5 percentage points improvement) whose observation can no longer be testified as significant. This is essentially caused by the fact that for such a small (but relevant) effect the sample set A0is simply no longer large enough to prove significance. Still, the probability for a randomly better performance is with 22% much higher than the CHAPTER 10. REAL-WORLD EXPERIMENTS 141 test-control-pair RROtest RROcontrol RROtest −RROcontrol signed rank 10.4215 0.6673 −0.2458 −28 20.5874 0.4758 0.1116 18 30.6572 0.4865 0.1708 25 40.5948 0.4773 0.1175 21 50.5491 0.4153 0.1338 24 60.5799 0.5117 0.0682 13 70.4833 0.5454 −0.0621 −12 80.4648 0.5124 −0.0476 −9 90.5051 0.4923 0.0128 2 10 0.4933 0.6094 −0.1162 −19 11 0.4926 0.4998 −0.0071 −1 12 0.4205 0.4706 −0.0501 −10 13 0.4352 0.3746 0.0607 11 14 0.7046 0.2860 0.4186 30 15 0.4547 0.5281 −0.0734 −14 16 0.5146 0.3846 0.1300 22 17 0.5285 0.4247 0.1038 17 18 0.4802 0.5081 −0.0279 −3 19 0.3562 0.4865 −0.1303 −23 20 0.4119 0.4496 −0.0377 −5 21 0.2195 0.2577 −0.0382 −6 22 0.4274 0.5437 −0.1163 −20 23 0.2262 0.6415 −0.4153 −29 24 0.4006 0.3252 0.0754 16 25 0.3779 0.4244 −0.0465 −8 26 0.4759 0.4008 0.0750 15 27 0.5926 0.3971 0.1955 26 28 0.4458 0.4116 0.0342 4 29 0.4540 0.4985 −0.0445 −7 30 0.5278 0.3050 0.2228 27 ∅0.4761 0.4604 0.0157 2.57 Table 10.5: RROs for the test-control-pairs – heavily cleaned-up data, 23 articles sample relative realized objective gross yield sales test 0.4761 0.6829 0.7951 control 0.4604 0.6744 0.8021 Table 10.6: Alternative performance metrics, heavily cleaned-up data. probability for a randomly worse performance given by P30(X≤271) = 78.6%. So far, we assessed the quality of the decisions of the various methods on the basis of our objective function that was carefully engineered together with our partner. Yet, it is interesting to see that the new two-stage method outperforms the old method in some very important criteria at the same time. In Table 10.5.4 we list average RRO per branch, relative gross yields, and relative sales for all test-control-pairs. For both revenue and gross yield we see improvements. In contrast to this, the number of sales is only minimally smaller for DISPO. Now, which decisions have been taken differently by the new method? Table 10.7 shows the differences in the lot-type designs of the new and the old method for the 23 remaining articles.3The most obvious effect is that the number of different lottypes used is usually smaller for the ISPO than for the LDP. Since the old method tries to approximate a fractional demand as closely as possible by a supply distribution on 3Since the lot-type design of the control branches had to be reconstructed from incomplete data – see Subsection 10.5.3, the multiplicities for the control branches do not always add up to 30. The lot-types are reliable, though. CHAPTER 10. REAL-WORLD EXPERIMENTS 142 -0,2 0 0,2 RRO_testRRO_control ordered test-control-pairs -0,5 0 0,5 RRO_testRRO_control ordered test-control-pairs Figure 10.5: RROtest −RROcontrol for ordered test-control-pairs – heavily cleaned up data, 23 articles the basis of suitable lot-types, it will usually use as many lot-types as possible, even if the improvements of a new lot-type are small. The goal of the new method is not to meet the demand as closely as possible but to earn as much money as possible. Obviously, an additional lot-type is not always justified by higher predicted profits in ISPO. Consequently, ISPO does not suggest to use such a new lot-type. In the table we clearly see that lot-type (1,...,1) is very often used. This is the result of the rule that each branch has to receive at least one piece in every size – a fact that reduces the potential for improvement and should be taken into account when the effect (1.5 to 2 percentage points improvement) of using the new method is assessed.