A Methodology for Consolidation Effects of Inventory Management with Serially Dependent Random Demand

Abstract

Most studies of inventory consolidation effects assume time-independent random demand. In this article, we consider time-dependence by incorporating an autoregressive moving average structure to model the demand for products. With this modeling approach, we analyze the effect of consolidation on inventory costs compared to a system without consolidation. We formulate an inventory setting based on continuous-review using allocation rules for regular transshipment and centralization, which establishes temporal structures of demand. Numerical simulations demonstrate that, under time-dependence, the demand conditional variance, based on past data, is less than the marginal variance. This finding favors dedicated locations for inventory replenishment. Additionally, temporal structures reduce the costs of maintaining safety stocks through regular transshipments when such temporal patterns exist. The obtained results are illustrated with an example using real-world data. Our investigation provides information for managing supply chains in the presence of time-patterned demands that can be of interest to decision-makers in the supply chain.

1. Introduction

This section provides the introduction, state-of-the-art, contributions, objectives, and plan of the sections used in this article.

1.1. Introduction and Bibliographical Review

In the supply chain, items are typically transferred from producers to distribution centers and then sold to retailers [1]. To create an efficient network, companies aim to minimize costs while maintaining a high level of service [2]. This can be achieved by making decisions regarding the facility location, number of distribution centers serving demand, and supply frequency [3].

One common approach is to set up a supply chain where each dedicated facility serves the demand per unit of time (DPUT), resulting in an independent system (IS) [4,5]. Thus, we define an IS as a framework where each decentralized location is served exclusively by a single centralized location. However, when demand for each dedicated facility in an IS exhibits high variability, costs associated with maintaining service levels, such as holding costs, tend to increase [6]. In such cases, inventory pooling can be implemented through centralized or regular supply using transshipment.

A centralized inventory is given when all decentralized locations are served by a single centralized location, while a regular transshipment (RT) occurs when all decentralized locations are served by all centralized locations. This approach can significantly reduce demand variability at individual points and often lead to cost savings compared to an IS [6,7,8].

The study conducted in [7] examined inventory pooling defined as a setting in supply chain to group goods minimizing the system’s cost in multiple locations. In [7], it was also investigated the cost reductions associated with inventory pooling and its relationship with demand variability assuming independent and identically distributed (IID) demand. The research [7] concluded that centralization is necessary for improving system operation when demand varies significantly across different areas. In a study given in [9], centralized and decentralized supply chains with uncertain supply were compared. The centralized network leveraged risk-pooling effects by combining demand variation from numerous sites, resulting in lower inventory costs and expected total costs (TC). In [10], the authors extended the concept of consolidated safety stock (${\mathrm{SS}}_{\mathrm{c}}$) to multiple storage points, while in [11], it was conducted a similar analysis for consolidated cycle stock (${\mathrm{CS}}_{\mathrm{c}}$). Recent techniques have explored cost reduction by identifying potential consolidation facilities that combine shipments with improved service levels [12].

The selection of appropriate demand allocation rules is crucial for considering the consolidation effect [6]. The Tyagi-Das (TD) rule is one of such approaches, where an equal demand is distributed to decentralized sites from a specific centralized location [5]. When using the TD rule and assuming IID DPUTs over time, the maximum consolidation impact is determined by the lead time (LT) between sites [6]. Consolidation is most effective when the consolidated TC (${\mathrm{TC}}_{\mathrm{c}}$) and inventory decrease. If the mean and standard deviation (SD) of the LT are the same at both centralized locations, the demand provided to decentralized sites should also be the same. Thus, when inventories are consolidated in a single place, the maximum consolidation impact is achieved by employing distinct means/SDs for LT and DPUT [6]. Another allocation approach is the Ballou-Burnetas (BB) rule, where a portion of the demand is met by a main location and the remaining demand by other secondary locations [4]. For IID DPUTs over time, RT is a suitable choice when DPUTs in each centralized site are directly correlated. This approach balances small/large mean/SD for LT and DPUT at various centralized locations and is particularly effective when the holding cost (HC) is relatively small  [6].

Inventory pooling primarily reduces costs by decreasing the inventory level required to meet a given service level. Nonetheless, inventory centralization (IC) does not only impact HC [6]. In general, the ${\mathrm{TC}}_{\mathrm{c}}$ and IS TC (${\mathrm{TC}}_{\mathrm{i}}$) in an inventory pooling system consist of four costs associated with SS and average inventory, including the ordering cost (OC), storage cost, and distribution cost (DC) over a given period [6,13]. For the ${\mathrm{TC}}_{\mathrm{c}}$, the four cost components are: (CC1) $\mathrm{HC}$ multiplied by ${\mathrm{SS}}_{\mathrm{c}}$; (CC2) $\mathrm{HC}$ multiplied by ${\mathrm{CS}}_{\mathrm{c}}$; (CC3) consolidated DC (${\mathrm{DC}}_{\mathrm{c}}$); and (CC4) consolidated OC (${\mathrm{OC}}_{\mathrm{c}}$) [6]; also refer to Section 2.3. The total SS under inventory pooling is crucially influenced by the randomness of DPUT and the LT period [13].

To reduce the ${\mathrm{TC}}_{\mathrm{c}}$, previous studies given in [6,11] used mathematical treatments of the consolidation effect by assuming random DPUT and the following: (a) all locations utilize a SS before and after consolidation, considering the uncertainty of LT demand (LTD); (b) all locations employ the same safety factor for the desired service level; (c) each location adopts a continuous-review model for lot size and reorder point in inventory control; (d) the total average demand of the system remains constant after consolidation; and (e) DPUT and LT are independent. Recent studies conducted in [14,15] extended these approaches by considering an asymmetric distribution for the DPUT and various degrees of kurtosis. These studies demonstrated that the ${\mathrm{TC}}_{\mathrm{c}}$ decreases under specific conditions, including the level of correlation among DPUTs at different locations, mean/SD of DPUTs and LT, safety factor for LTD, and costs associated with holding, ordering, and distribution [6].

1.2. Contribution, Objectives, and Organization of the Article

The main findings from prior research are: (a) small ${\mathrm{TC}}_{\mathrm{c}}$ values generated by centralized locations are linked to indirect correlations between the DPUTs and high values of DPUT SDs; (b) positive correlations between DPUTs and small values for DPUT SD are associated with the use of an IS supply; and (c) distinct SDs for the LT are managed in a better way by an RT due to the combination and balancing of the DPUTs that an RT serves. In [14,16], the authors confirmed the previous conclusions for non-Gaussian distributions. However, these works employed a continuous-review model with IID DPUT. Recent research has considered the dependence on DPUT [17]. Copulas methods have been the primary focus for modeling the joint DPUT distribution [11,18,19]. Nevertheless, it should be noted that DPUTs often exhibit time-dependence, as highlighted in [20]. In this dependence framework, an autoregressive moving average (ARMA) model can adequately describe time-dependent DPUT [21]. In this context, in [22], it was warned of a gap in the literature on inventory models between those that incorporate temporal dependence and those that do not. Also, the consideration of ARMA models makes it possible to better capture the behavior of customer demand in cases of temporal dependency. This dependency improves the accuracy of the forecasts. Thus, inventory management, which can be considered in customer relationship management (CRM) models to improve customer demand, enables us to plan seasonal marketing with product delivery on time.

One of the few works on consolidation effects in inventories under time-dependent DPUTs is attributed to [23]. These authors found the optimal SS by assuming homogeneity in the DPUT variability and considering the temporal correlation. In [24], optimal policies and approximations were proposed for inventories with temporally correlated demand in a similar setting, although not on the consolidation impact.

Given that ARMA models are commonly used to model time-dependent demand [25], its characteristic and how they can be incorporated into different allocation rules are of interest. One unique characteristic of the ARMA models is that the variance of the random variable (RV) of interest in a time-dependent framework is conditional on past data, and it is less than the marginal variance in an independent framework. Since decreasing demand variability results in cost savings under inventory pooling, this suggests that the precision of the inventory model increases when analyzing such demand time series, leading to a decrease in ${\mathrm{TC}}_{\mathrm{c}}$ through allocation rules. However, to the best of our knowledge, there is a gap in the literature because no studies have analyzed how the temporal dependence of demand affects the ${\mathrm{TC}}_{\mathrm{c}}$ of a pooled inventory for different allocation rules.

To address this gap, our investigation designs a procedure to analyze how the temporal dependence described by an ARMA model impacts the ${\mathrm{TC}}_{\mathrm{c}}$ of a pooled inventory for distinct allocation rules. Specifically, our study includes: (i) comparing the TC between a system with inventory pooling and IS under time-dependence; and (ii) analyzing changes in metrics of inventory TC, under the previously assumed scenarios but adding temporality.

The rest of the paper is as follows. Section 2 introduces our approach upon two pillars: (a) modeling assuming an IID DPUT and modeling DPUT based on an ARMA model for temporal dependence; and (b) the consolidation effects of inventory. In Section 3, we perform simulations to analyze how time-dependence affects ${\mathrm{TC}}_{\mathrm{c}}$ of inventories for distinct allocation rules. Section 4 analyzes an illustrative example with real data and guidelines on the theme are established. Section 5 discusses our investigation and states concluding remarks, limitations, and proposals to continue investigating the theme of the present study.

2. Proposed Methodology

