Data Science Solutions for Retail Strategy to Reduce Waste Keeping High Profit

Abstract

Waste reduction in retail is a fundamental problem for sustainability. Among waste reduction approaches such as recycling and donation, stock management based on demand estimation which leads to mitigate waste generation and maintain a high profit is expected to play an important role. However, demand estimation is generally difficult because fluctuations in sales are quite volatile, and stock-out leads to incomplete demand observation. Here, we propose data science solutions to estimate non-stationary demand with censored sales data including stock-outs and realize scientific stock management. Concretely, we extend a non-stationary time series analysis method based on Particle Filter to handle censored data, and combine it with the newsvendor problem formula to determine the optimal stock. Moreover, we provide a way of pricing waste reduction costs. A method to verify consistency between the statistical model and sales data is also proposed. Numerical analysis using actual Point-Of-Sales data in convenience stores shows food waste could be reduced several tenths percent keeping high profits in most cases. Specifically, in cases of foods disposed of frequently about 75% of working days, food waste decreases to about a quarter with the profit increases by about 140%. The way of pricing waste reduction costs tells new insights such as 27% waste reduction is achieved by 1% profit loss. Our method provides a practical solution for food waste reduction in the retail sector.

1. Introduction

1.1. Background

The United Nations predicts the world population reaches 8.6 billion in 2030 and 9.8 billion in 2050 [1]. Establishment of sustainable food systems for increasing food demand is essential. In this context, reduction of food loss and waste draws special attention since the one-third of foods are lost or wasted every year, which amounts to four times as much food required yearly for eliminating global hunger [2]. While food losses predominantly occur due to inefficient harvesting, storage facilities, infrastructural conditions and access to the market in low- and middle-income countries [3], food waste depends on such factors as cultural habits, meal planning, and shopping behavior [4]. The wasteful practices of the food industry and consumers in high-income countries are the predominant drivers for food waste [5].

Food waste reduction at retail is listed in the Sustainable Development Goals (SDGs) Target 12.3 [6]. In this study, we forcus on food waste reduction in the retail sector. We possess detailed sales data of a convenience store which contains records of disposal in addition to sales. Utilizing this data, we can consider specific measures for retail from the viewpoint of data science.

1.2. Literature Review on Food Waste Reduction

Component ratio analyses of total food waste [7,8,9] indicate that the household is the largest producer of food waste (53–71%), followed by food manufacturing and processing companies (17–30%), food service providers (9–12%) and retail (2–9%). Though the percentage of retail is relatively small, comprehensive initiatives including in retail are essential as the SDGs Target 12.3 covers the waste reduction in the retail sector. Studies on food waste constituent ratio in retail show it mainly consists of perishables, particularly, either bread or fruit and vegetables [10,11]. Concretely, a recent study shows that 30% originates from bread and 29% from fruit and vegetables in the total wasted mass [12]. Another recent research discussed the economic cost and the climate impact in addition to the mass of fruit and vegetable waste, and specified hotspot categories in Swedish food retail market [13].

Many strategies for food waste reduction have been proposed. For example, the concept of food recovery and food redistribution, specifically, food sharing including donating and bartering initiatives, have been widely spread in developed countries [14,15]. Recently, startups offer online platforms which help exchange surplus food [16,17] and cook and eat together [18,19]. Food waste valorization [20,21] is another primary strategy that converts food waste into chemicals, materials, and energy, or extracts high-value target compound from food waste. Specifically, techniques such as thermochemical, biochemical, fermentation, and extraction processes are used to obtain biofuels, drug, cosmetics, animal feed, fertilizer supply, compositing from food waste. Increasing demand for chemicals, materials, and fuels for sustainable future makes food waste valorization important approach. Customer awareness campaigns such as “Love Food Hate Waste” scheme introduced by WRAP in 2007 [22] and “I Value Food” project organized by Sustainable America [23], are also held for customer’s voluntary food waste reduction.

Actual practices taken for waste reduction at grocery retail stores in the United Kingdom has been studied recently [24]. Food donation to charities and food banks is conducted with making use of partnerships with surplus food distribution platforms such as FareShare [25] and local charities on food donations. Software and online platforms like “Food Cloud” [26] simplify the donating process, which helps retail managers grasp foods approaching their expiration date and communicate with registered charities. Other approaches include changing the labeling on products which shows just packaging date instead of “best before date” [27,28], investment in packaging and storage to improve food quality and shelf life, and price reduction for food items which do not meet company’s aesthetical standards or are approaching expiration date. Customer awareness campaigns, surplus foods for staff canteen, bakery waste sent for animal feeding, cooking oils diverted for biofuel production are also reported.

It is pointed out that waste management program which ignores social aspects of management and professional skills is prone to failure [29]. There exists a study identifying innovative food waste management practices and evaluating food service manager’s perspectives on these innovations [30]. In their study, incremental strategies in operational process and developments in current technology rather than radical innovations are reported to be dominant in the restaurants in Switzerland. Constraints to adoption of innovations are ranging from internal factors (lack of space, infrastructure, etc.), human capital factors (lack of time or staff, motivations, skills, etc.), external factors(uncollaborative customers, suppliers, etc.) and legal factors(local legislation, etc.). The primary driver for the adoption of food mitigation practices is financial savings and reputational gains [24,30]. Another research also reveals complex relationships in many drivers, namely, policy, private initiatives, consumer behavior, food waste impact, and food industry issues, and examine the effect of changing policies [31]. Technical proposal for waste reduction needs to take these managerial views into account for practical implementation.

1.3. Literature Review on Stock Management and Demand Estimataion

