8644
Ind. Eng. Chem. Res. 2010, 49, 8644–8656
Control Variance Amplification in Linear Time Invariant Decentralized Supply Chains: A Minimum Variance Control Perspective Hua Xu, Gang Rong,* Yiping Feng, and Yucheng Wu State Key Laboratory of Industrial Control Technology, Institute of Cyber-systems and Control, Zhejiang UniVersity, Hangzhou 310027, People’s Republic of China
This paper addresses demand uncertainty and its propagation in supply chains. The supply chain is considered as a linear time invariant (LTI) system driven by stochastic customer demand. Under general ARMA demand patterns and arbitrary lead times, a unified and structured framework based on the classical minimum variance control theory is proposed for decentralized supply chain management (SCM). Optimal forecasting, the traditional order-up-to policy and the generalized order-up-to policy are directly derived according to the minimum variance criterion. Given these strategies, stochastic properties of the supply chain are studied using LTI system theory in both the time and the frequency domain. Findings from previous literature are reinterpreted from a control-theory-oriented perspective, and new characteristics of the generalized order-up-to policy are deduced and analyzed. On the basis of the statistical analysis, an optimization model is constructed to minimize the variable operation costs which are related to the parameters of the SCM strategies. 1. Introduction A supply chain is a system of business enterprises that link together to move a product or a service from suppliers to customers. It is characterized by a forward flow of products and a backward flow of information. The objective of supply chain management (SCM) is to satisfy customer demand through the most efficient use of resources owned by supply chain partners (i.e., suppliers, manufacturers, distributors, and retailers). SCM requires an integration of resources and information from original suppliers to end customers. Herein integration dose not imply that one member controls the whole supply chain.1 All partners make local decisions independently. They cooperate with each other to bring down the total costs of the whole supply chain. In most cases, a supply chain is a decentralized system. The ultimate goal of a supply chain is customer satisfaction. In an increasingly competitive environment, customer demand changes rapidly and drives operations of a supply chain to vary with time. Enterprises in a supply chain need to respond promptly to customer demand fluctuations. Managers should make decisions to balance the economical potentials and the negative impacts of these rapid demand fluctuations. Therefore, it is necessary to model and analyze the dynamic characteristics of a supply chain in relation to customer demand. Meanwhile, consideration of dynamic characteristics of a supply chain is helpful to fully capture its stochastic nature, which mainly stems from uncertainty of demand fluctuations, lead times, and sales forecasting, etc. A negative impact of these uncertainties is imbalance between supply and demand in a supply chain, which is known as the “bullwhip effect”: downstream supply chain members propagate and amplify demand fluctuations upstream, thereby imposing high capacity and inventory costs. To describe, analyze, and find the solution to the “bullwhip effect”, statistical properties of a supply chain dynamic model should be carefully investigated. Recently dynamic and statistical analysis of supply chain systems has attracted a lot of attention. A volume of research has been devoted to this topic. Different researchers have different methodologies to obtain SCM strategies. * To whom correspondence should be addressed. Tel.: 86-57187953145. Fax: 86-571-8795227. E-mail:
[email protected].
One methodology is performing dynamic analysis directly on the simulation of complex supply chain systems. Perea-Lo´pez et al.2 presented a framework for dynamically modeled decentralized supply chains in continuous time domain. They evaluated performances of several heuristic control laws. Pitty et al.3 proposed a hybrid refinery supply chain dynamic model, which was able to consider various supply chain activities and stochastic variations. Mele et al. adopted a discrete event-driven system to model the supply chain.4 Their work was later extended by using an agent-based supply chain simulator that can take demand uncertainty into account.5 Further, given the replenishment rules, they provided a simulation-based optimization strategy to obtain appropriate parameters of the rules.4,5 Raj and Lakshminarayanan proposed a simulation-based troubleshooting framework for supply chain performance enhancement.6 When it comes to a supply chain with several business goals, they used multiobjective optimization to derive a desired solution on the basis of their framework.7 Recently they suggested an entropy-based optimization method to reduce the unpredictable nature of decentralized supply chain systems.8 The above simulation-based methodology can simulate the complex dynamical phenomenon of a supply chain. It provides decision tools for managers to simulate, evaluate, and improve various heuristic management rules. Another way to design dynamic SCM strategies is adopting the model predictive control (MPC) scheme. MPC is appreciated for its ability to deal with constraints and nonlinearities of multivariable dynamic systems. Therefore, it is natural to extend the methodology from process control to supply chain management for multiechelon supply chains.9-11 We refer readers to the work by Sarimveis et al. for a review of MPC application in SCM.12 MPC is also a flexible framework. Perea-Lo´pez et al.13 employed a MPC framework to deal with multiproduct, multiechelon supply chain networks. This framework contains a discrete mixed integer linear programming (MILP) model. Mestan et al.14 presented a MPC configuration with a mixed logical dynamical (MLD) system model. Puigjaner et al.15 developed a MPC scheme which permits integration of hierarchical decision levels (i.e., strategic, tactical, and operational). Li and Marlin16 proposed a MPC scheme involving a stochastic optimization formulation to cope with uncertainty in customer
10.1021/ie100167u 2010 American Chemical Society Published on Web 08/13/2010
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
demand and lead times. So far the MPC methodology has been prevalent in dynamic SCM. Both the simulation-based and MPC-based methods may belong to the “design via synthesis” methodology. Here “synthesis” is defined as structuring all factors that affect supply chain performances into a simulation or optimization model, so as to directly obtain optimal performance from the merged model. This methodology mainly addresses the optimization aspect in supply chain management. But it only offers numeric solutions which depend on specific data and settings. The implications of the optimal strategy and the impacts of its parameters can only be discovered through a large number of simulation trials. It may not contain enough analytical supports to investigate the characteristics of a supply chain and provide general guidelines to improve the supply chain performance. Motivated by quantifying the bullwhip effect, significant attention has been focused on more analytical investigation in supply chain dynamics and uncertainties, which provides a foundation for replenishment rules improvement. Such methodology may be referred to as “design via analysis”, which we will follow in this work. One approach of this methodology is based on time series analysis theory. Chen et al. addressed the role of demand forecasting on the bullwhip effect.17,18 They used an first-order autoregressive (AR(1)) process to describe the customer demand and studied the impact of moving average (MA) and exponential smoothing (ES) demand forecasting methods. Zhang derived a bullwhip effect measure (i.e., variance amplification ratio) using the minimum mean square error (MMSE) estimation method.19 Gilbert generalized the results of the previous research.20 They presented a multistage supply chain model that was based on autoregressive integrated moving average (ARIMA) models. They also derived ARIMA models of inventory and orders of the supply chain. All of the above analysis was performed on the assumption that the replenishment rule was the order-up-to (OUT) policy. Balakrishnan et al.21 modified the OUT policy and proposed MA and ES order smoothing policies to reduce order-size variability, thus reducing overall system costs. But their strategy was restricted to the assumption that customer demand was independently and identically distributed (iid) over time and the lead time was unity. Another analytical approach is based on the classical discrete control theory. Dejonckheere et al.22 used block diagram and transfer functions to model the replenishment rules and the forecasting method (moving average method, exponential smoothing method, and demand signal processing). On the basis of transfer function models, they adopted a frequency response method to link the replenishment rules and the bullwhip effect measure. Disney et al.23 summarized the relationships between the bullwhip effect measure, noise bandwidth, the sum of squares of the system’s impulse, and the system transfer function. For a single echelon supply chain, Disney et al.24 presented a general model of replenishment rules to reduce the bullwhip effect. On the basis of this model, they provided analytical expressions of the bullwhip effect measure and the inventory variance under an iid demand pattern. Also, they proposed “the golden ratio” to balance the inventory variance and the bullwhip effect. Lin et al.25 analyzed a supply chain model using the z transform technique. They used the stability analysis to explain the cause of bullwhip effect. They proposed a PI (proportional and integral) controller and cascade control structures to eliminate the bullwhip effect. Quyang and Daganzo26 presented a set of linear time invariant policies in multistage decentralized supply chains. They derived robust
8645
analytical conditions to predict the presence of the bullwhip effect. Their work was mainly built on the frequency domain. The above control theory based approach is well-suited for system analysis and replenishment rules design. But it may be more appropriate for deterministic signal analysis. It is a less straightforward approach to deal with stochastic demand, especially general ARMA demand patterns. As a sequel of the analytical works above, several researchers have adopted a combination of the time series analysis and the discrete control theory. They considered a replenishment rule as a controller to deal with stochastic signals and analyzed the performance of the controller using the time series analysis approach. Hosoda and Disney27 analyzed a three echelon supply chain with an OUT policy and MMSE forecasting. They presented closed-form expressions of the inventory variance and the orders variance under the AR(1) demand patterns. They exploited the contour integral approach28 to derive these variance expressions. Lin et al.29 presented a predictive control approach for a decentralized supply chain unit under the general ARIMA demand patterns. A replenishment rule with an explicit expression was directly derived from a given criterion. But their strategy was in essence a MPC configuration,7 and the forecaster might not be in a MMSE sense. Also, the expression of the replenishment rule was quite complex for further theoretical investigation. Hosoda and Disney30 proposed the “generalized OUT policy”, which could be considered as a proportional controller incorporated into the OUT policy. Assuming demand followed the AR(1) process, they combined the generalized OUT policy with MMSE forecasting and derived the exact formulas for the order and the inventory variance in a two-echelon supply chain. Chen and Disney31 integrated the ARMA(1,1) demand process, the MMSE forecasting, and the generalized OUT policy in a single-echelon supply chain model with a unit lead time. They pointed out that bullwhip effect avoidance always came at the cost of holding extra inventory. Boute et al.32 analyzed the impact of the generalized OUT policy on the performance of a two-echelon supply chain under an iid demand pattern. They also numerically studied the interaction between replenishment order sizes and lead times. Our work is an extension of the previous studies that are constructed on the discrete control theory and time series analysis. In previous studies, we can see that SCM strategies (the OUT policy, the generalized OUT policy, and the MMSE forecasting method) are determined in advance and investigated in specific settings (i.e., white noise, AR(1), or ARMA(1,1)). There is a lack of a more general and systematic framework that considers multiechelon supply chains, general ARMA demand patterns, and arbitrary lead times. To address the issue, we consider a supply chain as a linear time invariant (LTI) dynamic system whose inputs are stochastic processes. We propose a framework, which is on the basis of the welldeveloped minimum variance control theory, for decentralized SCM strategies design and analysis. The proposed framework consists of three parts (we follow the methodology presented by Åstro¨m33): (1) Optimal replenishment rules design: given a dynamic model and a criterion, find the control law which minimizes the criterion. (2) Analysis: what are the statistical properties of the system variables? (3) Parametric optimization: given a strategy with known structure but unknown parameters, how are the parameters to be adjusted?
8646
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
) 0 and Ot ) 0 are continuously maintained, the system remains in an equilibrium state such that the inventory keeps its target value T g 0. T is a prespecified value. We also assume that the system is in its equilibrium state for all t e 0 andI0 ) T. Define NSt ) It - T. We can obtain the following equation from eq 4: NSt ) NSt-1 + Ot-L - Dt
Figure 1. Supply chain model.
Moreover, the new framework may contribute a different understanding to previous findings. It also provides answers to some important questions: (1) Why should we choose the match between the OUT policy and the MMSE forecasting? Johnson and Thompson34 have proved that the match produced the minimum inventory variance for a single supply chain with a zero lead time. Is it still a myopic policy for arbitrary lead times and general ARMA demand patterns? More fundamentally, what are the relationships between the business forecasting and the replenishment rules? (2) So far the generalized OUT policy is only a heuristic strategy for bullwhip effect reduction. What are the theoretical characteristics of the policy under general ARMA demand patterns? (3) The generalized OUT policy offers more freedom to manipulate the dynamics and stochastic properties of the supply chain. How can we make use of this flexibility? 2. Model and Assumptions Consider a supply chain model26,27,29 depicted in Figure 1. It is a homogeneous chain. All suppliers are alike. In an arbitrary node of the chain, the material balance can be represented by the following difference equation: It ) It-1 + YSt - YCt
(1)
YSt
where is the amount of products delivered from its supplier in time period t, YCt is the amount of products delivered to its customer in time period t, and It is the inventory of period t. The relationship of information-material exchange is described as follows: YCt ) Qt
(2)
YSt ) Pt-L
(3)
where Pt is the order placed to its supplier (for a manufacturer, it represents production order quantity) in time period t, Qt is the customer demand (i.e., the order placed by its customer) in time period t, and L is the lead time (i.e., the period of time from the moment a customer places an order to the moment the customer receives the order). Here we assumeL g 1. The order placed by the downstream unit is the customer demand of the upstream unit (i.e., Pti-1 ) Qti). Substituting eqs 2 and 3 into eq 1, we obtain the following equation: It ) It-1 + Pt-L - Qt
(4)
There are obvious analogies between the above model and a process control system,10,11 which is depicted in Figure 2. Herein we assume the customer demand Qt and the orders Pt are both stationary Gaussian processes and E[Qt] ) E[Pt] ) µ.We can suppose that Dt and Ot are deviations from expected values of Qt and Pt, respectively. It means that if the conditionsDt
(5)
Note that Ot and Dt are all stationary Gaussian processes with expected values of zero. Meanwhile, Ot is selected to ensure that NSt is also a stationary Gaussian process with the expected value of zero. 3. Operation Cost Structure We divide the operation costs into two parts: fixed costs and variable costs. Fixed costs are related to constants T and µ, which are determined in advance. Here we focus on the variable costs. 3.1. Expected Inventory Cost. Since NSt is a stationary stochastic process, its expected value and variance do not change with time. Given its mean (i.e., 0) and the standard deviation (i.e., σNS), we can calculate the expected inventory cost per period. If NSt > 0 (i.e., It > T), a holding costh per unit is incurred. If NSt < 0 (i.e., It < T), a backlog costb per unit shortage is incurred. The expected inventory cost can be expressed as follows: e-x /2σNS (x) dx + b 0 σNS√2π (h + b)σNS
E[CINV] ) h )
2
∫
+∞
2
e-x /2σNS (-x) dx -∞ σNS√2π
∫
2
0
2
√2π
(6)
Note that if a system is in its equilibrium state (demand and order rates are constants), σNS ) zero; thus, we have no inventory variable costs. 3.2. Expected Production Cost. The manufacturer incurs a variable cost of V for every unit of additional order Ot (Ot > 0). If Ot exceeds the production capacity limit K, an external service capacity is needed and a unit is produced at cost g (g > V). If Ot < 0, no variable costs are incurred. The production cost structure can be seen in Figure 3. Given the standard deviation σO of the order series, we can calculate the expected production costs as follows: E[Corder] ) V )
∫
K
0
σO
e-x /2σO (x) dx + g σO√2π 2
2
∫
∞
K
e-x /2σO (x) dx σO√2π 2
2
[V + (g - V)e-K /2σO ] 2
2
√2π σO V (K . σO) ≈ √2π
(7)
From eqs 6 and 7, we can see that both variable costs are approximately proportional to the standard deviations. Thus, we can reduce the variable costs directly by controlling the variance of inventory and orders. 4. Minimum Variance Control Strategies 4.1. Preliminary Strategy. Herein we introduce a shift operatorz:
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
8647
Figure 2. Fluid representation of supply chain dynamics.
We can also rewrite A(z-1) and C(z-1) for convenience: A(z-1) ) 1 + a1z-1 + a2z-2 + ... + anz-n C(z-1) ) 1 + c1z-1 + c2z-2 + ... + cnz-n n ) max(p, q)
Suppose that we have observed all the customer demand up to period t (Ht ) {Dt,Dt-1,Dt-2,...}). We want to find an optimal k-step-ahead prediction for future demand. On the basis of eq 11, the future demand Dt+k can be expressed as follows:
Figure 3. Production cost structure.
z-1Xt ) Xt-1 zXt ) Xt+1
-1
-1
(1 - z )Xt ) Xt - Xt-1 1 Xt ) 1 - z-1
(8)
+∞
∑
Xt-h )
h)0
∑
) λFk(z-1)et+k + λ
Xh
h)-∞
∑
h)-∞
∑
(9)
Dh
h)-∞
t
∑O
t+L
h
h)-∞ t-1
) NSt +
∑
-
∑D
t-L
h
-(
h)-∞
∑O
t
h
-
h)-∞
∑ D) h
The k-step-ahead minimum mean squares error (MMSE) estimator is given by
h)-∞
t+L
Oh + Ot -
h)t-L+1
∑D
ˆ t+k/H ) λ D t
h
h)t+1
(10) t+L Define the lead-time demand DLt ) ∑h)t+1 Dh as the total demand in the future periods [t + 1, t + L]. We aim to divide it into two parts: predictable and unpredictable. Note that Dt is a general stationary stochastic process. It can be equivalently represented by an ARMA(p,q) model:
Dt ) λ
where the polynomials Fk(z-1) and Gk(z-1) satisfy the Diophantine equation: C(z-1) ) A(z-1) Fk(z-1) + z-kGk(z-1)
Then we can obtain the following equation: NSt+L ) NSt +
A(z-1)
(13) et
Gk(z-1) ) g0 + g1z-1 + g2z-2 + ... + gn-1z-(n-1)
t
Oh -
Gk(z-1)
Fk(z-1) ) 1 + f1z-1 + f2z-2 + ... + fk-1z-(k-1)
1 1 O Dt -1 t-L 1-z 1 - z-1 t-L
)
Dt+k ) λ C(z ) et+k A(z-1)
t
Thus, eq 5 can be reformulated as follows: NSt )
(12)
C(z-1) et ) λH0(z-1)et A(z-1)
(11)
A(z-1) ) 1 + a1z-1 + a2z-2 + ... + apz-p C(z-1) ) 1 + c1z-1 + c2z-2 + ... + cqz-q
where {et} is white noise and et ∈ N(0,1). λ is a constant.
Gk(z-1) A(z-1)
et )
Gk(z-1) C(z-1)
Dt
(14)
The forecast error is given by ˜ t+k/H ) λFk(z-1)et+k D t
(15)
ˆ t+k/H only depends on et, et-1, et-2, ..., while Note that D t ˜ Dt+k/Ht depends on et+k, et+k-1, ..., et+1. Thus, we have the following equation: ˆ t+k/H , D ˜ t+k/H ) ) 0 Cov(D t t
(16)
Equation 16 indicates that the optimal predictor is uncorrelated with its prediction error. We can divide the future demand Dt+k into two uncorrelated parts:
8648
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
ˆ t+k/H + D ˜ t+k/H Dt+k ) D t t
(17)
t+L ˆ tL ) ∑h)t+1 ˆ h/H . Define the lead-time MMSE estimator D D t L t+L ˜ h/H . ˜ Define the lead-time prediction error Dt ) ∑h)t+1 D t From eq 17, we can easily obtain the following equations:
connection between the optimal forecasting and the OUT policy. This conclusion is essentially a version of the well-known separation theorem in stochastic control theory. Substituting eq 24 into eq 10, we obtain the formula for inventory under the OUT policy:
ˆ Lt + D ˜ Lt DLt ) D
(18)
ˆ Lt - DLt ) -D ˜ Lt NSt+L ) D
ˆ Lt , D ˜ Lt ) ) 0 Cov(D
(19)
Combining eq 5 and eq 25, we can rewrite the OUT policy as follows:
4.2. Minimum Variance Control. Now we aim at controlling the inventory variance. Consider the following scenario. At the end of the period t, we have measured the inventory NSt. We have known all the past orders Ot-1, Ot-2, ... and customer demand Dt, Dt-1, Dt-2, .... We want to determine the order quantity Ot so that the variance of the inventory is as small as possible. We directly select the mean square value of the inventory as the performance index: J ) E(NSt+L2) ) σNS2
(20)
Substituting eq 10 into eq 20, we obtain the following equation:
Ot ) NSt+L - NSt+L-1 + Dt+L L L ˆ Lt - DLt - (D ˆ t-1 )D - Dt-1 ) + Dt+L L L ˆ ˆ ) Dt - Dt-1 + Dt
∑
J*1 ) min E(NSt+L2 + r02Ot2) Ot
ˆ Lt ) - D ˜ Lt ]2} Oh + Ot - D
(21)
h)t-L+1
ˆL Note that (NSt + ∑t-1 h)t-L+1 Oh + Ot - Dt ) only depends onet, ˜ Lt depends on et+L, et+L-1, ..., et+1. Thus, we et-1, et-2, ..., whileD have the decomposition: t-1
J ) E[(NSt +
∑
(22) According to the lemma proved by Åstro¨m33 (in the Appendix), we choose the optimal order quantityOt, which satisfies the following function:
t-1
∑
where r0 is a penalty factor and r0 g 0. Note that we can select a more general performance index to design more complex replenishment rules (i.e., a linear quadratic Gaussian regulator):35 J*2 ) min E{[P(z-1)NSt+L]2 + [R(z-1)Ot]2} Ot
ˆ Lt )2] + E[(D ˆ Lt )2] Oh + O t - D
h)t-L+1
(23) ˜ Lt )2]. The minimum value of eq 23 is J* ) E[(D The optimal order quantity can be represented by the following equation:
Here we choose the simplest one (i.e., eq 27). According to the criterion, we can obtain the following equation:
∑
Oh)
(24)
h)t-L+1
Here NSt is the on-hand stock (i.e., the stock immediately t-1 Oh available for satisfying the customer demand) and ∑h)t-L+1 is the on-order stock (i.e., stocks requisitioned but not yet t-1 received). Define the inventory position IPt ) NSt + ∑h)t-L+1 L ˆ Oh. Equation 24 becomes Ot ) Dt - IPt, which is the orderup-to (OUT) policy. It means that Ot is ordered in such a way that it raises the current inventory position IPt to the target level ˆ Lt , which is an estimate of the future customer demand. D From the above derivation, we can conclude that the OUT policy minimizes the inventory variance only if the MMSE estimator is adopted to predict the future demand. It provides a
(29)
where β ) 1/(1 + r02) and 0 < β e 1. The detail derivation can be seen in the Appendix. Equation 29 is also called the generalized order-up-to policy which was intuitively proposed by previous researchers for reducing the bullwhip effect. Here we derive it in a more rigorous way. Note that when β ) 1, eq 29 becomes the OUT policy. We can rewrite eq 29 as follows:
t-1
ˆ Lt - (NSt + Ot ) D
(28)
R(z-1) ) r0 + r1z-1 + r2z-2 + ... + rnrz-nr
ˆ Lt - IPt) Ot ) β(D
E(NSt+L2) ) E{min(NSt+L2)} J* ) min Ot Ot Ot
(27)
P(z-1) ) 1 + p1z-1 + p2z-2 + ... + pnpz-np
ˆ Lt )2] + E[(D ˜ Lt )2] Oh + Ot - D
h)t-L+1
) E[min(NSt +
(26)
4.3. Minimum Variance Control Where the Order Variance Is Restricted. The minimum variance controller (i.e., the OUT policy) previously discussed is designed to produce the minimum variance of the inventory. It is assumed there is no restriction in the variance of orders. It may be not acceptable by a manufacturer with high variable production costs. Thus, we choose the following objective to balance the inventory variance and the order variance.
t-1
J ) E{[(NSt +
(25)
Ot )
β L ˆ Lt - D ˆ t-1 (D + D t) 1 - (1 - β)z-1
(30)
The inventory under the generalized OUT policy can be expressed as follows: ˆ Lt - DLt NSt+L ) D
1-β L ˆ t-1 ˆ Lt - D + Dt) (D 1 - (1 - β)z-1 (31)
5. Statistical Properties of the Inventory and the Orders In this part we want to analyze the behavior of the supply chain when the replenishment rules are given. We consider the supply chain as a linear time invariant stochastic control system.
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
We introduce some basic concepts of the LTI system theory (i.e., impulse response, transfer function, and the signal spectrum, etc.) to describe its stochastic properties. From eqs 30 and 31, we can see that the variance of the inventory and the orders is closely related to three sources. The first term is the ˆ Lt - DLt over L periods. Note that the variance prediction error D of the prediction error is also the minimum inventory variance L ˆ t-1 ˆ Lt - D + (i.e., eq 25). The second term is the OUT policy D Dt, which is determined by the customer demand. The third is a “coordinator” (i.e., (1 - β)/(1 - (1 - β)z-1) or β/(1 - (1 β)z-1)), which is a smooth filter to adjust and balance the inventory and the order variances through different values of ˆ Lt - DLt depends on et+L, et+L-1, the parameter β. Also note that D ˆ Lt - D ˆ Lt-1 + Dt) depends ..., et+1 and [(1 - β)/(1 - (1 - β)z-1)](D on et, et-1, et-2, .... Thus, we can obtain the following equation: Var(NSt+L) )
+ 1-β L ˆ Lt - D ˆ t-1 Var (D + Dt) -1 1 - (1 - β)z ˆ Lt Var(D
[
DLt )
]
(32)
Equation 32 indicates that the inventory variance can be decomposed into two uncorrelated parts. One is the variance of the prediction error. The other is the variance of an additional signal introduced by the generalized OUT policy. 5.1. Prediction Error. According to eq 15, the prediction error can be represented in detail as follows: ˜ Lt ) λ[F1(z-1)et+1 + F2(z-1)et+2 + ... + FL(z-1)et+L] D ) λ[et+L + (1 + ψ1)et+L-1 + (1 + ψ1 + ψ2)et+L-2 + ... + (1 + ψ1 + ψ2 + ... + ψL-1)et+1]
(33) where
-1
8649
-1
C(z ), A(z ), and L. But we can still explore some general properties related to the bullwhip effect. ∞ can be From the viewpoint of LTI system theory, {ψi}i)0 considered as the impulse response of the dynamic system C(z-1)/A(z-1). Thus, the time series Dt can be written in an impulse response representation: ∞
Dt ) λ
C(z-1) et ) λ( ψiet-i) A(z-1) i)0
∑
(35)
where ψ0 ) 1. From eq 15 and eq 33 we can obtain the equation: ˜ t+k/H ) λFk(z-1)et+k ) λ(et+k + ψ1et+k-1 + D t ψ2et+k-2 + ... + ψk-1et+1)
(36)
According to eq 17 and eq 35, eq 14 can be rewritten as follows: ˆ t+k/H ) λ D t
Gk(z-1) A(z-1)
et ) λ(ψket + ψk+1et-1 + ... + ψk+jet-j + ...)
(37)
Then we can derive the impulse response formulation of the OUT policy, which was first presented by Gilbert,20 from the viewpoints of the minimum variance control theory: L ˆ Lt - D ˆ t-1 Ot ) D + Dt ) λ[(1 + ψ1 + ψ2 + ... + ψL)et + ψL+1et-1+ ψL+2et-2 + ψL+3et-3 + ...]
(38) Further, we derive the bullwhip effect measure M1:
-1
F1(z ) ) 1 -1
F2(z ) ) 1 + ψ1z
2 Var(Ot) - Var(Dt) ) λ (1 + ψ1 + ψ2 + ... + ψL) λ2(1 + ψ12 + ψ22 + ... + ψL2)
-1
F3(z-1) ) 1 + ψ1z-1 + ψ2z-2 l FL(z-1) ) 1 + ψ1z-1 + ψ2z-2 + ... + ψL-1z-(L-1)
L-1
) 2λ2(
∑
L-2
ψjψj+1 +
j)0 1
The variance of the prediction error can be calculated as follows:
∑ψψ j
∑ψψ j
j+2
+ ... +
j)0 0
j+L-1
+
j)0
∑ψψ j
j+L)
j)0
(39) ˜ Lt ) ) λ2[1 + (1 + ψ1)2 + ... + Var(D (1 + ψ1 + ψ2 + ... + ψL-1)2]
(34)
Equation 34 indicates that the longer lead time L will cause the larger prediction error, thus the larger inventory variance. 5.2. OUT Policy and the Bullwhip Effect: Time Domain Analysis. Herein we want to investigate the bullwhip effect of the OUT policy (eq 25). We use the difference between the order variance and the demand variance M1 ) Var(Ot) Var(Dt) to measure the bullwhip effect instead of the variance amplification ratio M2 ) Var(Ot)/Var(Dt), which is proposed by the previous researchers. The two measures are essentially the same with each other. But the measure M1 is more appropriate for analysis in the current work. If M1 > 0, the bullwhip effect exists. It means that the variance of the orders (retailer demand) is greater than the variance of the customer demand. The upstream supplier (a distributor or a manufacturer) observes more unstable demand and thus has to pay more in the variable costs. Note that the concrete expressions of Ot, Dt, and their variances depend on the specific structure and parameters of
From eq 39, we can see that if all ψi > 0 i ∈ [1,+∞), the bullwhip effect exists. The lead time is longer; the bullwhip effect is more obvious. We can also infer that if every member in a supply chain adopts the OUT policy, the variance amplification is unavoidable. However, if not all ψi are greater than zero, eq 39 may be less than zero and there is no bullwhip effect. Consider an AR(1) demand process Dt ) FDt-1 + et. If F e 0,M1 is equal or less than zero. In summary, the OUT policy does not necessarily lead to the emergence of the bullwhip effect. It depends on the specific patterns of the customer demand. Also note the definition of the autocovariance of a stationary stochastic process:36 ∞
γk ) Cov(Dt, Dt+k) ) λ2
∑ψψ j
j+k
(40)
j)0
Then eq 39 can be approximately written as follows: Var(Ot) - Var(Dt) ≈ 2λ2(γ1 + γ2 + ... + γL)
(41)
8650
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
Note that the autocovariance can be easily obtained from a realization of a stochastic process. Thus, eq 41 may be convenient for use if we cannot obtain the expressions of the polynomials C(z-1) and A(z-1) exactly. 5.3. OUT Policy and the Bullwhip Effect: Frequency Domain Analysis. Now we represent the OUT policy in a more compact form: L ˆ t-1 ˆ Lt - D + Dt ) Ot ) D
∑
Gk(z-1) + C(z-1)
L
et ) λH1(z-1)et
k)1
λ
A(z-1)
(42)
We characterize the stationary process Ot, Dt (eq 11 and eq 42) in the frequency domain and graphically display their second-order statistical properties using the signal spectrum. We refer readers to the work by Åstro¨m33 or Ljung36 for a detailed discussion about the signal spectrum analysis for LTI systems. Define the signal spectrum as follows: k)∞
Φ(ω) ) λ2
∑
γke-ikω ) 2λ2(γ0 + 2
k)-∞
∞
∑γ
k
cos(kω))
k)1
0eωeπ
(43)
We can see that the peak value of Φ(ω) is Φ(0) if all γi > 0. Such customer demand pattern will lead to the bullwhip effect. According to the theorem of transformation of spectra by linear systems, we can equivalently write eq 11 and eq 42 from the viewpoints of the signal spectrum: ΦD(ω) ) λ2 |H0(e-jω)| 2Φe(ω)
ΦO(ω) ) λ2 |H1(e-jω)| 2Φe(ω)
(44)
where ΦD(ω) and ΦO(ω) are the spectra of the time seriesDt and Ot, respectively. Φe(ω) is the spectrum of the white noise et and Φe(ω) ) 1. The variance of the orders Ot can be calculated through the following frequency domain formula: Var(Ot) )
λ2 π
∫
π
0
ΦO(ω) dω )
λ2 π
∫
π
0
|H1(e-jω)|2 dω
(45) We can also obtain the spectrum representation of eq 39: Var(Ot) - Var(Dt) )
λ2 π
∫
π
0
Here the lead time is set to 1. Figure 4 shows the spectra of the demand and the order process (under the OUT policy with MMSE forecasting). 5.4. Generalized OUT Policy. The order quantityOt under the generalized OUT policy can be expressed in the following equation: Ot )
L
(1 - z-1)
Dt - 0.8Dt-1 + 0.25Dt-2 ) et - 0.01et-2
(ΦO(ω) - ΦD(ω)) dω
[(1 - z-1) λ
β (1 - (1 - β)z-1)
-1
k
) + C(z-1)]
k)1
et ) A(z-1) λG1(z-1) H1(z-1)et (47)
where G1(z-1) ) β/(1 - (1 - β)z-1). Now let us draw the frequency domain plot of |G1(e-jω)|2with β ) 0.6. From Figure 5 we can see that G1(z-1) in the generalized OUT policy is a low-pass filter. It passes low-frequency signals but attenuates (reduces the amplitude of) signals with high frequencies. Smaller is the value β; more “noise” (i.e., highfrequency signals) is filtered out. For example, consider a singleechelon supply chain unit with a lead time 0. Assume the customer demand is a deterministic mixed sine wave Dt ) 2 sin(0.01πt) + sin(0.9πt). We use the generalized OUT policy with β ) 0.2. Note that eq 30 reduces to Ot ) [β/(1 - (1 β)z-1)]Dt. Figure 6 shows the impacts of the low-pass filterG1(z-1). We see that the filter almost eliminates the high-frequency signal sin(0.9πt), while retaining the low-frequency signal 2 sin(0.01πt). In the supply chain, such a phenomenon indicates that the generalized OUT policy makes the orders trace the slow stable customer demand changes which exhibits market trends. Meanwhile it rejects the rapid fluctuations and considers them as irregular. Note that when β ) 1 (i.e., the OUT policy), such rapid fluctuations are all transmitted to the upstream supplier, thus posing a heavy burden on the inventory and production of the upstream member. Although the generalized OUT policy can reduce the adverse of the bullwhip effect for the whole supply chain, it sacrifices the local interests. The high-frequency part of the customer demand is absorbed by the local inventory, thus amplifying the local inventory variance. According to eq 32, we can see that the inventory variance increases due to the additional signals introduced by the generalized OUT policy. Analogue to eq 47, the signal can be represented in the following equation:
(46) Equations 45 and 46 indicate that in the spectral plot the variance of Ot is proportional to the area under the curve ΦO(ω). The difference between the order variance and the customer demand variance is proportional to the area under the curve ΦO(ω) - ΦD(ω). Given the transfer function H0(z-1) and H1(z-1), Hosoda and Disney27 and Chen and Disney31 proposed a matrix method to obtain analytical expressions of the demand variance and the orders variance. But for a high-order transfer function, such formulas are not practical to handle. Thus, for an arbitrary highorder ARMA process, we adopt the signal spectrum plot to visually exhibit the statistical properties of the process. Meanwhile we use an efficient recursive numerical algorithm33 to compute its variance. For example, consider an ARMA(2,2) demand process:
∑ G (z
L
(β - 1)[(1 - z-1) ∆NSt+L ) λ
∑ G (z
-1
k
) + C(z-1)]
k)1
et ) (1 - (1 - β)z-1)A(z-1) λG2(z-1) H1(z-1)et (48)
where G2(z-1) )
β-1 1 - (1 - β)z-1
The variance of the additional signal can be calculated as follows: Var(∆NSt+L) )
λ2 π
∫
π
0
|G2(e-jω)|2 |H1(e-jω)|2 dω
(49)
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
min(cost(β)) ) (h + b)
σNS(β)
√2π
β
s.t. ratio(β) e c 0
