Managing inventory systems of slow-moving items

Author's Accepted Manuscript

Managing Inventory Systems of Slow-Moving Items G.J. Hahn, A. Leucht

www.elsevier.com/locate/ijpe

PII: DOI: Reference:

S0925-5273(15)00309-6 http://dx.doi.org/10.1016/j.ijpe.2015.08.014 PROECO6188

To appear in:

Int. J. Production Economics

Received date: 15 April 2014 Accepted date: 18 August 2015 Cite this article as: G.J. Hahn, A. Leucht, Managing Inventory Systems of SlowMoving Items, Int. J. Production Economics, http://dx.doi.org/10.1016/j. ijpe.2015.08.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.

Managing Inventory Systems of Slow-Moving Items G.J. Hahna , A. Leuchtb a

Professorship of Operations Management and Process Innovation, German Graduate School of Management and Law, Heilbronn/Germany b Institut f¨ ur Mathematische Stochastik, Professur f¨ ur Mathematische Stochastik mit Anwendungsbezug, Technische Universit¨ at Braunschweig/Germany

Abstract Slow-moving demand patterns frequently occur with spare parts as well as items in decentralized retail supply chains with large assortments. These patterns are commonly called lumpy since they exhibit comparably high demand variation and a high fraction of zero-demand events. In this paper, we examine two distribution-based approaches to model lumpy demand processes for inventory control: (i) a generalized hurdle negative binomial model, and (ii) a worst-case non-parametric model that is derived using a test-based approach. Considering a base stock inventory policy, we examine a set of lumpy time series from the industry to exemplify the suitability and benefit of the proposed approaches for managing inventory systems of slow-moving items. Keywords: Forecasting intermittent demand, inventory control, hurdle negative binomial distribution, Panjer recursion, non-parametric models

∗

Correspondence: Bildungscampus 2, 74076 Heilbronn/Germany Email addresses: [email protected] (G.J. Hahn), [email protected] (A. Leucht)

Preprint submitted to International Journal of Production Economics

August 20, 2015

1. Introduction and Research Background Spare parts and items in decentralized retail supply chains with large assortments frequently exhibit slow-moving demand patterns, i.e. rather low but highly fluctuating demand sizes and a high fraction of zero-demand periods. Extant metrics to characterize corresponding demand patterns include the average interval m between two consecutive positive demand events and the squared coefficient of variation of transaction sizes SCV, i.e., variance divided by squared mean both being calculated from positive demands only. Empirical evidence reports SCVs of around 0 to 1.7 and average demand intervals of around 3.8 to 24 periods on a monthly basis for the demand of spare parts at the Royal Air Force (Teunter et al., 2010). Snyder et al. (2012) analyze data from an US automotive company with an average SCV of 0.9 and an average demand interval of 2.5. According to the classification in Boylan and Syntetos (2008), the aforementioned time series are highly sporadic (m  1.34) and mostly show lumpy (SCV > 0.28) demand characteristics. Forecasting demand for inventory control of slow-moving items poses particular challenges to ensure service level compliance while keeping stock levels low (Ward, 1978; Snyder et al., 2012). Especially the inter-dependencies between demand forecasting and inventory control have received increasing attention in the literature more recently (e.g., see Teunter and Duncan, 2008; Tiacci and Saetta, 2009; Strijbosch et al., 2011; Rego and Mesquita, 2015). Three major approaches to model sporadic and lumpy demand processes are discussed in the literature (Willemain et al., 2004): (i) Croston’s extension of simple exponential smoothing (Croston, 1972), (ii) compound models with count or gamma-type distributions (e.g. Ward, 1978; Dunsmuir and Snyder, 1989), and (iii) various methods such as sampling and bootstrap that make direct use of the empirical data 2

