- Email: [email protected]
Sensitivity analysis in inventory models by means of ergodicity coefficients
Author’s Accepted Manuscript Sensitivity analysis in inventory models by means of ergodicity coefficients Boualem RABTA
www.elsevier.com/locate/ijpe
PII: DOI: Reference:
S0925-5273(17)30073-7 http://dx.doi.org/10.1016/j.ijpe.2017.03.014 PROECO6683
To appear in: Intern. Journal of Production Economics Received date: 23 September 2015 Revised date: 16 March 2017 Accepted date: 17 March 2017 Cite this article as: Boualem RABTA, Sensitivity analysis in inventory models by means of ergodicity coefficients, Intern. Journal of Production Economics, http://dx.doi.org/10.1016/j.ijpe.2017.03.014 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting galley proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Sensitivity analysis in inventory models by means of ergodicity coefficients Boualem RABTA ∗ Department of production management and logistics, Alpen–Adria–Universit¨ at Klagenfurt, Universitaetsstr. 65 - 67, 9020 Klagenfurt am W¨ orthersee - Austria. E-mail : [email protected]
Abstract We investigate the sensitivity of a class of (s, S) inventory models with respect to the perturbation of the demand distribution in terms of ergodicity coefficients. Ergodicity coefficients can be regarded as matrix norms that are useful to study qualitative as well as quantitative properties of some classes of stochastic systems. We obtain estimations of the absolute deviation of the stationary vector of the underlying Markov chain subject to perturbations. The particular structure of the transition matrix of this class of models allows us to derive simple closed form formulae for the computation of the ergodicity coefficients as well as the perturbation bounds. The perturbation bounds obtained can help us to decide whether or not the model remains an acceptable representation of the real system and thus decide whether or not it can be trusted for real life applications. Under perturbation, e.g. estimation errors or approximations, an optimal solution for the mathematical model may not be optimal for the real system and if implemented, the real performance measures may deteriorate considerably and deviate from the targeted values. Numerical examples are given to illustrate how the perturbation of the demand distribution may have a considerable impact on the optimal inventory policy and the performance
Preprint submitted to Elsevier
18 March 2017
measures in some cases. In a practical setting, the understanding of the sensitivity results for inventory models can help us to figure out which parameters have to be estimated with most attention and how to build models that are robust and as close as possible to the real systems they represent.
2
Key words: Inventory control, Sensitivity analysis, Markov chain, Ergodicity coefficient, Perturbation 1991 MSC: 60J10, 60J22, 90B05
1
Introduction
Many stochastic models including inventory models, are very complicated and analytical solutions are rarely available for them. In practical situations, the solution of such systems is obtained numerically and often after employing a series of approximations to obtain a mathematically tractable representation of the real system. Moreover, the model parameters are usually estimated from empirical data by means of statistical methods inducing further perturbations. Approximations must be used with care to avoid large deviations between the actual system and its mathematical representation (model). In this respect, the purpose of the sensitivity analysis is to determine whether or not the effect of these disturbances is limited, namely, a small perturbation of the model’s structure or parameters (inputs) generates only a small deviation of its characteristics (outputs), and to give an estimation of this effect. This is very important because only robust systems can be trusted and implemented for decision-making. A sensitive model can react unexpectedly to the perturbation of its input parameters and therefore cannot be considered reliable. The optimal solution of the mathematical model might only be sub–optimal for the real ∗ Department of production management and logistics, Alpen–Adria–Universit¨at Klagenfurt, Universitaetsstr. 65 - 67, 9020 Klagenfurt am W¨orthersee - Austria Email address: [email protected] (Boualem RABTA).
3
system. In many cases, this can result in over–investments in inventory and in unacceptable service levels [20]. Sensitivity analysis can help decision-makers to determine if the performance measures remain within acceptable bounds and to understand the importance of the various parameters in this context. If we obtain an affirmative response to the qualitative question (robustness), then we might be able to estimate the deviation of the characteristics of the system caused by the perturbations (quantitative question). In general, multiple parameters in the inventory system are unknown and must be approximated or estimated from empirical data. Customers’ demand, replenishment leadtimes and system costs (ordering, holding and shortage costs) are among the most important elements that are subject to perturbations. The classical Economic Order quantity model (EOQ) is well–known to be somewhat insensitive to the change in order quantity and against parameter estimates (see, e.g. [15]). However, it can be quite sensitive to forecast errors whenever the leadtime is non–zero [22]. Extensions of the classical EOQ model have been presented in numerous works and many of them also discuss the sensitivity analysis of the optimal order quantity (see, e.g. [5,7] and references therein). Wang [38] applied sensitivity analysis to the changes of optimal advertising expenditure and ordering decisions in a newsvendor model under advertising-dependent demand distribution. Khanra [19] reviewed the literature on sensitivity analysis of newsvendor inventory models and studied the sensitivity of the classical one. The author showed that it is more sensitive than the EOQ model to sub–optimal ordering decisions and identified conditions for symmetry/skewness of cost deviation found to be closely linked with symmetry/skewness of the demand density function. Inventory models were the first stochastic models for which monotonicity prop4
erties were proved [35]. Boylan [8] has proved that the solution of the optimal inventory equation depends continuously on its parameters. Sinha and Gupta [34] studied the effects of alternative forms of surplus and shortage cost functions and also investigated the effects of alternative forms of probability density function of demand in an inventory control situation. In [16], Glasserman and Tayur studied the sensitivity of base-stock levels in multi– echelon production–inventory systems, in particular with respect to demand parameters. Zheng [40] investigated the sensitivity of a standard single-item continuous-review (s, Q) inventory system with fixed leadtimes and found it to be more robust than its deterministic counterpart, the EOQ model. Chen and Zheng [12] established that the average cost in the (s, S) inventory model with exponentially distributed demand is quite insensitive to the changes in D = S − s. In [24], Rabta and A¨ıssani investigated the strong stability of an (s, S) inventory model with instantaneous deliveries. After proving the robustness of the underlying Markov chain with respect to the perturbation of the demand distribution, they obtained an upper bound for the deviation of the stationary distribution of the underlying Markov chain [24].
In many cases, sensitivity analysis is performed numerically (see, e.g. references in [14,39]). However, the numerical approach has proved to be insufficient. Chu and Chung [14] derived sensitivity formulas for an inventory model with backorders and showed that the numerical analysis might not lead to consistent sensitivity results. Hence, an analytical approach is required in order to obtain valid conclusions regarding the qualitative properties of the system. Local sensitivity analysis methods include the study of one parameter at a time and testing the effect when changing the value within a specific range (e.g. [27,26]) along with perturbation analysis (e.g. [4,3]) among others (see, 5
[36]). On the other hand, global methods test the sensitivity of the model in consideration of the uncertainty distribution reflecting the decision-maker state of knowledge in the parameters (see, e.g. [6] and the references therein). Borgonovo and Peccati [6] employed a global sensitivity analysis method based on Sobol functions and variance decomposition method to determine the most influential parameters on the model output. Borgonovo [5] introduced a new method to study of the effect of changes in one or more of the exogenous variables, defining sensitivity measurements that do not rest on differentiability and applied it to modified EOQ models. Sensitivity analysis can also be performed numerically (e.g. [14,39]) or by simulation, where results that can easily be analyzed statistically are generated (see, e.g. [10]). General techniques of stochastic analysis were also applied to stochastic inventory models. For instance, Stoyan [35] investigated their monotonicity properties while Rabta and A¨ıssani [24] analysed an (s, S) inventory model based on the results of the theory of the strong stability of Markov chains [2].
A large class of infinite horizon stochastic inventory systems can be represented by Markov chains (e.g. [9,24]). In this case, the focus is generally on the steady state of the system and the characteristics are expressed in terms of the stationary distribution. There is an abundant literature on the robustness of Markov chains with tools and results that can be applied to inventory models. For finite Markov chains, the problem of sensitivity to perturbations has been treated by many authors since Schweitzer (1968) [29]. Meyer singly and with co-authors used the group inverse notion [11] to derive several bounds (see, e.g. [21]). Hunter [17] used a wider class of generalized inverses. The coefficient of ergodicity has been used by Seneta [31,32]. Also, Cho and Meyer [13] obtained a perturbation bound expressed by means of mean first passage 6
times. The validity of most of these bounds extends easily to infinite discrete Markov chains [25] as well as Markov chains with general state space [2].
Among those results, only few have been applied to inventory models. The bound derived in [24] requires the computation of a part of the stationary distribution which is expressed in terms of the renewal function of the demand distribution.
In this paper, we investigate the sensitivity of an (s, S) inventory model as represented by Markov chains and subject to perturbations of the demand distribution. We use the so-called ergodicity coefficients with the aim of deriving bounds on the deviation of the stationary characteristics of the model that could be easily computed. The bound proposed in this paper, based on the ergodicity coefficient, does not involve the stationary distribution and is simple and easy to compute.
The generic (s, S) inventory model analyzed in the sequel assumes zero leadtime but the obtained results extend to a large set of similar models with only minor changes in the assumptions (see, remarks below). The particular structure of the transition matrices of this kind of model allows the computation of the corresponding ergodicity coefficient in a very simple and computationally efficient manner. Furthermore, we show that the obtained closed form formulas reduce to very simple forms for particular demand distributions. The obtained bounds are given in terms of the norm k.k1 of the stationary vector of the underlying Markov chain without the need to explicitly calculate this stationary distribution. We provide numerical examples to show the optimal inventory policy is affected by the change to the demand distribution and how this can affect the performance measures such as total cost and service level. 7
2
Preliminary and notations
In this article, we use the norm k.k1 for row vectors and the norm k.k∞ for column vectors and for matrices. Also, all vectors are column vectors (row vectors are transposed). Since it will always be clear if the norm corresponds to a matrix, a row vector or a column vector, all norms are noted k.k. Then, for a vector v = (vj )j∈E we have : kvk = kvk∞ = max |vi |, i∈E
kv T k = kv T k1 =
X
|vj |.
j∈E
For a matrix B = (bij )i,j∈E , we have : kBk = kBk∞ = max i∈E
X
|bij |.
j∈E
The following result (Paz’s inequality), given in [23, Ch. IIa] will be needed in the sequel : Lemma 1 Let c 6= 0 be a real vector such that cT e = 0. Then, for any complex vector d we have : |cT d| ≤ kcT k
maxi,j |di − dj | . 2
Here and in the sequel, eT = (1, 1, .., 1) is the vector of all ones. The ergodicity coefficient of a matrix B = (bij )i,j∈E with equal row sums is defined as ([30]): τ (B) =
maxT
uT e=0,
ku k<1
kuT Bk.
Seneta [30] gives an explicit expression of this quantity : τ (B) =
X 1 max |bik − bjk |. 2 i,j k
(1)
8
The ergodicity coefficient of a finite stochastic matrix P (with respect to the norm k.k1 ) satisfies the following properties (see, e.g. [33]) : • 0 ≤ τ (P ) ≤ 1. • τ (P ) = 0 if and only if P has rank 1 (i.e., P = ev T , for some probability vector v). • For every eigenvalue λ 6= 1 of P we have : |λ| ≤ τ (P ). Additionally, τ (P n ) ≤ (τ (P ))n . Consider a homogeneous irreducible Markov chain X with states in a finite space E. Denote by P the transition matrix of X and by π its stationary vector, i.e., the unique real vector satisfying π T P = π T , π T e = 1.
Suppose that P is perturbed to Q = P + ∆ such that Q corresponds to the transition matrix of another irreducible Markov chain Y with states in E. Denote by ν the stationary vector of the Markov chain Y . For finite Markov chains, Seneta [31] gives the following result. Theorem 2 For a finite irreducible Markov chain X satisfying τ (P ) < 1 we have : kν T − π T k ≤
kQ − P k . 1 − τ (P )
However, we observe that for two arbitrary probability vectors, the norm of their difference cannot exceed 2 since, kν T − π T k ≤ kν T k + kπ T k = 2. 9
Observe that for any two probability vectors π and ν, the vector (ν−π) satisfies (ν − π)T e = 0. Then, by taking d = ek (where ek is the k−th column of the identity matrix) and applying Paz’s inequality we obtain for every k ∈ E : 1 |νk − πk | ≤ kν T − π T k. 2 Hence, 1 max |νk − πk | ≤ kν T − π T k. k∈E 2
3
Perturbation bounds for the (s, S) inventory model
Consider the following single-item, single-echelon, inventory problem. We manage a non–perishable product in the face of a stationary stochastic demand. The inventory level is inspected every R time units. At the beginning of each period, the manager decides whether or not he should make a replenishment order and if yes, how much to order. Suppose that the outside supplier is perfectly reliable and that orders arrive immediately. We assume that the consecutive period demands are independent and identically distributed (i.i.d.) random variables. Remark 3 The period length parameter R does not play a direct role in the subsequent analysis. In fact, R might be random as in the case of the continuous review model (s, S) (e.g. time between consecutive arrivals of customers) or fixed as in the case of the periodic review model. The structure of the transition matrix given below remains the same and hence the results are valid for both cases. 10
For both continuous and periodic review models, it has been proved that the (s, S) inventory policy is optimal (see [28,18,37]). According to this control rule, the level is reviewed at dates tn = R1 + R2 + .. + Rn , n ≥ 0 where Ri can be deterministic (periodic review, Ri = R∀i) or stochastic (continuous review, Ri are i.i.d. random variables). If at the review moment, the inventory level Xn is below or equal to s (Order-Point) we order sufficient quantity to raise the inventory level to S (Order-Up-To level). Thus, the size of the order is equal to Zn = S − Xn . Many algorithms have been proposed for the computation of the optimal values of the parameters of this class of inventory models (see, for example [37]).
During period n, n ≥ 1, the total demand is a random variable ξn . Assume that ξn , n ≥ 1, are independent and identically distributed (i.i.d.) random variables with common probabilities given by
ak = P (ξ1 = k), k = 0, 1, 2, ...
So, ak is the probability to have a demand of k items during the period.
For the purpose of our analysis, it suffices to define the demand distribution by the finite sequence ak for k = 0, 1, .., S − 1 representing the first S demand probabilities. This observation allows us to deal with only finite sums while maintaining full accuracy and precision in the derived formulas. Therefore, we denote by αz the quantity:
αz = P (ξ1 ≥ z) = 1 −
z−1 X
ak
k=0
and assume that ak > 0, k = 0, .., S − 1 and αS > 0. 11
Additionally, we exclude the special case s = S − 1 (base stock policy) since for this model the ergodicity coefficient will be 0 and the Inventory Markov chain is already in its permanent regime where all rows of P are identical and the stationary vector is given by a closed form (any row of P ). Therefore, the deviation of the stationary vector can be computed exactly and the estimation of perturbation bounds is not necessary.
Remark 4 Considering the assumption of zero-leadtime, the model includes both lost sales and backorders cases (or a mixture policy). If the demand during a given period is greater than the available quantity, then the excess quantity is either backordered or the sales are lost. In both cases, we order sufficient quantity to raise the inventory level to the order-up-to point S and eventually satisfy the backordered demand. The transition matrix of the on-hand inventory process will not be affected.
Assume that the initial inventory level X0 > s (if not, we raise it to S by ordering sufficient quantity). The inventory level Xn+1 at the end of period n + 1 is given by : + (S − ξn+1 )
if Xn ≤ s,
Xn+1 =
+ (Xn − ξn+1 ) if Xn > s,
where (A)+ = max(A, 0).
The random variable Xn+1 depends only on Xn and ξn+1 , where ξn+1 is independent of n and the state of the system. Thus, X = {Xn , n ≥ 0} is a homogeneous Markov chain with values in E = {0, 1, ..., S} with transition 12
matrix P given by :
0
1
s
s+1
S
0
P∞
ak aS−1 · · · aS−s aS−s−1 · · · a0
1
P∞
ak aS−1 · · · aS−s aS−s−1 · · · a0
S
S
.. .
.. .
..
.
.. .
.. .
..
. . ..
P = P∞
s s+1
S
P∞
s+1
.. . S
ak aS−1 · · · aS−s aS−s−1 · · · a0
P∞ S
ak as · · · a1 .. .
...
.. .
a0
0
0
.. .
...
0
ak aS−1 · · · aS−s aS−s−1 · · · a0
The Markov chain X is finite and irreducible. It admits a unique stationary vector π.
3.1
The coefficient of ergodicity of P
(1) allows us to easily compute the coefficient of ergodicity of a finite Markov chain in general. Yet, the particular structure of the transition matrices of inventory models allows further simplifications. Observe that rows 0 to s of the transition matrix P are all identical and they are same as the S-th row. Therefore we might write: τ (P ) =
X 1 max |pik − pjk |. 2 s+1≤i
To compute |pik − pjk | we distinguish the following cases: 13
• k = 0, |pik − pjk | = |
P∞
l=i
al −
P∞
l=j
al | =
Pj−1 l=i
al .
• 1 ≤ k ≤ i, |pik − pjk | = |ai−k − aj−k |. • i + 1 ≤ k ≤ j, |pik − pjk | = aj−k . • j + 1 ≤ k ≤ S (for j < S), |pik − pjk | = 0. Thus,
j−1 j i X X X 1 τ (P ) = max al + |ai−k − aj−k | + aj−k 2 s+1≤i
Hence, Lemma 5 The coefficient of ergodicity of the transition matrix P is given by:
j−1 j−i−1 i−1 X X X 1 |ak − aj−i+k | + ak + ak . τ (P ) = max 2 s+1≤i
(2)
The coefficient of ergodicity of the transition matrix P could therefore be calculated exactly by means of a finite loop. Furthermore, we have : Lemma 6 The coefficient of ergodicity τ (P ) of the transition matrix P of the underlying Markov chain in the (s, S) inventory model satisfies : τ (P ) ≤ 1 − αS < 1.
PROOF.
j−1 j−i−1 i−1 X X X 1 τ (P ) = max |ak − aj−i+k | + ak + ak 2 s+1≤i
j−1 j−i−1 i−1 X X X 1 ≤ max (ak + aj−i+k ) + ak + ak 2 s+1≤i
j−1 j−1 j−i−1 X X X 1 ak + ak + ak ≤ max 2 s+1≤i
14
j−1 j−1 j−1 S−1 X X X X 1 = max ak + ak = max ak = ak s+1≤i
τ (P ) ≤ 1 − αS < 1.
3.1.1
Particular demand distributions
The expression of the ergodicity coefficient of the transition matrix P derived above could be further simplified for particular demand distributions. For instance, in the result below we impose conditions on the first S elements of the demand distribution. Lemma 7 If the sequence ak , k = 0, 1, .., S − 1 is non–increasing then: τ (P ) = 1 − αS−s−1 . Similarly, if the sequence ak , k = 0, 1, .., S − 1 is non–decreasing then: τ (P ) = αs+1 − αS .
PROOF. Suppose that ak , k = 0, 1, .., S − 1 is non–increasing then:
τ (P ) =
j−1 X
i X
j X
i X
1 max al + ai−k − aj−k + aj−k 2 s+1≤i
j−1 j−1 j−i−1 i−1 X X X X 1 = max al + ak − ak + ak 2 s+1≤i
j−1 j−1 j−1 X X X 1 = max ak − 2 ak + ak 2 s+1≤i
=
max
s+1≤i
j−1 X
ak −
k=0
j−1 X
ak
k=j−i
Thus, τ (P ) =
max
j−i−1 X
s+1≤i
k=0
15
ak =
S−s−2 X k=0
ak
Now, suppose that ak , k = 0, 1, .., S − 1 is non–decreasing then:
j−1 j i i X X X X 1 τ (P ) = max al − ai−k + aj−k + aj−k 2 s+1≤i
=
j−1 X
j−1 X
i−1 X
1 max al − ak + ak + 2 s+1≤i
j−1 X
j−1 X
j−1 X
j−i−1 X
ak
k=0
j−1 X 1 max 2 al − al + ak = max al = s+1≤i
Hence, τ (P ) =
S−1 X
al =
l=s+1
∞ X
al −
l=s+1
∞ X
al
l=S
The first condition of the previous lemma is satisfied by the following two discrete probability distributions.
3.1.1.1
Uniformly distributed demand In this case, the total demand
follows a (discrete) uniform distribution on the set {0, 1, .., M } with M ≥ S. ak =
1 , k = 0, 1, .., M. M +1
The coefficient of ergodicity τ (P ) of the transition matrix P writes: τ (P ) =
3.1.1.2
(S − s − 1) M +1
Geometrically distributed demand In this case, the demand
distribution is given by : ak = q k (1 − q), k = 0, 1, ... for 0 < q < 1. The coefficient of ergodicity τ (P ) of the transition matrix P writes: 16
τ (P ) = 1 − q S−s−1
3.2
Perturbation of the demand distribution
In this section, we replace the demand distribution in the previously described inventory model by another distribution close to the original one in a certain sense (we will define a measure of proximity (distance) of the two distributions). This case is encountered in many circumstances : (1) The demand distribution is fitted from empirical data. The perturbations result mainly from the estimation errors. (2) The demand distribution is replaced by another distribution for purposes such as simplification of the model and ease of the analysis. For instance, the set of phase-type distributions is dense in the set of positive continuous probability distributions [1]. Hence, a positive continuous probability distribution could be approximated (replaced) by a phase-type distribution at any desired precision level. (3) The form of the perturbed distribution is the same as the original one but differs in the value of one or more parameters. The difference comes for instance from errors of the estimation of the parameter by means of statistical methods. Consider another inventory model with the same structure but with period demands ξi0 , i ≥ 1 independent and identically distributed with common probabilities a0k = P (ξ10 = k), k = 0, 1, 2, ... 17
Denote by Y the inventory level Markov chain in the new system and by Q its transition matrix. First, we estimate the effect of the perturbation on the transition matrix of the original Markov chain. Namely, we compute a measure of the distance between the transition matrices of the two chains. Lemma 8 The deviation of the transition matrix P is given by : kQ − P k = |αS0 − αS | +
S−1 X
|a0k − ak |.
k=0
PROOF. kQ − P k = max i∈E
X
|qij − pij | = max(A, B)
j∈E
where, A = max
0≤i≤s
∞ X
=|
∞ X
a0k −
k=S
ak | +
k=S
S X
S X
|qij − pij |
j=0
|a0S−j − aS−j | = |αS0 − αS | +
j=1
S−1 X
|a0k − ak |
k=0
and B = max
s+1≤i≤S
"
≤ max
s+1≤i≤S
|
S X
|qij − pij | = max | s+1≤i≤S
j=0
∞ X k=S
a0k
−
∞ X k=S
ak | +
S−1 X
|a0k
∞ X
a0k −
k=i
− ak | +
k=i
∞ X k=i
i−1 X k=0
ak | +
i X
|a0i−j − ai−j |
j=1
#
|a0k
− ak | =
S−1 X 0 |αS −αS |+ |a0k −ak |. k=0
Similarly to the formula of the coefficient of ergodicity, the deviation of the transition matrix P subject to perturbation is given by a finite sum. Moreover, the computation involves only basic arithmetic and therefore the total computation time is linear with respect to the dimension of the matrix P . It 18
is not necessary to store the transition matrix or the complete demand distribution (only, the first S individual probabilities are required). Hence, the computations could be performed efficiently and scale well for large models. Using Theorem 2 we obtain the following result : Theorem 9 Let π and ν be the stationary distributions of chains X and Y respectively. Then : 0 |α0 − αS | + S−1 i=0 |ai − ai | . kν − π k ≤ S 1 − τ (P ) T
P
T
(3)
where τ (P ) is given by (2).
PROOF. This result holds by virtue of Theorem 2 and lemma 8.
3.3
Deviation of the mean inventory level
The previous perturbation bounds can be used to estimate the deviation of several quantities of interest. To illustrate this we compute the perturbation bounds for the average end period inventory level: ¯ (resp. Y¯ ) be the mean end of period inventory level in Theorem 10 Let X the base (resp. perturbed) system. Then : ¯ ≤ S kν T − π T k. |Y¯ − X| 2
PROOF. We have : ¯= X
S X j=0
19
jπj ,
Y¯ =
S X
jνj .
j=0
Then ¯ =| |Y¯ − X|
S X
jνj −
j=0
S X
jπj | = |(ν T − π T )f |,
j=0
where, f = (0, 1, .., S)T . Thus, by using Lemma 1 we get :
¯ ≤ |Y¯ − X|
3.4
S T kν − π T k. 2
Deviation of the inventory costs and service level
The individual components πi of the stationary vector π of the chain X give the probability that, at the equilibrium state, the on-hand inventory level is equal to i at the end of any period. The total cost is usually the sum of three components, namely, ordering (fixed and variable), holding and penalty costs. The latter being in general difficult to estimate, we might alternatively consider a service level in addition to an inventory cost function that is the sum of ordering and holding costs only.
Let’s assume lost-sales and the following cost structure : A fixed ordering cost K is incurred whenever a replenishment order is placed and is added to the variable acquisition cost of c per ordered unit of the product. The unit holding cost is h per unit of time. In addition to the replenishment and holding costs (C1 ) we consider a service level (ψ1 ). We illustrate the use of the above inequalities in the study of the sensitivity of the key performance measures of the considered inventory model. 20
3.4.1
Replenishment and holding cost
The overall cost function might have a different structure depending on the assumptions of the model. For the purpose of the following illustration we restrict the analysis to the case where R is constant and holding costs are computed from the average end period inventory level. Under the previous assumptions, we write the average total (period) replenishment and holding cost of the model :
C1 (π) =
S X
hiRπi +
i=1
s X
(K + c(S − i))πi = π T ζ,
i=0
where ζ is the vector for which components are given by K + c(S − i) + hiR if 0 ≤ i ≤ s
ζi =
hiR
if s < i ≤ S
Similarly, the average total replenishment and holding cost in the perturbed system is C1 (ν) = ν T ζ and by applying Paz’s inequality, we obtain |C1 (ν) − C1 (π)| = |(ν T − π T )ζ| ≤ kν T − π T k
maxi ζi − mini ζi . 2
Hence, |C1 (ν) − C1 (π)| ≤
kQ − P k (maxi ζi − mini ζi ) . 2(1 − τ (P ))
If we assume that c > hR (the unit purchase cost is greater than the unit period holding cost) then
max ζi = max {K + cS, K + cS + hR − c, hSR} = K + cS, i
21
and min ζi = min {K + cS, K + c(S − s) + hsR, h(s + 1)R} = h(s + 1)R i
Hence, Theorem 11 Under the previous assumptions on the cost function, we have :
|C1 (ν) − C1 (π)| ≤
3.4.2
|αS0 − αS | +
PS−1 i=0
|a0i − ai | (K + cS − h(s + 1)R) 2(1 − τ (P ))
.
Service level
We consider a service level ψ1 defined as follows: ψ1 = P {Inventory at the beginning of the period ≥ Period demand}. At the equilibrium, the inventory at the end of one period is equal to i with probability πi . Therefore, at the beginning of the next period the on hand inventory will be equal to S if i ≤ s and equal to i otherwise. Hence, ψ1 (π) =
s X
πi P (ξ1 ≤ S) +
s X i=0
πi (1 − αS+1 ) +
S X
πi P (ξ1 ≤ i)
i=s+1
i=0
=
S X
πi (1 − αi+1 ) = 1 −
s X i=0
i=s+1
πi αS+1 −
S X
πi αi+1
i=s+1
Thus s S s S X X X X 0 0 πi αS+1 − πi αi+1 − 1 − νi αS+1 − νi αi+1 |ψ1 (π) − ψ1 (ν)| = 1 − i=0 i=s+1 i=0 i=s+1 s S s S X X X X 0 0 = πi αS+1 + πi αi+1 − νi αS+1 + νi αi+1 i=0 i=s+1 i=0 i=s+1 X S s S X X X s = πi αS+1 + πi αi+1 − νi αS+1 + νi αi+1 i=0 i=s+1 i=0 i=s+1 s S s S X X X X 0 0 + νi αS+1 + νi αi+1 − νi αS+1 + νi αi+1 i=0 i=s+1 i=0 i=s+1
22
X S s S X X X s = πi αS+1 + πi αi+1 − νi αS+1 − νi αi+1 i=0 i=s+1 i=0 i=s+1 s X − νi α 0
S+1
+
i=0
S X
0 νi αi+1 −
i=s+1
s X
νi αS+1
i=0
− νi αi+1 i=s+1 S X
s S X X = (πi − νi )αS+1 − (πi − νi )αi+1 i=0 i=s+1 s X − νi (α0
S+1
S X
− αS+1 ) −
i=0
0 νi (αi+1
i=s+1
− αi+1 ) .
Whence,
|ψ1 (π) − ψ1 (ν)| = (π T − ν T )η − ν T (η 0 − η) , 0 ), where η (resp. η 0 ) is the vector having component i equals to αS+1 (resp. αS+1 0 ), if s + 1 ≤ i ≤ S. Consequently, if 0 ≤ i ≤ s and to αi+1 (resp. αi+1
|ψ1 (π) − ψ1 (ν)| ≤ (π T − ν T )η + ν T kη 0 − ηk .
Applying Lemma 1, we obtain :
|ψ1 (π) − ψ1 (ν)| ≤ kπ − νk
maxi,j |ηi − ηj |
T
0 + ν kη − ηk , 2
0 where maxi,j |ηi − ηj | = αs+1 −αS+1 , ν T = 1 and kη 0 − ηk = maxi=s+1,..,S |αi+1 −
αi+1 |. Finally : Theorem 12 Under the previous assumptions, the deviation of the service level ψ1 satisfies: 0 |α0 − αS | + S−1 i=0 |ai − ai | αs+1 − αS+1 0 −αi+1 |. |ψ1 (π) − ψ1 (ν)| ≤ S + max |αi+1 i=s+1,..,S 1 − τ (P ) 2
P
23
4
Effect of the perturbations on the optimal inventory policy: Numerical examples
Consider the following base scenario: under the model’s assumptions above, we assume that the inventory is reviewed periodically with period R = 7. Customers arrive according to a stationary Poisson distributed process with rate λ = 6 (per unit of time) and all customers request exactly one unit of the product. The cost function elements are : Fixed ordering cost K = 25 per order, variable ordering cost c = 3 per unit ordered and holding cost h = 0.1 per unit in stock at the end of the period and per unit of time. Additionally, we either consider a shortage cost p = 5 per unit of unsatisfied demand or impose a service level constraint (ψ1 = 0.95). In the following examples, we illustrate the effect of the perturbations on the optimal inventory policy as well as the performance measures of the system. We implemented a linear search algorithm that determines the optimal inventory policy in the following two cases:
• Minimization of the sum of ordering and holding costs under service level constraint. • Minimization of the total cost (sum of ordering, holding and shortage costs).
Obviously, the solution might be different in the two cases and we expect a different impact on the performance measures. In the first case, the performance of the inventory system is more affected in terms of costs than in service level which is limited by the constraint.
In each of the above mentioned cases, we introduce a small perturbation on 24
the demand distribution by varying the customers’ arrival rate. Two issues are investigated. We first show the impact of the perturbation of the demand distribution on the optimal parameters s∗ andS ∗ of the inventory policy and we compare the performance measures of the original and perturbed system (optimal policy in the model vs. optimal policy in the real system). However, the situation is different in practice. The decision is made based on the solution of the mathematical model and the obtained optimal solution is implemented in the real system. Since the parameters of the real system are different from those of the mathematical model, the implemented policy is not optimal anymore and another error component is added to the approximation error. We illustrate this situation by comparing the outcome of the same policy in both the original and perturbed model. Therefore, in this case we compare the optimal solution of the mathematical model with the sub–optimal policy in the real system. In each (sub–)case, we proceed as follows: • Calculate the optimal policy, total cost and service level under the base scenario’s assumptions. • Perturb the demand distribution and calculate key performance measures in the perturbed model. • Compare the original and the perturbed models in terms of performance measures. • Repeat these steps for different perturbation amplitudes. To perturb the demand distribution, we act on the customers’ arrival rate and vary it by a small amount each time. The perturbation amplitude will range between -10% and +10% of the value of the arrival rate. This is sufficient to 25
show the impact of the perturbation while remaining in a realistic situation. In each, we report the values of the optimal inventory policy parameters s∗ , S ∗ , the total cost C ∗ and the service level ψ1∗ for both the original and perturbed ∗
models as well as the relative deviation of these performance measures ( ∆C C∗ and
4.1
∆ψ1∗ ψ1∗
in percentage). Discussion on the outcome of the experiments follows.
Example 1 : Sensitivity of the optimal inventory policy under minimization of ordering and holding costs with service level constraint
In practice, shortage costs are very difficult to estimate. Along with the profit loss, shortages can lead to shutdowns of the production or influence customer satisfaction and retention which may have a negative impact on future demand and revenue. Therefore, in many cases a target service level is fixed and the system is optimized under this constraint rather than with respect to the total cost which includes an estimated (imprecise) shortage cost component. Let us set the service level to 95%. The algorithm will search for the optimal inventory policy by minimizing the sum of ordering and holding costs under the service level constraint. The results are reported in Figure 1. Figure 1 goes here. Fig. 1. Deviation of the optimal policy parameters with respect to the perturbation of the arrival rate.
The effect of the perturbation in this case remains somewhat limited. For the maximum amplitude (10%) of the perturbation of the customers’ arrival rate, we observe only 8.5% deviation of the sum of ordering and holding costs while the service level is kept near to the target level by the constraint. The change in the cost follows naturally the change in the optimal inventory policy 26
parameters s∗ and S ∗ . Next, we inject the solution of the mathematical model into the perturbed system and compare the performance measures. Results are represented in Figure 2. Figure 2 goes here. Fig. 2. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
The effect of the perturbation on the cost remains limited in this case as well but we observe a considerable change in the service level. A 10% increase in the customers’ arrival rate induced a drop of the service level from 0.95 to 0.86 (10%). In many real situations, the impact of the deterioration of the service level by this magnitude can be considerable.
4.2
Example 2 : Sensitivity of the optimal inventory policy under total cost minimization
In this example, we aim to show how the perturbation of the demand distribution (say, estimation errors) can affect the optimal values of the (s, S) inventory policy when the latter is obtained by minimization of the total cost function, i.e., the sum of ordering, holding and shortage costs. The results of the procedure are shown in Figures 3 and 4. Figure 3 goes here. Fig. 3. Deviation of the optimal inventory parameters with respect to the perturbation of the arrival rate.
The perturbation of the demand distribution has a direct impact on the op27
Figure 4 goes here. Fig. 4. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
timal inventory policy (optimal s and S). In the example above, the optimal inventory policy in the base scenario is (25, 43). When varying the arrival rate between -10% and +10% of its nominal value, the optimal inventory policy changes between (21, 39) to (29, 48). Naturally, this will lead to the deviation of the performance measures of the (s, S) inventory such as the total inventory cost and the service level. For example, a 10% error in the customers’ arrival rate generates around 8.2% deviation of the total inventory cost and up to 20.7% on the service level. In the following example, we implement the optimal inventory policy obtained in the base scenario in the perturbed model and report the results in Figure 5. Figure 5 goes here. Fig. 5. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
Notice that in this case, the deviation in the stationary performance measures might be higher than that in previous case. For example, a 10% error in the customers’ arrival rate can lead to a deviation of the service level as high as 46.4% (compared to 20.7% previously). Observe as well that the service level has dropped to 47% which means a high level of unsatisfied demand that can lead to serious consequences (e.g. loss of customers, production shutdown,...). Here, the deviation resulting from the approximation of the demand distribution is added to the deviation resulting from the use of a non–optimal policy in the real situation. 28
Finally, observe that without constraint, the service level performance measure might be at a low level even with optimal inventory policy. This could be explained by the fact that the shortage costs are poorly estimated which in turn increases the tolerance of the system to shortages allowing more unsatisfied demand. In the first case, i.e., minimization of the cost under service level constraint, the (s, S) inventory system looks more robust. If one can ensure that the perturbation is kept small then the constraint will help to keep the service level near to the target value. It seems also that costs are least affected in this case. The examples above show that perturbations of the demand distribution can affect the output of the (s, S) inventory model in various ways. In these examples, the optimal inventory policy (s∗ , S ∗ ) is affected and the values of its optimal parameters change with visible consequences for the total cost and the service level. If the perturbation is large enough, the difference between the real system and the model in terms of performance measures might be significant, especially when the optimal policy obtained as a solution of the model is implemented in the real system. When modeling inventory systems, attention must be paid to the estimation of the demand distribution but other parameters are also very important (e.g. shortage costs) and errors in the estimation of their values can greatly impact the performance of the system. Note however, that the performance of an (s, S) inventory model can be influenced by numerous parameters including but not limited to, the demand distribution and the cost structure. Furthermore, interactions between the various factors might exist, i.e., the effect of the changes in one factor can be modified by the variation of another parameter. Hence, numerical examples can demonstrate the existence of particular relationships but remain insuffi29
cient to get a complete understanding and draw conclusions for the general context or to quantify the effects of the perturbation of various parameters. Chu and Chung [14] have already pointed out the limitations of the numerical approach and the importance of analytical results in this context.
5
Conclusion
Approximations are usually employed when modeling real systems in order to obtain a mathematically tractable model which is a fair representation of the reality. Additionally, model parameters are usually estimated from empirical data by means of statistical methods. The analysis of the sensitivity of the mathematical model is necessary in order to ensure that the characteristics of the model are within an acceptable distance from those of the represented real system and that perturbations have only a bounded effect on the performance measures. In this paper, we study the sensitivity of an (s, S) inventory system represented by Markov chains. Firstly, we provide closed form expressions for the ergodicity coefficient of its transition matrix. Then, we estimate the deviation of the stationary characteristics of the model under the perturbation of the demand distribution. In the common practical situations decision-makers deal with model results obtained by approximation methods and with perturbed parameters (e.g. estimated from empirical data). Sensitivity analysis is very important because it allows us to know in which situation the constructed models are trustable and the obtained results are reliable. Sensitive models cannot be trusted as their output might change unexpectedly when the input is perturbed. The numerical examples show that the impact of the perturbations on the perfor30
mance of the system can be very high and lead to serious consequences. The perturbation bounds obtained in this paper are given in closed form and are simple to calculate. They provide a way for the quick and precise assessment of the sensitivity of the considered inventory model and help in decision-making. The results can be easily extended to large class of similar inventory models.
References
[1] Asmussen, S., 2003. Applied Probability and Queues. 2nd ed. Springer, New York. [2] A¨ıssani, D., Kartashov, N.V., 1983. Ergodicity and stability of Markov chains with respect to operator topology in the space of transition kernels. Dokl. Akad. Nauk. Ukr. SSR, ser. A 11, 3–5. [3] Bogataj, M., Bogataj, L., 2004. On the compact presentation of the lead times perturbations in distribution networks. International Journal of Production Economics, 88(2), 145-155. ˘ [4] Bogataj, L., Cibej, J.A., 1994. Perturbations in living stock and similar biological inventory systems. International Journal of Production Economics, 35(13), 233–239. [5] Borgonovo, E., 2010. Sensitivity analysis with finite changes: An application to modified EOQ models. European Journal of Operational Research, 200(1), 127–138. [6] Borgonovo, E., Peccati, L., 2007. Global sensitivity analysis in inventory management. International Journal of Production Economics, 108(1-2), 302– 313.
31
[7] Borgonovo, E., Peccati, L. (2011). Finite Change Comparative Statics for Risk Coherent Inventories. International Journal of Production Economics, 131(1), 52–62. [8] Boylan, E.S., 1969. Stability theorems for solutions to the optimal inventory equation. Journal of Applied Probability, 6, 211–217. [9] Browne, S., Zipkin, P., 1991. Inventory Models with Continuous, Stochastic Demands. The Annals of Applied Probability, 1(3), 419–435. [10] Butler,J., Jia, J., Dyer, J., 1997. Simulation techniques for the sensitivity analysis of multi-criteria decision models. European Journal of Operational Research, 103, 531-546, [11] Campbell,
S.L.,
Meyer,
C.D.,
1991.
Generalized inverses of linear transformations. Dover Publications Inc., New York. [12] Chen, F., Zheng, Y.S., 1997. Sensitivity analysis of an (s, S) inventory model. Operations Research Letters, 21(1), 19–23. [13] Cho, G.E., Meyer, C.D., 2000. Markov chain sensitivity measured by mean first passage times. Linear Algebra and its Applications, 316, 21-28. [14] Chu, P., Chung, K.-J., 2004. The sensitivity of the inventory model with partial backorders. European Journal of Operational Research, 152, 289-295. [15] Dobson, G., 1988. Sensitivity of the EOQ Model to Parameter Estimates. Operations Research,36(4) , 570–574. [16] Glasserman, P., Tayur, S.R., 1995. Sensitivity Analysis for Base-Stock Levels in Multi-Echelon Production-Inventory Systems. Management Science, 41, 263– 281.
32
[17] Hunter, J.J., 1986. Stationary distributions of perturbed Markov chains. Linear Algebra and its Applications, 82, 201–214. [18] Iglehart, D., 1963. Dynamic programming and stationary analysis of inventory problems, multistage inventory models and techniques. In : Scarf, H., Gilford, D., Shelly, M., editors, Multistage Inventory Models and Techniques, Stanford University Press, Stanford (Calif.), 1–31. [19] Kharna, A., Soman, C., Bandyopadhyay, T., 2014. Sensitivity analysis of the newsvendor model. European Journal of Operational Research, 239(2): 403-412. [20] Matheus, P., Gelders, L., 2000. The (R, Q) inventory policy subject to a compound Poisson demand pattern. International Journal of Production Economics, 68, 307–317. [21] Meyer, C.D., 1980. The condition of a finite Markov chain and perturbation bounds for the limiting probabilities. SIAM J. Algebraic Discrete Methods, 1(3): 273–283. [22] Mykytka, E.F, Ramberg, J.S., 1984. On the Sensitivity of the EOQ to Errors in the Forecast of Demand. IIE Transactions, 16(2), 144–151. [23] A. Paz, 1971. Introduction to Probabilistic Automata. Academic, New York. [24] Rabta, B., A¨ıssani, D., 2005. Strong Stability in an (R, s, S) Inventory Model. International Journal of Production Economics, 97(2), 159–171. [25] Rabta, B., A¨ıssani, D., 2008. Strong stability and perturbation bounds for discrete Markov chains. Linear Algebra and its Applications, 428(89), 19211927. [26] Ray, J., Chaudhuri, K.S., 1997. An EOQ model with stock–dependent demand, shortage, inflation and time discounting. International Journal of Production Economics, 53(2), 171–180.
33
[27] Ray, P.K., Sahu, S., 1992. Productivity measurement in multi-product manufacturing firms: Evaluation and control through sensitivity analysis. International Journal of Production Economics, 28(1), 71–84. [28] Scarf, H., 1960. The optimality of (S, s) policies in the dynamic inventory problem. In : Arrow,K.J., Karlin,S., Suppes, P., editors. Mathematical methods in the social sciences. Stanford University Press, Stanford (Calif.), 196–202. [29] Schweitzer, P., 1968. Perturbation theory and finite Markov chains. Journal of Applied Probability, 5, 401–413. [30] E. Seneta, 1984. Explicit forms for ergodicity coefficients and spectrum localization, Linear Algebra and its Applications, 60, 187–197. [31] E. Seneta, 1988. Perturbation of the stationary distribution measured by ergodicity coefficients. Advances in Applied Probability, 20, 288–230. [32] E. Seneta, 1991. Sensitivity analysis, ergodicity coefficients and rank-one updates for finite Markov chains. In : W.J. Stewart, ed., Numerical Solution of Markov Chains (Dekker, New York) pp. 121–129. [33] E. Seneta, 1993. Explicit forms for ergodicity coefficients of stochastic matrices. Linear Algebra and its Applications, 191, 245–252. [34] Sinha, B.K., Gupta, S.K., 1969. Sensitivity of Cost Functions in Inventory Models. AIIE Transactions, 1(3), 252–256. [35] Stoyan, D., 1983, Comparison Methods for Queues and Other Stochastic Models. English translation, D.J. Daley, Editor ( J. Wiley and Sons, New York). [36] Turanyi, T., Rabitz, H., 2000. Local Methods. In: Saltelli, A., Chan, K., Scott, E.M. (Eds.), Sensitivity Analysis, John Wiley & Sons, Chichester. [37] Veinott, A., Wagner, H., 1965, Computing optimal (s, S) policies, Management Sciences 11, 525–552.
34
[38] Wang, T., 2011. Sensitivity analysis for newsvendor model with general advertisement-sensitive demand. Procedia Engineering, 15, 1209-1213. [39] Yang, G. K., 2007, Note on sensitivity analysis of inventory model with partial backorders. European Journal of Operational Research, 177(2), 865-871. [40] Zheng, Y.S., 1992. On Properties of Stochastic Inventory Systems. Management Science, 38(1), 87–103.
35
Fig. 1. Deviation of the optimal policy parameters with respect to the perturbation of the arrival rate.
Fig. 2. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
1
Fig. 3. Deviation of the optimal inventory parameters with respect to the perturbation of the arrival rate.
Fig. 4. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
2
Fig. 5. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
3
www.elsevier.com/locate/ijpe
PII: DOI: Reference:
S0925-5273(17)30073-7 http://dx.doi.org/10.1016/j.ijpe.2017.03.014 PROECO6683
To appear in: Intern. Journal of Production Economics Received date: 23 September 2015 Revised date: 16 March 2017 Accepted date: 17 March 2017 Cite this article as: Boualem RABTA, Sensitivity analysis in inventory models by means of ergodicity coefficients, Intern. Journal of Production Economics, http://dx.doi.org/10.1016/j.ijpe.2017.03.014 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting galley proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Sensitivity analysis in inventory models by means of ergodicity coefficients Boualem RABTA ∗ Department of production management and logistics, Alpen–Adria–Universit¨ at Klagenfurt, Universitaetsstr. 65 - 67, 9020 Klagenfurt am W¨ orthersee - Austria. E-mail : [email protected]
Abstract We investigate the sensitivity of a class of (s, S) inventory models with respect to the perturbation of the demand distribution in terms of ergodicity coefficients. Ergodicity coefficients can be regarded as matrix norms that are useful to study qualitative as well as quantitative properties of some classes of stochastic systems. We obtain estimations of the absolute deviation of the stationary vector of the underlying Markov chain subject to perturbations. The particular structure of the transition matrix of this class of models allows us to derive simple closed form formulae for the computation of the ergodicity coefficients as well as the perturbation bounds. The perturbation bounds obtained can help us to decide whether or not the model remains an acceptable representation of the real system and thus decide whether or not it can be trusted for real life applications. Under perturbation, e.g. estimation errors or approximations, an optimal solution for the mathematical model may not be optimal for the real system and if implemented, the real performance measures may deteriorate considerably and deviate from the targeted values. Numerical examples are given to illustrate how the perturbation of the demand distribution may have a considerable impact on the optimal inventory policy and the performance
Preprint submitted to Elsevier
18 March 2017
measures in some cases. In a practical setting, the understanding of the sensitivity results for inventory models can help us to figure out which parameters have to be estimated with most attention and how to build models that are robust and as close as possible to the real systems they represent.
2
Key words: Inventory control, Sensitivity analysis, Markov chain, Ergodicity coefficient, Perturbation 1991 MSC: 60J10, 60J22, 90B05
1
Introduction
Many stochastic models including inventory models, are very complicated and analytical solutions are rarely available for them. In practical situations, the solution of such systems is obtained numerically and often after employing a series of approximations to obtain a mathematically tractable representation of the real system. Moreover, the model parameters are usually estimated from empirical data by means of statistical methods inducing further perturbations. Approximations must be used with care to avoid large deviations between the actual system and its mathematical representation (model). In this respect, the purpose of the sensitivity analysis is to determine whether or not the effect of these disturbances is limited, namely, a small perturbation of the model’s structure or parameters (inputs) generates only a small deviation of its characteristics (outputs), and to give an estimation of this effect. This is very important because only robust systems can be trusted and implemented for decision-making. A sensitive model can react unexpectedly to the perturbation of its input parameters and therefore cannot be considered reliable. The optimal solution of the mathematical model might only be sub–optimal for the real ∗ Department of production management and logistics, Alpen–Adria–Universit¨at Klagenfurt, Universitaetsstr. 65 - 67, 9020 Klagenfurt am W¨orthersee - Austria Email address: [email protected] (Boualem RABTA).
3
system. In many cases, this can result in over–investments in inventory and in unacceptable service levels [20]. Sensitivity analysis can help decision-makers to determine if the performance measures remain within acceptable bounds and to understand the importance of the various parameters in this context. If we obtain an affirmative response to the qualitative question (robustness), then we might be able to estimate the deviation of the characteristics of the system caused by the perturbations (quantitative question). In general, multiple parameters in the inventory system are unknown and must be approximated or estimated from empirical data. Customers’ demand, replenishment leadtimes and system costs (ordering, holding and shortage costs) are among the most important elements that are subject to perturbations. The classical Economic Order quantity model (EOQ) is well–known to be somewhat insensitive to the change in order quantity and against parameter estimates (see, e.g. [15]). However, it can be quite sensitive to forecast errors whenever the leadtime is non–zero [22]. Extensions of the classical EOQ model have been presented in numerous works and many of them also discuss the sensitivity analysis of the optimal order quantity (see, e.g. [5,7] and references therein). Wang [38] applied sensitivity analysis to the changes of optimal advertising expenditure and ordering decisions in a newsvendor model under advertising-dependent demand distribution. Khanra [19] reviewed the literature on sensitivity analysis of newsvendor inventory models and studied the sensitivity of the classical one. The author showed that it is more sensitive than the EOQ model to sub–optimal ordering decisions and identified conditions for symmetry/skewness of cost deviation found to be closely linked with symmetry/skewness of the demand density function. Inventory models were the first stochastic models for which monotonicity prop4
erties were proved [35]. Boylan [8] has proved that the solution of the optimal inventory equation depends continuously on its parameters. Sinha and Gupta [34] studied the effects of alternative forms of surplus and shortage cost functions and also investigated the effects of alternative forms of probability density function of demand in an inventory control situation. In [16], Glasserman and Tayur studied the sensitivity of base-stock levels in multi– echelon production–inventory systems, in particular with respect to demand parameters. Zheng [40] investigated the sensitivity of a standard single-item continuous-review (s, Q) inventory system with fixed leadtimes and found it to be more robust than its deterministic counterpart, the EOQ model. Chen and Zheng [12] established that the average cost in the (s, S) inventory model with exponentially distributed demand is quite insensitive to the changes in D = S − s. In [24], Rabta and A¨ıssani investigated the strong stability of an (s, S) inventory model with instantaneous deliveries. After proving the robustness of the underlying Markov chain with respect to the perturbation of the demand distribution, they obtained an upper bound for the deviation of the stationary distribution of the underlying Markov chain [24].
In many cases, sensitivity analysis is performed numerically (see, e.g. references in [14,39]). However, the numerical approach has proved to be insufficient. Chu and Chung [14] derived sensitivity formulas for an inventory model with backorders and showed that the numerical analysis might not lead to consistent sensitivity results. Hence, an analytical approach is required in order to obtain valid conclusions regarding the qualitative properties of the system. Local sensitivity analysis methods include the study of one parameter at a time and testing the effect when changing the value within a specific range (e.g. [27,26]) along with perturbation analysis (e.g. [4,3]) among others (see, 5
[36]). On the other hand, global methods test the sensitivity of the model in consideration of the uncertainty distribution reflecting the decision-maker state of knowledge in the parameters (see, e.g. [6] and the references therein). Borgonovo and Peccati [6] employed a global sensitivity analysis method based on Sobol functions and variance decomposition method to determine the most influential parameters on the model output. Borgonovo [5] introduced a new method to study of the effect of changes in one or more of the exogenous variables, defining sensitivity measurements that do not rest on differentiability and applied it to modified EOQ models. Sensitivity analysis can also be performed numerically (e.g. [14,39]) or by simulation, where results that can easily be analyzed statistically are generated (see, e.g. [10]). General techniques of stochastic analysis were also applied to stochastic inventory models. For instance, Stoyan [35] investigated their monotonicity properties while Rabta and A¨ıssani [24] analysed an (s, S) inventory model based on the results of the theory of the strong stability of Markov chains [2].
A large class of infinite horizon stochastic inventory systems can be represented by Markov chains (e.g. [9,24]). In this case, the focus is generally on the steady state of the system and the characteristics are expressed in terms of the stationary distribution. There is an abundant literature on the robustness of Markov chains with tools and results that can be applied to inventory models. For finite Markov chains, the problem of sensitivity to perturbations has been treated by many authors since Schweitzer (1968) [29]. Meyer singly and with co-authors used the group inverse notion [11] to derive several bounds (see, e.g. [21]). Hunter [17] used a wider class of generalized inverses. The coefficient of ergodicity has been used by Seneta [31,32]. Also, Cho and Meyer [13] obtained a perturbation bound expressed by means of mean first passage 6
times. The validity of most of these bounds extends easily to infinite discrete Markov chains [25] as well as Markov chains with general state space [2].
Among those results, only few have been applied to inventory models. The bound derived in [24] requires the computation of a part of the stationary distribution which is expressed in terms of the renewal function of the demand distribution.
In this paper, we investigate the sensitivity of an (s, S) inventory model as represented by Markov chains and subject to perturbations of the demand distribution. We use the so-called ergodicity coefficients with the aim of deriving bounds on the deviation of the stationary characteristics of the model that could be easily computed. The bound proposed in this paper, based on the ergodicity coefficient, does not involve the stationary distribution and is simple and easy to compute.
The generic (s, S) inventory model analyzed in the sequel assumes zero leadtime but the obtained results extend to a large set of similar models with only minor changes in the assumptions (see, remarks below). The particular structure of the transition matrices of this kind of model allows the computation of the corresponding ergodicity coefficient in a very simple and computationally efficient manner. Furthermore, we show that the obtained closed form formulas reduce to very simple forms for particular demand distributions. The obtained bounds are given in terms of the norm k.k1 of the stationary vector of the underlying Markov chain without the need to explicitly calculate this stationary distribution. We provide numerical examples to show the optimal inventory policy is affected by the change to the demand distribution and how this can affect the performance measures such as total cost and service level. 7
2
Preliminary and notations
In this article, we use the norm k.k1 for row vectors and the norm k.k∞ for column vectors and for matrices. Also, all vectors are column vectors (row vectors are transposed). Since it will always be clear if the norm corresponds to a matrix, a row vector or a column vector, all norms are noted k.k. Then, for a vector v = (vj )j∈E we have : kvk = kvk∞ = max |vi |, i∈E
kv T k = kv T k1 =
X
|vj |.
j∈E
For a matrix B = (bij )i,j∈E , we have : kBk = kBk∞ = max i∈E
X
|bij |.
j∈E
The following result (Paz’s inequality), given in [23, Ch. IIa] will be needed in the sequel : Lemma 1 Let c 6= 0 be a real vector such that cT e = 0. Then, for any complex vector d we have : |cT d| ≤ kcT k
maxi,j |di − dj | . 2
Here and in the sequel, eT = (1, 1, .., 1) is the vector of all ones. The ergodicity coefficient of a matrix B = (bij )i,j∈E with equal row sums is defined as ([30]): τ (B) =
maxT
uT e=0,
ku k<1
kuT Bk.
Seneta [30] gives an explicit expression of this quantity : τ (B) =
X 1 max |bik − bjk |. 2 i,j k
(1)
8
The ergodicity coefficient of a finite stochastic matrix P (with respect to the norm k.k1 ) satisfies the following properties (see, e.g. [33]) : • 0 ≤ τ (P ) ≤ 1. • τ (P ) = 0 if and only if P has rank 1 (i.e., P = ev T , for some probability vector v). • For every eigenvalue λ 6= 1 of P we have : |λ| ≤ τ (P ). Additionally, τ (P n ) ≤ (τ (P ))n . Consider a homogeneous irreducible Markov chain X with states in a finite space E. Denote by P the transition matrix of X and by π its stationary vector, i.e., the unique real vector satisfying π T P = π T , π T e = 1.
Suppose that P is perturbed to Q = P + ∆ such that Q corresponds to the transition matrix of another irreducible Markov chain Y with states in E. Denote by ν the stationary vector of the Markov chain Y . For finite Markov chains, Seneta [31] gives the following result. Theorem 2 For a finite irreducible Markov chain X satisfying τ (P ) < 1 we have : kν T − π T k ≤
kQ − P k . 1 − τ (P )
However, we observe that for two arbitrary probability vectors, the norm of their difference cannot exceed 2 since, kν T − π T k ≤ kν T k + kπ T k = 2. 9
Observe that for any two probability vectors π and ν, the vector (ν−π) satisfies (ν − π)T e = 0. Then, by taking d = ek (where ek is the k−th column of the identity matrix) and applying Paz’s inequality we obtain for every k ∈ E : 1 |νk − πk | ≤ kν T − π T k. 2 Hence, 1 max |νk − πk | ≤ kν T − π T k. k∈E 2
3
Perturbation bounds for the (s, S) inventory model
Consider the following single-item, single-echelon, inventory problem. We manage a non–perishable product in the face of a stationary stochastic demand. The inventory level is inspected every R time units. At the beginning of each period, the manager decides whether or not he should make a replenishment order and if yes, how much to order. Suppose that the outside supplier is perfectly reliable and that orders arrive immediately. We assume that the consecutive period demands are independent and identically distributed (i.i.d.) random variables. Remark 3 The period length parameter R does not play a direct role in the subsequent analysis. In fact, R might be random as in the case of the continuous review model (s, S) (e.g. time between consecutive arrivals of customers) or fixed as in the case of the periodic review model. The structure of the transition matrix given below remains the same and hence the results are valid for both cases. 10
For both continuous and periodic review models, it has been proved that the (s, S) inventory policy is optimal (see [28,18,37]). According to this control rule, the level is reviewed at dates tn = R1 + R2 + .. + Rn , n ≥ 0 where Ri can be deterministic (periodic review, Ri = R∀i) or stochastic (continuous review, Ri are i.i.d. random variables). If at the review moment, the inventory level Xn is below or equal to s (Order-Point) we order sufficient quantity to raise the inventory level to S (Order-Up-To level). Thus, the size of the order is equal to Zn = S − Xn . Many algorithms have been proposed for the computation of the optimal values of the parameters of this class of inventory models (see, for example [37]).
During period n, n ≥ 1, the total demand is a random variable ξn . Assume that ξn , n ≥ 1, are independent and identically distributed (i.i.d.) random variables with common probabilities given by
ak = P (ξ1 = k), k = 0, 1, 2, ...
So, ak is the probability to have a demand of k items during the period.
For the purpose of our analysis, it suffices to define the demand distribution by the finite sequence ak for k = 0, 1, .., S − 1 representing the first S demand probabilities. This observation allows us to deal with only finite sums while maintaining full accuracy and precision in the derived formulas. Therefore, we denote by αz the quantity:
αz = P (ξ1 ≥ z) = 1 −
z−1 X
ak
k=0
and assume that ak > 0, k = 0, .., S − 1 and αS > 0. 11
Additionally, we exclude the special case s = S − 1 (base stock policy) since for this model the ergodicity coefficient will be 0 and the Inventory Markov chain is already in its permanent regime where all rows of P are identical and the stationary vector is given by a closed form (any row of P ). Therefore, the deviation of the stationary vector can be computed exactly and the estimation of perturbation bounds is not necessary.
Remark 4 Considering the assumption of zero-leadtime, the model includes both lost sales and backorders cases (or a mixture policy). If the demand during a given period is greater than the available quantity, then the excess quantity is either backordered or the sales are lost. In both cases, we order sufficient quantity to raise the inventory level to the order-up-to point S and eventually satisfy the backordered demand. The transition matrix of the on-hand inventory process will not be affected.
Assume that the initial inventory level X0 > s (if not, we raise it to S by ordering sufficient quantity). The inventory level Xn+1 at the end of period n + 1 is given by : + (S − ξn+1 )
if Xn ≤ s,
Xn+1 =
+ (Xn − ξn+1 ) if Xn > s,
where (A)+ = max(A, 0).
The random variable Xn+1 depends only on Xn and ξn+1 , where ξn+1 is independent of n and the state of the system. Thus, X = {Xn , n ≥ 0} is a homogeneous Markov chain with values in E = {0, 1, ..., S} with transition 12
matrix P given by :
0
1
s
s+1
S
0
P∞
ak aS−1 · · · aS−s aS−s−1 · · · a0
1
P∞
ak aS−1 · · · aS−s aS−s−1 · · · a0
S
S
.. .
.. .
..
.
.. .
.. .
..
. . ..
P = P∞
s s+1
S
P∞
s+1
.. . S
ak aS−1 · · · aS−s aS−s−1 · · · a0
P∞ S
ak as · · · a1 .. .
...
.. .
a0
0
0
.. .
...
0
ak aS−1 · · · aS−s aS−s−1 · · · a0
The Markov chain X is finite and irreducible. It admits a unique stationary vector π.
3.1
The coefficient of ergodicity of P
(1) allows us to easily compute the coefficient of ergodicity of a finite Markov chain in general. Yet, the particular structure of the transition matrices of inventory models allows further simplifications. Observe that rows 0 to s of the transition matrix P are all identical and they are same as the S-th row. Therefore we might write: τ (P ) =
X 1 max |pik − pjk |. 2 s+1≤i
To compute |pik − pjk | we distinguish the following cases: 13
• k = 0, |pik − pjk | = |
P∞
l=i
al −
P∞
l=j
al | =
Pj−1 l=i
al .
• 1 ≤ k ≤ i, |pik − pjk | = |ai−k − aj−k |. • i + 1 ≤ k ≤ j, |pik − pjk | = aj−k . • j + 1 ≤ k ≤ S (for j < S), |pik − pjk | = 0. Thus,
j−1 j i X X X 1 τ (P ) = max al + |ai−k − aj−k | + aj−k 2 s+1≤i
Hence, Lemma 5 The coefficient of ergodicity of the transition matrix P is given by:
j−1 j−i−1 i−1 X X X 1 |ak − aj−i+k | + ak + ak . τ (P ) = max 2 s+1≤i
(2)
The coefficient of ergodicity of the transition matrix P could therefore be calculated exactly by means of a finite loop. Furthermore, we have : Lemma 6 The coefficient of ergodicity τ (P ) of the transition matrix P of the underlying Markov chain in the (s, S) inventory model satisfies : τ (P ) ≤ 1 − αS < 1.
PROOF.
j−1 j−i−1 i−1 X X X 1 τ (P ) = max |ak − aj−i+k | + ak + ak 2 s+1≤i
j−1 j−i−1 i−1 X X X 1 ≤ max (ak + aj−i+k ) + ak + ak 2 s+1≤i
j−1 j−1 j−i−1 X X X 1 ak + ak + ak ≤ max 2 s+1≤i
14
j−1 j−1 j−1 S−1 X X X X 1 = max ak + ak = max ak = ak s+1≤i
τ (P ) ≤ 1 − αS < 1.
3.1.1
Particular demand distributions
The expression of the ergodicity coefficient of the transition matrix P derived above could be further simplified for particular demand distributions. For instance, in the result below we impose conditions on the first S elements of the demand distribution. Lemma 7 If the sequence ak , k = 0, 1, .., S − 1 is non–increasing then: τ (P ) = 1 − αS−s−1 . Similarly, if the sequence ak , k = 0, 1, .., S − 1 is non–decreasing then: τ (P ) = αs+1 − αS .
PROOF. Suppose that ak , k = 0, 1, .., S − 1 is non–increasing then:
τ (P ) =
j−1 X
i X
j X
i X
1 max al + ai−k − aj−k + aj−k 2 s+1≤i
j−1 j−1 j−i−1 i−1 X X X X 1 = max al + ak − ak + ak 2 s+1≤i
j−1 j−1 j−1 X X X 1 = max ak − 2 ak + ak 2 s+1≤i
=
max
s+1≤i
j−1 X
ak −
k=0
j−1 X
ak
k=j−i
Thus, τ (P ) =
max
j−i−1 X
s+1≤i
k=0
15
ak =
S−s−2 X k=0
ak
Now, suppose that ak , k = 0, 1, .., S − 1 is non–decreasing then:
j−1 j i i X X X X 1 τ (P ) = max al − ai−k + aj−k + aj−k 2 s+1≤i
=
j−1 X
j−1 X
i−1 X
1 max al − ak + ak + 2 s+1≤i
j−1 X
j−1 X
j−1 X
j−i−1 X
ak
k=0
j−1 X 1 max 2 al − al + ak = max al = s+1≤i
Hence, τ (P ) =
S−1 X
al =
l=s+1
∞ X
al −
l=s+1
∞ X
al
l=S
The first condition of the previous lemma is satisfied by the following two discrete probability distributions.
3.1.1.1
Uniformly distributed demand In this case, the total demand
follows a (discrete) uniform distribution on the set {0, 1, .., M } with M ≥ S. ak =
1 , k = 0, 1, .., M. M +1
The coefficient of ergodicity τ (P ) of the transition matrix P writes: τ (P ) =
3.1.1.2
(S − s − 1) M +1
Geometrically distributed demand In this case, the demand
distribution is given by : ak = q k (1 − q), k = 0, 1, ... for 0 < q < 1. The coefficient of ergodicity τ (P ) of the transition matrix P writes: 16
τ (P ) = 1 − q S−s−1
3.2
Perturbation of the demand distribution
In this section, we replace the demand distribution in the previously described inventory model by another distribution close to the original one in a certain sense (we will define a measure of proximity (distance) of the two distributions). This case is encountered in many circumstances : (1) The demand distribution is fitted from empirical data. The perturbations result mainly from the estimation errors. (2) The demand distribution is replaced by another distribution for purposes such as simplification of the model and ease of the analysis. For instance, the set of phase-type distributions is dense in the set of positive continuous probability distributions [1]. Hence, a positive continuous probability distribution could be approximated (replaced) by a phase-type distribution at any desired precision level. (3) The form of the perturbed distribution is the same as the original one but differs in the value of one or more parameters. The difference comes for instance from errors of the estimation of the parameter by means of statistical methods. Consider another inventory model with the same structure but with period demands ξi0 , i ≥ 1 independent and identically distributed with common probabilities a0k = P (ξ10 = k), k = 0, 1, 2, ... 17
Denote by Y the inventory level Markov chain in the new system and by Q its transition matrix. First, we estimate the effect of the perturbation on the transition matrix of the original Markov chain. Namely, we compute a measure of the distance between the transition matrices of the two chains. Lemma 8 The deviation of the transition matrix P is given by : kQ − P k = |αS0 − αS | +
S−1 X
|a0k − ak |.
k=0
PROOF. kQ − P k = max i∈E
X
|qij − pij | = max(A, B)
j∈E
where, A = max
0≤i≤s
∞ X
=|
∞ X
a0k −
k=S
ak | +
k=S
S X
S X
|qij − pij |
j=0
|a0S−j − aS−j | = |αS0 − αS | +
j=1
S−1 X
|a0k − ak |
k=0
and B = max
s+1≤i≤S
"
≤ max
s+1≤i≤S
|
S X
|qij − pij | = max | s+1≤i≤S
j=0
∞ X k=S
a0k
−
∞ X k=S
ak | +
S−1 X
|a0k
∞ X
a0k −
k=i
− ak | +
k=i
∞ X k=i
i−1 X k=0
ak | +
i X
|a0i−j − ai−j |
j=1
#
|a0k
− ak | =
S−1 X 0 |αS −αS |+ |a0k −ak |. k=0
Similarly to the formula of the coefficient of ergodicity, the deviation of the transition matrix P subject to perturbation is given by a finite sum. Moreover, the computation involves only basic arithmetic and therefore the total computation time is linear with respect to the dimension of the matrix P . It 18
is not necessary to store the transition matrix or the complete demand distribution (only, the first S individual probabilities are required). Hence, the computations could be performed efficiently and scale well for large models. Using Theorem 2 we obtain the following result : Theorem 9 Let π and ν be the stationary distributions of chains X and Y respectively. Then : 0 |α0 − αS | + S−1 i=0 |ai − ai | . kν − π k ≤ S 1 − τ (P ) T
P
T
(3)
where τ (P ) is given by (2).
PROOF. This result holds by virtue of Theorem 2 and lemma 8.
3.3
Deviation of the mean inventory level
The previous perturbation bounds can be used to estimate the deviation of several quantities of interest. To illustrate this we compute the perturbation bounds for the average end period inventory level: ¯ (resp. Y¯ ) be the mean end of period inventory level in Theorem 10 Let X the base (resp. perturbed) system. Then : ¯ ≤ S kν T − π T k. |Y¯ − X| 2
PROOF. We have : ¯= X
S X j=0
19
jπj ,
Y¯ =
S X
jνj .
j=0
Then ¯ =| |Y¯ − X|
S X
jνj −
j=0
S X
jπj | = |(ν T − π T )f |,
j=0
where, f = (0, 1, .., S)T . Thus, by using Lemma 1 we get :
¯ ≤ |Y¯ − X|
3.4
S T kν − π T k. 2
Deviation of the inventory costs and service level
The individual components πi of the stationary vector π of the chain X give the probability that, at the equilibrium state, the on-hand inventory level is equal to i at the end of any period. The total cost is usually the sum of three components, namely, ordering (fixed and variable), holding and penalty costs. The latter being in general difficult to estimate, we might alternatively consider a service level in addition to an inventory cost function that is the sum of ordering and holding costs only.
Let’s assume lost-sales and the following cost structure : A fixed ordering cost K is incurred whenever a replenishment order is placed and is added to the variable acquisition cost of c per ordered unit of the product. The unit holding cost is h per unit of time. In addition to the replenishment and holding costs (C1 ) we consider a service level (ψ1 ). We illustrate the use of the above inequalities in the study of the sensitivity of the key performance measures of the considered inventory model. 20
3.4.1
Replenishment and holding cost
The overall cost function might have a different structure depending on the assumptions of the model. For the purpose of the following illustration we restrict the analysis to the case where R is constant and holding costs are computed from the average end period inventory level. Under the previous assumptions, we write the average total (period) replenishment and holding cost of the model :
C1 (π) =
S X
hiRπi +
i=1
s X
(K + c(S − i))πi = π T ζ,
i=0
where ζ is the vector for which components are given by K + c(S − i) + hiR if 0 ≤ i ≤ s
ζi =
hiR
if s < i ≤ S
Similarly, the average total replenishment and holding cost in the perturbed system is C1 (ν) = ν T ζ and by applying Paz’s inequality, we obtain |C1 (ν) − C1 (π)| = |(ν T − π T )ζ| ≤ kν T − π T k
maxi ζi − mini ζi . 2
Hence, |C1 (ν) − C1 (π)| ≤
kQ − P k (maxi ζi − mini ζi ) . 2(1 − τ (P ))
If we assume that c > hR (the unit purchase cost is greater than the unit period holding cost) then
max ζi = max {K + cS, K + cS + hR − c, hSR} = K + cS, i
21
and min ζi = min {K + cS, K + c(S − s) + hsR, h(s + 1)R} = h(s + 1)R i
Hence, Theorem 11 Under the previous assumptions on the cost function, we have :
|C1 (ν) − C1 (π)| ≤
3.4.2
|αS0 − αS | +
PS−1 i=0
|a0i − ai | (K + cS − h(s + 1)R) 2(1 − τ (P ))
.
Service level
We consider a service level ψ1 defined as follows: ψ1 = P {Inventory at the beginning of the period ≥ Period demand}. At the equilibrium, the inventory at the end of one period is equal to i with probability πi . Therefore, at the beginning of the next period the on hand inventory will be equal to S if i ≤ s and equal to i otherwise. Hence, ψ1 (π) =
s X
πi P (ξ1 ≤ S) +
s X i=0
πi (1 − αS+1 ) +
S X
πi P (ξ1 ≤ i)
i=s+1
i=0
=
S X
πi (1 − αi+1 ) = 1 −
s X i=0
i=s+1
πi αS+1 −
S X
πi αi+1
i=s+1
Thus s S s S X X X X 0 0 πi αS+1 − πi αi+1 − 1 − νi αS+1 − νi αi+1 |ψ1 (π) − ψ1 (ν)| = 1 − i=0 i=s+1 i=0 i=s+1 s S s S X X X X 0 0 = πi αS+1 + πi αi+1 − νi αS+1 + νi αi+1 i=0 i=s+1 i=0 i=s+1 X S s S X X X s = πi αS+1 + πi αi+1 − νi αS+1 + νi αi+1 i=0 i=s+1 i=0 i=s+1 s S s S X X X X 0 0 + νi αS+1 + νi αi+1 − νi αS+1 + νi αi+1 i=0 i=s+1 i=0 i=s+1
22
X S s S X X X s = πi αS+1 + πi αi+1 − νi αS+1 − νi αi+1 i=0 i=s+1 i=0 i=s+1 s X − νi α 0
S+1
+
i=0
S X
0 νi αi+1 −
i=s+1
s X
νi αS+1
i=0
− νi αi+1 i=s+1 S X
s S X X = (πi − νi )αS+1 − (πi − νi )αi+1 i=0 i=s+1 s X − νi (α0
S+1
S X
− αS+1 ) −
i=0
0 νi (αi+1
i=s+1
− αi+1 ) .
Whence,
|ψ1 (π) − ψ1 (ν)| = (π T − ν T )η − ν T (η 0 − η) , 0 ), where η (resp. η 0 ) is the vector having component i equals to αS+1 (resp. αS+1 0 ), if s + 1 ≤ i ≤ S. Consequently, if 0 ≤ i ≤ s and to αi+1 (resp. αi+1
|ψ1 (π) − ψ1 (ν)| ≤ (π T − ν T )η + ν T kη 0 − ηk .
Applying Lemma 1, we obtain :
|ψ1 (π) − ψ1 (ν)| ≤ kπ − νk
maxi,j |ηi − ηj |
T
0 + ν kη − ηk , 2
0 where maxi,j |ηi − ηj | = αs+1 −αS+1 , ν T = 1 and kη 0 − ηk = maxi=s+1,..,S |αi+1 −
αi+1 |. Finally : Theorem 12 Under the previous assumptions, the deviation of the service level ψ1 satisfies: 0 |α0 − αS | + S−1 i=0 |ai − ai | αs+1 − αS+1 0 −αi+1 |. |ψ1 (π) − ψ1 (ν)| ≤ S + max |αi+1 i=s+1,..,S 1 − τ (P ) 2
P
23
4
Effect of the perturbations on the optimal inventory policy: Numerical examples
Consider the following base scenario: under the model’s assumptions above, we assume that the inventory is reviewed periodically with period R = 7. Customers arrive according to a stationary Poisson distributed process with rate λ = 6 (per unit of time) and all customers request exactly one unit of the product. The cost function elements are : Fixed ordering cost K = 25 per order, variable ordering cost c = 3 per unit ordered and holding cost h = 0.1 per unit in stock at the end of the period and per unit of time. Additionally, we either consider a shortage cost p = 5 per unit of unsatisfied demand or impose a service level constraint (ψ1 = 0.95). In the following examples, we illustrate the effect of the perturbations on the optimal inventory policy as well as the performance measures of the system. We implemented a linear search algorithm that determines the optimal inventory policy in the following two cases:
• Minimization of the sum of ordering and holding costs under service level constraint. • Minimization of the total cost (sum of ordering, holding and shortage costs).
Obviously, the solution might be different in the two cases and we expect a different impact on the performance measures. In the first case, the performance of the inventory system is more affected in terms of costs than in service level which is limited by the constraint.
In each of the above mentioned cases, we introduce a small perturbation on 24
the demand distribution by varying the customers’ arrival rate. Two issues are investigated. We first show the impact of the perturbation of the demand distribution on the optimal parameters s∗ andS ∗ of the inventory policy and we compare the performance measures of the original and perturbed system (optimal policy in the model vs. optimal policy in the real system). However, the situation is different in practice. The decision is made based on the solution of the mathematical model and the obtained optimal solution is implemented in the real system. Since the parameters of the real system are different from those of the mathematical model, the implemented policy is not optimal anymore and another error component is added to the approximation error. We illustrate this situation by comparing the outcome of the same policy in both the original and perturbed model. Therefore, in this case we compare the optimal solution of the mathematical model with the sub–optimal policy in the real system. In each (sub–)case, we proceed as follows: • Calculate the optimal policy, total cost and service level under the base scenario’s assumptions. • Perturb the demand distribution and calculate key performance measures in the perturbed model. • Compare the original and the perturbed models in terms of performance measures. • Repeat these steps for different perturbation amplitudes. To perturb the demand distribution, we act on the customers’ arrival rate and vary it by a small amount each time. The perturbation amplitude will range between -10% and +10% of the value of the arrival rate. This is sufficient to 25
show the impact of the perturbation while remaining in a realistic situation. In each, we report the values of the optimal inventory policy parameters s∗ , S ∗ , the total cost C ∗ and the service level ψ1∗ for both the original and perturbed ∗
models as well as the relative deviation of these performance measures ( ∆C C∗ and
4.1
∆ψ1∗ ψ1∗
in percentage). Discussion on the outcome of the experiments follows.
Example 1 : Sensitivity of the optimal inventory policy under minimization of ordering and holding costs with service level constraint
In practice, shortage costs are very difficult to estimate. Along with the profit loss, shortages can lead to shutdowns of the production or influence customer satisfaction and retention which may have a negative impact on future demand and revenue. Therefore, in many cases a target service level is fixed and the system is optimized under this constraint rather than with respect to the total cost which includes an estimated (imprecise) shortage cost component. Let us set the service level to 95%. The algorithm will search for the optimal inventory policy by minimizing the sum of ordering and holding costs under the service level constraint. The results are reported in Figure 1. Figure 1 goes here. Fig. 1. Deviation of the optimal policy parameters with respect to the perturbation of the arrival rate.
The effect of the perturbation in this case remains somewhat limited. For the maximum amplitude (10%) of the perturbation of the customers’ arrival rate, we observe only 8.5% deviation of the sum of ordering and holding costs while the service level is kept near to the target level by the constraint. The change in the cost follows naturally the change in the optimal inventory policy 26
parameters s∗ and S ∗ . Next, we inject the solution of the mathematical model into the perturbed system and compare the performance measures. Results are represented in Figure 2. Figure 2 goes here. Fig. 2. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
The effect of the perturbation on the cost remains limited in this case as well but we observe a considerable change in the service level. A 10% increase in the customers’ arrival rate induced a drop of the service level from 0.95 to 0.86 (10%). In many real situations, the impact of the deterioration of the service level by this magnitude can be considerable.
4.2
Example 2 : Sensitivity of the optimal inventory policy under total cost minimization
In this example, we aim to show how the perturbation of the demand distribution (say, estimation errors) can affect the optimal values of the (s, S) inventory policy when the latter is obtained by minimization of the total cost function, i.e., the sum of ordering, holding and shortage costs. The results of the procedure are shown in Figures 3 and 4. Figure 3 goes here. Fig. 3. Deviation of the optimal inventory parameters with respect to the perturbation of the arrival rate.
The perturbation of the demand distribution has a direct impact on the op27
Figure 4 goes here. Fig. 4. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
timal inventory policy (optimal s and S). In the example above, the optimal inventory policy in the base scenario is (25, 43). When varying the arrival rate between -10% and +10% of its nominal value, the optimal inventory policy changes between (21, 39) to (29, 48). Naturally, this will lead to the deviation of the performance measures of the (s, S) inventory such as the total inventory cost and the service level. For example, a 10% error in the customers’ arrival rate generates around 8.2% deviation of the total inventory cost and up to 20.7% on the service level. In the following example, we implement the optimal inventory policy obtained in the base scenario in the perturbed model and report the results in Figure 5. Figure 5 goes here. Fig. 5. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
Notice that in this case, the deviation in the stationary performance measures might be higher than that in previous case. For example, a 10% error in the customers’ arrival rate can lead to a deviation of the service level as high as 46.4% (compared to 20.7% previously). Observe as well that the service level has dropped to 47% which means a high level of unsatisfied demand that can lead to serious consequences (e.g. loss of customers, production shutdown,...). Here, the deviation resulting from the approximation of the demand distribution is added to the deviation resulting from the use of a non–optimal policy in the real situation. 28
Finally, observe that without constraint, the service level performance measure might be at a low level even with optimal inventory policy. This could be explained by the fact that the shortage costs are poorly estimated which in turn increases the tolerance of the system to shortages allowing more unsatisfied demand. In the first case, i.e., minimization of the cost under service level constraint, the (s, S) inventory system looks more robust. If one can ensure that the perturbation is kept small then the constraint will help to keep the service level near to the target value. It seems also that costs are least affected in this case. The examples above show that perturbations of the demand distribution can affect the output of the (s, S) inventory model in various ways. In these examples, the optimal inventory policy (s∗ , S ∗ ) is affected and the values of its optimal parameters change with visible consequences for the total cost and the service level. If the perturbation is large enough, the difference between the real system and the model in terms of performance measures might be significant, especially when the optimal policy obtained as a solution of the model is implemented in the real system. When modeling inventory systems, attention must be paid to the estimation of the demand distribution but other parameters are also very important (e.g. shortage costs) and errors in the estimation of their values can greatly impact the performance of the system. Note however, that the performance of an (s, S) inventory model can be influenced by numerous parameters including but not limited to, the demand distribution and the cost structure. Furthermore, interactions between the various factors might exist, i.e., the effect of the changes in one factor can be modified by the variation of another parameter. Hence, numerical examples can demonstrate the existence of particular relationships but remain insuffi29
cient to get a complete understanding and draw conclusions for the general context or to quantify the effects of the perturbation of various parameters. Chu and Chung [14] have already pointed out the limitations of the numerical approach and the importance of analytical results in this context.
5
Conclusion
Approximations are usually employed when modeling real systems in order to obtain a mathematically tractable model which is a fair representation of the reality. Additionally, model parameters are usually estimated from empirical data by means of statistical methods. The analysis of the sensitivity of the mathematical model is necessary in order to ensure that the characteristics of the model are within an acceptable distance from those of the represented real system and that perturbations have only a bounded effect on the performance measures. In this paper, we study the sensitivity of an (s, S) inventory system represented by Markov chains. Firstly, we provide closed form expressions for the ergodicity coefficient of its transition matrix. Then, we estimate the deviation of the stationary characteristics of the model under the perturbation of the demand distribution. In the common practical situations decision-makers deal with model results obtained by approximation methods and with perturbed parameters (e.g. estimated from empirical data). Sensitivity analysis is very important because it allows us to know in which situation the constructed models are trustable and the obtained results are reliable. Sensitive models cannot be trusted as their output might change unexpectedly when the input is perturbed. The numerical examples show that the impact of the perturbations on the perfor30
mance of the system can be very high and lead to serious consequences. The perturbation bounds obtained in this paper are given in closed form and are simple to calculate. They provide a way for the quick and precise assessment of the sensitivity of the considered inventory model and help in decision-making. The results can be easily extended to large class of similar inventory models.
References
[1] Asmussen, S., 2003. Applied Probability and Queues. 2nd ed. Springer, New York. [2] A¨ıssani, D., Kartashov, N.V., 1983. Ergodicity and stability of Markov chains with respect to operator topology in the space of transition kernels. Dokl. Akad. Nauk. Ukr. SSR, ser. A 11, 3–5. [3] Bogataj, M., Bogataj, L., 2004. On the compact presentation of the lead times perturbations in distribution networks. International Journal of Production Economics, 88(2), 145-155. ˘ [4] Bogataj, L., Cibej, J.A., 1994. Perturbations in living stock and similar biological inventory systems. International Journal of Production Economics, 35(13), 233–239. [5] Borgonovo, E., 2010. Sensitivity analysis with finite changes: An application to modified EOQ models. European Journal of Operational Research, 200(1), 127–138. [6] Borgonovo, E., Peccati, L., 2007. Global sensitivity analysis in inventory management. International Journal of Production Economics, 108(1-2), 302– 313.
31
[7] Borgonovo, E., Peccati, L. (2011). Finite Change Comparative Statics for Risk Coherent Inventories. International Journal of Production Economics, 131(1), 52–62. [8] Boylan, E.S., 1969. Stability theorems for solutions to the optimal inventory equation. Journal of Applied Probability, 6, 211–217. [9] Browne, S., Zipkin, P., 1991. Inventory Models with Continuous, Stochastic Demands. The Annals of Applied Probability, 1(3), 419–435. [10] Butler,J., Jia, J., Dyer, J., 1997. Simulation techniques for the sensitivity analysis of multi-criteria decision models. European Journal of Operational Research, 103, 531-546, [11] Campbell,
S.L.,
Meyer,
C.D.,
1991.
Generalized inverses of linear transformations. Dover Publications Inc., New York. [12] Chen, F., Zheng, Y.S., 1997. Sensitivity analysis of an (s, S) inventory model. Operations Research Letters, 21(1), 19–23. [13] Cho, G.E., Meyer, C.D., 2000. Markov chain sensitivity measured by mean first passage times. Linear Algebra and its Applications, 316, 21-28. [14] Chu, P., Chung, K.-J., 2004. The sensitivity of the inventory model with partial backorders. European Journal of Operational Research, 152, 289-295. [15] Dobson, G., 1988. Sensitivity of the EOQ Model to Parameter Estimates. Operations Research,36(4) , 570–574. [16] Glasserman, P., Tayur, S.R., 1995. Sensitivity Analysis for Base-Stock Levels in Multi-Echelon Production-Inventory Systems. Management Science, 41, 263– 281.
32
[17] Hunter, J.J., 1986. Stationary distributions of perturbed Markov chains. Linear Algebra and its Applications, 82, 201–214. [18] Iglehart, D., 1963. Dynamic programming and stationary analysis of inventory problems, multistage inventory models and techniques. In : Scarf, H., Gilford, D., Shelly, M., editors, Multistage Inventory Models and Techniques, Stanford University Press, Stanford (Calif.), 1–31. [19] Kharna, A., Soman, C., Bandyopadhyay, T., 2014. Sensitivity analysis of the newsvendor model. European Journal of Operational Research, 239(2): 403-412. [20] Matheus, P., Gelders, L., 2000. The (R, Q) inventory policy subject to a compound Poisson demand pattern. International Journal of Production Economics, 68, 307–317. [21] Meyer, C.D., 1980. The condition of a finite Markov chain and perturbation bounds for the limiting probabilities. SIAM J. Algebraic Discrete Methods, 1(3): 273–283. [22] Mykytka, E.F, Ramberg, J.S., 1984. On the Sensitivity of the EOQ to Errors in the Forecast of Demand. IIE Transactions, 16(2), 144–151. [23] A. Paz, 1971. Introduction to Probabilistic Automata. Academic, New York. [24] Rabta, B., A¨ıssani, D., 2005. Strong Stability in an (R, s, S) Inventory Model. International Journal of Production Economics, 97(2), 159–171. [25] Rabta, B., A¨ıssani, D., 2008. Strong stability and perturbation bounds for discrete Markov chains. Linear Algebra and its Applications, 428(89), 19211927. [26] Ray, J., Chaudhuri, K.S., 1997. An EOQ model with stock–dependent demand, shortage, inflation and time discounting. International Journal of Production Economics, 53(2), 171–180.
33
[27] Ray, P.K., Sahu, S., 1992. Productivity measurement in multi-product manufacturing firms: Evaluation and control through sensitivity analysis. International Journal of Production Economics, 28(1), 71–84. [28] Scarf, H., 1960. The optimality of (S, s) policies in the dynamic inventory problem. In : Arrow,K.J., Karlin,S., Suppes, P., editors. Mathematical methods in the social sciences. Stanford University Press, Stanford (Calif.), 196–202. [29] Schweitzer, P., 1968. Perturbation theory and finite Markov chains. Journal of Applied Probability, 5, 401–413. [30] E. Seneta, 1984. Explicit forms for ergodicity coefficients and spectrum localization, Linear Algebra and its Applications, 60, 187–197. [31] E. Seneta, 1988. Perturbation of the stationary distribution measured by ergodicity coefficients. Advances in Applied Probability, 20, 288–230. [32] E. Seneta, 1991. Sensitivity analysis, ergodicity coefficients and rank-one updates for finite Markov chains. In : W.J. Stewart, ed., Numerical Solution of Markov Chains (Dekker, New York) pp. 121–129. [33] E. Seneta, 1993. Explicit forms for ergodicity coefficients of stochastic matrices. Linear Algebra and its Applications, 191, 245–252. [34] Sinha, B.K., Gupta, S.K., 1969. Sensitivity of Cost Functions in Inventory Models. AIIE Transactions, 1(3), 252–256. [35] Stoyan, D., 1983, Comparison Methods for Queues and Other Stochastic Models. English translation, D.J. Daley, Editor ( J. Wiley and Sons, New York). [36] Turanyi, T., Rabitz, H., 2000. Local Methods. In: Saltelli, A., Chan, K., Scott, E.M. (Eds.), Sensitivity Analysis, John Wiley & Sons, Chichester. [37] Veinott, A., Wagner, H., 1965, Computing optimal (s, S) policies, Management Sciences 11, 525–552.
34
[38] Wang, T., 2011. Sensitivity analysis for newsvendor model with general advertisement-sensitive demand. Procedia Engineering, 15, 1209-1213. [39] Yang, G. K., 2007, Note on sensitivity analysis of inventory model with partial backorders. European Journal of Operational Research, 177(2), 865-871. [40] Zheng, Y.S., 1992. On Properties of Stochastic Inventory Systems. Management Science, 38(1), 87–103.
35
Fig. 1. Deviation of the optimal policy parameters with respect to the perturbation of the arrival rate.
Fig. 2. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
1
Fig. 3. Deviation of the optimal inventory parameters with respect to the perturbation of the arrival rate.
Fig. 4. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
2
Fig. 5. Relative cost and service level deviation with respect to the perturbation of the arrival rate.
3