We explore specific solutions for food waste reduction in retail using the sales data. Scientific stock management based on precise demand estimation with sales data is a candidate for waste prevention. Related researches have been widely conducted. The newsvendor problem is a representative example. Perishable products such as newspapers and foods are disposed of when remaining unsold at the end of the day, while stock-out of those result in opportunity loss. The newsvendor problem tells the optimal stock to maximize profit under the fluctuating demand. There exist a wide variety of extension [32,33,34]. For demand uncertainty, Bayesian approach, specifically Bayesian Dynamic Programming (BDP), is proposed to update the uncertain demand distribution with sales data [35,36,37]. After the adverse effect of stock-outs in sales data on demand estimation is shown [38,39], BDP for censored sales data because of stock-outs is proposed [40,41]. Since BDP with censored data is tractable only when demand distribution is assumed to be exponential and Weibull distribution, more general demand distribution such as Poisson distribution cannot be assumed. Instead of solving BDP rigorously, heuristics such as the myopic policy which optimize the stock without considering future demand, are also developed [42,43]. Though the myopic policy is proved to order more than the optimal stock [44,45], it is shown that the difference between the myopic and the optimal stock is generally bellow 1% in most cases [46]. While BDP needs to assume the functional form of the demand distribution, some non-parametric approaches are also proposed [47,48,49,50]. Although vast researches have been conducted, these BDP and non-parametric approaches for censored demand assume the time dependency of the demand is stationary [51].

For non-stationary demand, there exist a variety of studies. For example, Markov Modulated Demand Process (MMD) provides non-stationary demand estimation [52,53,54,55]. BDP for sales time series with a single change point is also studied [56]. For non-stationary demand with auto-correlation, Auto Regressive Integrated Moving Average (ARIMA) is adapted [57,58,59,60]. Machine learning approaches for sales time series with underlying patterns is also studied [61,62,63]. However, these non-stationary demand analyses are not ready to handle censored sales data.

As for the conventional studies on waste reduction using sales data, there exists research on quantifying the amount of waste reduction, that reveals the amount of packaging waste of the shampoo and hair conditioner is reduced by introducing refill packages [64]. Other researches using food sales data, such as milk [65], bread [66,67], foods in convenience stores [68,69] and grocery retailing [70], mainly aim at estimating demand and a way for waste reduction is not proposed.

1.4. The Aim and Contribution to the Literature

Our aim and contribution are to provide a practical solution for demand estimation and stock management under non-stationary and censored demand, which realizes food waste reduction keeping high profit. We recently proposed a non-stationary time series analysis method [71] which is suitable for sales data analysis since it incorporates a recently discovered nature of the fluctuation of the sales data, namely, Taylor’s law [72]. Specifically, the sales process is basically the Poisson process of which fluctuation width is known to be the square root of the mean, however, the fluctuation width is proportional to the mean when the mean is large, which is caused by the effect of fluctuation in the population (customer) [73]. Our study [71] shows that the parameter (demand) estimation of non-stationary sales data without considering Taylor’s law leads to a misreading of the abnormal fluctuation in sales data as non-stationary changes of the demand, that results in poor accuracy of demand estimation.

In this study, we extend this non-stationary time series analysis method to be applicable to the demand estimation with censored sales data. We combine this method with the newsvendor problem formula to determine the optimal stock for each time. Moreover, we extend the newsvendor problem formula to determine the stock realizing the desired waste reduction, with pricing the cost of the waste reduction. We verify the effectiveness of our method for food waste reduction keeping high profit with the actual sales and disposal data.

2. Methods

Our method is based on the two methods, namely, the newsvendor problem [32,33,34] and the parameter estimation method for non-stationary time series [71]. In this section, we review these methods briefly and explain the detail of our method.

2.1. Newsvendor Problem

In the retail sector, managing the stock of perishable items such as newspaper and daily food is a challenging task, since those remaining unsold at the end of the day are to be disposed of, while sold-out of those results in opportunity loss. The difficulty in determining the optimal stock lies in the fluctuation of the demand. Namely, the demand k is not supposed to be constant but follow some Probability Density Function (PDF) $m(k)$. The newsvendor problem tells us how much stock s should we prepare a day for the fluctuating demand $m(k)$ to obtain the maximum profit. Assuming the cost c and price p, the expected profit $R(s)$ is written as follows:

$$\begin{array}{c}\hfill R(s)=p\{\sum _{k=0}^{s}km(k)+\sum _{k=s+1}^{\infty }sm(k)\}-cs\end{array}$$

(1)

Here, we assume the stock s is discrete and the unsold products at the end of the day are disposed of. The stock $s^{*}$ for the maximum profit is obtained by calculating the differential value of the Equation (1).

$$\begin{array}{c}\hfill R(s+1)-R(s)=p\sum _{k=s+1}^{\infty }m(k)-c=0\end{array}$$

(2)

Thus, $s^{*}$ is written with the Inverse of Cumulative Distribution Function(CDF) $M^{-1}(k)$, where $M(k)=Prob(K>k)$,44]

$$\begin{array}{c}\hfill s^{*}=M^{-1}(\frac{c}{p})\end{array}$$

(3)

Equation (3) means the stock $s^{*}$ for the maximum profit can be calculated with the cost ratio $c/p$, and the demand distribution $m(k)$ which determines $M^{-1}(k)$.

2.2. Assumption of the Demand Distribution

To solve the Equation (3), the functional form of the demand distribution needs to be assumed. In conventional researches, distributions such as the Poisson distribution [38], the Normal distribution [39], and the negative binomial distribution [74] are assumed. Meanwhile, our recent study [73] on sales data of a convenience store reveals it follows Taylor’s law [72]. That is, underlying fluctuation in the population of the observing objects (customers) affects the fluctuation of the observation of the Poisson process. Sales is basically a Poisson process of which fluctuation is known to be the square root of the mean value $\lambda$, but the fluctuation of the number of the customer leads to the fluctuation of the demand following the demand mean value $\lambda$ linearly, when the demand mean is large. The standard deviation $\sigma$ of Taylor’s fluctuation scaling is generally written as follows:

$$\begin{array}{c}\hfill \sigma =\sqrt{\lambda +{(\gamma ·\lambda )}^2}\end{array}$$

(4)

where $\gamma$ is the proportionality constant of Taylor’s fluctuation scaling. We incorporate this knowledge into the assumption of the demand distribution $m(k|\lambda )$ for given demand mean value $\lambda$. Namely, we adopt the Poisson distribution $Po(\lambda )$ for small $\lambda$, otherwise the Normal distribution $N(\lambda ,\sigma )$ with the fluctuation $\sigma$ corresponding to the Equation (4).

