Optimization of Contribution Margins in Food Services by Modeling Independent Component Demand

We propose a methodology useful for food services, allowing contribution margins to be optimized. This is based on statistical tools, inventory models and ﬁnancial indicators. To reduce the gap between theory and practice, we apply this methodology to the case study of a Chilean company to show its potential. We conduct a real-world demand data analysis for perishable and non-perishable products in the company’s inventory assortment. Then, we use suitable inventory models to optimize the associated costs. We compare the proposed optimized system with the non-optimized system currently employed by the company, using ﬁnancial indicators

Table 1: Acronyms used through the paper.

Introduction
Supply systems and inventory policies affect company logistics positively, minimizing the costs involved, and reducing inefficiencies in their management.It is known that total inventory cost is a function of purchasing (PC), ordering (OC) and storing (SC) costs; see Hillier & Lieberman (2005).Several authors have discussed the importance of having optimal supply and inventory policies in a company together with efficient logistics management; see Blankley, Khouja & Wiggins (2008) and Kogan & Tell (2009).These aspects of logistics are also present in collective food service companies; see Ramirez (2013).Such services prepare menu food portions according to diverse specifications, including nutritional and sanitary issues, based on the different types of clients who consume this menu; see Marambio, Parker & Benavides (2005).The increase in food services is generating an important source of employment in countries and providing multiple market opportunities.This is attributed to people needing to eat out due to activities related to businesses, factories, hospitals, schools and universities.Because of service diversity, the complexity of this food industry has grown considerably, requiring professional management and regulation by government agencies; see MINSAL (2004).
In Chile food services are deemed as small and medium enterprises.Many of these Chilean services are not optimizing their supply of raw materials.Such materials make up their inventory assortment that is divided in perishable (fruits, meats, vegetables) and non-perishable products with greater storing capacity subject to shortage; see Grant, Karagianni & Li (2006).Logistics of these raw materials is based on monthly planning of the menus guided by nutritional considerations.The management of logistics can be improved by using inventory policies that allow contribution margins (CMs) of company to be increased; see Soman (2006) and Nicolau (2005) for a case study of hotels.CMs are the gross profits of a company and summarize the movements of income and costs, which may be direct (variable costing) and indirect (absorption costing).These margins vary depending on product units sold, unit costs, ratio between them, and the total costs and fixed costs involved; see Ramanathan (2006).
An optimal inventory policy can be attained choosing the most adequate inventory model, which involves several aspects; see Botter & Fortuin (2000), Braglia, Eaves & Kingsman (2004), Braglia, Grassi & Montanari (2004), Wanke (2011) and Wanke (2012).When non-perishable (multi-period) products are considered, inventory models are classified in two types: pull and push, ranging from the economic order quantity (EOQ) to the just in time (JIT) supply; see Wanke (2009).The EOQ model is the cornerstone of several software packages for inventory control and is widely used in practice; see Nahmias (2001).The JIT method is useful for raw materials that can be supplied as timely as they are required, however it imposes constraints on logistics limiting its use for certain types of products in food services; see Carter, Carter, Monczka, Slaight & Swan (2000) and Wanke, Arkader & Rodrigues (2008).Chiu (2010) discussed models for multi-period products where shortage is not permitted, seeking to find the EOQ and reorder point (ROP), appropriate for groceries often used by food services.Considering lead time (LT) in the modeling renders its assumptions to be more adherent to real world settings; see Ben-Daya & Raouf (1994).The EOQ model is used altogether with the ROP in inventory control to determine safety stocks (SS) under both random LT and demand, which randomness directly affects the operation of a logistics system; see Speh & Wagenheim (1978) and Wanke (2009).Perishable (single-period) products can only be stored during a limited period.These products usually are fruits, meats and vegetables, which are essential raw materials in food services.When contemplating this type of products, the model based on the critical ratio (CR) or service level is often employed; see Hillier & Lieberman (2005, pp. 961-975).
Multiple and single period models must consider that the demanded quantity of a product cannot accurately be predicted owing to several factors, making it to be a random variable (RV) and, therefore, its behavior must be described by a statistical distribution (or probabilistic model); see Johnson, Kotz & Balakrishnan (1994).The Gaussian (or normal) distribution is often used for describing data of three RVs involved in inventory models: demand, LT, and lead-time demand.It is known that this distribution is validly used for RVs that take negative and positive values.Therefore, quantities less than zero can be admitted in the modeling, although, in practice, this is not possible for the three mentioned RVs; see Keaton (1995) and Nahmias (2001).Mentzer & Krishnan (1988) studied the nonnormality effect on inventory control, indicating that demand for products presenting a normal distribution is found in few practical cases.This is because demand data often follow asymmetric distributions; see Moors & Strijbosch (2002).In any case, the normality assumption must be checked by goodness-of-fit methods; see Castro-Kuriss, Kelmansky, Leiva & Martinez (2009), Castro-Kuriss, Kelmansky, Leiva & Martinez (2010), Barros, Leiva, Ospina & Tsuyuguchi (2014) and Castro-Kuriss, Leiva & Athayde (2014).Thus, using the normal distribution to model the demand and LT, and to determine the ROP and SS, may provoke wrong results, leading to stock shortage or excess.Some non-normal distributions used for describing demand or LT in inventory models are the gamma or Erlang, inverse Gaussian (IG), lognormal (LN), uniform and Weibull; see Burgin (1975), Tadikamalla (1981), Lau (1989), Wanke (2008) and Cobb, Rumí & Salmerón (2013).
A good statistical modeling of demand data and a scientific management of inventories for food service companies can maximize their CMs, resulting in better competitiveness, efficiency and profitability of these companies.This can be helpful in making optimal decisions.
The main objective of this paper is to propose a methodology useful for food services allowing CMs to be optimized.This methodology is based on statistical tools, inventory management models and financial indicators.Specifically, the methodology uses probabilistic models that describe the behavior of demand data for raw materials employed by food services that prepare a daily menu.Hereafter, we refer to these raw materials as components (or products) forming part of a food menu.Then, the logistics process is optimized by using inventory models that depend on the type of product from the corresponding assortment.Hence, the CMs of the company are measured by using absorption costing and improved by means of logistics management.Such an improvement is evidenced when comparing the financial results obtained from the optimized system with respect to the nonoptimized system used by the food service.Since certain authors have stressed the need for conducting case studies to reduce the gap between theory and practice and enable researchers to increase their background (Wagner & Lindemann 2008), we apply this methodology to the case study of a food company that serves the staff of a Chilean hospital.This paper is organized as follows.In Section 3, we propose our methodology.In Section 4, we conduct a case study for a Chilean food service.In Section 5, we provide an illustration for one product from the inventory assortment.Finally, in Section 6, we discuss the conclusions of this study.