observations (Snyder, 2002; Willemain et al., 2004). Our research contribution is twofold and adds to the second and third stream of literature from a methodological viewpoint as explained below. For this purpose, we adopt the approach of evaluating demand forecasting models from an inventory control perspective using a base stock inventory policy and a set of time series from the industry. As the first contribution, we investigate count distributions and apply a compound Bernoulli model in order to accommodate for additional zero-demand events. Extant compound Bernoulli approaches as described in Dunsmuir and Snyder (1989), Janssen et al. (1998), and Strijbosch et al. (2000) make use of the gamma and Erlang distributions. Count data models are frequently used in actuarial science (Klugman et al., 2012), but have rarely been used for modeling demand size distributions in inventory control even though they are most natural. Although Johnston and Boylan (1996) and Strijbosch et al. (2000) advocate the applicability of count models for slow-moving items, Snyder et al. (2012) is the only paper to the best of our knowledge that uses count models for this purpose. In their paper, they consider the Poisson, negative binomial, and hurdle Poisson distribution which corresponds to a compound Bernoulli model with positive demand following a Poisson distribution. Snyder et al. (2012) conclude that there is no superior model to be used. This actually supports our idea of deriving a generalized hurdle negative binomial model which encompasses the aforementioned distributions and allows modeling of even more over-dispersed demand processes since further distribution models are included in this approach. Second, instead of fitting a distribution model to the demand data first, the empirical distribution is frequently used to estimate the inventory control parameters directly. This approach can be problematic for lumpy demand

3

patterns with only comparably few positive demand observations and potentially lacking demand data for the right tail when data samples are small. Inspired by a recent paper of Klabjan et al. (2013), we propose a nonparametric approach to finding a worst-case demand distribution from the empirical data that a goodness-of-fit test would not reject. This worst-case distribution can be used to evaluate the risk of a downward deviation from the target service level and to determine more robust control parameters. By this, we pursue an alternative way in contrast to Ahmed et al. (2007) who take the demand distribution as given and resort to coherent risk measures to obtain robust results considering the risk aversion of the decision-maker. The remainder of this paper is organized as follows: in Section 2, we discuss the generalized hurdle negative binomial model and the non-parametric worst-case approach. The numerical study is conducted in Section 3 using exemplary time series from the industry. We conclude the paper in Section 4 and give an outlook for further research. Mathematical proofs and corresponding derivations are deferred to the Appendix. 2. Modeling Demand Processes of Slow-moving Items Throughout the paper, we assume demand D to be a sequence of independent and identically distributed (iid) random variables, i.e., we have a stationary process with no dependencies between the demand observations. 2.1. Generalized Hurdle Negative Binomial Model In the following, we derive a generalization of the hurdle Poisson model that is discussed in Snyder et al. (2012). For this purpose, we use a general class of count distributions that was first described in Panjer (1981) and is also referred to as the (a,b,0) class (Klugman

4

et al., 2012). The (a,b,0) class contains the binomial, Poisson, and negative binomial (NB) distributions with the pdf being defined as   b · pk−1 pk = a + k = 1, 2, . . . , (1) k  obviously requiring p0 = 1 − k≥1 pk . The parameters a and b are presumed to satisfy a < 0 and b = −q · a for some q ∈ N+ for the binomial distribution, a = 0 and b > 0 for the Poisson distribution, and 0 < a < 1 and b > −a for the NB distribution (Hess et al., 2002). In what follows, we allow for a zero-modification to independently model p0 using the parameter π. More precisely, we assume a compound Bernoulli model such that   b · pk−1 and pk = a + k = 2, 3, . . . , (2) p0 = π k  with p1 = 1 − π − k≥2 pk which is referred to as the (a,b,1) class in the literature (Klugman et al., 2012). In accordance with Snyder et al. (2012), we exclude the binomial model (and the degenerate case of all mass being accumulated in 0) from our considerations and thus restrict the parameter θ = (π, a, b) to the parameter space Θ := { θ = (π, a, b) | 0 ≤ π < 1, 0 ≤ a < 1, b > −2 · a} .

(3)