Given the presence of correlated demand over time, we propose using the following methodology to consolidate inventories. To carry out our study, we have followed the approach given in [6] which assumed IID demands. However, to derive a consolidation effect on time-dependent demand, we relax the IID assumption by assuming an ARMA model for DPUT. Next, this section states the new approach introduced in the present investigation.

2.1. Modeling DPUT without Temporal Dependence

Suppose that an IID RV Z models the randomness of the DPUT with a normal distribution. Then, Z has a probability density function (PDF) given by

$$f_Z\left(z\right)=\frac{1}{\sqrt{2\pi {\sigma }_Z^2}}exp\left(-\frac{{(z-{\mu }_Z)}^2}{2{\sigma }_Z^2}\right),z\in \mathbb{R},{\mu }_Z\in \mathbb{R},{\sigma }_Z>0,$$

with ${\mu }_Z$ and ${\sigma }_Z^2$ denoting the statistical parameters related to the marginal mean and variance of Z, respectively. In the usual regression model, ${\mu }_Z$ can be described by the values of r explanatory variables (also called independent variables or covariates) summarized in the $(r+1)\times 1$ vector of values $\mathit{x}={(1,x_1,\cdots ,x_r)}^{\top }$ by means of

$${\mu }_Z=\mathrm{E}\left(Z\right)={\mathit{x}}^{\top }\mathit{\beta },$$

(1)

where $\mathit{\beta }={({\beta }_0,{\beta }_1,\cdots ,{\beta }_r)}^{\top }$ is the coefficient vector related to each element of $\mathit{x}$.

Let $Z_1$ and $Z_2$ denote IID RVs which have a cumulative distribution function (CDF) $F_{Z_1,Z_2}$. Hence, we get a relation given by

$$C(c,d;\varrho )=F_{Z_1,Z_2}(F_{Z_1}^{-1}\left(c\right),F_{Z_2}^{-1}\left(d\right)),c,d\in [0,1],$$

(2)

where

$$\varrho ={\varrho }_{Z_1,Z_2}=\frac{\mathrm{Cov}(Z_1,Z_2)}{\sqrt{{\sigma }_{Z_1}^2{\sigma }_{Z_2}^2}},$$

(3)

with $\varrho$ denoting a dependence parameter between $Z_1,Z_2$, and ${\sigma }_{Z_1}^2=\mathrm{V}\left(Z_1\right)$ and ${\sigma }_{Z_2}^2=\mathrm{V}\left(Z_2\right)$. Hence, from the formula stated in (2) and in [16], we reach that the PDF related to $C(u,v,\varrho )$ is established as

$$c(F_{Z_2}\left(z_2\right),F_{Z_1}\left(z_1\right),{\varrho }_{Z_1,Z_2})=\frac{{\partial }^{2}C}{\partial F_{Z_2}\partial F_{Z_1}},$$

where ∂ denotes the corresponding partial derivatives. For more details about the use of copulas to model the correlated DPUTs in consolidated effects, see [11].

2.2. Modeling DPUT with Temporal Dependence

To describe the DPUT when it is time-dependent, we can employ an RV $Z_t$ measured at time t with $t\in \{1,\cdots ,N\}$. For this purpose, we may consider the conditional normal distribution of $Z_t$ depending on a past data set ${\mathit{H}}_t\equiv \{{\mathit{x}}_1,\cdots ,{\mathit{x}}_t,z_1,\cdots ,z_{t-1}\}$ by an ARMA model of p (autoregressive) and q (moving average) orders, that is, ARMA($p,q$), with ${\mathit{x}}_1,\cdots ,{\mathit{x}}_t$ and $z_1,\cdots ,z_{t-1}$ being the observed values of $\mathit{x}$ and Z, respectively, at the corresponding instants of time. Through this model, the mean of $Z_t$ can be formulated as

$${\mu }_{Z_t}=\mathrm{E}\left(Z_t\right)={\mathit{x}}_t^{\top }\mathit{\beta }+\sum _{l=1}^p{\varphi }_l(z_{t-l}-{\mathit{x}}_{t-l}^{\top }\mathit{\beta })+\sum _{m=1}^q{\vartheta }_m(z_{t-m}-{\mu }_{t-m}),$$

(4)

with ${\varphi }_l,{\vartheta }_m$ being the l-th and m-th elements of an ARMA structure, respectively, of $p,q$ orders, and $\mathit{\beta }$ is defined as in (1), but now related to r covariates depending over time, that is, ${\mathit{x}}_t={(1,x_{1t},\cdots ,x_{rt})}^{\top }$.

If we define a martingale residual as $U_t=Z_t-{\mu }_{Z_t}$ for the model error and state $V_t=Z_t-{\mathit{x}}_t^{\top }\mathit{\beta }$, from (4), we can take the expression given by

$$V_t=\sum _{h=1}^p{\varphi }_h{v}_{t-h}+\sum _{j=1}^q{\vartheta }_j{u}_{t-j}+U_t.$$

(5)

Now, if we consider the lag operators presented as

$$\begin{array}{ccc}\hfill \mathsf{\Phi }\left(M\right)& =& 1-{\varphi }_1{L}^{\left(1\right)}-\cdots -{\varphi }_p{L}^{\left(p\right)},\hfill \end{array}$$

(6)

$$\begin{array}{ccc}\hfill \vartheta \left(M\right)& =& 1+{\vartheta }_1{M}^{\left(1\right)}+\cdots +{\vartheta }_q{M}^{\left(q\right)},\hfill \end{array}$$

(7)

then, we may rewrite the model defined in (5) as

$$V_t=\frac{\vartheta \left(M\right)}{\mathsf{\Phi }\left(M\right)}{U}_t=\mathsf{\Psi }\left(M\right)U_t,$$

with

$$\mathsf{\Psi }\left(M\right)=\frac{\vartheta \left(M\right)}{\mathsf{\Phi }\left(M\right)}=1+{\psi }_1{M}^{\left(1\right)}+{\psi }_2{M}^{\left(2\right)}+\cdots ,$$

(8)

and $\mathsf{\Phi }\left(M\right)$ being an invertible function. In [26], it was proved that, for the RV Z, its marginal mean and variance are stated as

$$\begin{array}{ccc}\hfill {\mu }_{Z_t}& =& \mathrm{E}\left(Z_t\right)=\mathrm{E}({\mathit{x}}_t^{\top }\mathit{\beta }+V_t)=\mathrm{E}\left({\mathit{x}}_t^{\top }\mathit{\beta }\right)+\mathrm{E}\left(V_t\right)={\mathit{x}}_t^{\top }\mathit{\beta },\forall t,\hfill \\ \hfill {\sigma }_{Z_t}^2& =& \mathrm{V}\left(Z_t\right)=\mathrm{V}\left(V_t\right)=\mathrm{E}\left(V_t^2\right)={\mathsf{\Psi }}^{\left(2\right)}\left(M\right)\mathrm{V}\left(Z_t\right|{\mathit{H}}_t),\forall t,\hfill \end{array}$$

(9)

where ${\mathsf{\Psi }}^{\left(2\right)}\left(M\right)=1+{\psi }_1^2{L}^{\left(1\right)}+{\psi }_2^2{L}^{\left(2\right)}+\cdots$, for stationary and invertible time series.

Remark 1. 

Note that ${\mathsf{\Psi }}^{\left(2\right)}\left(M\right)\ge 1$. Thus, the variance conditional on past data is less than the marginal variance, that is, when no temporal disposition of the data is considered. In the modeling, consider the assumption that ${\sigma }_{Z_t}^2$ is invariant over time.

Now, if we define an RV T as the random LT, with T being independent of $Z_t$, then the random sum S of $Z_t$ for an item under inventory until reaching LT is expressed as $S={\sum }_{t=1}^T{Z}_t$. We denote the PDF of S as $f_S$, which is defined on $[0,\infty ]$ (non-negative support), and its CDF is given by $F_S\left(r\right)={\int }_0^r{f}_S\left(x\right)\mathrm{d}x$, with its quantile function being stated as $F_S^{-1}\left(x\right)$, for $x\in [0,1]$. The conditional mean and variance of S on past data, ${\mathit{H}}_t$ namely, are defined, respectively, as

$$\begin{array}{ccc}\hfill {\mu }_{S_t}& =& \mathrm{E}\left(S\right|{\mathit{H}}_t)={\mu }_T\mathrm{E}\left(Z_t\right|{\mathit{H}}_t)={\mu }_T{\mu }_{Z_t},\hfill \end{array}$$

(10)

$$\begin{array}{ccc}\hfill {\sigma }_{S_t}^2& =& \mathrm{V}\left(S\right|{\mathit{H}}_t)={\sigma }_T^2{\mathrm{E}}^2\left(Z_t\right|{\mathit{H}}_t)+{\mu }_T\mathrm{V}\left(Z_t\right|{\mathit{H}}_t)={\sigma }_T^2{\mu }_{Z_t}^2+{\mu }_T{\sigma }_{Z_t}^2,\hfill \end{array}$$

(11)

where ${\mu }_T=\mathrm{E}\left(T\right)$ and ${\sigma }_T^2=\mathrm{V}\left(T\right)$. Note that for IID DPUTs, the expressions given in (10) and (11) reduce to

$$\begin{array}{ccc}\hfill {\mu }_S& =& \mathrm{E}\left(S\right)=\mathrm{E}\left(T\right)\mathrm{E}\left(Y\right)={\mu }_T{\mu }_Y,\hfill \\ \hfill {\sigma }_S^2& =& \mathrm{V}\left(S\right)=\mathrm{V}\left(T\right){\mathrm{E}}^2\left(Y\right)+\mathrm{E}\left(T\right)\mathrm{V}\left(Y\right)={\sigma }_T^2{\mu }_Y^2+{\mu }_T{\sigma }_Y^2.\hfill \end{array}$$