Methodology
In this section, we provide a methodology for food services that allows CMs to be optimized.First, we discuss the assumptions and the limitations of our methodology.Second, we mention how demand data for components of a food menu should be collected.Third, we present the statistical tools needed to fit a demand data set to a suitable distribution.Fourth, we detail inventory management models to be used for optimizing the supply system based on the selected distributions.Fifth, we describe the financial indicators of our methodology.An algorithm that summarizes this methodology is provided.

Assumptions and Limitations
The main assumptions of our methodology are (i) random demand, (ii) demand time series free of seasonality and trend, (iii) independent component demand, (iv) constant LT and (v) the need to ascertain managerial costing calculations.Some limitations of our methodology are related to (i) additional research needed to improve the results, especially by introducing aspects to better reflect real world settings, such as issues related to seasonality, trend and independence, and (ii) relevant costs of operation characteristics.
Note that shortage costs for non-perishable products might be unavailable and therefore not incorporated in the analysis.Unlike situation proposals by Silver, Pyke & Peterson (1998) and Zipkin (2010), in our methodology, there are no CMs or penalties imposed to a product with unsatisfied demand.In practice, for food companies, a product is replaced by another when unavailable and customers continue consuming their meals.We set a target level of service based on a safety factor (SF), instead of the simultaneous optimization of EOQ and SF.This is due to the eminently practical nature of our study, which aim is, among others, the transfer of knowledge and management of inventory policy over time to the studied company.Using the SF in an inventory policy necessarily requires the manager to think in terms of service level and inventory segmentation by levels of criticality with respect to shortage of items.When the simultaneous optimization of EOQ and SF is carried out, these issues are less explicit for the manager.

Recording the Data
We recommend designing and implementing a system of record for all products comprising the inventory assortment of raw materials of the food service company.This system must be based on identification codes, unit PC, demanded quantity, price, date and time of entry and exit of products used in the preparation of food portions; see Harvey (2002) and Yajiong (2005).The system of record must be based on individual identification using bar codes and developed for demand profiles of the products in the inventory assortment, throughout the period planned for the study.We recommend a period of six months (26 weeks).Initially, we considered 26 weeks (half a year) as a sample of convenience based on the project's budget for data collection.However, we were able to collect data for one week more.This additional week was considered to increase the sample size.