The parameter space in (a, b) is visualized in Figure 1 highlighting the different distributions that are included in this group of count distribution models (see Klugman et al., 2012): Poisson (a = 0, b > 0), negative binomial (NB: 0 < a < 1, b > −a), logarithmic (0 < a < 1, b = −a), and extended truncated negative binomial (ETNB: 0 < a < 1, −2 · a < b < −a). Including the ETNB model allows for additional modeling flexibility concerning overdispersed demand data compared to extant approaches as described in Snyder 5

Poisson

Negative Binomial

Geometric

Logarithmic

Figure 1: Parameter space in (a, b) of the generalized hurdle negative binomial model

et al. (2012). In summary, we obtain a generalized hurdle negative binomial (hNB) distribution. For the expectation, it holds that E [X] =

p1 + (a + b) · (1 − π) 1−a

where p1 has the closed form representation ⎧ ⎪ (1 − π) · ebb−1 if a = 0, b > 0, ⎪ ⎪ ⎨ b p1 = (1 − π) · log(1+b) if 0 < a < 1, b = −a, ⎪ ⎪ ⎪ a+b ⎩ (1 − π) · if 0 < a < 1, b > −2 · a, b = −a. (1−a)−(1+b/a) −1

(4)

(5)

Equation (5) reinforces the idea of a generalized hNB distribution that involves the hurdle Poisson (case 1) and the hurdle logarithmic model (case 2) as the two limiting cases of the hurdle (ET)NB model which 6

is characterized by case 3 (Klugman et al., 2012). The second moment is given by E [X] · (a + b + 1) E X2 = 1−a

(6)

(see Lemma 1 in the Appendix). Using Steiner’s formula leads to the variance V ar [X] =

E [X] · (a + b + 1) − (E[X])2 . 1−a

(7)

We can fit generalized hNB models to data adapting the so-called minimum quadratic distance estimator approach (see Luong and Garrido, 1993) for recursively defined distributions with finite support. The basic idea is comparable to least square estimation in regression. Even though it can be verified that this approach yields consistent parameter estimators, it shows poor performance in finite samples. Consequently, we propose generalized method-of-moment estimators or more precisely plug-in estimators as an alternative approach. Our method is based on empirical analogues of (4) and (6) as well as of p0 and p1 . Suppose that the true parameter lies in Θ and let pk denote the empirical relative frequency of k. Further, we denote the empirical counterparts of E[X] and   E[X 2 ] by M1 = n−1 nt=1 Xt and M2 = n−1 nt=1 Xt2 . Considering only practical relevant samples with at least two observations of different positive demands, we obtain the following plug-in estimators π  = p0 ,  a = max

M2 · (1 − π ) − M1 · (M1 + 1 − π  − p1 ) , 0 , M2 · (1 − π ) − M12 b = M2 · (1 −  a) − 1 −  a. M1

(8) (9)

(10)

The justification of these equations as well as the proof of consistency of the parameters is deferred to Lemma 2 in the Appendix. 7

2.2. Non-Parametric Robust Worst-Case Approach Modeling (intermittent) demand processes using comparably few observation from short time series bears the risk of misspecifying the demand size distribution that is used to derive inventory control parameters. The consequences of inappropriate control parameters are either excess inventories or poor service level compliance. Focusing on the issue of poor service level compliance which is more prevalent in spare parts management, we are interested in deriving a worst-case demand size distribution given incomplete information in order to evaluate the potential downward impact on inventory service and to determine more robust inventory control parameters. In the field of information theory, maximum entropy approaches aim at finding a probability distribution that maximizes uncertainty in the sense of Shannon (1948) subject to case-dependent constraints (see Golan et al., 1996, for an overview). Inspired by Klabjan et al. (2013), we pursue a related approach. Subject to certain inequality conditions, we derive a non-parametric worst-case model given the empirical frequencies. However, worst-case does not mean high uncertainty in the sense of Shannon (1948), but much probability mass in the right tail of the distribution which goes along with a high risk for poor service level compliance. For this purpose, we apply the Kolmogorov-Smirnov (KS) goodness-of-fit test to find a hypothetical model with the cumulative distribution function F0 that is not rejected for a given level of significance δ. max! F0−1 (0.5) K max 