The data generated from a bivariate ARMA model conditional on past data are presented in [27]. Generalized volatile multivariate time series models by using the Cholesky decomposition were presented in [28]. With this generalization, it is possible to parameterize the variance-covariance matrix of a bivariate ARMA model for $Z_{1_t}={\mu }_{Z_{1_t}}+U_{1_t}$ and $Z_{2_t}={\mu }_{Z_{2_t}}+U_{2_t}$, which can be related by using the regression model given by

$$U_{2_t}=\beta u_{1_t}+{\epsilon }_{2_t},$$

(12)

where $\beta$ is a coefficient to be determined, for example, with the ordinary least squares method; $u_{1_t}$ is the value of the covariate $U_{1_t}$; and ${\epsilon }_{2_t}$ is the model error with mean zero.

For the model defined in (12), it is possible to prove that the corresponding conditional variance over time is $\mathrm{V}\left({\epsilon }_{2_t}\right|{\mathit{H}}_t)=\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_t)(1-{\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}^2)$, where ${\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}$ is the correlation coefficient between both variables conditional on past data defined as

$${\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}=\frac{\mathrm{Cov}(Z_{1_t},Z_{2_t}|{\mathit{H}}_t)}{\sqrt{\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_t)\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_t)}},$$

(13)

with $\mathrm{Cov}(Z_{1_t},Z_{2_t}|{\mathit{H}}_t)=\widehat{\beta }\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_t)$ being the conditional covariance on past data obtained from the Cholesky decomposition. Based on (9), we have that

$$\begin{array}{ccc}\hfill {\sigma }_{Z_{1_t}|{\mathit{H}}_t}^2& =& \mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_t)=\mathrm{V}\left(Z_{1_t}\right)/{\mathsf{\Psi }}_{Z_{1_t}}^{\left(2\right)}\left(M\right),\hfill \\ \hfill {\sigma }_{Z_{2_t}|{\mathit{H}}_t}^2& =& \mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_t)=\mathrm{V}\left(Z_{2_t}\right)/{\mathsf{\Psi }}_{Z_{2_t}}^{\left(2\right)}\left(M\right)=\mathrm{V}\left({\epsilon }_{2_t}\right)/((1-{\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}^2){\mathsf{\Psi }}_{Z_{2_t}}^{\left(2\right)}\left(M\right)),\hfill \end{array}$$

(14)

where ${\mathsf{\Psi }}_{Z_{j_t}}^{\left(2\right)}\left(M\right)$, for $j\in \{1,2\}$, are as given in (8), whereas the expressions stated in (14) are the conditional variances on past data for $Z_{1_t}$ and $Z_{2_t}$, which have a correlation coefficient given by (13). Note that ${\mathsf{\Psi }}_{Z_{1_t}}^{\left(2\right)}\left(M\right)$ and ${\mathsf{\Psi }}_{Z_{2_t}}^{\left(2\right)}\left(M\right)$ are autoregressive parametric functions ($\mathsf{\Phi }$), with moving average ($\vartheta$), describing the temporal dependence of the DPUT. If there is no temporal dependence (that is, $\mathsf{\Phi }=0$ and $\vartheta =0$), the variances of $Z_{1_t}$ and $Z_{2_t}$ are the marginal variances of each variable. Therefore, their conditional correlation on past data is ${\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}=0$ based on (3).

Remark 2. 

The results given in (14), in addition to Remark 1, are fundamental for improving the precision of the estimation of conditional variances on past data and their respective correlations. Note that these variances and correlations are involved in the inventory TC for both pooled and supply IS. Thus, we should have more precision when comparing both systems if temporal dependence is considered.

2.3. Consolidation Effect

The components that form the ${\mathrm{TC}}_{\mathrm{c}}$ are defined as

$$\begin{array}{ccc}\hfill {\mathrm{SS}}_{\mathrm{c}}& =& k\sum _{f=1}^m\sqrt{{\left(\mathrm{E}\left(\sum _{j=1}^n{\varsigma }_{j,f}{Z}_{j_t}\right)\right)}^2\mathrm{V}\left(T_f\right)+\mathrm{V}\left(\sum _{j=1}^n{\varsigma }_{j,f}{Z}_{j_t}\right)\mathrm{E}\left(T_f\right)},\hfill \end{array}$$

(15)

$$\begin{array}{ccc}\hfill {\mathrm{CS}}_{\mathrm{c}}& =& \frac{1}{\sqrt{2\mathrm{HC}}}\sum _{f=1}^m\sqrt{{\mathrm{OC}}_f\mathrm{E}\left(\sum _{j=1}^n{\varsigma }_{j,f}{Z}_{j_t}\right)},\hfill \end{array}$$

(16)

$$\begin{array}{ccc}\hfill {\mathrm{OC}}_{\mathrm{c}}& =& \frac{1}{2{\mathrm{CS}}_{\mathrm{c}}}\sum _{f=1}^m{\mathrm{OC}}_f\mathrm{E}\left(\sum _{j=1}^n{\varsigma }_{j,f}{Z}_{j_t}\right),\hfill \end{array}$$

(17)

$$\begin{array}{ccc}\hfill {\mathrm{DC}}_{\mathrm{c}}& =& \sum _{f=1}^m\sum _{j=1}^n{\mathrm{DC}}_{{\mathrm{T}}_{f,j}}{\varsigma }_{j,f}{Z}_{j_t},\hfill \end{array}$$

(18)

where ${\mathrm{OC}}_f$ is the OC measured in $ per order for the centralized site/location f; ${\mathrm{HC}}_f=\mathrm{HC}$ denotes the HC (in $ per time unit); $k_f=k$ represents the safety factor; ${\varsigma }_{j,f}$ is the part of the DPUT mean from the site j, which is decentralized, linked to the site f, which is centralized; $\mathrm{E}\left(T_f\right)$ and $\mathrm{V}\left(T_f\right)$ represent the LT mean and variance, respectively; and ${\mathrm{DC}}_{{\mathrm{T}}_{f,j}}$ denotes the unit cost to move an article from site f to decentralized site j. Observe that ${\varsigma }_{j,f}$ belongs to $[0,1]$, for all j and f from 1 to n and m, respectively, ${\sum }_{f=1}^m{\varsigma }_{j,f}=1$, where $m,n$ stand for the numbers of centralized/decentralized sites, respectively. Also note that $Z_{j_t}$ is the DPUT at decentralized site j over time $t\in \{1,\dots ,N\}$; otherwise if it is independent over time, we simply denote it by $Z_j$.

We compare the DPUT described by an ARMA model through the conditional DPUT on past data (${\mathit{H}}_t$) versus a DPUT assuming an IID framework. From Remark 2, note that $\mathrm{V}\left({\sum }_{j=1}^n{\varsigma }_{j,f}{Z}_{t_j}\right)\mathrm{E}\left(T_f\right)$ given in (15) changes whether the DPUT has temporal dependence or not, altering the optimal results for each allocation rule. Then, expressions given in (15)–(18) have different functional forms according to the allocation employed for inventory pooling.

We show the case of $m=n=2$ to provide a basis to compare performance. However, it may be generalized to any m and n, where $m\le n$, following the works presented in [4,5,6,11]. If $n=m=2$, the TD rule gives ${\varsigma }_{1,f}={\varsigma }_{2,f}={\varsigma }_f$, for $0\le {\varsigma }_f\le 1$, with ${\sum }_{f=1}^m{\varsigma }_f=1$ and ${\varsigma }_2=1-{\varsigma }_1$. Nonetheless, according to the BB rule, a centralized facility f provides equal demand to its main decentralized site j (that is, ${\varsigma }_{1,1}={\varsigma }_{2,2}={\varsigma }_r$; ${\varsigma }_{1,2}={\varsigma }_{2,1}=1-{\varsigma }_r$), where r is the site with the largest portion of the demand.

The optimal solution using the BB rule is not only different from the optimal solution obtained with the TD rule, but it also leads to different combinations of inventories [6]. Under the BB rule, if the optimal solution of ${\varsigma }_r$ is zero or one, the supply chain reduces to an IS with both decentralized centers attended by a dedicated facility for each one. Nevertheless, if ${\varsigma }_r$ is between zero and one, RT occurs implying that all centralized sites serve both dispersed places. Under the TD rule, if ${\varsigma }_1\in \{0,1\}$, decentralized sites have a single facility serving all decentralized sites. Figure 1 exemplifies the aforementioned.

When the optimal solution of ${\varsigma }_r$ is zero or one (supply IS), both sites being decentralized are attended by a facility under the BB rule. Nevertheless, if ${\varsigma }_1\in \{0,1\}$, decentralized sites have a single facility employing the TD rule. Moreover, if the optimum value ${\varsigma }_r$ is between zero and one, RT occurs when the BB rule is utilized. This implies that all centralized sites serve both places.