Demand Statistical Distributions
Based on the system of record mentioned in Subsection 3.2, demand data needed to model the distribution of the demand for each component must be collected and then the demand distribution fitted.Until very recently, one of the problems for using a demand distribution different from the normal model was the limitation of statistical software.However, today this is not a problem, first because currently we have a number of statistical software that has implemented several statistical distributions and, second, the scientific community has at its disposal a non-commercial and open source software for statistics and graphs, named R, which can be obtained at no cost from www.r-project.org.The statistical software R is nowadays very popular in the international scientific community.Then, to perform a statistical analysis of demand data, we use the R software and also employ some of its packages to carry out more specific statistical analysis.As mentioned, the gamma, IG, LN, uniform and Weibull distributions have been used for modeling the demand or the LT in inventory problems and they are implemented in R software packages named gamlss and ig; see Stasinopoulos & Rigby (2007) and Leiva, Hernandez & Sanhueza (2008).Statistical analysis based on BS distributions, including a version known as the BS-Student-t (BS-t) distribution, which has been proven to provide robust estimates of its parameters against outliers (Paula et al. 2012), can be conducted by means of an R software package named gbs; see Barros, Paula & Leiva (2009).Next, we provide some useful results for all of these distributions; see details in Johnson et al. (1994).
The BS-t Distribution.A RV D following the BS-t distribution with shape α > 0, ν > 0 and scale β > 0 parameters is denoted by D ∼ BS-t(α, β, ν).In this case, the PDF and CDF of D are where I a (b, c) is the incomplete beta function ratio.The corresponding QF is again d(q) = F −1 D (q) = β[αz(q)/2 + (αz(q)/2) 2 + 1] 2 , for 0 < q < 1, but now z(q) is the QF of the Student-t distribution with ν degrees of freedom.Note that β is once again the median or 50th percentile of the distribution.The mean and variance of In this case, W = Z 2 follows a Fisher distribution with one degree of freedom in the numerator and ν degrees of freedom in the denominator, which also is useful for goodness of fit purposes.Some of its properties are: c D ∼ BS-t(α, c β, ν), with c > 0, and 1/D ∼ BS-t(α, 1/β, ν).
The Gamma Distribution.A RV D following the gamma distribution with shape α > 0 and scale β > 0 parameters is denoted by D ∼ Gamma(α, β).In this case, the PDF and CDF of D are where Γ(•) and γ(•, •) denote the usual and incomplete gamma functions, respectively.The corresponding QF given by d(q) = F −1 D (q), for 0 < q < 1, must be obtained by solving this equation with an iterative numerical method.The mean and variance of D are E[D] = β and Var[T ] = α 2 β 2 , respectively.The gamma distribution also shares the property c D ∼ Gamma(α, c β), with c > 0.
The Inverse Gaussian Distribution.A RV D following the IG distribution with mean λ > 0 and scale β > 0 parameters is denoted by D ∼ IG(λ, β).In this case, the PDF and CDF of D are 2dλ 2 and and once again the corresponding QF given by d(q) = F −1 D (q), for 0 < q < 1, must be obtained by solving this equation with an iterative numerical method.The The Weibull Distribution.A RV D following the Weibull distribution with shape α > 0 and scale β > 0 parameters is denoted by D ∼ Wei(α, β).In this case, the PDF and CDF of D are The QF of D is d(q) = β[− log(1 − q)] 1/α , for 0 < q < 1, and its mean and variance are E Data Analysis, Parameter Estimation and Goodness-of-Fit of Distributions.As mentioned, R is a free software environment for statistical computing and graphics.Using this software (i) exploratory data analysis (EDA) can be conducted for diagnosing the statistical features present in the demand data; (ii) estimation of the parameters of the BS, BS-t, gamma, IG, LN and Weibull distributions can be carried out by the popular maximum likelihood (ML) method, and (iii) goodness-of-fit of a distribution to a demand data set can be performed by Anderson-Darling (AD) and Kolmogorov-Smirnov (KS) tests and probability plots.Next, we describe the R commands of the gbs, ig and basics packages and briefly illustrate their use.
First, the R software must be downloaded from CRAN.r-project.org and installed as any other software.Second, this software can be used in a simple interactive form with the R commander by installing the Rcmdr package.Third, the gbs and ig packages must be also installed.Data analyses based on the BS and BS-t distributions can be performed with the gbs package, for the IG distribution with the ig package and for the gamma, LN and Weibull distributions with the basics or gamlss packages.Thus, once these packages are installed, they must be loaded into the R software, for example, by the command library(gbs) typing them at the R prompt of the R commander, or with any editor program that the user is considering.Once all these instructions are completed, the data, for example, "component1", must be loaded as data(component1).The data can also be directly typed by the R commander such as an Excel sheet or imported from text files, from other statistical software or from Excel.Table 2 provides examples of some commands that allow us to work with the BS distribution, similar instructions may be used for other distributions; for more details about how to use the gbs package, see Barros et al. (2009).