s.t.

(11)

pk = 1

(12)

k=0

8

1 ≤ pk γ·n

k = 1, . . . , Kmax

1 − KS ( p0 , . . . , pKmax , F0 ) ≤ δ

(13)

(14)

To fit the non-parametric distribution F0 , we use the non-linear optimization model in (11) to (14). The probabilities pk are the decision variables on the domain {0, . . . , Kmax } for some Kmax specified below. We maximize the median F0−1 (0.5) in (11) to shift the probability mass towards the right tail. Obviously probabilities must sum up to one in (12). To ensure sufficient probability mass on empty bins and at the same time allowing for a sufficiently large support, we define a minimum probability depending on the sample size n and some scaling factor γ in (13). A higher scaling factor γ leads to a larger solution space for the decision model and thus corresponds to a setting with incomplete information resulting in different possible distribution models. In (14), we require the KS test to comply with the significance level δ; the function KS determines the corresponding pvalue given the empirical relative frequencies pk and the hypothetical model F0 . Initially, we set Kmax = max{k | pk > 0} and then gradually increase it as long as the solver finds a feasible solution to finally obtain the worst-case distribution model. 3. Numerical Study 3.1. General Approach and Datasets The numerical study is structured along the two contributions of the paper: first, we compare the generalized hNB model proposed above against the hurdle Poisson (hP) and the negative binomial (NB) model as the extant count models in the literature for inventory control when demand is 9

intermittent (see Snyder et al., 2012). For this purpose, we use ex-post inventory performance measures (achieved service level, physical stock) to compare suitability and accuracy of the different models. The hP and NB models as well as their corresponding method of moments estimators for the parameters are implemented according to Klugman et al. (2012). As a reference point, we include the ordinary empirical distribution (OE) in our investigations. Second, the robust worst-case (rWC) approach is illustrated by calculating expected service level (ex-ante measure) to quantify the actual risk of a downward deviation in service level when information quality of demand data is low. For both analyses, a simple base stock policy (review cycle r = 1) with order-up-to level (OUL) S, fill rate β, and a deterministic lead time L is assumed. Let Y and Z denote random variables of the demand during the lead time and during the entire risk period (lead time plus review cycle). S ∗ = min{S | E [Z − S]+ − E [Y − S]+ ≤ (1 − β) · E[D]}

(15)

An optimal OUL S ∗ can be found by minimizing S in (15) such that the expected shortage does not exceed 1 − β times the expected demand per period (Brown, 1959). For this purpose, we need to calculate the expected shortage both at the end and the beginning of the replenishment cycle. +

E[Z − S] =



(k − S) ·

pZk

= E [Z] − S −

k≥S+1

S 

(k − S) · pZk

(16)

k=0

The expression E[·]+ represents the first-order loss function which can be calculated for discrete distributions according to (16) and pZk denotes the (r+L)-fold convolution of the demand distribution for the risk period. Recall that for two random variables D1 , D2 taking their values in N0 , the 10

convolution can be calculated as P (D1 +D2 = k) =

k

l=0

P (D1 = l)·P (D2 =

k − l). Using (15), one can easily derive the expected fill rate for a given optimal OUL. Given the optimal OUL, we need to conduct a simulation of the inventory system to obtain ex-post results for the average achieved service level and the physical stock. The following sequence of events is assumed in each period t of the simulation: 1. Replenishment order is issued 2. Order placed in period t − L arrives 3. Backlogged demand is (partially) satisfied 4. Demand of period t is served from stock 5. Shortage quantity is backlogged 6. Service level for period t is determined 7. Physical inventory at the end of period t is determined In the numerical analysis, we examine three exemplary sets of industry time series with the characteristics summarized in Table 1. In general, the time series cover a broad range of fractions of zero-demand events (0.21 to 0.73) and also include higher levels of demand variation (SCVs up to 1.05). We calculated the empirical auto-correlation function with respect to time lag 1 which shows only minor time-dependency in the data (ranges from -0.22 to 0.26) and thus does not contradict the iid assumption. The first half of the observations of the time series are used to determine the order-up-to levels and the remaining observations are the input for the inventory simulation. We consider lead times of L = 1 and L = 3 periods and assume fill rates β in the range of 0.90 to 0.99 separated in steps of 0.03. In the rWC approach, we apply a scaling factor of γ = 1.5 which has been determined in advance using simulation and allows 11