$$\begin{array}{c}\hfill m(k|\lambda )=\left\{\begin{array}{cc}Po(\lambda ),\hfill & \mathrm{if}\lambda <20\hfill \\ N(\lambda ,\sqrt{\lambda +{(\gamma ·\lambda )}^2}),\hfill & \mathrm{otherwise}\hfill \end{array}\right.\end{array}$$

(5)

Here, we approximate the Poisson distribution to the Normal distribution for $\lambda \ge 20$. Although the approximation can be performed for $\lambda \ge 10,$ [75] we adopt $\lambda \ge 20$ since the skewness of the Poisson distribution of $\lambda =20$ is $30\%$ smaller than that of $\lambda =10$, therefore, closer to 0, that is the value of skewness of the normal distribution, whereas the change ratio of the standard deviation by Taylor’s fluctuation scaling at $\lambda <20$ is below $10\%$, where $\gamma =0.1$. The specific $\gamma$ value for the sales time series used in this research is 0.12 [73].

The Equations (3) and (5) allow calculating the optimal stock $s^{*}$ for a demand mean value $\lambda$ for each cost ratio. For example, assuming $\lambda =10$, the $M^{-1}(k)$ in the Equation (3) is calculated using the Poisson distribution with $\lambda =10$. The optimal stocks $s^{*}$ for the cost ratio 0.9, 0.7, 0.5 are determined to be 6, 8, 10, respectively. Assuming $\lambda =50$ and the proportionality constant of Taylor’s fluctuation scaling $\gamma =0.1$, the Normal distribution in the Equation (5) is used to determine the optimal stock $s^{*}$. In this case, the optimal stock $s^{*}$ for the cost ratio 0.9, 0.7, 0.5 are determined to be 39, 45, 50, respectively. The optimal stock depends on the cost ratio and the demand distribution in this way. The detailed information is provided in the Appendix A (Figure A8).

2.3. Parameter Estimation of the Demand Distribution

To determine $m(k|\lambda )$ in the Equation (5), we need to estimate the demand parameter $\lambda$ for each time t. We have to consider that the demand is generally non-stationary and the fluctuation of the sales data follows Taylor’s law. We recently proposed parameter estimation method for non-stationary time series with Taylor’s fluctuation scaling [71].

In this method, the State Space Model (SSM) [76] is adopted for statistical modeling of sales data. The SSM assumes an observation data $y_t$ at time t is obtained stochastically by underlying invisible state $x_t$, and $x_t$ is also stochastically determined by $x_{t-1}$. Non-stationarity of a time series is modeled by a transition function of $x_{t-1}$ to $x_t$, which is called the system model. The system model we proposed for tracking non-stationary parameter is as follows:

$$\begin{array}{c}\hfill x_t=\delta (x_{t-1}+{\nu }_t)\end{array}$$

(6)

$$\begin{array}{c}\hfill \delta (x)=\left\{\begin{array}{cc}x,\hfill & \mathrm{if}x>0\hfill \\ 0,\hfill & \mathrm{otherwise}\hfill \end{array}\right.\end{array}$$

(7)

$$\begin{array}{cc}\hfill {\nu }_t\sim & (1-m)·N(0,\alpha ·x_{t-1})+m·U(-\beta ,\beta )\end{array}$$

(8)

where ${\nu }_t$ in the Equation (8) is a stochastic variable which governs the transition from $x_{t-1}$ to $x_t$ in the Equation (6). The value of $x_t$ is limited to be non-negative by the delta function in the Equation (7) since the Poisson parameter $\lambda$ (or $x_t$) needs to be positive. The tilde (∼) in the Equation (8) denotes that ${\nu }_t$ follows the distribution in the right side of the formula, namely, the superposition of Normal distribution $N(0,\alpha ·x_{t-1})$ with the mean value zero and the standard deviation $\alpha ·x_{t-1}$ ($0<\alpha$), and Uniform distribution $U(-\beta ,\beta )$ which takes values between $-\beta$ and $\beta$ ($0<\beta$). m characterizes the ratio of superposition of the Normal and the Uniform distribution ($0\le m\le 1$). The superposition of two distributions considers the case where $x_t$ generally follows the Normal distribution, but occasionally changes up to $\beta$ with a small probability m. The value of $\beta$ is determined so that the range of the Uniform distribution is much larger than the fluctuation of the Normal distribution ($\alpha ·x_{t-1}$), that allows a large change of $x_t$. m, $\alpha$, $\beta$ are hyper parameters. In this research, we adopted $m=0.05$, $\alpha =0.005$ and $\beta =4x_{t-1}$ by examining Root Mean Squared Error of $\lambda$ estimation of some artificial time series. The detailed information of the determination of the hyper parameters is provided in the Appendix of the reference [71].

We specify the observation model in the SSM, which describes the observation of data $y_t$ from the state $x_t$ as follows:

$$\begin{array}{c}\hfill y_t\sim \left\{\begin{array}{cc}Po(x_t),\hfill & \mathrm{if}\lambda <20\hfill \\ N(x_t,\sqrt{x_t+{(\gamma ·x_t)}^2}),\hfill & \mathrm{otherwise}\hfill \end{array}\right.\end{array}$$

(9)

Here, Taylor’s fluctuation scaling at large $x_t$ is considered as in the Equation (5). Wide fluctuation caused by Taylor’s fluctuation scaling is correctly treated as fluctuation, instead of a non-stationary change of $x_t$, by this term.

This SSM can be solved with the Particle Filter [77,78] to determine the state $x_t$ at time t. The Particle Filter is a kind of Monte Carlo simulation which approximates an arbitrary PDF of the state $x_t$ as the distribution of Monte Carlo sample values, which is called the particle distribution. The particle distribution is updated with the observation data $y_t$ at time t in accordance with the likelihood of each particle value.

The likelihood function for our observation model is PDF of the Poisson and Normal distribution in the Equation (9). We extend our model for the case of a censored observation. We adopt the following likelihood function $F(k_t|\theta )$ for a censored observation $k_t$ [79] for updating the particle distribution.

$$\begin{array}{c}\hfill F(k_t|\theta )={\int }_{s_t}^{\infty }f(k|\theta )dk\end{array}$$

(10)

Thus, the likelihood function for our SSM is described as follows:

$$\begin{array}{c}\hfill p(y_t|x_t)=\left\{\begin{array}{cc}\frac{1}{\sqrt{2\pi {\sigma }_t^2}}exp(\frac{-{(y_t-x_t)}^2}{2{\sigma }_t^2}),\hfill & \mathrm{if}{x}_t\ge 20\mathrm{and}\mathrm{non-censored}\hfill \\ \frac{1}{\sqrt{2\pi {\sigma }_t^2}}{\int }_{m=y_t}^{\infty }exp(\frac{-{(m-x_t)}^2}{2{\sigma }_t^2})\hfill & \mathrm{else}\mathrm{if}{x}_t\ge 20\mathrm{and}\mathrm{censored}\hfill \\ \frac{x_t^{y_t}·exp(-x_t)}{y_t!},\hfill & \mathrm{else}\mathrm{if}{x}_t<20\mathrm{and}\mathrm{non-censored}\hfill \\ \sum _{m=y_t}^{\infty }\frac{x_t^m·exp(-x_t)}{m!},\hfill & \mathrm{otherwise}\hfill \end{array}\right.\end{array}$$

(11)

where ${\sigma }_t=\sqrt{x_t+{(\gamma ·x_t)}^2}$.

The number of Monte Carlo samples is set to N = 10,000 and the initial particles are generated by the Equations (6)–(8) assuming $x_0=y_1$, namely sales $k_1$ for simplicity. In the case $k_1=0$, the initial particles are generated assuming $x_0=1$. The $\lambda$ value for each time t is obtained with calculating the median, instead of the mean, of the Monte Carlo samples of $x_t$ since its robustness against variations of particle values. While the SSM can be solved with Markov Chain Monte Carlo (MCMC) [80], the Particle Filter is selected since it can stably solve the SSM even with a small amount of data.

2.4. Pricing the Cost of Waste Reduction

We formulate the expected disposal amount $d_t^{*}$ as the following equation.

$$\begin{array}{c}\hfill d_t^{*}=\sum _{k=0}^{s_t^{*}-1}(s_t^{*}-k)m(k|{\lambda }_t)\end{array}$$

(12)

where $s_t^{*}$ is the stock for the maximum profit at time t which is determined with the Equation (3), and $m(k|{\lambda }_t)$ is the demand distribution shown in the Equation (5). Thus, the stock $s{(\alpha )}_t$ for reducing disposal $\alpha (0\le \alpha \le 1)$ times compared to $d_t^{*}$ is expressed as follows:

$$\begin{array}{c}\hfill s{(\alpha )}_t=\underset{0\le s\le s^{*}}{arg\; min}(|\alpha d_t^{*}-\sum _{k=0}^{s-1}(s-k)m(k|{\lambda }_t)|)\end{array}$$

(13)

The profit $R(s{(\alpha )}_t)$ is estimated as follows:

$$\begin{array}{c}\hfill R(s{(\alpha )}_t)=p(s{(\alpha )}_t-\alpha d_t^{*})-cs{(\alpha )}_t\end{array}$$

(14)

Thus, we can estimate the cost of the waste reduction by examining $R(s{(\alpha )}_t)$ / $R(s_t^{*})$.

2.5. Stock Determination for Fine Control of Waste Reduction

The demand distribution $m(k|{\lambda }_t)$ in the Equation (13) is assumed to be the Poisson distribution at small $\lambda$ as shown in the Equation (5). Since the Poisson distribution is discrete, the stock $s{(\alpha )}_t$ obtained with the Equation (13) is inevitably discrete, that leads to poor controllability of the waste reduction. For example, assuming the demand mean $\lambda =10$ and the target disposal ratio $0.5\le \alpha \le 1.0$, the stock value $s{(\alpha )}_t$ can only be two values, namely, 7 or 8. We modify the Poisson distribution with the $\mathsf{\Gamma }$ function to handle a real number stock value as follows:

$$\begin{array}{c}\hfill m(k|\lambda )=\left\{\begin{array}{cc}\frac{{\lambda }^k·exp(-\lambda )}{\mathsf{\Gamma }(k+1)},\hfill & \mathrm{if}\lambda <20\hfill \\ N(\lambda ,\sqrt{\lambda +{(\gamma ·\lambda )}^2}),\hfill & \mathrm{otherwise}\hfill \end{array}\right.\end{array}$$

(15)

There is one more thing to consider for fine control of waste reduction. The actual stock values for products in the retail, such as riceballs and meat buns, need to be discrete. We propose a probabilistic selection of the two stock value around $s{(\alpha )}_t$ to approximate the real number of $s{(\alpha )}_t$ under many trials. For example, if the target stock $s{(\alpha )}_t$ is $10.3$, the stock value 10 with probability 0.7 and the stock value 11 with probability 0.3 is adopted to reproduce approximately 10.3 under many trials. Namely, the probabilistic determination of the stock for a real number s is written with the following formula.

$$\begin{array}{c}\hfill \left\{\begin{array}{c}P(X=\lfloor s\rfloor )=1-(s-\lfloor s\rfloor )\hfill \\ P(X=\lfloor s\rfloor +1)=(s-\lfloor s\rfloor )\hfill \end{array}\right.\end{array}$$

(16)

where $\lfloor s\rfloor$ is the integer part of a real number s.

The effectiveness of the stock determination using the modified Poisson distribution in the Equation (15) and the probabilistic determination of the stock in the Equation (16), along with the method of pricing the cost of waste reduction are shown in the Appendix A (Figure A1).

2.6. Summary of the Procedure to Determine the Stock

We summarize the procedure to determine the stock $s{(\alpha )}_t$ with non-stationary and censored sales time series. Schematic diagram is shown in Figure 1.

• Determine the target disposal ratio $\alpha$ with considering the cost of the waste reduction by the Equations (13) and (14). The proportionality constant $\gamma$ of Taylor’s fluctuation scaling needs to be estimated empirically with the procedure described in the reference [71].
• Set time $t=1$. The initial stock $s_1$ needs to be determined empirically. In this research, $s_1$ is calculated with the Equation (13) using sales value $k_1$ for simplicity.
• Observe the sales $k_t$ and disposal $d_t$ of a product at time t. In the case of sold-out, the observed sales $k_t$ is the same value as the stock $s_t$.
• Estimate the demand mean value $\lambda$ at time t with the SSM and the Particle Filter, namely Equations (6)–(11). In this estimation, the sales and disposal of the sales time series up to the day t is used.
• Specify the demand distribution $m(k_t|{\lambda }_t)$ at the day t with the Equation (15).
• Determine the stock $s{(\alpha )}_{t+1}$ with the Equations (13) and (16). Here, the demand mean value $\lambda$ at time $t+1$ is assumed to be equal to the $\lambda$ at time t.
• Inclement $t=t+1$ and go back to step 3.

The verification of the method for estimating demand with censored sales time series and its effectiveness for determining the stock for each time is discussed in the Appendix A (Figure A2, Figure A3, Figure A4 and Figure A5).

2.7. Point-Of-Sales (POS) Data

The POS data are taken from 326 chain stores of a Japanese leading convenience store company, Seven-Eleven Japan Co., Ltd., during the 153 days from 1 June to 31 October 2010. The 326 stores are located in Kanagawa and Yamaguchi Prefecture in Japan. Individual owners mainly franchise the stores and sell extensive ranges of products, such as processed food, fast food, daily food, and non-food items including magazines, newspapers, stationeries, batteries, kitchen utensils, pet foods, liquors, and cigarettes. Sales of food items consist of 67.4% of total store sales and revenue from operations par shop a year was about 42 millions of Yen on average in 2010. The average floor space par shop is about 123 square meters [81].

The Pos data contains the record of each purchase at cash registers, and the daily amount of delivery and disposal for each shop and product. We obtained daily sales, delivery and disposal time series for each product in each store for the simulation of our method. We choose sales data of the food products which are disposed of more than 75% of sales days. There are two reasons for this selection. The first one is that our research aims to reduce food waste and verify the method can be applied to reduce an actual large amount of food waste. The other one is a technical reason for the verification of our method. The sales data with many sold-out days do not represent the underlying demand, which means the result of the demand estimation with our method cannot be verified, especially in the case when a higher demand than sales data is estimated with our method. The selected sales data are mainly on food items which are processed in the shops, such as fried chickens, frankfurter sausages and french fries. We limit the data to these items processed in the shops. Though these items are produced and disposed of several times a day, we aggregate production and disposal for each day for simplicity. The data with lost records in 153 days or unusual records because of promotions are excluded. We use 215 sales time series which are summarized in

Table 1. POS data used in this research.

Product	Number of Shops	Total Disposal	Total Sales	Price (yen)
frankfurter sausages A	11	6030	22,200	150
frankfurter sausages B	8	2670	4291	150
frankfurter sausages C	2	635	1970	150
french fries A	22	8532	21,653	155
french fries B	25	11,939	30,350	155
french fries C	2	921	2615	165
fried chicken A	72	40,140	171,193	165
fried chicken B	33	17,098	152,790	105
fried chicken C	10	4464	14,869	105
grilled chicken A	3	1164	3075	105
grilled chicken B	12	4253	8651	105
grilled chicken C	13	4873	19,272	105
potato croquette	1	452	759	80
skewered beef	1	297	331	120
(Total number)	215	103,468	454,019	-

for the simulation. The cost ratio is assumed to be 0.7, which is the typical value for the convenience store that is shown in the company’s investor relations [81].

The actual POS data with the record of disposal is rare in academic study. Though the POS data we possess is not the latest one, the statistical properties of sales data, namely, following the Poisson process and Taylor’s law, are universal and not time bounded. With difficulty in estimating fluctuating demand precisely, wastes from inappropriate stock management are assumed to occur all the same as formerly. Our POS data is qualified for the consideration of the method for waste reduction at present.

3. Results

3.1. Optimal Stock Determination Using POS Data

Sales time series is generally non-stationary because of many factors such as weather condition and social trends. Figure 2a,b shows examples of daily sales number for fried chickens (fried chicken A in Table 1) in two shops. The green line indicates the sales. Non-stationary changes on demand are observed, such as the continuous rise after day 100 in Figure 2a, and the abrupt up-and-downs in Figure 2b. The orange line shows the stock and the blue line denotes the disposal for each day. In Figure 2a, the disposal are observed in 78% of 153 days. Specifically, 327 fried chickens in total are disposed of in 153 days. As for the profit, 1656 fried chickens in total are sold at 165 yen, so the amount of sales is 273,240 yen. The production cost of 1983 fried chickens is 229,037 yen (assuming cost ratio 0.7). Thus, the profit in the period is 44,203 yen. In Figure 2b, fried chickens are disposed of in 76% of 153 days, and the disposal amounts to 438 in total. The total profit is 117,068 yen since 3387 fried chickens are sold at 165 yen. Our method aims at the reduction of these large amounts of disposal with keeping a high profit by estimating non-stationary demand and determining the optimal stock for each day.

Figure 2c,d are the simulation results of our method using the same sales data as in Figure 2a,b. In this simulation, as we mentioned in the Methods Section 2.6, the stock at day t is determined with the sales and disposal of the sales time series up to the day $t-1$. Specifically, we determine the mean of the demand at day $t-1$ with the time series analysis method described in the Section 2.3, namely Equations (6)–(11). In the case of sold-out, the potential sales number without sold-out is expected to be the same or more than the observed sales number. We estimate the likelihood of the demand when sold-out with the likelihood function for incomplete (censored) observation [79], namely Equations (10) and (11) in the Section 2.3, that take the expectation of more sales than stock into account. Assuming the demand at day t equals to the demand at day $t-1$, the demand distribution at day t is determined with the Equation (15) in the Section 2.5. The optimal stock at day t is determined with the Equation (3) in the Section 2.1, using the demand distribution and the cost ratio. In Figure 2c,d, the gray line indicates the underlying demand. As the disposal are observed almost every day, we assume the sales can be regarded as the demand. The green line shows the observed sales, the orange line is the stock, the magenta line is the estimated mean of the demand and the blue line denotes the disposal. The disposal and the profit in Figure 2c are 69 and 54,599 yen, respectively. The disposal is reduced to 21.1% (69/327) and the profit is increased to 124% (54,599/44,204) with our method. In Figure 2d, the disposal is reduced to 43.9% (192/438) and the profit is decreased to 93.7% (109,692/117,068).

We applied our method to the 215 sales time series mentioned above. In each case daily stock is determined with the Equation (3) using sales and disposal time series up to the day. Figure 2e shows the profit and the disposal of shops for each 215 product, which is normalized as follows. We simulated the daily sales time series with our method for each 215 product, and obtained the total profit and disposal for each product. For comparison, we divide the profit and the disposal of shops by those of our method, respectively. The horizontal axis denotes the total disposal of a sales time series of a shop normalized by that of our method, and the vertical axis shows the total profit of a sales time series of a shop normalized by that of our method. The plots whose value exceed 1 in the horizontal axis, which is shown in the dashed gray line, correspond to the result that the waste amount simulated by our method is smaller than that of actual shops. The blue plots show the cases that the profit of our method is higher than that of real shops, and the red plots show the cases that our method gives lower profits. The results of Figure 2c,d are marked with black circles. We confirmed our method improved profit in 83.7% (180/215) of the plots. The negative value in the vertical axis corresponds to the negative profit of a shop, that results in the negative ratio of the profit. By aggregating the profit and the waste in Figure 2e, the waste of the shops (total amount of waste 103,468) is reduced to 23.2% (total waste 23,967), and the profit of the shops (total profit 8,255,002 yen) is increased to 140% (total profit 11,558,387 yen).

3.2. Further Waste Reduction

As we mentioned in the Section 2.4 and Section 2.5, our method provides a way to estimate the cost of waste reduction, which helps reduce waste further than the waste level shown in Figure 2e. Figure 3a shows the result of controlling waste reduction with estimating the cost. Each black plot represents the median and interquartile range of the profit and the disposal for each target disposal ratio from 0.5 to 0.9 assuming the waste level in Figure 2e is 1.0. Here, the Equation (13) is used to determine the daily stock. The dashed blue line represents the theoretical cost, namely profit loss, for the reduction of disposal, which is calculated with Equations (13) and (14). We confirmed that the food waste reduction by 27% is realized with 1% profit loss, and that of 45% is realized with 3.3% profit loss with this actual data.

In Figure 3a, the deviation between the target disposal ratio and the obtained disposal ratio is observed, such as the obtained disposal ratio 0.73 for target disposal ratio 0.7 and obtained disposal ratio 0.55 for target disposal ratio 0.5. Though the deviation is relatively small especially for a large target disposal ratio, we can empirically compensate for this deviation by regression analysis with the obtained and target disposal ratio. Figure 3b is the result of the regression analysis. The regression function is estimated to be $g(\alpha )=1.09\alpha -0.09$, where $g(\alpha )$ is the target disposal ratio and $\alpha$ is the obtained disposal ratio. This function tells the disposal ratio $\alpha$ in the sales data analysis is obtained with target disposal ratio $g(\alpha )$, so it can be used to compensate the target disposal ratio to obtain the desired disposal ratio in the sales data analysis.

The magenta plots in Figure 3a are the results using the compensated target disposal ratio $g(\alpha )$, specifically from 0.455 to 1.0, to determine the optimal stock. We can see the deviation between the obtained and target disposal ratio is suppressed.

In the real shop operation, shop managers can estimate the regression function occasionally considering the deviation of the target and obtained disposal ratio. Besides, the regression function can be estimated in advance of the actual sales operation by using some artificial sales time series generated by random number simulation. The example of this estimation is shown in the Appendix A (Figure A6 and Figure A7).

3.3. Detailed Examination on Our Method

We assume that the cost ratio is 0.7 which is the typical value of the convenience store [81]. Meanwhile, the cost ratio is supposed to be different for each food industry and product. As a small cost ratio allows a large amount of disposal for obtaining the same profit, the effectiveness of our method for various cost ratio should be examined. We examined the dependency of the profit and the disposal on the cost ratio with the 215 sales time series as a case study. The stock for each day is determined with the procedure described in the Section 2.6, and with the Equation (3) to determine the optimal stock as in Figure 2. Figure 4a is the estimated dependency of the total profit on the cost ratio. The green line denotes the total profit of the shops, and the magenta line shows that of the simulated result with our method. The total profit of our method is higher than that of the shops in the case that the cost ratio is above 0.55. Figure 4b is the total disposal of the shops (the green line) and that of our method (the magenta line). The disposal of our method is smaller than that of the shops above the cost ratio 0.25. At the cost ratio 0.55, the same profit as the shops and the waste reduced to 40.6% (42,016/103,468) compared to the shops are obtained with our method.

Figure 4a shows the profit of the shops surpasses that of our method below the cost ratio 0.55, whereas our method provides the optimal stock. One reason is that the sales data contains stock-out up to 25% of 153 sales days as mentioned in the Section 2.7. The increase in the profit with the optimal stock is supposed to be underestimated. We discuss this effect in detail later. Another reason is that there exist fluctuation that is not considered in our assumption. Figure 5a,b shows the examples of actual sales time series. The green line indicates sales and the magenta line denotes the estimated demand mean, assuming the cost ratio 0.1. We examine the probability of each sales data realized by the assumed demand distribution, namely the Equation (15). The red plots indicate that the value of Cumulative Distribution Function (CDF) of the sales value realized by the demand distribution is greater than or equal to 0.95. These peaks with CDF greater than or equal to 0.95 appear more than theoretically expected. Figure 5c shows the histogram of the CDF values for each sales in the 215 sales time series. Though the frequencies of CDF values are to be uniform if the assumed distribution fits the sales data, the frequency of CDF greater than or equal to 0.95 is more than two times larger than that of other bins of the histogram. Figure 5d is the histogram of the count of CDF greater than or equal to 0.95 in each 215 sales time series. The theoretical realization probability of CDF greater than or equal to 0.95, for 15 days in 153 days is 1.0%. But 84.2% of the sales time series show the CDF values more than or equal to 15 days. One of the sales time series contains sales peak over $100\sigma$. These peaks are supposed to arise from some irregular reasons such as bulk orders for events. Since shops cover these abrupt sales peaks which occur with very low probability, some prior knowledge or notification in advance of the sales peaks might have been provided somehow to the shop.

To discuss the potential of our method, we remove the abnormal sales time series whose count of CDF greater than or equal to 0.95, are observed more than or equal to 15 days. A sales time series with a peak over $100\sigma$ is also removed. Figure 6a,b are the estimated dependency of the profit and the disposal on the cost ratio with the 33 sales time series obtained. The green line is the result of shops, and the magenta line denotes that of our method. The navy line illustrates the result of our method with compensation of the sales underestimation because of the stock-out. Here we calculate the expected sales in stock-out days using the Probability Density Function (PDF) of the demand distribution, namely the Equation (15). Specifically, we get the expected sales above the stock for each sold-out day by calculating the sum product of the PDF and sales count greater than or equal to the stock, and devide the obtained value by the summation of the PDF whose sales count greater than or equal to the stock. Then, we compensate the sales for each sold-out days up to disposal count at the day. In Figure 6a, we can see the compensated profit of our method (the navy line) approaches asymptotically to that of shops (the green line) as the cost ratio decreases. In Figure 6b, the compensated disposal of our method (the navy line) decreased compared to the original disposal of our method (the magenta line) since the additional sales by the compensation decrease the disposal.

Figure 6c shows the histogram of CDF with the 33 sales time series. The excessive count of CDF greater than or equal to 0.95 is suppressed. Since the data size is relatively small, the histogram appears not to be uniform. To see the consistency with the result and the theoretically expected frequency, we simulated the theoretical histogram with random numbers whose distribution follows Equation (15) and each demand mean is the estimated values for each 33 sales time series. Figure 6d denotes an example of the simulated theoretical histogram, which shows a similar tendency with the result of actual sales data. We confirmed the equivalency of the theoretical and the actual distribution with Kolmogorov–Smirnov test. Namely, the p-value comparing the CDF distribution obtained with the actual sales data, and that of random number simulation, repeated 200 times with changing random number seed, is 0.07 on average, which shows the difference of the theoretical and the actual distribution of the CDF is not statistically significant at the p-value 0.05.

Figure 6a tells the profit and the disposal with our method compared to shops are equivalent at the cost ratio 0.2, and the disposal 59.7% (8247/13,822) with the same profit compared to shops are obtained at the cost ratio 0.4. At the cost ratio 0.7 which is typical for the convenience store, the profit is 208.0% (1,260,691/606,070 yen) and the disposal is 22.4% (3097/13,822) compared to the shops.

4. Discussion

As we mentioned in the Introduction section, some methods are already proposed for demand estimation. Here, we examine the conventional methods and clarify the advantages of our method.

BDP and some related heuristics [35,36,37,40,41,42,43] are developed for demand estimation and stock management under demand uncertainty and censoring. Meanwhile, BDP for non-stationary demand is not proposed [51], except for the recent research which proposes BDP for demand (without censoring) with a single change point [56]. Non-parametric methods proposed for demand estimation and stock management [47,48,49,50] also assume stationary demand. Developing these methods for non-stationary demand is one of the open research areas [51].

MDD, ARIMA and NN are applicable to non-stationary demand. In these methods, MDD proposed in the literature [52,53,54,55] assume the transition probability matrix and the demand distribution within each Markov state is known and fixed. This assumption limits the practicability since non-stationarity is usually unpredictable. ARIMA [57,58,59,60] and NN [62] can be used to non-stationary demand with underlying patterns such as auto-correlation and seasonality. However, the underlying patterns may change with such factors as the social trend and weather condition. Moreover, MDD, ARIMA and NN are not ready for handling censored sales data.

One may think of a more simple non-stationary time series analysis method, Moving Average and their family, such as Exponentially Weighted Moving Average (EWMA) [82]. These smoothing methods are also not applicable to censored sales data.

We examine above mentioned methods which are applicable to a non-stationary sales time series, namely ARIMA, NN and EWMA. In these methods, ARIMA and NN are executed with the statistics software R, specifically the software packages $forecast$ and $nnet$, respectively. EWMA is implemented as $D_t=\alpha Y_{t-1}+(1-\alpha )D_{t-1}$, where $D_t$ is estimated demand for day t, $Y_{t-1}$ is sales data at day $t-1$ and $\alpha$ is the smoothing constant [82].

Here, we explain the procedure of the examination. The POS data shown in Figure 2c is used for the examination. As ARIMA and NN require training data, the demand data within the initial 90 days are used for initial training and no censoring is considered in this training data. For each day from day 91, all previous data up to the day is used for retraining the ARIMA and NN model. The demand for the next day is predicted with the model. For simplicity, the stock for each day is determined as the same value as the predicted demand, that is equivalent to the case where the cost ratio is $0.5$ for a symmetric demand distribution.

Specific ARIMA model for each day is determined to be $(p,d,q)=(0,1,1)$ with Akaike Information Criterion (AIC) [83], where p is the order of AR, d is the degree of differencing and q is the order of MA [57]. As for the NN model, the number of neurons in the hidden layer is set to be 10. As a large number of neurons results in overfitting and instability of the estimation, while small number leads to underfitting, we select one of the balanced condition that can generally follow the trend with stability.

As EWMA needs no training, daily demand is estimated from day 2 with the EWMA equation. $D_1$ is assumed to be $Y_1$ and the value of $\alpha$ is set to be 0.12 with try and error to follow the trend.

Figure 7a–c is the result of demand estimation with ARIMA, NN and EWMA, respectively. The gray line indicates underlying demand. The blue line shows the estimated demand with sales data without censoring, and the magenta line shows that with censoring. The estimated demand without censoring generally follows the trend. Meanwhile, the estimated demand with censoring fails to follow the trend. Using censored sales data without any compensation results in an underestimate of the underlying demand.

Since these conventional time series analysis and machine learning are not ready for censored demand, some exogenous compensation of censoring is needed for demand estimation with censored sales data, such as compensation with intuition, experience or considering some patterns in sales. Our method compensates the censored demand with the statistical likelihood, namely Equation (11), and needs no exogenous compensation. To our knowledge, our approach firstly provides practical solutions for non-stationary and censored demand estimation.

As the examination within our method, we compare the model with and without considering Taylor’s fluctuation scaling (TFS). Figure 7d is the result of demand estimation with and without considering TFS. Note that demand $\lambda$ in Figure 7d is relatively large (i.e., $\lambda >20$), so Taylor’s fluctuation scaling is assumed in the sales data. The blue line denotes the demand estimation with TFS, specifically the Equation (9) is used for the estimation. The magenta line indicates the demand estimation without TFS, namely only Poisson distribution in the Equation (9) is assumed. The estimated demand without TFS is smaller than that with TFS. We confirmed the profit obtained without TFS is $3.3\%$ smaller than that with TFS. As shown in the reference [71], large fluctuations caused by TFS are regarded as non-stationary changes if TFS is not considered in the conventional Poisson model. Unlike large values, small values of the demand distribution are observed under censoring. These small values are considered as non-stationary changes in the Poisson model, which lowers the estimation.

5. Conclusions

The non-stationary parameter estimation method based on Particle Filter which considers Taylor’s fluctuation scaling is extended for demand estimation with censored sales data. By combining this method with the newsvendor problem formula to determine the optimal stock, reduction of food waste keeping high profit is realized. The extension of the newsvendor problem formula to examine the relationship between the amount of waste and profit provides new insights into the cost of waste reduction.

The retail shops handle many products and have to determine the stock for each product for each day. For example, the convenience stores generally deal with over 500 food products which require disposal management considering the expiration date. Our method helps retailers estimate the non-stationary demand and determine the optimal stock for each product for each day, which enables reducing food waste with keeping a high profit in the real activity. The stock determination with intuition [66] or constant inventory level that result in a large amount of waste can be avoided by our method.

Meanwhile, multilateral managerial viewpoints, especially constraints such as internal factors (lack of space, infrastructure, etc.), human capital factors(lack of time or staff, motivations, skills, etc.), external factors(uncollaborative customers, suppliers, etc.) and legal factors(local legislation, etc.) are needed to be considered for practical implementation [30]. As for the internal, external, and legal factors, our method is accessible to ordinary retail stores. The method is based on data analysis on a personal computer. More concretely, the data used in our method is sales data of the shop itself. Furthermore, the data analysis requires no sizeable computational resource such as a large server or cloud computing. Namely, our method can estimate the demand for 153 days of a product in a few seconds with a personal computer for consumer use (i.e., Intel Core i5-6400 2.7GHz CPU and 16GB Memory). Short computational time is important since ordinary shops are not necessarily accessible to huge computational resources and have to manage stocks for many products.

In the factors mentioned above, human capital factors, specifically, the requirement for some Information Technology (IT) and data science knowledge, can be a barrier to the introduction of our method. Human resources who are knowledgeable of IT and data science are not sufficient worldwide. Providing some computer software or online platform that is easy to use without programming skills, along with holding some hands-on seminars, could be a realistic solution for implementation. Some other factors, such as perception in retail managers that stock-out leads to lost not only sales but also the loyalty of the customer [66] could be an obstacle since our method could increase stock-outs by eliminating over-stocking. Awareness campaigns for managers as well as consumers to change the perception of stock-out might be required. Low recognition of technologies and reluctance of retail managers to adopt IT tools [30] should also be addressed. We need to accumulate successful cases and present the merit properly. Such cultural barriers as a lack of consumer interest and awareness as well as a hesitant company culture are also considered to be the primary obstacles in the case of promoting the circular economy [84]. As a long time measure, policies to increase the experts on IT and data science who can make full use of the company’s data is also essential since scientific consideration based on data, which is not limited to the proposed method, should realize effective measures for food sustainability.

Financial issues for implementing green innovations should be carefully examined. In the case of innovations which require a substantial amount of investment, financial actors such as banks, venture capitals, business angels, and private investors play an important role. Recent researches reveal relationships between the types of funding and innovations for renewable energy [85], the pressures of various financial actors in the greening of the financial system in Italy [86], and factors that prevent investment decisions of green companies in biomass production sector [87]. Financial sustainability in cases of waste management of municipalities is also discussed [88]. Since retailers now work on green innovations mainly in their operational processes and technologies available, taking advantage of investments from green finance to overcome financial issues is not reported in surveys on the actual waste reduction practices in retail [24,30]. Our method needs the initial cost of introducing software that is not supposed to require a substantial amount of investment as the cases of biomass and renewable energy. Meanwhile, a small amount of investment could be an obstacle for small shops. Promoting green finance for small-scale investment needs, which is reported to be scarce [87], could help spread our method and other green innovations in retail.

As for the technical aspects, we demonstrated an approach to verify the consistency of the assumed statistical model and the actual sales data using CDF in the Section 3.3. Specifically, we obtained a CDF histogram which corresponds to the realization probability of sales value by the assumed demand distribution. This CDF histogram tells mismatch between the assumed model and actual data. Besides its information can be used to remove the abnormal data, it can be used to modify the statistical model.

There exist limitations and rooms for future research in our method. Periodicities in sales time series such as weekly cycle and seasonality that are specific for each type of industry and product are not included in our model. If the systematic deviations because of periodicities are substantial, taking such periodicities into account in our model may yield better accuracy of the demand estimation. We assume a constant price for each product. Retailers reduce the price of foods which are approaching their expiry date. Incorporating pricing decisions into the model makes the model more efficient. Demand substitution which expects to mitigate the stock-out effect is not considered in our model. Understanding the interaction for each sales amount of product considering stock-out, and incorporating this knowledge into the model may lead to further waste reduction.

The adjustment of Taylor’s scaling factor $\gamma$ in the Equation (4), which governs fluctuation of the demand distribution, may contribute to the estimation accuracy. Statistically, the correlation of the underlying population (customers) causes Taylor’s scaling exponent takes a value between 0.5 and 1.0 whereas the Equation (4) assumes the scaling exponent to be 1 (linear). Some examples of the scaling exponent below 1 are found in non-food products on the POS data, such as stamps which are often purchased in bulk and magazines that are specifically bought on the launch day [73]. Influence of a TV program or social trend may cause the correlation of customers for food products. In the future work, the modification of Equation (4) corresponding to the scaling exponent below 1 will provide a more general form of our method.

Our method provides practical solutions for estimating non-stationary and censored demand. The method has the potential to apply to sales data of other sectors than retail, specifically, food manufacturing and supplying industry, and food service sectors such as restaurants, cafeterias, and catering services since the property of the sales data, that follows Poisson process and Taylor’s law, is assumed to be universal.