Inventory Management Models
Once the most suitable demand distribution is chosen from the pool presented in Subsection 3.3, the appropriate inventory management model must be applied taking into account if the product (component) is (i) perishable -under a single period-, (ii) non-perishable -under multiple periods-, or (iii) supplied under the JIT method.Thus, depending on the type of product, we obtain the quantity to be replenished that minimizes the PCs, OCs and SCs according to one of the following inventory models: (M1) Model for non-perishable products or (Q, r): considers that the quantity needed to optimize the OCs and SCs is based on the EOQ model given by where λ is the demand rate in units of the product per time unit, calculated as the mean (expected value) according to the distribution that adequately fits the data.In the model M1, we must also consider the ROP, which is the level that an inventory must have in stock when a purchase order is placed, calculated as r = l λ, where λ is defined in (1) and l is the constant LT.However, given that random consumption occurs causing demand uncertainties, to protect against such uncertainties, it is necessary to include a SS and the ROP becomes where µ D l = E[D l ] = l λ is the mean of the demand during the LT (D l ) and SS = k q σ D l , with k q being the SF associated with a service level q × 100%, or amount of standard deviations (SD) of the demanded quantity during the LT given by σ As noted in (2), it is necessary to know the demand distribution during the LT to determine the SF; see Keaton (1995).This factor can be established by using some percentile of the demand during the LT.To be protected against any unexpected logistics situation, the 95th percentile is usually considered, that is, q = 0.95.Thus, k 0.95 must be obtained from the statistical distribution that adequately fits the demand data during LT.
Note that in the model M1 is not considered a shortage cost of the product, because, in the event of a shortage, it is possible to produce an emergency menu, to avoid any unsatisfied demand of the final product (menu); see details about this model and it assumptions in Hillier & Lieberman (2005, pp. 956-961).Also, we recall that no simultaneous optimization of Q and r is carried out due to the practical nature of our methodology; see details in Subsection 3.2.
(M2) Model for perishable products: it considers the quantity needed to optimize the cost of ordering one unit less (generating temporary shortage), in contrast to ordering one unit more (generating temporary overstock), based on the CR in this case given by CR = [UC − PC]/[UC + HC], where UC is the unsatisfied demand shortage cost per unit, that includes lost revenue and loss cost of customer goodwill, PC is expressed as a purchasing cost per unit of the product, and HC is the holding cost per unit, per day, that includes the SC minus a salvage value of a product unit.The numerator UC − PC results in a decrease in profit, due to not ordering a unit that could have been sold during such period, whereas the denominator UC + HC results also in decrease in profit, but due to ordering a unit that could not be sold during such period.Thus, the single period model for perishable products allows us to obtain the optimum stored quantity from the optimum service level given by where F D (•) is the CDF of the demanded quantity and d 0 the optimum quantity of ordered units; see details in Hillier & Lieberman (2005, pp. 961-975).
(M3) JIT model: it is the just quantity for production, it does not take storage into account and is used for specific products requested for completing the daily menu of the food service company.A Kanban type information system can be used in this case, which allows the availability of the product to be harmonically coordinated; see Carter et al. (2000).