No. of time series No. of observations

A

B

C

9

12

12

128

120

96

Fraction of

Min

0.21

0.47

0.26

zero-demands

Median

0.26

0.64

0.41

Max

0.38

0.73

0.61

Mean positive

Min

2.50

1.94

2.24

demand

Median

3.03

3.26

3.49

Max

5.34

4.06

6.24

Min

0.22

0.28

0.07

Median

0.26

0.54

0.29

SCV

ACF[1]

Max

0.89

1.05

0.44

Min

-0.21

-0.17

-0.22

Median

-0.02

-0.06

-0.02

0.10

0.18

0.26

Max

Table 1: Characteristics of the sets of time series

for a sufficiently large solution space. Furthermore, we consider two different levels of significance of δ = 0.01 and δ = 0.05. Calculations are conducted using the statistics package R in version 3.0.2 (R Foundation, 2014). For the non-parametric robust approach, we resort to the package dgof (Arnold and Emerson, 2011) that provides a variant of the KS test for testing with hypothesized discrete distributions and the non-linear solver of the package solnp by Ghalanos and Theussl (2012) that implements an augmented Lagrange multiplier method. 3.2. Results and Discussion Generalized hNB model. First, we examine achieved vs. target fill rate for the three datasets and the two lead time configurations (see Figures 2 to 4). Across all time series and parameter configurations, we observe that achieved fill rates on average overshoot target fill rates in a range of 3.3 to 6.3 percentage points (pp) for a target fill rate β = 0.90. The overshoot decreases for higher target fill rates to a range of 0.4 to 1.0 pp (β = 0.99) 12

13 Figure 2: Achieved vs. target fill rate for L = 1 and L = 3 (Dataset A)

14 Figure 3: Achieved vs. target fill rate for L = 1 and L = 3 (Dataset B)

15 Figure 4: Achieved vs. target fill rate for L = 1 and L = 3 (Dataset C)

and is generally lower for longer lead times. The overshoot can be explained by discrete and mostly single-digit OULs in our case which lead to higher overshoots when target fill rates are comparably low. Lower achieved fill rates for the case of longer lead times are the result of decreasing model accuracy when working with convolutions of single-period demand data. Comparing the different forecasting methods, we find that the NB model mostly shows the highest overshoot (on avg. 0.9 to 4.9 pp) while the OE and hNB models on the contrary exhibit lower overshoots (on avg. 0.7 to 4.4 pp). Next, we investigate the trade-off between achieved fill rate vs. average physical stock comparing the forecasting methods along the four different target fill rates and for the two different lead time configurations as defined above (see Figure 5). Average physical stock ranges between 4.9 to 10.5 and 7.0 to 14.0 items (dataset A), 6.6 to 15.0 and 8.6 to 18.0 items (B), and 5.2 to 11.5 and 7.5 to 15.2 items (C) for lead times equal to 1 and 3 periods, respectively. Recalling the aforementioned overshoot of the NB model, we also observe disproportionally higher average physical stock especially for higher target fill rates (on average plus 1.1 to 2.1 items). In most of the parameter configurations, the trade-off curves for the hP, hNB, and OE model closely match with a better performance of the hNB model and (in fewer cases) of OE model especially for lower target fill rates. Having observed similar demand process characteristics across the three datasets, we join them in order to compare subsets containing the lower and upper quartiles of the time series with respect to the SCV of demand (see Figure 6) and the fraction of zero-demand events (see Figure 7). Considering the two cases with low (SCV < 0.26) vs. high (SCV > 0.45) demand variation, we observe lower average physical stock for a low SCV since less variation induces lower safety stock and vice versa. Excluding the NB model,