Under the TD rule, intermediate values of ${\varsigma }_1$ yield the same pattern, but they are non-optimal solutions [11]. In what follows, we state how to compute ${\mathrm{CS}}_{\mathrm{c}}$, ${\mathrm{SS}}_{\mathrm{c}}$, ${\mathrm{OC}}_{\mathrm{c}}$ and ${\mathrm{DC}}_{\mathrm{c}}$, employing both allocation rules, considering the temporal dependence of the DPUTs. However, the same expressions are valid when facing IID DPUTs, but they must be adapted to the marginal variance of the DPUTs, and, therefore, of the LTDs. Then, under the TD rule, assuming the temporal dependence of $Z_{j_t}$ on ${\mathit{H}}_{j_t}$ at decentralized site j, formulas defined in (15)–(18) have structures stated as

$$\begin{array}{ccc}\hfill {\mathrm{SS}}_{\mathrm{c}}& =& k\left(\sqrt{{\varsigma }_1^2{\left(\sum _{j=1}^2{\mu }_{j_t}\right)}^2\mathrm{V}\left(T_1\right)+{\varsigma }_1^2\mathrm{V}\left(\sum _{j=1}^2{Z}_{j_t}|{\mathit{H}}_t\right)\mathrm{E}\left(T_1\right)}\right.\hfill \end{array}$$

$$\begin{array}{ccc}& & \left.+\sqrt{{(1-{\varsigma }_1)}^2{\left({\displaystyle \sum _{j=1}^2}{\mu }_{j_t}\right)}^2\mathrm{V}\left(T_2\right)+{(1-{\varsigma }_1)}^2\mathrm{V}\left({\displaystyle \sum _{j=1}^2}{Z}_{j_t}|{\mathit{H}}_t\right)\mathrm{E}\left(T_2\right)}\right),\hfill \end{array}$$

(19)

$$\begin{array}{ccc}\hfill {\mathrm{CS}}_{\mathrm{c}}& =& \sqrt{\frac{{\sum }_{j=1}^2{\mu }_{j_t}}{2\mathrm{HC}}}\left(\sqrt{{\mathrm{OC}}_1{\varsigma }_1}+\sqrt{{\mathrm{OC}}_2(1-{\varsigma }_1)}\right),\hfill \end{array}$$

(20)

$$\begin{array}{ccc}\hfill {\mathrm{OC}}_{\mathrm{c}}& =& \frac{{\mathrm{OC}}_1{\varsigma }_1({\mu }_{1_t}+{\mu }_{2_t})}{2{\mathrm{CS}}_{\mathrm{c}}}+\frac{{\mathrm{OC}}_2(1-{\varsigma }_1)({\mu }_{1_t}+{\mu }_{2_t})}{2{\mathrm{CS}}_{\mathrm{c}}},\hfill \end{array}$$

(21)

$$\begin{array}{ccc}\hfill {\mathrm{DC}}_{\mathrm{c}}& =& {\mathrm{DC}}_{{\mathrm{T}}_{1,1}}{\mu }_{1_t}+{\mathrm{DC}}_{{\mathrm{T}}_{1,2}}{\mu }_{2_t},\hfill \end{array}$$

(22)

with ${\mu }_{1_t},{\mu }_{2_t}$ being formulated as in (4), and 

$$\mathrm{V}\left(Z_{1_t}+Z_{2_t}|{\mathit{H}}_t\right)=\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_{1_t})+\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_{2_t})+2{\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}\sqrt{\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_{1_t})\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_{2_t})}.$$

We assume that the inventory is consolidated at site 1 to compute ${\mathrm{DC}}_{\mathrm{c}}$. Under the BB rule, considering the temporal dependence of $Z_{j_t}$ on ${\mathit{H}}_{j_t}$ at decentralized site j, (15)–(18) take, respectively, the forms given by

$$\begin{array}{ccc}\hfill {\mathrm{SS}}_{\mathrm{c}}& =& k\left(\right({({\varsigma }_r{\mu }_{1_t}+(1-{\varsigma }_r){\mu }_{2_t})}^2\mathrm{V}\left(T_1\right)+({\varsigma }_r^2\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_{1_t})+{(1-{\varsigma }_r)}^2\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_{2_t}).\hfill \\ & & +2{\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}\sqrt{\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_{1_t})\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_{2_t})}{\varsigma }_r(1-{\varsigma }_r)\left)\mathrm{E}\left(T_1\right)\right){)}^{1/2}\hfill \\ & & +k\left(\right({((1-{\varsigma }_r){\mu }_{1_t}+{\varsigma }_r{\mu }_{2_t})}^2\mathrm{V}\left(T_1\right)+({(1-{\varsigma }_r)}^2\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_{1_t}).\hfill \\ & & +{\varsigma }_r^2\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_{2_t})2{\varrho }_{Z_{2_t},Z_{1_t}|{\mathit{H}}_t}\sqrt{\mathrm{V}\left(Z_{1_t}\right|{\mathit{H}}_{1_t})\mathrm{V}\left(Z_{2_t}\right|{\mathit{H}}_{2_t})}{\varsigma }_r(1-{\varsigma }_r)\left)\mathrm{E}\left(T_2\right)\right){)}^{1/2},\hfill \end{array}$$

(23)

$$\begin{array}{ccc}\hfill {\mathrm{CS}}_{\mathrm{c}}& =& \frac{1}{\sqrt{2\mathrm{HC}}}\left(\sqrt{{\mathrm{OC}}_1({\varsigma }_r{\mu }_{1_t}+(1-{\varsigma }_r){\mu }_{2_t})}+\sqrt{{\mathrm{OC}}_2\left((1-{\varsigma }_r){\mu }_{1_t}+{\varsigma }_r{\mu }_{2_t}\right)}\right),\hfill \end{array}$$

(24)

$$\begin{array}{ccc}\hfill {\mathrm{OC}}_{\mathrm{c}}& =& \frac{{\mathrm{OC}}_1\left({\mu }_{1_t}{\varsigma }_r-{\mu }_{2_t}({\varsigma }_r-1)\right)}{2{\mathrm{CS}}_{\mathrm{c}}}+\frac{{\mathrm{OC}}_2\left({\varsigma }_r{\mu }_{2_t}-({\varsigma }_r-1){\mu }_{1_t}\right)}{2{\mathrm{CS}}_{\mathrm{c}}},\hfill \end{array}$$

(25)

$$\begin{array}{ccc}\hfill {\mathrm{DC}}_{\mathrm{c}}& =& \left({\varsigma }_r{\mathrm{DC}}_{{\mathrm{T}}_{1,1}}+(1-{\varsigma }_r){\mathrm{DC}}_{{\mathrm{T}}_{2,1}}\right){\mu }_{1_t}+\left({\varsigma }_r{\mathrm{DC}}_{{\mathrm{T}}_{2,2}}+(1-{\varsigma }_r){\mathrm{DC}}_{{\mathrm{T}}_{1,2}}\right){\mu }_{2_t}.\hfill \end{array}$$

(26)