Determination of Financial Indicators
Once an appropriate inventory management model is chosen from M1, M2 or M3 (see Subsection 3.4), the CM for each of p products (components) of the inventory assortment used in the preparation of a food menu portion must be calculated, based on the income obtained during w weeks for the company corresponding to this menu.The quantities of each component used in the preparation of the menu (ingredients) are determined with its respective consumption measuring unit; see Table 12 in Appendix for an example on the equivalence among these units for the products of the case study included in Section 4.
The prorated demand of the product i in the jth week can be obtained by means of the proportion that each product of the food portion holds weekly in the menu calculated according to where DQ i,j is the demanded quantity of the product i in the jth week and DQ j is the demanded quantity for all the products during that week.The income of the company for all the portions of the food menu sold during the jth week is where N j is the number of menus sold and S j is the price of the food menu portion, both of them in the jth week.Thus, the prorated income derived from the product i during the jth week is obtained as PI i,j = I j PD i,j , i = 1, . . ., p, j = 1, . . ., w, where PD i,j and I j are defined in (4) and ( 5), respectively.The PC for the product i in the jth week is PC i,j = NC i,j PQ i,j , i = 1, . . ., p, j = 1, . . ., w, where NC i,j and PQ i,j are the unit net cost and the purchased quantity of the product i during the jth week, respectively.Note that, for the optimized system with the inventory model for non-perishable products, PQ i,j must be estimated from Q given in (1), whereas, in the case of perishable products, PQ i,j must be estimated from d 0 − L j given in (3), with L j being the stock level at the beginning of the jth week.For the non-optimized system, this value can be empirically calculated.Once financial indicators PI i,j and PC j defined in ( 6) and ( 7) are obtained, the variable contribution margin (VCM) of the product i during the jth week must be computed as The OC for the product i during the jth week can be obtained as OC i,j = OC i /52, i = 1, . . ., p, j = 1, . . ., w, where OC i is the annual OC of the product i given by OC i = 3 h=1 OC h i OQ i , with OQ i being the annual order quantity and OC h i the cost of type h given in Table 3, both for the product i.Note that, for the optimized system with the inventory model for non-perishable products, OQ i must be estimated from Q given in (1) using the expression λ/Q for each product (with λ being expressed as a demand rate per year), whereas in the case of perishable products OQ i = 52, for all i = 1, . . ., p.For the non-optimized system, this value can be empirically calculated.Administrative costs associated with order movements (input and general service costs with respect to order generation).OC 2 Inspection and receiving costs (social security contributions and warehouseman wages) of movements associated with an order.OC 3 Transportation costs related solely to order generation.Source: generated by the authors based on Hernández-González (2011).
The SC for the product i during the jth week is where SC i is the annual SC of the product i given by SC i = 5 k=1 SC k i /SQ i , with SC k i being the annual SC of type k defined in Table 4 and SQ i = 52 j=1 SQ i,j the annual stored quantity, both for the product i, and SQ i,j is the stored quantity of the product i in the jth week.Note that, for the optimized system with the inventory management model for non-perishable products, SQ i,j must be estimated from SQ = Q/2 + SS, where Q and SS are given in (1) and (2), respectively, whereas, in the case of perishable products, SQ i,j must be estimated from the expected inventory level by single period.For the non-optimized system, this value can be empirically calculated.Annual cost of energy spent on the storehouse, including battery charging necessary for handling, data processing equipment and lighting.SC 5 Annual cost of rental of equipment and facilities, during insurance, storage and communications, and taxes.Source: generated by the authors based on Morillo (2009).
Revista Colombiana de Estadística 38 (2015) 1-30 We consider CMs as absorbable by sales with respect to indirect costs, which are subtracted from the VCM given in (8) to obtain the total CM of the product i during the jth week as where VCM i,j , OC i,j and SC i,j are given in ( 8), ( 9) and ( 10), respectively.Thus, we collect a series of CMs for p products (one for each of them).Hence, the CM of all the products of the inventory assortment during the jth week is CM j = p i=1 CM i,j , for j = 1 . . ., w, where CM i,j is given in (11).Therefore, the total CM of the inventory system is Note that the objective function to be maximized is the sum of CM i,j for the product i in the jth week, during the entire period of study totalizing w weeks, for the menu composed by p products with independent demand.Here, the margins VCM i,j and costs OC j and SC i,j depend on the inventory model of the product i.This function is expressed as Since our approach to calculating (i) CMs from the differential revenues and (ii) costs from the movements in and out of the inventory assortment is based on independent components (products) and not on the menu, absorbable costs for ordering and storing are also calculated using the same criteria of independence and considering the spread of demand from the proportion of components used in the menu.This approach turns out to be more streamlined, because it does not consider the correlations that might exist between components of the menu, which is a source for future work; see Section 6.

Summary of the Methodology
Algorithm 1 summarizes our methodology in six main steps divided into 13 substeps based on the aspects detailed in Subsections 3.2 to 3.5, from the collection of data until the establishment of the CMs to evaluate the optimized system in relation to the current (non-optimized) system.We recall this algorithm considers the demand for independent components, but once all the components are considered, the total contribution of the components used in the service are maximized.
Algorithm 1 Main methodological steps 1: Collect demand data for the product i in each day of the w weeks (i = 1, . . ., p). 2: For the statistical analysis: 2.1 Carry out a correlation study for data collected in Step 1 to detect possible seasonality, trend or dependence.If neither autocorrelation nor correlation between components are detected, then an EDA for independent data must be conducted.Otherwise, these seasonality and/or trend must be removed using suitable techniques.