16

17 Figure 5: Achieved fill rate vs. average physical stock for the three datasets

18

Figure 6: Achieved fill rate vs. average physical stock for low vs. high SCVs of demand

19

Figure 7: Achieved fill rate vs. average physical stock for low vs. high fractions of zero-demand events

we find that the trade-off curves closely match for demand processes with a low fraction of zero-demand events and low levels of demand variation. In contrast, we observe poor performance of the hP model for a high fraction of zero-demand events and higher levels of demand variation which is obviously due to limited modeling flexibility. On the other hand, the NB model shows reasonable performance in these settings. Summarizing the previous findings, we conclude the following: the OE approach can serve as a starting point for modeling slow-moving demand processes and shows good performance across a diverse set of parameter constellations as well as time series characteristics. The same is true for the hNB model which represents a generalization of the hP and NB models and thus resolves the aforementioned limitations and shortcomings of the two specific count distributions. This finding confirms our initial idea of using a generalized model that provides additional modeling flexibility and does not require selecting a specific count distribution model upfront. Non-parametric rWC approach. We apply the rWC approach in order to quantify the risk potential of a downward deviation in inventory service when the quality of demand data is low (see Figure 8). The optimal OUL of the hNB approach (and the equivalent expected fill rate) serve as the reference point to calculate the downward deviation. We find that the downward risk of the fill rate ranges between -10.2 to -1.5 and -22.2 to -3.2 pp (dataset A), -22.2 to -9.9 and -27.1 to -13.6 pp (B), and -14.0 to -7.3 and -17.2 to -9.8 pp (C) for lead times equal to 1 and 3 periods, respectively. Across the three datasets, we observe that the downward risk is higher when lead times are longer and target fill rates are lower. The reasoning behind these findings is similar to the case of achieved vs. target fill rate as discussed above. Concerning the different levels of significance, we make two observations:

20

21 Figure 8: Downward risk in inventory service for the three datasets

first, we find higher levels of downward deviations (except for one case) when the level of significance is high (δ = 0.05) and thus less strict when testing worst-case demand patterns. Consequently, the level of significance can be used to consider different levels of information quality in the rWC approach. Second, we observe different magnitudes for the offset between the two levels with a distinct range of 1.6 to 11.4 pp for dataset A in contrast to moderate ranges of 0.1 to 2.9 and 0.2 to 3.5 pp for datasets B and C. These results confirm that varying levels of information quality can have an impact on downward risk and thus should be considered. In general, the rWC approach can provide structured decision support to especially risk-averse decisionmakers that aim for more robust inventory control parameters when demand data is suspected to involve low information quality. 4. Conclusion and Outlook for Further Research In this paper, we examined methods to model demand processes of slowmoving items for inventory control. Basically, two novel distribution-based models were derived and evaluated with respect to their impact on inventory performance metrics. A comparative numerical study using exemplary industry data was conducted. The following insights were derived: the proposed generalized hurdle negative binomial model provides a unified approach that proved to be beneficial for modeling demand process across a wide range of parameter settings and time series characteristics. This is especially true for time series with higher levels of demand variation and higher fractions of zero-demand events. The ordinary empirical distribution also proved to be a eligible for inventory control of slow-moving items even for settings with shorter time series that are often found in practice. 22

The non-parametric worst-case model introduces a robust decision-making view on determining inventory control parameters. This novel approach allows addressing the issue of incomplete information and its impact on inventory performance which is of high relevance in practical settings especially when time series are short. The test-based procedure provides a structured way to identify the hidden risk potential in the right tail of the demand distribution which can be translated into additional required safety stock. With respect to further research, next steps would be to conduct further numerical analyses especially using larger sets of industry data and investigating other inventory policies. A more rigorous approach could be to examine inventory systems with stochastic lead times. Moreover, time-dependent approaches with rolling updates of control parameters could be considered to relax the restrictive assumption of iid demand processes. Appendix Lemma 1. Suppose that X admits a generalized hurdle negative binomial distribution, i.e. p0 = π, pk = (a + b/k) · pk−1 ,

