Analytical Solution of Safety Stock Determination in Case of Uncertain Unimodal Lead-Time Demand

As companies state that a delivery service is important to their customers, an out-of-stock is considered harmful and therefore they keep safety stock in case of uncertain demand. For decision making on the level of safety stock a complete formulation of the distributional form of the demand during lead time is required. In practice, this information may not be available. In such a case, only partial information on the distribution might be available, such as the range, the mode, the mean or the variance. Given a value for a service performance measure, the decision maker, in this case, is not confronted with a single value for the safety stock but rather with an interval. The present research shows how upper and lower bounds of the safety stock are obtained in an analytical way, given a pre-specified service level using a service performance measure, called ‘expected number of units short’. The technique is also illustrated and compared within the framework of the research.


I. INTRODUCTION
In the business world, investment in inventories might be very high.Furthermore, companies are confronted with fluctuations in inventories in time and with uncertainties both in demand and supply, which directly influence decisions on inventories.Logistics managers use sophisticated systems for inventory management that enable them to take correct and timely decisions.Software implements decision models, which concentrate on the determination of inventory replenishment quantities based on relevant costs, such as order costs or storage costs.Most of these models assume deterministic demand patterns.But in real life uncertainties may appear both in demand and in supply or even in the quality of the delivered goods.The present research investigates in detail a specific case of inventory management decisions in the case of demand uncertainty.
The decisions are made on the basis of optimisation models taking a performance measure into consideration, which might be cost-oriented or service-oriented.Performance measures of the service-oriented type may be expressed relatively as a probability of a stock-out during a certain replenishment period, or may be expressed absolutely in terms of the number of units short, which is a direct indication for lost sales.These types of measures do not explicitly include a cost for additional paperwork or transport cost or loss of goodwill as, in practice, experience shows that these costs are hard to determine.Both in the cost-oriented and in the service-oriented performance evaluation, a specific integral plays an important role in the decision-making process.This integral is defined by Silver et al. [1, p. 258] as the expected shortage per replenishment cycle (ESPRC): in which it is assumed that the demand x in a replenishment lead time has a probability density function ) (x f and an order is placed at some time when the inventory position is at level t .If ordered per quantity Q , the fraction backordered is equal to Q ESPRC/ and a performance measure, indicated as 2 P , is defined as (see [1]): (2) In classical textbooks, the decision on the level of safety stock is based on the expected value and the standard deviation of the lead-time demand.Furthermore, it is assumed that the demand follows a Normal distribution, so that, given a service level and the knowledge about the first two moments, the required safety stock can be determined.However, normality of the distribution is sometimes violated.Assuming normality, when there is not true distribution behind the data, it might lead to either a degraded service level or an increase in the storage cost.
Regarding the violation of the normality assumption, the scientific literature is contradictory.Naddor [2] finds that false assumptions about the distribution may lead to higher cost in the case of extreme distributions, but, with realistic distributions, only the first and second moments are essential.On the other hand, Bartezzaghi et al. [3] show a significant impact of the shape of the demand distribution on the service level, based on a large set of experiments.Their analysis demonstrates that the shape of the distribution is a primary factor in the determination of inventories.Lau and Zaki [4] note that mean and variance are not sufficient for safety stock calculation, since skewness and kurtosis should also be accounted for.Furthermore, Käki et al. [5] show the impact of the demand distribution shape on replenishment, based on experiments with qualitative shape characteristics (normal, positively skewed, negatively skewed, and bimodal).In Janssens and Ramaekers' research [6], an approach has been developed to obtain the reorder point based on the knowledge of the range, mean and variance of the demand distribution only, which is the same information as required for the use of the normal distribution (as many times used in commercial software).
Given a mean and a standard deviation, and without knowledge of the real probability distribution behind the data, many distributions are candidate distributions.Each of them would lead to a different safety stock, given a preset service Information Technology and Management Science _______________________________________________________________________________________________ 2018/21 76 level.A conservative approach would be to determine the safety stock in the worst case distribution.In case only partial information on the probability density function is available (such as range, mean or variance), it is possible to determine the worst-case distribution.The mathematical development, however, shows that the worst-case distributions (in case mean and variance are known) are discrete distributions, in which the mass is concentrated in two points.This property has already been described by Scarf [7].
It might be argued that this approach is too conservative for practical purposes.Let us think of convincing a logistics manager to determine the safety of stock of a specific product based on the idea that the demand distribution has mass only in two points.A manager certainly will agree that some probability distribution, including uncertainty, has to be taken into account, but the distribution does not have extreme shapes.A manager can believe in asymmetry of the distribution, but it is rather uncommon that the manager will believe that the distribution is not unimodal.This means that the 'conservative' approach should be modified to a less conservative one, in which the 'worst case' distributions have a unimodal shape.
The present paper aims at investigating how such a distribution can be determined and which techniques are available to achieve it.It would be outstanding if a set of analytical formulas were available to determine the safety stock, given the knowledge of range, mean and variance of the demand during lead time, with the constraint that the resulting distribution should be unimodal.However, up to now, this simple solution is not available and other methods need to be explored.In the following section, solution methods first are explored in the case of knowledge of range, mean and variance but without constraints on unimodality.Afterwards, it is investigated which of these methods can be used with the constraint of unimodality.Illustrations, opportunities and limitations of the use of the methods are provided for the unimodal case.

II. SOLUTIONS FOR THE CASE WITHOUT CONSTRAINT ON UNIMODALITY
In this section, the incomplete information on the demand during lead time includes the finite range of the distribution, and the first and second moments.Let the size of the demand X for a specific product in a finite period have a distribution F with first two moments 1 = ( ) From a mathematical point of view, the problem is to find the following bounds: where φ is the class of all distribution functions F which have moments 1 µ and 2 µ , and which have support in + ℜ .
Further the following abbreviations are used: It has been shown that the upper bound corresponds to a 2atomic distribution, i.e., the probability mass is concentrated in two points (see Heijnen and Goovaerts [8] and Janssen, Haezendonck and Goovaerts [9]).Given a reorder point, the upper bound on the expected number of units short is given in Table I.I has all ingredients to determine the reorder point (and by this the safety stock), the decision problem is different from what the table offers.The table provides the upper bound on the expected number of units short given a value of t, but the decision problem is looking for a value of t, given a maximum value of the expected number of units short.There are several solutions proposed to solve the decision problem: -if an analytical equation exists, the equation for t should be solved in an analytical or numerical way (Method 1); -if an analytical solution exists, the equation should be solved by a one-dimensional search procedure (Method 2); -if no analytical procedure exists, a linear program should be constructed to find an approximative solution and use a one-dimensional search procedure (Method 3).
As can be noticed from Table I, the equations (listed in the second column called 'Upper bounds') are non-linear but not hard non-linear functions.Thus, for this illustrative example, the first of the three solution methods is used.This means, as the conditions (in the first column) cannot be checked in Information Technology and Management Science _______________________________________________________________________________________________ 2018/21 77 advance, all four equations need to be solved and the solutions need to be validated with the conditions.The conditions are mutually exclusive, so no confusion can exist.The solutions can be found by any non-linear equation solver, including spreadsheets like Microsoft Excel, using the Solver.Let us illustrate this way of working by the following example.The parameters of the problem are: a = 0, b = 50, 1 µ =30, 2 µ = 1200, so c = 25, 2 σ = 300 and 2 t µ σ = 325.Let us put the upper bound of the expected number of units short equal to 12.
The solutions of the four equations are listed in Table II.The function, describing the upper bound, is a decreasing function in t.That makes it easy to apply Method 2, given that t = 0 offers the expected value (30 in the example) and t = 50 offers zero.The procedure of cutting the relevant range into halves offers an efficient way of finding the solution.The pseudocode of the procedure is shown below.Target refers to the required value of the expected number of units short.Required_precision refers to the precision required on the value of t.Running_t refers to the running value of t as it changes during the procedure.Precision refers to the running value of the precision, calculated as the difference between the previous and current value of upper_bound.Upper_bound refers to the value of the upper bound as calculated from Table I.Interval refers to the size of the interval, which is used for adapting the value of t.

Parameters target, precision_required, a, b interval
The example is worked out in Table III.A number of steps in the procedural computations are shown.The procedure stops whenever the required precision is obtained.The third method will be illustrated in the unimodal case (in section III).

III. SOLUTIONS FOR THE CASE WITH A CONSTRAINT ON UNIMODALITY
Sometimes it is considered difficult to observe or to estimate both first and second moments, but experts have an opinion on the unimodality of a demand distribution, either with or without additional knowledge on the expected value.This idea has been common in expert opinion on the duration of activities in a project, where experts estimate three values: an optimistic, a most likely (or modal) and a pessimistic duration [10].
Let I be an interval on the real line ℜ .A fixed point I m ∈ is called the mode.A m-unimodal density function is a density function increasing (not necessarily strict) at the left of m and decreasing (not necessarily strict) at the right of m .Let us consider first the case in which the range, expected value and mode are known.The domain of the parameters is: ( + ).In this case, analytical upper bounds exist [11,Chapter 5].The bounds are shown in Table IV.The type of distributions that correspond to the extreme distributions leading to the upper bound are distributions, which consist of rectangular parts.The probability density in the left part (between 0 and the mode on the x-axis) is smaller than the density on the right part (between the mode and the upper bound of the range).In the example, the probability density is given by () = 0 if 0 ≤  ≤ 10 and () = 1/40 if 10 <  ≤ 50.
Method 2 is also applicable in the case of unimodal distributions with a known mode.In Table VI, a number of steps in the procedural computations are shown.The procedure is similar to the pseudocode as presented in Section III: the 'case' component has only two choices, i.e., related to the conditions mentioned in Table IV.In case, more constraints are involved, it is hardly possible to find analytical solutions and, therefore, to find the way to obtain the best value of safety stock or reorder point.In this situation, the authors of the present paper propose using a search procedure based on a linear programming approach.This corresponds to Method 3 mentioned in Section I.The method is elaborated and illlustrated for the case of a unimodal distribution with mode known but, additionally to the knowledge of the expected value, also the constraint of a known variance is included.A separate section is dedicated to this topic.In our application, the constraints are moment constraints, i.e., the first and second moment equalities and the obvious constraint because any member of φ is a probability distribution.which means that:

IV. LINEAR PROGRAMMING SOLUTION FOR THE CASE WITH
First, let us introduce the method for the case of the problem with two moments known.The extreme distributions for the supremum problem has been shown to be a two-point distribution.The integral (4) may be approximated by a sum making use of finite masses i p in a large number of points i x .
Its formulation looks like this: ( ) ., subject to

1, =
and However, any refinement in granularity (more values of xi) leads to an increased number of variables both in the objective function (5) and in the constraints ( 6)-( 8).This phenomenon does not guarantee any convergence towards the exact upper bound.
The optimisation problem (4) has a dual program of the type: Mostly the set J is infinite, so the number of linear constraints on y is infinite.The idea is to replace J with a large finite subset of J and then to solve the obtained linear program.( ) A numerical example demonstrates the use of bounds on the expected number of units short.The following information on demand during the replenishment period is known: the first moment 1 µ = 25, the second moment 2 µ = 725 and the range of demand is [0, 50].The upper bounds on the number of stockout units are presented in  In case the constraint of a unimodal distribution with a known mode m is introduced, the functions  ̂(. ) take a different form: With the same data, and additional values for the mode m, the problem (9) can be formulated.The upper bounds on the expected number of units short are presented in Table VIII.Note that the number of constraints is linear in terms of the granularity k.This keeps the solution effort of the linear programs tractable.From the tables, it can be seen that with a reasonable value of k a very good approximation can be obtained.The effort is higher than the one required for Methods 1 and 2. But for these methods an analytical solution of the integral, given the value of t, is required.

V. CONCLUSION
The research aims at proposing methods for obtaining a safety stock, given a customer service level, in a situation where the demand during lead time is not completely known.The incomplete information may include the range, first and second moments, unique mode, or tail probability.In some special cases, the expected number of units short can be expressed in an analytical way.In these cases, the incomplete information can be used to obtain an upper bound on the safety stock given a pre-described service level.Two methods are explained in the case of two moments known and in the case with a unique mode and the mean known.In most other cases, no analytical formulae exist on which our methods are based.For these situations, a linear programming approach is presented.The situation is more complex, but also leads to a good approximation of the required solution.

TABLE I UPPER
BOUNDS FOR THE EXPECTED NUMBER OF UNITS SHORT also applicable in this case.With the same data (but no variance known), let us add the additional information of the mode m = 10.The result of the application of Method 1 is shown in TableV.
Through a careful selection of the θ -values, convergence to the supremum is guaranteed as exemplified by Janssens & A Q : Table VII for different values of t and varying sizes k of the evaluation point set.The exact values for the upper bounds are 16.37931 for t = 10, 5.0 for t = 25 and 1.37931 for t = 40.

TABLE VIII SOLUTIONS
OF THE LP APPROXIMATION METHOD: UNIMODAL CASE