Case Study
In this section, given need to conduct case studies focusing on their applicability in firms to reduce the gap between theory and practice, we apply the methodology summarized in Algorithm 1 to an anonymous Chilean food company, which serves the staff of a hospital in the city of Valparaiso.This case study enables researchers to increase their practical knowledge of the aspects involved and understand the complexity of this environment and the managerial efforts made by firms become evident.
This study was leaded by Fernando Rojas and Victor Leiva in the University of Valparaiso-Chile (www.uv.cl) by means of the project grant DIUV 14/2009, during w = 27 weeks covering the period since 20-Nov-2011 to 26-May-2012 (189 days).Details of p = 89 products of the inventory assortment considered in this study are provided in Table 12 with their respective equivalence units.
We bring to mind that unsatisfied demand shortage costs for non-perishable products are unavailable for their incorporation into the analysis, on account that there are no CMs or penalties imposed for a product with unsatisfied demand.If a product is missing, it is replaced by a similar item.Moreover, owing to the practical nature of this study, for non-perishable products, we set a service level based on a SF instead of simultaneously optimizing Q and r.
As mentioned, the data were collected during the period indicated above following the system of record in Subsection 3.2.Note that, in the type of data that we analyze (food services for hospitals), seasonality or trend factors usually are not present; see Step 2 of Algorithm 1.We have also explored the correlation between some products and only a small correlation but marginally not significant was detected, hence we discarded this aspect.In any case, some comments in this line are provided in the conclusions of this study.Moreover, demand data are usually observed over time.Then, one must check whether these data are time dependent or not.An autocorrelation graphical analysis detected that the corresponding autocorrelations are very small, so that dependence on time can be discarded too.This graphical analysis can be corroborated by the Durbin-Watson test and its bootstrapped p-value to examine independence in these data.
Second, we carry out a statistical analysis of these data from the EDA until the selection of the most appropriate distribution for the demand data of each product under study is defined, following Subsection 3.3.Table 5 display a summary of the statistical results for 89 products of the inventory assortment of the Chilean food service.This summary indicates, among other aspects, the statistical distribution that fits the demand data best for each product.
Third, once we have selected the most appropriate distribution to model the demand data, we then use an adequate inventory management model to determine the optimum level in stock to place an order of products, and the optimum ordered quantity to minimize total inventory costs, following Subsection 3.4.Table 5 also show the optimal quantity of replenishment and the ROP obtained by applying the appropriate inventory management model.In Table 5, note that "P" and "EOQ" are the "inventory model" for the corresponding perishable and non-perishable product, respectively; "statistical distribution" corresponds to the fitted demand distribution for the indicated product according to the ID detailed in Table 12; λ and σ are the estimated demand mean rate and SD given in (1) and (2); k 0.95 is the SF for a service level of 95% given below (2); "SS" is the safety stock given in (2); Q is the EOQ given in (1); CR and Revista Colombiana de Estadística 38 (2015) 1-30 d 0 are given in (3); r is the ROP given in (2); and "JIT" is considered when this method is used for such a product.Note also that the symbol "-" is used when the corresponding value must not be calculated; "average" indicates that any distribution can be fitted for such a product and then the normal distribution is used.From Table 5, we describe 38 of the 89 products (components) from the inventory assortment with the EOQ model given in (1), 47 with the perishable model given in (3), and 4 with the JIT method, indicating that the total inventory is made up of mostly perishable type products.BS distributions were adequate for several of the demand data sets for those products allowing a distribution to be fitted (not JIT).
Fourth, once we have calculated the elements of the inventory models by using equations ( 1) and ( 3), following the financial approach detailed in Subsection 3.5, we compute the differences between direct and indirect costs (unit and annual), with weekly and annual ordering, and obtain the differences between the CMs with respect to the entire product inventory assortment of non-optimized and optimized systems, by using equations ( 4) to ( 12).Table 6 shows the annual and weekly OCs in both systems, where OCs for the optimized system increase 54.64%.Morillo (2009).
Table 7 shows the annual SCs of the non-optimized and optimized systems, which are diminished by 84.05%, translating our proposal into significant savings due to the improvements obtained by using the inventory management models.the non-optimized and optimized systems, for each of these financial indicators.A positive value of the differential indicates savings detected for the indicated product, using the optimized system.A negative value indicates the attained optimization is unfavorable for the indicated product.In Table 11, the set of critical products, which account for close to 80% of the optimized values obtained in these financial indicators have been delimited by a line.This is established as a cumulative percentage (CP) of optimized values regarding the total optimization attained in the differential profit or saving of the financial indicator, using the classification ABC; see details in Ramanathan (2006).

Illustration
In this section, we illustrate the optimized analysis for one of the 89 products in the inventory assortment of the case study presented in Section 4. We select this product due to its statistical features, so that a practitioner can better understand how the analysis is produced for a component, and thus replicated for other components.This analysis is divided into three parts following Steps 2, 3 and 4 of Algorithm 1.