k = 2, 3, . . . for some

θ = (π, a, b) ∈ Θ. Then equations (4) and (6) hold. Proof. Invoking the definition of expectation for discrete random variables and the Panjer recursion, we get E [X] = p1 + a ·

∞ 

(k + 1) · pk + b ·

k=1

∞ 

pk

k=1

= p1 + (a + b) · (1 − π) + a · E [X] =

p1 + (a + b) · (1 − π) , 1−a

23

(17)



E X

2



 b · pk (k + 1) · a + = p1 + k+1 k=1

= p1 + a · E X 2 + (2 · a + b) · E [X] + (a + b) · (1 − π) ∞ 

=



2

(18)

E [X] · (a + b + 1) . 1−a

Lemma 2. Suppose that X1 , . . . , Xn are independent, identically distributed observation admitting a generalized hurdle negative binomial distribution, i.e. p0 = π, pk = (a + b/k) · pk−1 ,

k = 2, 3, . . . for some θ = (π, a, b) ∈ Θ.

The empirical distribution is assumed to assign mass to at least two positive bins. Then, the plug-in estimators given in (8) to (10) are well-defined and consistent for θ, i.e. θ = ( π,  a, b) converges stochastically to θ. Remark: Note that (by the Glivento-Cantelli-Theorem) the assumption of the empirical distribution to assign positive mass to at least two possible bins is satisfied with probability tending to one and hence not restrictive. Proof. The parameter estimators are obtained by the plug-in principle and suitable expansions of p0 , p1 , E[X1 ], and E[X12 ]. As π < 1 implies E[X1 ] > 0, equation (18) yields E[X 2 ] · (1 − a) − 1 − a. b= E[X]

(19)

Plugging this relation into equation (17), straight-forward calculations imply a=

E[X 2 ] · (1 − π) − E[X] · {E[X] + 1 − π − p1 } . E[X 2 ] · (1 − π) − (E[X])2

(20)

To assure that the right-hand side is well-defined, we require that the denominator is strictly positive. Since 1 − π > 0, this holds true if E[X 2 ] 1 > 2 (E[X]) 1−π 24

(21)

which is equivalent to a+b+1 1 > p1 + (a + b) · (1 − π) 1−π

(22)