From (24), we must assume that the largest parts ${\varsigma }_r$ of ${\mu }_{j_t}$ are provided with the smallest ${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$ and ${\mathrm{DC}}_{{\mathrm{T}}_{2,2}}$.

In the case of decentralized site j, considering a supply IS, OC, CS, DC, and SS are expressed as  

$$\begin{array}{ccc}\hfill {\mathrm{SS}}_j& =& k\sum _{f=1}^m\sqrt{\mathrm{E}\left(T_f\right)\mathrm{V}\left(\sum _{j=1}^n{Z}_{j_t}|{\mathit{H}}_t\right)+\sum _{j=1}^n{\mu }_{j_t}^2\mathrm{V}\left(T_f\right)},\hfill \end{array}$$

(27)

$$\begin{array}{ccc}\hfill {\mathrm{CS}}_j& =& \sqrt{\frac{{\mathrm{OC}}_{{}_j}{\mu }_{j_t}}{2\mathrm{HC}}},{\mathrm{OC}}_j=\sum _{f=1}^m\frac{{\mathrm{OC}}_f{\mu }_{j_t}}{2{\mathrm{CS}}_j},{\mathrm{DC}}_j={\mu }_j{\mathrm{DC}}_{{\mathrm{T}}_{f,j}}.\hfill \end{array}$$

(28)

As in the centralized case, formulas stated in (28) and (27) consider temporal dependence for the DPUTs. Nevertheless, the same expressions are valid when facing IID DPUTs, but they must be adapted to the marginal variance of the DPUTs, and therefore of the LTDs.

Thus, the objective function for a centralized system to be minimized is stated as

$${\mathrm{TC}}_{\mathrm{c}}=\mathrm{HC}({\mathrm{CS}}_{\mathrm{c}}+{\mathrm{SS}}_{\mathrm{c}})+{\mathrm{DC}}_{\mathrm{c}}+{\mathrm{OC}}_{\mathrm{c}},$$

(29)

under the TD rule, that is, for an IC. In contrast, the objective function under a supply IS for site j is given by

$${\mathrm{TC}}_j=\mathrm{HC}({\mathrm{CS}}_j+{\mathrm{SS}}_j)+{\mathrm{DC}}_j+{\mathrm{OC}}_j,j\in \{1,2\},$$

(30)

where the cost coefficients of a centralized system are given in (19)–(22) for the TD rule, in (24)–(26) for the BB rule, and in (29) and (30) for decentralized site j. As the objective is to determine which IC, IS or RT reduces the costs, the decision variables are related to the proportion of demand to be assigned at each site, that is, to find the optimal values of ${\varsigma }_{j,f}$ in the TD rule, and the optimal values of ${\varsigma }_r$ for the BB rule.

Following the approach given in [13], we state a consolidation effect metric associated with “the percentage reduction in aggregate SS made possible by a consolidation effect of inventory from multiple locations into one location”. This inventory metric for site j is presented as

$${\mathrm{CM}}_j=1-\frac{{\mathrm{CS}}_{\mathrm{c}}+{\mathrm{SS}}_{\mathrm{c}}}{{\mathrm{CS}}_j+{\mathrm{SS}}_j},j\in \{1,\dots ,n\}.$$

(31)

Figure 2 shows a flow diagram of the methodology proposed in Section 2.

3. Simulation and Analysis of Results

The new methodology stated in the current investigation is evaluated by employing simulations. In this section, we detail our computational setting and different scenarios.

3.1. Computational Setting and Simulation Framework

We use the R software to implement the proposed methodology (www.R-project.org, accessed on 29 June 2020). For the employment of the R language in inventory structures, see [14,29,30,31]. We employ the copula package of R to construct the joint distribution of the DPUTs [32]. This package provides us with the most popular copulas like Archimedian, elliptical, and Gaussian, permitting us to simulate correlated random numbers and fit a copula by goodness-of-fit methods. A generator of bivariate time series described by an ARMA model was programmed by the authors in R code and is available upon request. In addition, the generation of random numbers from bivariate ARMA models also was carried out with R codes, implemented by the authors, and available under request. We use the nloptr package of R to optimize the inventory TC of pooling models. Again, we simulate the scenario $m=n=2$, which may be generalized to other values such as $m\le n$.

The simulation considers inventory policies related to allocation and supply; modeling of the DPUT under independent and dependent cases; and optimization of costs using mathematical programming [33,34,35,36]. First, for inventory policies, we utilize three types of supply (IS/RT/IC) and two allocation rules (TD/BB). Second, the statistical models consider IID DPUTs and temporal dependence by an ARMA structure. We assume several ARMA temporal structures with bivariate normal distributions for the DPUT, whereas the IID framework is studied under marginal normal distributions. In addition, diverse simulated scenarios for mean, SDs, and correlations are considered to generate both types of DPUTs (independent and dependent) in this statistical modeling. The parameters are uniformly distributed into $[a,b]$, denoted as U$[a,b]$, which are chosen from values obtained from the literature; see, for example, Table 2 and Appendix C of [11].

Following the analysis done in [6], we contrast IC and supply IS under the TD rule based on the expressions given in (29) and (30), respectively, but assuming temporal dependence for the DPUT, which we called Case A. Moreover, the simulation examines the supply IS, IC, and RT as shown in [6], using the inventory TC defined as

$${\mathrm{TC}}_{\mathrm{c}}=\mathrm{HC}({\mathrm{CS}}_{\mathrm{c}}+{\mathrm{SS}}_{\mathrm{c}})+{\mathrm{DC}}_{\mathrm{c}},$$

(32)

under both the TD and BB rules. In this case, we consider the formula stated in (32) as the objective function to be optimized in a mathematical program and compare it to the optimization of the function defined in (30), which we name Case B.

Note that the ${\mathrm{TC}}_{\mathrm{c}}$ given in (29) is conformed by CC1 to CC4, whereas that the ${\mathrm{TC}}_{\mathrm{c}}$ stated in (32) considers CC1, CC2, and CC3. Similarly for ${\mathrm{TC}}_j$ formulated in (30), under both IID DPUTs and temporal structure, the components of costs considered are the same. We consider 10,000 random scenarios for $m=n=2$ minimizing the corresponding TC. Then, we obtain the optimal values of CC1 to CC4. The metrics of the simulation are detailed in

Table 1. Type of metric of the simulation study.

Statistical Parameters	Costs
${\mu }_{Z_{1_t}},{\mu }_{Z_{2_t}}\sim \mathrm{U}(80,120)$	${\mathrm{DC}}_{{\mathrm{T}}_{1,1}},{\mathrm{DC}}_{{\mathrm{T}}_{2,2}}\sim \mathrm{U}(0.20,0.25)$
${\sigma }_{Z_{1_t}|{\mathit{H}}_t}\sim \mathrm{U}(3,30)$	${\mathrm{DC}}_{{\mathrm{T}}_{1,2}},{\mathrm{DC}}_{{\mathrm{T}}_{2,1}}\sim \mathrm{U}(0.25,0.30)$
${\sigma }_{Z_{2_t}|{\mathit{H}}_t}\sim \mathrm{U}(3,30)\times {(1-{\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}^2)}^{1/2}$	$\mathrm{HC}\sim \mathrm{U}(1,1000)$
${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}\sim \mathrm{U}(-1,1)$	${\mathrm{OC}}_1\sim \mathrm{U}(17,67)$
${\mu }_{T_1},{\mu }_{T_2}\sim \mathrm{U}(1,5)$	${\mathrm{OC}}_2\sim \mathrm{U}(20,140)$
${\sigma }_{T_1},{\sigma }_{T_2}\sim \mathrm{U}(0.5,2)$
ARMA parameters	Safety factor
${\varphi }_1,{\varphi }_2,{\vartheta }_1,{\vartheta }_2\sim \mathrm{U}(-1,1)$	$k\sim \mathrm{U}(1,3)$

.

Remark 3. 

Similar metrics are used for the IID framework, but discarding the ARMA parameters and past data. Also, ${\sigma }_{Z_{1_t}}$, ${\sigma }_{Z_{2_t}}\sim \mathrm{U}(3,30)$, whose values were taken from [6,11].

3.2. Simulation Study I

We examine which supply system yields the smallest costs in the 10,000 scenarios for Cases A and B, as well as independent and dependent frameworks of DPUT. The TC defined in (29) and (32) for Cases A and B, respectively, are the smallest ones under an IC. We obtain 8621 optimal scenarios recommended an IC, and 1379 scenarios with a supply IS as the best solution for Case A under an IID DPUT. Also, we attain 8628 scenarios deciding an IC as optimal supply, and 1372 scenarios with a supply IS for Case A with temporal dependence on DPUT. In Case B with IID DPUT, we reach 8912, 444, and 644 optimal values for IC, RT, and supply IS, respectively. Under temporal dependence on DPUT, the results are 8889, 697, and 414, respectively.

For the 10,000 simulations performed, we compared the number of occasions in which cost minimization was given by an IC, RT, or IS, both for the model under the assumption of an IID DPUT and for the case where it is taken into account a time-dependent DPUT. For this experiment, for each of the possible decision-making (IC, RT, IS), we applied the ${\chi }^2$ test for the difference in proportions between the model considering IID DPUT and DPUT with time-dependence. The results showed a slight difference only in the case of the proportion of occasions in which the respective models opted for the choice of an IS. For this particular analysis, we detected a difference with a statistical significance of 13% (which can be considered non-significant). The above results indicate that the temporal dependence of the DPUTs reduces the TC with the supply IS according to Case B. This finding is explained by the fact that the variances of the DPUTs conditional on past data are less than for IID DPUTs, as indicated in Remarks 1 and 2.

For the 10,000 scenarios under both IID and dependent frameworks, we statistically compare the parts of CC1 to CC4 of the TC for IC, RT, and supply IS. As mentioned, note that Case B analyses the four cost components CC1-CC4, while Case A only considers CC1-CC3. To graphically compare these cost components between the combinations of IC, RT or IS supply, and IID or ARMA DPUT, we show the distribution of the corresponding cost values in Figure 3. From this figure, we can see that component CC1 has the largest median cost. For the dependent framework, we analyzed the number of scenarios with CC1 as the higher cost than the other components, with 6347, 7721, and 8858 cases for RT, IC, and supply IS in Case B. We found similar results under IID DPUTs. Furthermore, the cost values of CC1-CC4 were again tested using a ${\chi }^2$ statistic for differences on proportions. Significant differences were detected in favor of CC1 under IC, RT, and supply IS, for DPUTs using both with or without time-dependence frameworks. Therefore, we detected CC1 as the crucial component of the TC. Similar results were found for cases A by comparing CC1-CC3. Then, we used the Wilcoxon statistic to test the medians of the costs associated with CC1 obtained by IC, RT, and supply IS for DPUTs considering frameworks with or without time-dependence. We conclude that there are no statistically significant differences at 10% between medians for components of TC under both frameworks.

As in [14], we used the Kruskal-Wallis (KW) statistic to test medians of metrics associated with statistical parameters, costs, and safety factors regarding IC, RT, and supply IS. This test is performed if a minimal TC is attained in Cases A (IC/supply IS) and B (IC/RT/supply IS), verifying how these metrics change between DPUTs with and without temporal structure.

Table 2. Median and KW p-value for the mentioned model, metric, and supply, employing the TC in Case A.

IID Model				ARMA Model
Metric	IC	IS	KW p-Value	Metric	IC	IS	KW p-Value
${\mu }_{Z_1}$	99.61	102.69	<0.001	${\mu }_{Z_{1_t}}$	99.58	101.61	<0.001
${\mu }_{Z_2}$	100.26	100.42	0.71	${\mu }_{Z_{2_t}}$	99.90	99.68	0.37
${\sigma }_{Z_1}$	16.41	16.09	0.44	${\sigma }_{Z_{1_t}|{\mathit{H}}_t}$	16.44	16.31	0.64
${\sigma }_{Z_2}$	16.65	16.56	0.52	${\sigma }_{Z_{2_t}|{\mathit{H}}_t}$	11.98	11.97	0.25
${\varrho }_{Z_1,Z_2}$	−0.014	0.04	0.04	${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$	−0.008	0.034	0.2381
${\mu }_{T_1}$	2.97	2.92	0.21	${\mu }_{T_1}$	2.99	2.96	0.61
${\mu }_{T_2}$	3.01	2.89	0.03	${\mu }_{T_2}$	2.98	2.96	0.75
${\sigma }_{T_1}$	1.22	1.33	<$0.001$	${\sigma }_{T_1}$	1.23	1.33	<$0.001$
${\sigma }_{T_2}$	1.26	1.13	<$0.001$	${\sigma }_{T_2}$	1.27	1.17	<$0.001$
${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$	0.125	0.123	0.02	${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$	0.125	0.124	0.33
${\mathrm{DC}}_{{\mathrm{T}}_{2,2}}$	0.125	0.124	0.86	${\mathrm{DC}}_{{\mathrm{T}}_{2,2}}$	0.125	0.125	0.47
${\mathrm{DC}}_{{\mathrm{T}}_{1,2}}$	0.375	0.376	0.13	${\mathrm{DC}}_{{\mathrm{T}}_{1,2}}$	0.375	0.376	0.34
${\mathrm{DC}}_{{\mathrm{T}}_{2,1}}$	0.375	0.376	0.02	${\mathrm{DC}}_{{\mathrm{T}}_{2,1}}$	0.374	0.376	0.011
k	1.99	1.90	<$0.001$	k	2.01	1.89	<$0.001$
$\mathrm{HC}$	0.38	0.048	<$0.001$	$\mathrm{HC}$	0.38	0.05	<$0.001$
${\mathrm{OC}}_1$	42.37	38.51	<$0.001$	${\mathrm{OC}}_1$	42.76	38.47	<$0.001$
${\mathrm{OC}}_2$	80.75	72.62	<$0.001$	${\mathrm{OC}}_2$	80.87	71.54	<$0.001$
-	-	-	-	${\varphi }_1$	−0.004	−0.010	0.87
-	-	-	-	${\varphi }_2$	−0.0001	−0.0024	0.69
-	-	-	-	${\vartheta }_1$	0.001	0.0191	0.55
-	-	-	-	${\vartheta }_2$	−0.003	0.0005	0.92

and

Table 3. Median and KW p-value for the mentioned model, metric, and supply, employing the TC in Case B.

IID Model					ARMA Model
Metric	IC	IS	RT	KW p-Value	Metric	IC	IS	RT	KW p-Value
${\mu }_{Z_1}$	99.96	100.58	100.14	0.22	${\mu }_{Z_{1_t}}$	99.90	99.65	100.76	0.78
${\mu }_{Z_2}$	100.08	101.62	102.08	0.06	${\mu }_{Z_{2_t}}$	99.55	103.45	101.42	<$0.001$
${\sigma }_{Z_1}$	16.69	12.01	16.61	<$0.001$	${\sigma }_{Z_{1_t}|{\mathit{H}}_t}$	16.83	10.78	17.89	<$0.001$
${\sigma }_{Z_2}$	16.97	12.05	16.85	0.001	${\sigma }_{Z_{2_t}|{\mathit{H}}_t}$	12.33	8.41	11.91	<$0.001$
${\varrho }_{Z_1,Z_2}$	−0.05	0.53	0.06	<$0.001$	${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$	−0.03	0.42	0.06	<$0.001$
${\mu }_{T_1}$	2.98	2.81	2.88	0.034	${\mu }_{T_1}$	3.01	2.74	2.87	<$0.001$
${\mu }_{T_2}$	3.02	2.62	2.94	<$0.001$	${\mu }_{T_2}$	3.01	2.64	3.03	<$0.001$
${\sigma }_{T_1}$	1.22	1.36	1.37	<$0.001$	${\sigma }_{T_1}$	1.23	1.41	1.30	<$0.001$
${\sigma }_{T_2}$	1.25	1.22	1.14	0.002	${\sigma }_{T_2}$	1.25	1.27	1.19	0.16
${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$	0.124	0.124	0.123	0.86	${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$	0.124	0.123	0.125	0.38
${\mathrm{DC}}_{{\mathrm{T}}_{2,2}}$	0.124	0.124	0.124	0.99	${\mathrm{DC}}_{{\mathrm{T}}_{2,2}}$	0.125	0.126	0.124	0.44
${\mathrm{DC}}_{{\mathrm{T}}_{1,2}}$	0.375	0.376	0.376	0.63	${\mathrm{DC}}_{{\mathrm{T}}_{1,2}}$	0.375	0.374	0.378	0.012
${\mathrm{DC}}_{{\mathrm{T}}_{2,1}}$	0.375	0.374	0.377	0.11	${\mathrm{DC}}_{{\mathrm{T}}_{2,1}}$	0.376	0.374	0.376	0.13
k	2.00	1.80	1.88	<$0.001$	k	2.01	1.80	1.81	<$0.001$
$\mathrm{HC}$	0.37	0.07	0.04	<$0.001$	$\mathrm{HC}$	0.37	0.07	0.044	<$0.001$
${\mathrm{OC}}_{{}_1}$	42.07	39.77	42.25	0.035	${\mathrm{OC}}_{{}_1}$	42.27	40.73	41.64	0.21
${\mathrm{OC}}_{{}_2}$	80.17	75.47	73.08	0.025	${\mathrm{OC}}_{{}_2}$	80.55	73.31	69.54	<$0.001$
-	-	-	-	-	${\varphi }_1$	−0.008	0.003	0.007	0.36
-	-	-	-	-	${\varphi }_2$	−0.00028	−0.0098	0.024	0.95
-	-	-	-	-	${\vartheta }_1$	0.006	−0.021	−0.032	0.50
-	-	-	-	-	${\vartheta }_2$	0.0017	−0.02	−0.047	0.21

report this comparison for Cases A and B.

We observe that the medians of ${\varrho }_{Z_1,Z_2}$, ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$, ${\sigma }_{Z_1}$, ${\sigma }_{Z_2}$, ${\sigma }_{Z_{1_t}|{\mathit{H}}_t}$ and ${\sigma }_{Z_{2_t}|{\mathit{H}}_t}$ are related to a decrease in the TC for Case B; see Table 3. Positive correlations and low variability of the demand are linked to the supply IS, validating what was reported in [6], but now with temporal dependence for DPUTs. However, as pointed out, this dependence causes a decrease in the medians of correlations and variances of the DPUTs, when they are conditional on past data.

Table 2 and Table 3 show that it is not possible to observe significant statistical associations for medians of the ARMA parameters ($\varphi ,\vartheta$) in Cases A and B. As stated in Section 2, such parameters affect the correlations and conditional variances of DPUTs on past data, which are components of ${\mathrm{SS}}_{\mathrm{c}}$ and ${\mathrm{SS}}_{\mathrm{j}}$.

Table 4. Median and KW p-value for the mentioned metric, and supply, employing an ARMA model for DPUT based on CC1.

Metric	IC	IS	RT	KW p-Value
${\sigma }_{Z_{1_t}|{\mathit{H}}_t}$	16.46	9.41	17.01	<$0.001$
${\sigma }_{Z_{2_t}|{\mathit{H}}_t}$	11.97	8.17	13.39	<$0.001$
${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$	−0.001	0.50	−0.058	<$0.001$
${\varphi }_1$	−0.007	−0.06	0.021	0.570
${\varphi }_2$	−0.005	0.21	0.090	0.007
${\vartheta }_1$	0.002	0.31	−0.024	0.030
${\vartheta }_2$	−0.001	0.01	−0.019	0.680

reports the medians of some parameters that favor supply IS, RT, and IC for $\mathrm{CC}1=\mathrm{HC}\times \mathrm{SS}$. We think the temporal dependence favors lower $\mathrm{CC}1$ with RT systems concerning the TC stated under the expression stated in (32), which better balances the conditional variances of the DPUTs. Note that ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$ cannot be highly negative for obtaining a low $\mathrm{CC}1$ in IC. In addition, ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$ increases to high positive values favoring a low $\mathrm{CC}1$ in the supply IS, and slightly negative values to reach lower $\mathrm{CC}1$ in the RT. Slightly negative values of the AR parameter $\varphi$ favor an IC, moderately positive values favor a supply IS, and values close to zero favor the RT. The MA coefficients $\vartheta$ close to zero favor an IC, moderately positive values favor a supply IS, and slightly negative values favor an RT.

3.3. Simulation Study II

According to [6], the consolidated effect metric formulated in (31) may be a response variable modeled by a linear function of $R=\mathrm{V}\left(T_2\right)/\mathrm{V}\left(T_1\right)$. From the scatter-plots of Figure 4 with (left) independent and (right) dependent frameworks for the DPUTs, note that approximately linear forms are identified. Then, we consider a linear regression model to analyze the behavior of the consolidated effect metric. The estimated regression coefficients, standard errors, and t-test p-value, with temporal dependence for both $R>1$ and $R<1$, related to the consolidation effect, are reported in

Table 5. Indicators of the t-test for the listed covariate parameter on the consolidated effect with an ARMA model utilizing DPUT and the mentioned R.

Covariate	$\mathit{R}>$ 1			$\mathit{R}<$1
	Estimate	Standard Error	p-Value	Estimate	Standard Error	p-Value
Constant	$3.256\times {10}^{-1}$	$1.333\times {10}^{-2}$	<$0.001$	$1.998\times {10}^{-1}$	$1.423\times {10}^{-2}$	<$0.001$
${\mu }_{Z_{1_t}}$	$-6.482\times {10}^{-4}$	$2.876\times {10}^{-5}$	<$0.001$	$3.449\times {10}^{-4}$	$2.995\times {10}^{-5}$	<$0.001$
${\mu }_{Z_{2_t}}$	$4.408\times {10}^{-4}$	$2.878\times {10}^{-5}$	<$0.001$	$-4.365\times {10}^{-4}$	$2.983\times {10}^{-5}$	<$0.001$
${\sigma }_{Z_{1_t}|{\mathit{H}}_t}$	$3.127\times {10}^{-4}$	$4.209\times {10}^{-5}$	<$0.001$	$7.747\times {10}^{-5}$	$4.424\times {10}^{-5}$	0.080
${\sigma }_{Z_{2_t}|{\mathit{H}}_t}$	$8.543\times {10}^{-5}$	$4.222\times {10}^{-5}$	0.0431	$3.036\times {10}^{-4}$	$4.451\times {10}^{-5}$	<$0.001$
${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$	$-7.013\times {10}^{-3}$	$5.685\times {10}^{-4}$	<$0.001$	$-6.852\times {10}^{-3}$	$5.915\times {10}^{-4}$	<$0.001$
${\mu }_{T_1}$	$-6.603\times {10}^{-4}$	$2.842\times {10}^{-4}$	<$0.001$	$2.576\times {10}^{-3}$	$3.006\times {10}^{-4}$	<$0.001$
${\mu }_{T_1}$	$1.570\times {10}^{-3}$	$2.859\times {10}^{-4}$	<$0.001$	$-1.249\times {10}^{-3}$	$3.034\times {10}^{-4}$	<$0.001$
${\sigma }_{T_1}$	$-3.185\times {10}^{-1}$	$1.267\times {10}^{-3}$	<$0.001$	$2.091\times {10}^{-1}$	$1.123\times {10}^{-3}$	<$0.001$
${\sigma }_{T_2}$	$2.009\times {10}^{-1}$	$1.161\times {10}^{-3}$	<$0.001$	$-2.981\times {10}^{-1}$	$1.284\times {10}^{-3}$	<$0.001$
${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$	$4.184\times {10}^{-2}$	$2.341\times {10}^{-2}$	<$0.001$	$6.839\times {10}^{-3}$	$2.403\times {10}^{-2}$	0.776
${\mathrm{DC}}_{{\mathrm{T}}_{2,2}}$	$6.057\times {10}^{-3}$	$3.189\times {10}^{-2}$	0.8494	$4.562\times {10}^{-3}$	$2.397\times {10}^{-2}$	0.849
${\mathrm{DC}}_{{\mathrm{T}}_{1,2}}$	$2.424\times {10}^{-2}$	$2.331\times {10}^{-2}$	0.2985	$4.596\times {10}^{-4}$	$2.425\times {10}^{-2}$	0.985
${\mathrm{DC}}_{{\mathrm{T}}_{2,1}}$	$2.810\times {10}^{-2}$	$2.322\times {10}^{-2}$	0.2263	$-1.952\times {10}^{-3}$	$2.385\times {10}^{-2}$	0.935
k	$-2.020\times {10}^{-2}$	$6.098\times {10}^{-4}$	<$0.001$	$-8.566\times {10}^{-5}$	$6.119\times {10}^{-4}$	0.889
$\mathrm{HC}$	$-7.082\times {10}^{-2}$	$1.907\times {10}^{-3}$	<$0.001$	$2.374\times {10}^{-3}$	$1.831\times {10}^{-3}$	0.195
${\mathrm{OC}}_{{}_1}$	$-1.015\times {10}^{-3}$	$2.345\times {10}^{-5}$	<$0.001$	$1.151\times {10}^{-3}$	$2.346\times {10}^{-5}$	<$0.001$
${\mathrm{OC}}_{{}_2}$	$8.449\times {10}^{-4}$	$9.813\times {10}^{-6}$	<$0.001$	$-6.272\times {10}^{-4}$	$1.054\times {10}^{-5}$	<$0.001$
${\varphi }_1$	$-4.756\times {10}^{-4}$	$5.899\times {10}^{-4}$	<$0.001$	$-5.691\times {10}^{-5}$	$6.008\times {10}^{-4}$	0.925
${\varphi }_2$	$-3.858\times {10}^{-4}$	$5.592\times {10}^{-4}$	<$0.001$	$-2.969\times {10}^{-4}$	$5.989\times {10}^{-4}$	0.620
${\vartheta }_1$	$-4.507\times {10}^{-4}$	$5.770\times {10}^{-4}$	<$0.001$	$-4.469\times {10}^{-4}$	$5.926\times {10}^{-4}$	0.451
${\vartheta }_2$	$2.567\times {10}^{-4}$	$5.945\times {10}^{-4}$	<$0.001$	$4.247\times {10}^{-4}$	$5.972\times {10}^{-4}$	0.477

.

Comparing the regressions obtained with $R>1$ and $R<1$ in Table 5, note that both coefficients of the intercept, ${\mu }_{2_t},{\mathrm{OC}}_{{}_2}$ and ${\sigma }_{T_2}$ decrease and maintain their significance, while the coefficients of ${\mu }_{1_t},{\mathrm{OC}}_{{}_1}$ and ${\sigma }_{T_1}$ increase keeping their significance. First, the regression coefficients of ${\sigma }_{Z_{1_t}|{\mathit{H}}_t}$ and ${\mu }_{T_2}$ decrease in magnitude and significance, when $R>1$ and $R<1$ are compared, whereas an opposite result occurs with ${\sigma }_{Z_{2_t}|{\mathit{H}}_t}$ and ${\mu }_{T_1}$. Second, the regression coefficients of the distribution costs (${\mathrm{DC}}_{{\mathrm{T}}_{1,1}}$${\mathrm{DC}}_{{\mathrm{T}}_{2,2}},{\mathrm{DC}}_{{\mathrm{T}}_{1,2}},{\mathrm{DC}}_{{\mathrm{T}}_{2,1}}$), as well as the AR (${\varphi }_1,{\varphi }_2$) and MA (${\vartheta }_1,{\vartheta }_2$) coefficients, are not significant in all cases. Third, the coefficient of ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$ maintains its magnitude and significance with $R>1$ and $R<1$, whereas its R-squared coefficients are similar. An analogous analysis without time-dependence can be consulted in Table 7 of [6].

4. Empirical Application and Implications

This section applies our methodology to real data and provides management implications in decision making.

4.1. Illustrative Example

To illustrate our methodology, we apply it to real data taken from kaggle.com/c/competitive-data-science-predict-future-sales (accessed on 29 June 2020).

The data correspond to sales of one of the biggest software firms in Russia named 1C Company. This data set is free to utilize for any purpose, including commercial use. We employed the total monthly sales over 24 months for the product with the highest demand in Moscow. The identifier of this product is ID 20949.

We considered the two demand points that contain sufficient history to fit an ARMA model: the demand points #27 and #31. We call these demand points the decentralized stores. The data considered are shown in

Table 6. Data set corresponding to sales of the 1C Company for the indicated month and year.

Month	Store #27			Store #31
	2013	2014	2015	2013	2014	2015
January	N/A	557	552	N/A	899	743
February	N/A	337	300	N/A	941	180
March	N/A	360	237	N/A	776	481
April	N/A	218	280	N/A	597	460
May	N/A	272	197	N/A	602	434
June	N/A	300	126	N/A	625	482
July	271	268	N/A	799	528	N/A
August	270	348	N/A	820	591	N/A
September	260	335	N/A	950	639	N/A
October	283	312	N/A	978	634	N/A
November	348	426	N/A	989	772	N/A
December	598	792	N/A	1305	1209	N/A

Where N/A is not applicable, not available, not assessed, or no answer.

and graphically represented in Figure 5. Descriptive statistics for both demand points are detailed in

Table 7. Descriptive statistics for the data set corresponding to sales of the 1C Company for the indicated store.

Store	Minimum	Q1	Median/Q2	Q3	Maximum	Mean	SD
#27	128.00	269.50	300.00	351.00	792.00	343.71	147.16
#31	180.00	575.25	691.00	909.50	1305.00	726.42	258.42

Where Q1, Q2, and Q3 are the first, second, and third quartiles, respectively.

. From this table, graphically supported by Figure 5, note that demand point #31 has demand mean and SD greater than stock point #27. Plots of the autocorrelation function (ACF) and partial ACF (PACF) are simple techniques to find temporal DPUT patterns. The reader is suggested to see [28] for an explanation of these graphical plots. To check the temporal dependence of sales, we used the ACF and PACF for both demand points #27 and #31. Their graphical plots are shown in Figure 6, which suggest time-dependence in both time series. Given this result, we utilized the forecast library of the R software with its auto.arima function to fit an ARMA model and estimate its order as well as its parameters. Also, we used the Akaike information criteria (AIC) to select the better model, where a less value of AIC indicates a better fit of the model. The results coincide with the recommendation provided in Figure 6, since the fitted ARMA models are AR($\varphi =0.59$) and MA($\vartheta =0.69$) for demand points #27 and #31, respectively. Furthermore, the estimated SDs conditioned on past data are 118.88 and 211.12 units, respectively.

To evaluate the behavior of our model, we considered low, medium, and high values for the cost parameters of the proposed model, employing the metrics detailed in Section 3.1. Following [6], we stated the safety factor parameter: $k\in \{1,2,3\}$; OC: OC${}_1\in \{20,40,60\}$ and OC${}_2\in \{20,80,140\}$; mean LTs: ${\mathrm{LT}}_1={\mathrm{LT}}_2\in \{1,3,5\}$; and SDs: $S_{{\mathrm{LT}}_1}=S_{{\mathrm{LT}}_2}\in \{0.5,1,2\}$. The proposed model was implemented for all possible combinations. The TCs were evaluated under the TD and BB rules and compared with TC under an IS. In this experiment, overall, 56.06% of the occasions have lower TC under the TD rule, while in 43.9% of the cases, the best decision was the BB rule. There are no lower TC under the IS. One of the most interesting findings resides in the choice of k. It is expected that, as k increases, the TC increases as well, regardless of the evaluated system. However, when $k=1$, most decisions minimizing costs were made under the BB rule. In other cases, when $k\in \{2,3\}$, the TD rule minimized the costs in 59.12% and 69% of the cases, respectively. These results are shown in Figure 7.

Without loss of generality, when the choice is to centralize, the inventories are consolidated in location 1. Unsurprisingly, such a choice tends to be made when OCs for that centralized location are lower, while OCs for the second centralized location are high. No significant differences were observed concerning the average LT.

Nevertheless, it is possible to notice that, when the variability of the LT in centralized location 1 is low compared to the variability of the LT of location 2, the best decision also tends to make centralization in one location. Then, we applied a decision tree using the R software to detect deep patterns in the decision framework. Figure 8 shows the corresponding results. From this plot, we can observe that most decisions are related to LT and variability of the SS. Furthermore, we can observe a particular case when the OC${}_1$ is high with low OC${}_2$. Then, the decision is to apply a centralization under the TD rule. This is expected because this allocation rule centralizes the inventories mainly in one dedicated facility. In this specific case, the solution is centralized in the facility with lower OC, that is, location 1, except for the presence of a high LT variance for facility 1 and lower LT variance for facility 2. In the case of the IID framework, the results are very similar. Note that examination of Figure 8 and Figure 9 shows that the safety factor is more critical. We can notice that there are no different patterns.

4.2. Implications in Management

Following the analysis done in [14], we explain the consequences of the obtained findings for logistics management when decision making. First, when inventory consolidation is desired, it is critical to understand the temporal patterns that DPUT data might have. We conclude from our results that these temporal patterns directly affect the metrics associated with the correlation and variance of the demand when historical data are included. These metrics directly influence component CC1 of the inventory TCs, which is the essential component in all pooled inventory and supply IS strategies; see Figure 3. Next, some of the main implications are provided:

• When facing DPUTs with a temporal structure, the supply IS is favored as it is related to decreasing variances and correlations conditioned on past data.
• IID DPUTs favor IC systems with an increase of the HC of $\mathrm{CS}$.
• Temporal dependence favors lower $\mathrm{CC}1$ with RT systems.
• The metrics with a high effect on percentual consolidate effect for $R>1$ are ordered as: ${\sigma }_T^2>\mathrm{HC}>k>{\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}>{\mathrm{OC}}_{{}_f}>{\mu }_{Z_t}>{\mu }_T>\mathrm{V}\left(Z_t\right|{\mathit{H}}_t)$.
• The metrics with a high effect on the percentual consolidate effect for $R<1$ are: ${\sigma }_T^2>{\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}>{\mathrm{OC}}_{{}_f}>{\mu }_T>{\mu }_{Z_t}>\mathrm{HC}>\mathrm{V}\left(Z_t\right|{\mathit{H}}_t)>k$.

As explained at the end of Section 3.2, one of our findings is that the medians of ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$ associated with IC and supply IS is less than ${\varrho }_{Z_1,Z_2}$. This finding can be explained because ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$ should be regulated by the temporal structure described by the parameters $\varphi$ and $\vartheta$. Then, we propose regressions, for $R>1$ and $R<1$, between the response variable related to the consolidate effect and ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$, as well as with other covariates conditional on given values of $\mathsf{\Phi }=\varphi$ and $\vartheta =\vartheta$, to establish when this regulation is greater or less. These results are shown qualitatively for better comprehension in

Table 8. Regulation effect of ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}$ given the temporal structure values of $\mathsf{\Phi }$ and $\vartheta$ on the consolidate effect.

$R>1$	$\mathsf{\Phi }$	$[-1.0-0.5]$	$[-0.5,0.0]$	$[0.0,0.5]$	$[0.5,1.0]$
$\vartheta$	$[-1,-0.5]$	−	+	−	+
	$[-0.5,0.0]$	+	−	−	$++$
	$[0.0,0.5]$	$++$	$++$	−	$+++$
	$[0.5,1.0]$	$+++$	$+++$	−	−
$R<1$	$\mathsf{\Phi }$	$[-1.0,-0.5]$	$[-0.5,0.0]$	$[0.0,0.5]$	$[0.5,1.0]$
$\vartheta$	$[-1.0,-0.5]$	$+++$	$+++$	$+++$	−
	$[-0.5,0.0]$	$++$	$++$	$++$	−
	$[0.0,0.5]$	+	+	+	−
	$[0.5,1.0]$	−	−	−	−

The regulation effect goes from more to less according to the notation $+++>++>+>-$.

, in which the regulation effect of the temporal structures is given in intervals for parameters $\mathsf{\Phi }$ and $\vartheta$, indicating that this effect goes from more to less according to the notation $+++>++>+>-$. Recall that the IC is associated with ${\varrho }_{Z_{1_t},Z_{2_t}|{\mathit{H}}_t}<0$, and its regulations are considered more significant to favor the IC, while the corresponding regression coefficient is more negative.

From the analysis provided in Table 4, it is essential to highlight some differences that occur when comparing the medians of the ARMA metrics. These metrics are associated with the different supply systems that decrease CC1. There, models with pure positive AR coefficients (without considering the MA part) or positive coefficients of pure MA (without considering the AR part) would be more associated with an IS supply. This generates visible “wave” patterns, which would be related to the IS supply. In addition, coefficients of positive MA are related to upward trends of the time series, which would be associated with an IS supply as well. In the case of positive AR models, the consecutive residuals share similar values and signs, showing large areas of overestimation and underestimation.

5. Discussion, Concluding Remarks, Limitations, and Further Work

In this section, we discuss our investigation and provide conclusions, limitations, and ideas for further research.

5.1. Discussion

Our findings suggest new decision-making patterns to be considered in comparison with previous works. We model the temporal disposition of demand data, due to its effect on the variance conditional on past data. This effect suggests that the rule with the best behavior is a pooled inventory in the safety stock that was consolidated, which is obtained by multiplying with the holding cost and then by the total cost.

Specifically, given that safety stocks increase as the variance is greater, the existence of temporary structures for demand data under an ARMA model implies lower levels of safety stocks, due to the conditional demand overtime on the data is less than the marginal variances typically assumed under an approach of independent and identically distributed random variables. Regarding the total cost, we suggest selecting dedicated facilities for supply. Regular transshipment outperformed the independent supply system statistically for the component associated with the safety stock that was consolidated, which is multiplied by the holding cost.

We detect that when the demand per unit of time of different decentralized sites has low variance and positive correlation among them, an independent system is usually the best option to consider in the supply chain. This result is consistent with the conclusions given in [6] for an independent and identically distributed random demand, but now studied in a time-dependent framework.

5.2. Conclusions

In this study, we have analyzed the effect of consolidation on inventory costs compared to a system without consolidation. We have formulated an inventory structure based on continuous-review using allocation rules for regular transshipment and centralization, which has established temporal structures of demand. Our methodology proposed in this work was implemented in R language developed by the authors, which is available to the readers upon request. This allowed us to conduct simulations employing synthetic data with and without autoregressive moving average structures and to prove our model under a real-world data scenario. Numerical simulations have demonstrated that, under time-dependence, the conditional demand variance, based on past data, is less than the marginal variance. Our finding favors dedicated locations for inventory replenishment. In addition, temporal structures have reduced the costs of maintaining safety stocks through regular transshipments when such temporal patterns exist. These results have been also illustrated in an example employing real-world data.

In summary, in the previous research [11,14,16,17,18,19,21] to our study, such as was warned in [22], there was a gap on inventory models under consolidation effects between those that incorporate temporal dependence and those that do not. In this context, one of the few works on consolidation effects in inventory models is [23], while in [24], optimal policies were proposed for inventories with temporally correlated demand in a similar setting, although not on the consolidation impact. In the present study, we have covered such a gap. Our investigation has provided essential information for managing supply chains in the presence of time-patterned demands. These findings are of interest to decision-makers in the supply chain.

5.3. Limitations and Future Research

Although this work considers a critical aspect, such as the time-dependence on the demand, one of its limitations is the assumption of a Gaussian distribution, which is usually violated in practice. To address this problem, this work may be expandable to more general models, for example, to the largely known generalized linear modeling. This type of model is widely flexible, permitting functional forms that may be linear and non-linear [37]. Furthermore, probability distributions different from the Gaussian or normal model can be considered this kind of statistical model assumed for the demand, enabling other shapes of the demand statistical distribution to be assumed [31]. Moreover, multivariate distributions may also be considered [38]. It is essential to state these new statistical methods when studying a cost reduction through the use of locations for facilities of consolidation, which combine diverse shipments enhancing the service level [12]. These ideas are being explored by the authors and their findings will be documented in new publications. Other directions and applications on our proposal are presented in [39,40,41,42,43].