Statistical Analysis
The data correspond to the demanded amount (D) of the ground beef product (in kg) with ID = P42, which was collected during the period under study.Table 8 displays a descriptive summary of the demand data that includes the sample median (50th percentile), mean ( d), SD, coefficients of variation (CV), skewness or asymmetry (CS) and kurtosis (CK), and sample size (n), among other statistics.This summary is obtained by the command descriptiveSummary() of the gbs package.From Table 8, we note that the CS and CK for P42 data show a distribution with positive skewness and moderate kurtosis.Figure 1 shows the histogram, boxplot and graph of the empirical CDF (ECDF) for P42 data.These graphs are built with the command histgbs() of the gbs package and boxplot() and ecdf() of the base R package.Note that: (i) the histogram shows a PDF with positive skewness and moderately heavy tails; see also Table 8; and (ii) the boxplot displays some outliers.Based on the EDA results, BS distributions seem to be good options for modeling P42 data, because they can accommodate their outliers and degrees of variability, skewness and kurtosis.BS, BS-t, gamma, IG, LN and Weibull parameters can be estimated using the ML method; see Barros et al. (2009).For this purpose, commands mlegbs() and gamlss() of the gbs and gamlss packages, respectively, can be used.The goodness of fit of the model to P42 data can be checked using the AD and KS tests, which compare the ECDF and the theoretical CDF assumed for the data (within BS, BS-t, gamma, IG, LN and Weibull models).The command used for obtaining these results is ksgbs() of the gbs package, and its corresponding adaptations to the gamma, IG, LN and Weibull distributions.Table 9 provides the p-values of the AD and KS tests for P42 data, from which we note that almost all of these distributions seem to be reasonable models for these data.However, based on the AD test results, which is more powerful than the KS test (Barros et al. 2014), only the BS, BS-t and LN models fit the data well at a significance level of 1%.The fit of the model to P42 data is visually illustrated in Figure 1, from where the ECDF (gray line) and the theoretical BS-t CDF (black dots) are compared on the right, whereas the histogram with the estimated BS-t PDF is plotted on the left.Probability plots with envelopes are shown in Figure 1, where "envelopes" are bands constructed by a simulation process facilitating the display setting.In the case of BS distributions, these envelopes are built using the formulas given in Subsection 3.3; see details in Leiva, Athayde, Azevedo & Marchant (2011) and references therein.The commands used to obtain these graphs are envelopebs() and envelopegbs() and their corresponding adaptations to the gamma, IG, LN and Weibull distributions.From these graphs, we note the appropriate fitting provided by BS, BS-t and LN distributions proposed for modeling P42 data, but also the gamma and Weibull models (omitted here) are suitable, where all the points are inside of their envelopes, therefore corroborating the results provided in Table 9.However, as can be seen from the boxplot given Figure 1(center), there are some outliers that can introduce an adverse effect on the ML estimates of the parameters of the distributions detected as suitable by the goodness-of-fit methods.Nevertheless, as mentioned, only the BS-t distribution has been proven to provide estimates robust to these outliers.Thus, we choose the BS-t distribution as the most suitable within the distributions proposed for describing P42 data.

Inventory Analysis
Once we have selected the BS-t distribution as the most suitable to describe the demanded quantity P42, we use the perishable product model for single period to determine the optimum quantity to be ordered for minimizing the total cost of inventory.First, we estimate the demand rate from the BS-t distribution as λ = 23.19 kg/day.Then, with this value, we determine the optimum replenishment quantity as d 0 = 15.45 kg, by using the formula given in (3), whose value must be applied as refueling.We consider a constant LT of l = 3 days, which is the same for all the products of the inventory assortment.Thus, at the beginning of each week, the stock level must be checked, and then a quantity of (15.45 − L j ) kg, for j = 1, . . ., 27, of the product must be ordered.
Note that, for the case of non-perishable products, once again we first estimate the demand rate and, then, with this estimate, we calculate the optimum Q and r using the formulas given in (1) and (2), respectively.Thus, when the stock level is in r units, we generate an order of Q units of this type of product.

Financial Analysis
Once we have chosen the appropriate inventory model for the P42 product, we determine its CMs according to the expression given in (11).First, we obtain the VCMs based on (8) following the definitions and the sequence of equations given in ( 4)-( 7).Second, we calculate the corresponding OC and SC using the formulas displayed in ( 9) and (10), obtaining OC 42,j = US$0.88per order and SC 42,j = US$0.085per stored unit of the P42 product.With this, we obtain CM 42 = −US$5242.72(optimized value) in comparison to −US$5613.85(nonoptimized value), reaching a reduction of 6.61% for this product.Table 10 provides the weekly values of SQs, VCMs, OCs, SCs and CMs for the non-optimized and optimized systems in the P42 product.