by (4) and (6). The latter finally reduces to p1 < 1 − π which holds true in view of (5) since we excluded the Bernoulli distribution. Hence, we obtain parameter estimators by substituting π, p1 , E[X], and E[X 2 ] by their empirical counterparts, p0 , p1 , M1 , and M2 . It remains to check that the resulting estimators are well-defined and lie in the parameter space Θ. Since the empirical distribution is assumed to assign mass to at least two positive bins, π  < 1. For the same reasoning  a < 1 if additionally M2 · (1 − π ) − (M1 )2 > 0. The latter follows from ) − (M1 )2 = M2 · (1 − π

∞ 

pk · pl · (k 2 − k · l) =

k,l=1



pk · pl · (k − l)2 > 0. (23)

k
Also by definition  a is non-negative. It remains to investigate the behavior of b. First, note that in view of M2 − M1 > 0, b > 0 if  a = 0. To verify b > −2 ·  a we invoke its equivalence to M2 · (1 −  a) > 1 −  a. M1

(24)

The latter is satisfied if M2 > M1 which trivially holds true in the present situation. Note that π ,  a, and b are continuous functions of p0 , p1 , M1 , and M2 which in turn converge almost surely to p0 , p1 , E[X], and E[X 2 ] by the strong law of large numbers. Combining these two facts finally yields consistency. 5. References Ahmed, S., C ¸ akmak, U., Shapiro, A., 2007. Coherent risk measures in inventory problems. European Journal of Operational Research 182 (1), 226– 238. 25

Arnold, T. B., Emerson, J. W., 2011. Nonparametric goodness-of-fit tests for discrete null distributions. The R Journal 3 (2), 34–39. Boylan, J. E., Syntetos, A. A., 2008. Forecasting for inventory management of service parts. In: Complex System Maintenance Handbook. Springer, pp. 479–506. Brown, R. G., 1959. Statistical forecasting for inventory control. McGrawHill. Croston, J. D., 1972. Forecasting and stock control for intermittent demands. Operational Research Quarterly 23, 289–303. Dunsmuir, W. T. M., Snyder, R. D., 1989. Control of inventories with intermittent demand. European Journal of Operational Research 40, 16–21. Ghalanos, A., Theussl, S., 2012. Rsolnp: General non-linear optimization using augmented Lagrange multiplier method. R package version 1.14, http://CRAN.R–project.org/package=Rsolnp. Golan, A., Judge, G. G., Miller, D., 1996. Maximum entropy econometrics: Robust estimation with limited data. Series in financial economics and quantitative analysis. Wiley, New York. Hess, K., Liewald, A., Schmidt, K., 2002. An extension of Panjer’s recursion. ASTIN Bulletin 32 (2), 283–297. Janssen, F., Heuts, R., de Kok, T., 1998. On the (R, s, Q) inventory model when demand is modelled as a compound Bernoulli process. European Journal of Operational Research 104, 423–436. Johnston, F. R., Boylan, J. E., 1996. Forecasting for items with intermittent demand. Journal of the Operational Research Society 47 (1), 113–121. 26

Klabjan, D., Simchi-Levi, D., Song, M., 2013. Robust stochastic lot-sizing by means of histograms. Production and Operations Management 22 (3), 691–710. Klugman, S. A., Panjer, H. H., Willmot, G. E., 2012. Loss models: From data to decisions, 4th Edition. Wiley, New York. Luong, A., Garrido, J., 1993. Minimum quadratic distance estimation for a parametric family of discrete distributions defined recursively. Australian Journal of Statistics 35, 59–67. Panjer, H. H., 1981. Recursive evaluation of a family of compound distributions. ASTIN Bulletin 12, 22–26. R Foundation, 2014. The R Project for Statistical Computing. URL: http://www.r-project.org/, last accessed: 15/04/2014. Rego, J. R. d., Mesquita, M. A. d., 2015. Demand forecasting and inventory control: A simulation study on automotive spare parts. International Journal of Production Economics 161, 1–16. Shannon, C. E., 1948. A mathematical theory of communication. Bell System Technical Journal 27, 379–423. Snyder, R. D., 2002. Forecasting sales of slow and fast moving inventories. European Journal of Operational Research 140, 684–699. Snyder, R. D., Ord, J. K., Beaumont, A., 2012. Forecasting the intermittent demand for slow-moving inventories: A modelling approach. International Journal of Forecasting 28, 485–496.

27

Strijbosch, L. W. G., Heuts, R. M. J., van der Schoot, E. H. M., 2000. A combined forecast-inventory control procedure for spare parts. Journal of the Operational Research Society 51 (10), 1184–1192. Strijbosch, L. W. G., Syntetos, A. A., Boylan, J. E., Janssen, E., 2011. On the interaction between forecasting and stock control: The case of non-stationary demand. International Journal of Production Economics 133 (1), 470–480. Teunter, R. H., Duncan, L., 2008. Forecasting intermittent demand: A comparative study. Journal of the Operational Research Society 60 (3), 321– 329. Teunter, R. H., Syntetos, A. A., Babai, M. Z., 2010. Determining orderup-to levels under periodic review for compound binomial (intermittent) demand. European Journal of Operational Research 203, 619–624. Tiacci, L., Saetta, S., 2009. An approach to evaluate the impact of interaction between demand forecasting method and stock control policy on the inventory system performances. International Journal of Production Economics 118 (1), 63–71. Ward, J. B., 1978. Determining reorder points when demand is lumpy. Management Science 24 (6), 623–632. Willemain, T. R., Smart, C. N., Schwarz, H. F., 2004. A new approach to forecasting intermittent demand for service parts inventories. International Journal of Forecasting 20, 375–387.

28