Concluding Remarks
We proposed a methodology useful for food service companies that allows their contribution margins to be optimized.The methodology was based on statistical tools, inventory management models and financial indicators.Its main steps were synthesized in Algorithm 1.Because there is a need to conduct case studies focusing on their applicability in firms to reduce the gap between theory and practice, and to transfer knowledge to the industry, we applied this methodology.Specifically, the case study was conducted with a Chilean food service company, which showed the importance of considering inventory models and statistical aspects for improving its supply and inventory policies, increasing its contribution margins.Inventory management models for perishable products showed to adjust adequately the demand for fruits and vegetables, which have the greatest unit contribution margins, so that such products can be considered as critical in the inventory assortment.With respect to the non-perishable products, it is noteworthy that the EOQ model fitted them considerably well.Such products can be stored indefinitely, without losses occurring until expiry, which has been used to reduce ordering costs.Products fitted by the model for perishables presented an optimized quantity similar to the demand rate, because such products have a shelf life and cannot be stored for a long period.A small amount of products (about 5%) used a JIT method.In summary, we validated improvements in logistics management applying an appropriate inventory model, increasing the contribution margins of a Chilean company.This result agrees with that reported by Ramanathan (2006), who linked inventory cost minimization and contribution margin maximization with type A products in the ABC classification.These products correspond to around 20% the total of products in the inventory assortment and are responsible for a proportion close to 80% of the total contribution margin.It is noteworthy that, although we attained an improvement of 10.47% in contribution margin for the studied Chilean company, 1.8% in total variable contribution margin, 54.64% in ordering costs, and 84.05% in storing costs, using the proposed methodology, still some aspects can be improved.For example, it is possible to explore the statistical dependence among products.Seasonality and trend factors, as well as dependence on the time of the demand, can be considered in the modeling by using time series models.Statistical dependence among products can be analyzed by means of multivariate structures for the models considered in this work.In fact, the authors of the paper are planning to collect real-world data of this type for a future study.In addition, from the practical standpoint, another future study considering demand data for a menu instead of its components is being considered by the authors.Studies of this type have been presented in the literature under names such as: assemble to order systems, inventory models with correlated demand, inventory models with multivariate demand, inventory models with multi-item demand, inventory models with multi-component demand and inventory models with component commonality, among others; see Agrawal & Cohen (2001) and Lu & Song (2005).In any case, incorporating all of these elements in the modeling can improve the precision of results, but can also increase its statistical complexity, rendering its use less attractive to a practitioner.This last aspect provides relevance to the present work.

Appendix: Tables with Results and ID Products
mean and variance of D are E[D] = λ and Var[D] = λ 3 /β, respectively.The IG distribution also shares the scale property, that is, c D ∼ IG(c λ, c β), with c > 0. The Lognormal Distribution.If Y = log(D) has a normal distribution with mean µ and variance α 2 , that is, Y = log(D) ∼ N(µ, α 2 ), then the RV D follows the LN distribution with shape α > 0 and scale β = exp(µ) > 0 parameters.The notation D ∼ LN(α, β) is used in this case.Thus, the PDF and CDF of D are

Figure 1 :
Figure 1: [first panel] histogram with estimated BS-t PDF (left), boxplot (center) and ECDF with estimated BS-t CDF (right) and [second panel] plots of probability with envelopes for the indicated distribution using P42 data.

Table 2 :
Basic functions of the gbs package.

Table 3 :
Costs involved in generating a purchase order (OC h ).

Table 4 :
Annual costs involved in the storage of a product (SC k ).
Compute the VCM for the product i in the jth week of the optimal policy obtained in Step 3.2.4.2 Determine the corresponding OC for the product i in the jth week.4.3 Calculate the SC for the product i in the jth week.4.4 Obtain the CM for the product i in the jth week.5: Repeat steps 1 to 4 until completing p products.6: Establish the optimized total CM and compare it with the non-optimized total CM.

Table 5 :
Summary of statistical and inventory models for the indicated product.

Table 6 :
Annual and weekly OCs (in US$) for the indicated system.

Table 7 :
Morillo (2009) US$) for the indicated system.SC(a) is the annual SC in US$.In non-optimized and optimized systems, 919962.6 and 192342.9unit/yeararestored, respectively; (2) SC(w) is the weekly SC in US$; (3) SC(au) is the annual SC in US$.Source: generated by the authors based onMorillo (2009).

Table 11 (
see Appendix) show the differential of VCM, OC, SC and CM values for all the products in a descendent order, obtained by subtracting the results from Fernando Rojas, Víctor Leiva, Peter Wanke & Carolina Marchant

Table 9 :
p-values of the indicated method and distribution for P42 data.

Table 11 :
Differential of optimized values and CP of VCM, OC, SC and CM for the indicated ID.